本期分享一篇2017年12月21号发表在 Molecular Cell 上的偏分析的文章，文章的题目是 A Massively Parallel Reporter Assay of 30 UTR Sequences Identifies In Vivo Rules for mRNA Degradation

理解本文大致需要的知识有：

• 简单的微分方程以及线性回归概念。
• 对胚胎发育过程有一定了解，对简单的分子生物学实验有一定了解，如PCR、显微注射等。

本文的主要结论：

本文研究的关注点在胚胎发育的早期，即卵细胞受精一段时间以后，mRNA的降解速率与3’UTR序列中可能存在的顺式元件的关系，通过分析发现了一些motif可能会影响mRNA的稳定性，如促进mRNA稳定的poly-U，导致不稳定的ARE、miR-430、PUM motif等。

首先我简单阐述一下生物学问题：3’UTR的某些序列或者说是motif会影响mRNA的稳定性，比如poly-U、miR-430(miRNA结合位点，可以引起mRNA降解)，这篇文章的主要亮点是把UTR序列对mRNA稳定性的影响同发育过程联系了起来。

以斑马鱼的胚胎发育为例简单阐述一下胚胎发育过程中的转录情况，为了方便之后的描述，我们定义受精时的时间为坐标原点，即 0hr。在0-3hr(卵细胞受精后的3个小时以内)，受精卵细胞核中的DNA是无法转录形成新的mRNA的，此时细胞质中存在的主要是先前母亲的卵细胞留下来的mRNA，而在3hr之后，即可以转录形成新的mRNA，也就意味着3hr之后一些新转录形成的miRNA可以参与到调控mRNA稳定性的过程中去。

这篇文章的机制做的不是很深，但分析的方法可圈可点，因此这一期的解读主要集中在算法层面。

实验样品与流程

首先基于斑马鱼内源基因的3’UTR序列合成了90000条长度为110nt的人工设计的3’UTR序列，同时两端各预留出20bp进行同源重组和PCR扩增(这种合成方法的极限是150bp)，将这些人工合成的3’UTR接到GFP荧光蛋白基因的3’端，最终构建成如图中所示的载体。将得到的载体通过体外转录的方法得到mRNA，即图中所示的mRNA reporters，同时为了检验poly-A尾对降解速率的影响，有分为A+(在末尾加19bp的A)与A-(不添加poly-A)两类，之后将得到的mRNA通过显微注射的方法注射到斑马鱼的受精后的胚胎中，在不同的时间点如(0hr,1hr,2hr….)取样，通过RNA-seq的方法测量此时胚胎中mRNA的含量，即可得到Fig1.B中所示的点图，之后通过不同时间点测量所得的mRNA的含量估计出对应mRNA的降解速率。

通过不同时间点的mRNA含量估计mRNA的降解速率

作者在这一块的分析的亮点在于首先考虑到了可能存在两种降解类型的mRNA，一类是不受基因组转录的miRNA影响的mRNA，另一类是收到影响的mRNA(即3hr之后降解速率发生变化的mRNA)，所以在估计降解率时要进行区分，如Fig2.A B所示。利用mRNA在不同时间段的测量值进行降解速率的估计要用到微分方程的一些相关知识，不过不是很难，我着重讨论一下。

首先假设mRNA的降解是一级化学反应，则可以列出如下微分方程 $$\frac{dX}{dt}=-b(t) \cdot X(t)$$ 其中$X(t)$表示t时刻测量得到的mRNA的含量，$b(t)$表示在t时刻的降解速率。对于Fig2A中所示的降解类型的mRNA，$b(t)$是一个定值，我们用$\beta$来表示,则新的微分方程可以写作: $$\frac{dX}{dt}=-\beta \cdot X(t)\Rightarrow X(t) = X_0 \cdot e^{-\beta t} \Rightarrow log X(t) =log(X_0) - \beta \cdot t$$ 因此对于第一类mRNA(如Fig2A所示)可以通过对log X(t)和t做回归，此时斜率的绝对值即是降解率，截距的指数即是初始值。

