这次是分享李婧翌团队的一篇综述《Modeling and analysis of RNA-seq data: a review from a statistical perspective》，从统计学的角度理解RNA-seq的分析

分析的方向

目前正对RNA-seq的数据主流的有4个方向（当然事实上不止这些，可以辛苦读者慢慢收集整理，欢迎与我讨论）

1. 基因sample-level，这里主要是看生物学处理间，基因表达模式的相似性，通常用Pearson或Spearman相关系数进行表示
2. Gene-level，这里涉及到基因表达的定量
3. Transcript-level，这里涉及到对不同转录本的定量
4. Exon-level，这里涉及到差异可变剪切的检测

接下来作者主要围绕这四块内容进行在统计学上的理解

1). Sample-level

基于sample的分析，目的是检测不同sample的表达模式的相似性，通常可以利用Pearson and Spearman correlation coefficients来衡量。如果是利用全部基因来计算相关系数，管家基因的存在势必会 "夸大" 相关系数，因此比较好的方法是利用相关基因而不是全部的基因来计算，而R包TROM就是用来解决这类问题的，TROM通过计算TROM分数来选择出相关基因后，进行sample间相关系数的计算

除了计算相关系数，我们可以利用非线性的方法t-SNE或UMAP来进行降维聚类，以观测样本间的相似性

2). Gene-level

Gene层面的研究主要是对基因进行定量，并且进行差异表达分析，差异表达分析基本统计学模型的假设为，某个基因的count（表达量）在各个sample中的分布服从泊松分布或者负二项分布（如果是log以后的值一般认为服从正态分布）：

其中：

1. Y_k,ij 代表的是 condition k 中第 j 个sample gene i 的表达量
2. S_kj 代表 condition k 中第 j 个sample 的size factor
3. θ_ki 代表 condition k 中 gene i 的真实表达水平（可理解为在 condition k 的条件下， gene i 在各个 sample 中的平均表达水平）
4. Φ_i 表示 gene i 的dispersion

其基本假设为：

image.png

上图表示某个基因A在所有sample中表达量的分布（但由于生物学的sample较少，所以统计学家往往直接利用负二项分布去拟合），均值为S_kjθ_ki
经过统计学检验两个分布的差异，显然该基因在condition 2的表达量要小于condition 1中的，p值的计算可以考虑用置换检验来从两个分布中抽样计算p值

另外一种就是基于的共表达分析：

其中：

1. A_ij 代表gene i 与 gene j 的相关性矩阵
2. k 代表 gene k
3. d_ij = 1 - T_ij，用于表征基因之间的相似性距离

2). Transcript-level

一个基因可能有不同的转录本，基于转录本水平的分析主要是对一个基因的不同转录本进行定量

而对转录本定量往往存在一个问题，那就是对于同一个基因来说，一部分转录本的序列有overlap，那么reads在比对回去的时候，很难区分这些reads到底来自哪一个转录本，因此统计学家往往采用EM算法进行转录本的定量

并作出如下定义：

θ_j 表示 reads 来自于 isoform j 的概率
定义isoform的集合为：{ 1，2，3，.....，J }

Region based：
假设 X={ X_s | s∈S }，X_s代表map到region s上总的reads数，假设map到region s上总的reads数服从λ_s的泊松分布：

这里假设参数λ_s满足线性关系：
假设如下例子：

一共有三个isoform，这里的 X_s 特指map到外显子上的reads，而该例子中一共有4个外显子，X_s = X_s1 + X_s2 + X_s3 + X_s4
对于每一个转录本来说，如果该转录本缺乏某一个外显子，那么这个外显子上的reads数为0，似然函数：

相应的外显子区域的多项式值为 1（相当于没有贡献），利用极大似然估计的思想，我们的目的是确定似然函数 L() 取得最大值的时候参数 λ_s 的值，而 λ_s与θ_j满足线性关系**，即确定 λ_s 后利用EM算法对 θ_i 进行分配，原理参见：《用简单的EM算法模型理解RSEM算法》

经过计算后，我们可以得到：例如 θ₁=0.37，θ₂=0.33，θ₃=0.3，相当于一共有100条reads分配到该区域（该基因），isoform 1 表达37条，isoform 2 表达33条，isoform 3 表达30条

Reads based：
基本模型如下：

这种模型的基本思想是计算 reads i 来自于 isoform j 的概率，根据条件概率公式，

表征同时选中 reads i 和 isoform j 的概率，也就是定量结果

Regression-based：

回归的方法和 Region based 的方法理解相似，只不过 Region based 利用极大似然的方法估计参数；而 Regression-based 基于最小二乘的思想求解参数

4). Exon-level

这一块主要分析的是可变剪切事件，那么可变剪切事件的PSI定义为：

其中：

1. C_I denotes the number of reads supporting the inclusion isoform
2. C_E denotes the number of reads supporting the exclusion isoform
3. L_I and L_E denote the lengths or the adjusted lengths

而可变剪切的统计学模型是：

例如 inclusion 事件的reads的分布满足于总reads数为 n = C_I + C_E，reads属于 inclusion 的概率为 ψ（PSI）的二项分布（均值μ = n×p），而判断差异可变剪切事件为：
构建不同condition的二项分布，对于某个基于来说
经过统计学检验两个分布是有差异的（C_Ik的分布是有差异的），因而判断为差异可变剪切事件
经过统计学检验两个分布是没有差异的（C_Ik的分布是没有差异的），因而判断为非差异可变剪切事件