edgeR基因表达差异分析
正式文件的总结
==注意：:==- edgeR实际上是为了读取计数而设计的。 不建议将预测的转录丰度输入edgeR而不是实际计数。
引数
可以使用read.table ()、read.delim (如果文件很多，则为readdge )和columns=)一次读取多个文件。 readDGE ) )必须指定两列：计数列和基因标识符列
x
DGEList对象，构建组
DGEList是可以包含各种内容和统计的列表。 DGEList至少有必要的元素。 counts、samples (包括group组信息和lib.size文库的大小)、counts用于存储表示矩阵，samples用于标记样本信息和库的大小，group声明组
创建包含#组标记的DGEList
y
组
y
在DGEList中需要分组，但也可以将其他信息，例如lane通道、基因长度、基因注释等信息追加到1个注释化的DGEList中
过滤，删除低表达基因
低表达基因不仅不能使用，而且会干扰结果，要清除任何样品都没有足够数量序列片段的基因，应从下游分析中排除。 1 .低表达没有生物学意义，2 .通过去除低表达数据可以更准确地推测数据中平均-方差关系，3 .观察差异表达的下游分析中的运算量减少
edgeR包的filterByExpr函数提供了一种自动过滤基因的方法，可以尽可能多地保留表达数量足够的基因。 该函数默认选择最小组内的样本数量作为最小样本数量，并至少在该数量的样本中保存有计数10个以上序列片段的基因。 过滤基准以最小对内的样本数为基准(在本例中，最组内的样本为3 )，在所有样本中表达数(count )小于10的基因的个数超过最组内的样本数的情况下，过滤的基准为其遗传基因换算成cpm为cut.off.cpm=10/
CPM比例变换
常用的标度变换有CPM(counts per million )、log-CPM、FPKM、RPKM. CPM，将counts变换为CPM指数countspermillion，消除测序深度的影响log-CPM是CPM的逻辑二进制化。 cpm函数在进行log2转换之前将补全值与cpm值相加。 缺省完成值为2/L，2为“预计数”，l为采样总序列数(百万)的平均值，因此log-CPM值通过log2(CPM2/L )从CPM值计算。 对于一个基因，CPM值为1，相当于序列约2千万个样本中有20个计数，或序列约7.6千万个计数76个。
cpm
lcpm
手动过滤
经验以cpm=1位为cutoff点。 但是，这不是最准确的。 因为随着测序深度的增加，例如20million(2千万)，cpm=1意味着counts=20。 门槛可能有点高。 如果测序深度较低，例如2million(2百万)，cpm=1意味着counts=1。 阈值可能过低。 在这种情况下，可以使用自动过滤，也可以通过cut.off.cpm=10/进行计算。
CPM(y )是
keep_cpm 1)=3 #在此示例中，将组中的采样数设置为3
x
自动过滤
使用edgeR:filterByExpr ()自动过滤
keep.exprs
x
正规化
虽然归一化不是绝对必要的，但建议进行归一化。 在重复样品中，本不应具有生物学意义的外部因素会影响单个样品的表达。 例如，中第一批制备的样品总体表达高于第二批制备的样品，假设所有样品的表达值范围和分布应该相似，需要进行归一化以确保整个实验中各样品的表达分布相似。 可以使用edgeR:calcNormFctors ()函数用TMM方法归一化，在对象为DGElist的情况下，通过计算归一化系数来确定文库大小的缩放系数为x$samples$norm.factors
测序深度
文库的大小不同，测序深度也不同。 这是基本建模过程的一部分，它会自动进入倍数变化或p值计算。 始终存在，不需要用户干预。
有效库大小
在某些情况下，如果某些高表达基因消耗了特定样本文库的大部分大小，则该样本的剩余基因采样不足。 除非调整该效应，否则该样品中的剩余基因可能错误地表现为下调。 (由于高表达基因多，测序资源被大量消耗，非高表达基因相对较低。 )函数edgeR:calcNormFactors ()，缺省情况下使用TMM规格化，规格化后将为采样指定缩放因子。 原始库大小与比例因子的乘积称为有效库大小。 有效库大小将替换所有下游分析中的原始库大小。
y
==注意：==规范化不会直接修改为counts值，规范化因子会自动存在于x$samples$norm.factors中。
GC含量
基因的长度
为了~
GC含量和基因长度在每个基因间不会改变，因此这些都是相对的，因此对差异分析影响很小。MDS图形展示 样本无监督聚类
这种图表使用无监督聚类方法展示出了样品间的相似性和不相似性，能让我们在进行正式的检验之前对于能检测到多少差异表达基因有个大致概念。
# 图形展示
library(RColorBrewer)
lcpm
col.group
levels(col.group)
col.group
plotMDS(lcpm,labels=group,col=col.group)
#或者直接使用
plotMDS(y)
负二项式模型
计算生物变异系数
edgeR中对于离散度的估计主要有两种模型，其中比较经典的是==qCML==。但qCML方法仅适用于==单因素==的数据集。例如，比较癌和正常样本间的基因表达差异，而不考虑其他因素(如年龄、性别等)的影响。qCML后续会使用与fisher精确检验比较类似的 ==exact test== 进行差异表达分析。 如果使用qCML估计离散度，则不需要预先设定实验矩阵！(因为就只看一个因素，不必大费周折) 而对于更复杂的实验设计，则推荐使用GLMs模型，此时就需要提供实验矩阵design。 qCML方法：仅适用于具有单因素设计的数据集，但是有着更可靠的性能，尤其在小样本下有很好的表现，比较切合NGS测序。 通过estimateDisp()来估计公共离散度和 tagwise离散度(一个命令来运行)
y < -estimateDisp(y)
或者先计算common离散度，再计算tagwise离散度
y
y
计算差异基因
et
topTags(et)
广义线性模型(Glm)
对于更复杂的实验设计(有多个因素)，可以通过广义线性模型来及性能你和。
计算离散度
通过下面来估计common离散度、trended离散度、tagwise离散度。
y
或者分开依次进行
y
y
y
关于design的构建(?)：
#group是区别分组的factor对应不同样本，lane是不同样本的的lane道
design
colnames(design)
计算DE基因
edgrR涉及到差异表达分析的函数有很多： exactTest、glmFit、glmLRT、glmQLFit、glmQLFTest。 qCML估计离散度需要搭配 exact test 进行差异表达分析，对应 exactTest 函数。 而其他四个glm*都是与GLM模型搭配使用的函数。其中，glmFit 和 glmLRT 函数是配对使用的，用于 likelihood ratio test (似然比检验)，而 glmQLFit和 glmQLFTest则配对用于 quasi-likelihood F test (拟极大似然F检验)。 有两个方法glmQLFit()和glmQLFTest(),在两个检验方法中，首选QLFit，因为它反映了估计每个基因的离散度时的不确定性。当重复次数较少时，它提供了更强大和可靠的错误率控制
# 例子
group
design
fit
如果没有重复样本？
首先确定一个良好的BCV，人类数据一般设定为0.4 一个例子：
bcv
counts
y
et
如果是人类数据：
y_bcv
bcv
et
注： 如果差异基因过少，可以尝试下调BCV
输出结果
result = topTags(et, n = nrow(et$table))$table
查看统计
summary(de
因为我的笔记分布CSDN、简书、知乎专栏等比较零散，管理起来比较麻烦，思考再三申请了一个 微信公众号，会更加方便地发布更多有关生信息、统计方面内容，如果你觉得有需要可以尝试性和我缔结契约 ^ . ^ 。公众号如下：
对我唯一的精神鼓励可能就是下方的点赞了吧 ^^