发布时间 : 星期日 文章R语言DESeq 包介绍更新完毕开始阅读be5e159ab9f3f90f76c61ba8
R语言DESeq包介绍
分析RNA序列数据的一个主要任务是探测基因的差异表达，DESeq包提供了测试差异表达的方法，应用负二项分布和收缩的分布方程估计。
1. 包的安装
输入如下命令，DESeq和相关的包就可以自动下载和安装。 > source(\> biocLite(\
相关的包会自动下载安装，安装的包如下：
中间会有个选择需要更新相关的包，选择更新全部，更新的包有：
另外还要安装一个数据包，供下面介绍包中的方法时使用，包名为pasilla.
2. 输入数据和准备
2.1 计数表
数据表的第i行第j列元素表示第j个样本的第i个基因有多少个reads。本文使用的数据来
自于pasilla数据包，函数system.file告诉我们数据文件保存的路径。 > datafile = system.file( \> datafile
[1] \在R中读取这个文件，使用read.table函数。
> pasillaCountTable = read.table( datafile, header=TRUE, row.names=1 ) > head( pasillaCountTable )
2.2 元数据
没有元数据的数据是没有用的，元数据可以分为三组，分别是样本(行)，特征(列)和整个实验的信息。
首先需要样本的描述信息，data.frame的列表示各种信息，行表示7个样本。 > pasillaDesign = data.frame(
+ row.names = colnames( pasillaCountTable ),
+ condition = c( \+ \
+ libType = c( \+ \> pasillaDesign
这边简单地使用R代码进行设定，通常情况是从扩展表中读取这些数据。为了分析这些样本，我们需要解释single-end和paired-end的方法，这边首先简单的分析paired-end样本。 > pairedSamples = pasillaDesign$libType == \> countTable = pasillaCountTable[ , pairedSamples ] > condition = pasillaDesign$condition[ pairedSamples ] 现在，我们有下面的数据输入 > head(countTable)
> condition
对于自己的数据，可以简单的创建因子。 > #not run
> condition = factor( c( \我们现在举例CountDataSet，DESeq包的核心数据结构 > library( \
> cds = newCountDataSet( countTable, condition )
2.3 规范化
函数estimateSizeFactors估计统计数据的大小因子 > cds = estimateSizeFactors( cds ) > sizeFactors( cds )
如果统计数据的每列除以这列的大小因子，这样统计值就变成同一规模，使它们具有可比性。函数counts可以做这个计算。
> head( counts( cds, normalized=TRUE ) )
3. 方差估计
DESeq推断依靠估计典型的数据间的方差和平均关系，或者等效的离差和均值。离差可以理解为生物变异系数的平方。估计离差可以使用以下命令。 > cds = estimateDispersions( cds )
函数estimateDispersions做了三步，首先估计每条基因的离差，然后，通过估计匹配一条曲线，最后，每个基因分配一个离差，从每条基因估计值和匹配值选一个。为了让用户知道中间过程，fitInfo对象被储存下来。 > str( fitInfo(cds) )
函数 plotDispEsts可以画出每条基因的估计值和平均正常统计值的关系。 > plotDispEsts( cds )
图1 经验的(黑点)和匹配的(红线)离差值与平均正常统计值的关系图
在任何情况下，可以被子序列测试使用的离差值被存储在cds的特征数据集里。 > head( fData(cds) )
4. 推断：称之为差异表达
4.1 两个实验条件下的标准比较
为了看在条件”untreated\和”treated\下是否有差异表达，我们简单的使用nbinomTest函数， > res = nbinomTest( cds, \> head(res)