对于第二类mRNA，即假设b(t)是随着时间变化的分段函数，即为： $$ b(t)=\left \{ \begin{array}{c} \beta_0 \ \ \ t \leq t_0\\ \beta \ \ \ t \gt t_0 \end{array} \right. $$ 其中$\beta_0为t_0$之前的降解速率，此时需要拟合的数据变为3个，即$\beta_0, \ \beta,\ t_0$，可以采用分段函数线性回归的方法，我第一时间想到的是MATLAB的fminsearch函数，然而实际操作后已被日翻，看来还是要用作者手动选取$t_0$的方法用来估计，顺便搜了一下发现python有个包pwlf可以做这个事，并不是很难，就不展开叙述了，同时把我的MATLAB代码放在下面，如果有读者可以改进，请直接pull request，与此同时，再一次被Python社区活跃的氛围所震撼。

对于一个mRNA report，通过似然比检验的方法选择一个最符合当前mRNA降解模式的模型，进而得到这一段110nt序列所对应的降解速率。在这里，我必须要吐槽一下，为什么不提供代码啊。

Fig.2 主要结论

Fig.2 是全文比较重要的一张图，首先Fig.2 确实发现了预期中存在的两类mRNA，即图中所示的early-onset、late-onset，并且late-onset的估计所得的$t_0$主要集中在3hr-5hr(Fig.2D),并且含有poly-A尾的mRNA的半衰期确实要长(Fig.2E)。但我不能不能理解为什么不用虚线和实线分别表示early-onset、late-onset这个颜色看的真是要急死人，打成黑白的之后更是看不懂了。

%参数拟合代码，效果奇差，欢迎改进
tdata_1 = 0:0.1:5;
ydata_1 = 0.5*tdata_1 + 1;
tdata_2 = 5.1:0.1:10;
ydata_2 = 1.5*tdata_2 - 4.15;
tdata=[tdata_1,tdata_2]
ydata=[ydata_1,ydata_2]
%函数部分另外保存为一个新文件，MATLAB不能再脚本内定义函数
function sse = sseval(x,tdata,ydata)
A1 = x(1);
A2 = x(2);
lambda_1 = x(3);
lambda_2 =x(4);
t0=x(4);
sse = sum((ydata - ((tdata <= t0).*(tdata*A1+lambda)+(tdata <= t0).*(tdata*A2+lambda))).^2);
%根据数据进行参数拟合
fun = @(x)sseval(x,tdata,ydata);
x0 = rand(5,1);
bestx = fminsearch(fun,x0)

Fig.3

Fig.3 作者把这些不同的mRNA分成了三类：

• 第一类是在A+(含有poly-A)与A-(不含有poly-A)中都属于early-onset类型的，即受精后即开始降解的。
• 第二类是在A+(含有poly-A)与A-(不含有poly-A)分属于不同类型的。
• 第三类是在A+(含有poly-A)与A-(不含有poly-A)都属于late-onset类型的。

分为这三类的目的是为了之后找这三类的motif，也就是Fig.4的结果，此处结果处理的不太好，不如之后通过回归模型的结果可信，不继续探讨了。

作者提出的假说

因为发现了存在两种降解模型的mRNA，作者提出在胚胎发育过程中存在两种降解的信号，一是受精后0-3hr(核内DNA未发生转录)存在的maternal degradation singal，另一类是3hr之后的存在的zygotic singal。为了验证这一假说，作者把前文提到的体外合成的mRNA report library再一次注射到细胞中，只不过这一次是没有受精的卵细胞(没有了zygotic singal)，最终的实验的结果表明,对于只受到maternal degradation singal影响的early-onset的mRNA，其在卵细胞(Fig5.A纵坐标)和受精卵(Fig5.A横坐标)变化并不大，相关性很好，r有0.81，但在late-onset类型的mRNA中，其降解速率在卵细胞中明显低于胚胎，相关性也较差0.21。算是间接的验证了作者假说的正确性，那么接下来就是找机制了，但由于这篇文章找的机制都是之前已经报道过的，所以机制也没有发现新的东西，我猜也可能是这篇文章没有发到 Cell 上的原因，但后续的分析真的做的很出彩。我当时读完的感觉是，线性回归还能这么玩，强行找到了motif，大写的服气。

通过mRNA report序列特征与降解速率的回归分析找motif

这一块不细说了，涉及细节较多，讲不清楚，如有同好感兴趣，可以私下讨论。简单捋一下，首先把110nt的序列(字符串)切成3-7 mer的字符串，共计21824个自变量，与此序列对应的降解速率做回归。首先通过Lasso回归的方法选择权重不为0的自变量，再用选取得到的变量再做一次回归，之后再通过线性回归的模型找motif，最后将找到的motif做了实验验证，感觉做的很扎实，但对有有两个Y值的late-onset的回归模型没有详细说明究竟怎么确定的因变量，感觉是藏了点什么东西？

通过线性回归找到了late-onset(核基因转录开始之后才开始降解的那一类mRNA)中的motif： miR-430的靶定位点，可以看出确实符合预期即转录开始后产生miR-430，导致了这一类mRNA的降解，但是尴尬的是，这个机制就是本文的第二个通讯作者，也就是他老板找到的。

总的来说，我很喜欢这篇文章的分析部分，首先是能想到把mRNA降解分为两类，之后通过回归的方式居然能想到办法找motif，我之前也用过文中的回归方法，但我真的是只做了个回归，这一点学到了，但前段时间看的一篇用deep learning中CNN网络找motif的方法也很好，与本文相比，计算量要小很多，这一块的知识现在不是很会，等学会了也可以试一下。最后再吐槽一些作者神一般让人看不懂的写作与令人抓狂的图片展示方式。