基因芯片主要分为双通道cDNA芯片和高密度寡核苷酸芯片。
双通道cDNA芯片:每个微阵列产生两个探针水平的数据集(红色通道和绿色通道)。
高密度寡核苷酸芯片:每个微阵列产生一个探针水平的数据集。一些探针是匹配探针(Perfect match,PM),一些探针是错配探针(Mismatch,MM。不过有些芯片无MM探针)。
高密度寡核苷酸芯片,以Affymetrix芯片为例。这些芯片有多种型号,如常见的Affymetrix 3’ -biased Arrays 有HG_U95Av2、HG-U133A_2、HG-U133_Plus_2等型号,Affymetrix Exon ST Arrays则有HuGene-1_0-st-v1等型号。通常Affymetrix 3’ -biased Arrays有PM、MM探针,而Affymetrix Exon ST Arrays只有PM探针。很多时候可以根据芯片型号判断出该芯片所属物种,如:HG、HuGene是人类,Mouse是鼠,Zebrafish是斑马鱼。常见的芯片数据库有GEO、ArrayExpress。若想在GEO获得指定型号芯片的有关信息,访问http://www.ncbi.nlm.nih.gov/gds?term=affymetrix%20AND%20GPL%5BETYP%5D。如点击[HG-U133_Plus_2]Affymetrix Human U133 Plus 2.0 Array,可在Web link找到该型号芯片的Affymetrix产品页面链接(该页面可下载产品说明,注释文件等,不过下载需要先注册一个账号)。
所有的芯片用Affymetrix扫描仪进行扫描,用Affymetrix软件初始量化特征,该软件涉及到一些文件格式,它们是:
[EXP]包含实验的基本信息
[DAT]芯片的扫描图像
[CEL]特征的初始量化(每个探针的荧光强度)
[CDF]探针在芯片中的定位信息,探针到探针组的映射
[CHP]包含基因表达水平(用affy软件评估)
通常由DAT扫描图像得到CEL文件,在后续的计算中,我们只需要CEL文件和CDF文件。在这里不讨论EXP和CHP文件。
基因芯片技术的特点是使用寡聚核苷酸探针检测基因。使用ReadAffy函数读取CEL文件获得的数据是探针水平的(probe level),即杂交信号,而芯片数据预处理的目的是将杂交信号转成表达数据(即表达水平数据,expression level data)。存储探针水平数据的是AffyBatch类对象,而表达水平数据为ExpressionSet类对象。基因芯片探针水平数据处理的R软件包有affy, affyPLM, affycomp, gcrma等.
Affy芯片数据的预处理一般有三个步骤:
- 背景处理(background adjustment)
- 归一化处理(normalization,或称为“标准化处理”)
- 汇总(summarization)。
最后一步获取表达水平数据。需要说明的是,每个步骤都有很多不同的处理方法(算法),选择不同的处理方法对最终结果有非常大的影响。
Affy芯片基础知识
Affy基因芯片的探针长度为25个碱基,每个mRNA用11~20个探针去检测,检测同一个mRNA的一组探针称为probe sets。由于探针长度较短,为保证杂交的特异性,affy公司为每个基因设计了两类探针,一类探针的序列与基因完全匹配,称为perfect match(PM)probes,另一类为不匹配的探针,称为mismatch (MM)probes。PM和MM探针序列除第13个碱基外完全一样,在MM中把PM的第13个碱基换成了互补碱基。PM和MM探针成对出现。
背景处理
R软件包affy用于芯片背景噪声消减的函数是bg.correct(),而MAS和RMA方法是最常用的两种方法。
MAS方法将芯片分为k(默认值为16)个网格区域,用每个区域使用信号强度最低的2%探针去计算背景值和噪声。RMA方法的原理比较复杂,可以参看文献:
R. A. Irizarry, B. Hobbs, F. Collin, et al. Exploration, normalization, and summaries of high density oligonucleotide array probe level data. Biostatistics, 4:249–64, 2003b. 11, 18, 27, 232, 241, 432, 443。
data.rma <- bg.correct(data.raw, method="rma")
data.mas <- bg.correct(data.raw, method="mas")
class(data.rma)
## [1] "AffyBatch"
## attr(,"package")
## [1] "affy"
class(data.mas)
## [1] "AffyBatch"
## attr(,"package")
## [1] "affy"
class(data.raw)
## [1] "AffyBatch"
## attr(,"package")
## [1] "affy"
可以看到:ReadAffy()读入的CEL芯片数据以AffyBatch类数据形式存储,而背景消减后得到的依然是AffyBatch类数据。
MAS方法应用后PM和MM的信号强度都被重新计算。RMA方法仅使用PM探针数据,背景调整后MM的信号值不变。
head(pm(data.raw)-pm(data.mas),2)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501131 79.34 62.93 63.23 72.87 62.48 64.43 62.97 60.86
## 251604 77.57 63.06 63.73 80.69 63.07 66.53 63.33 60.34
head(pm(data.raw)-pm(data.rma),2)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501131 111.2 102.21 93.23 116.36 93.75 76.15 102.82 85.21
## 251604 103.9 88.69 72.57 90.75 82.23 72.64 87.75 85.32
head(mm(data.raw)-mm(data.mas),2)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501843 79.34 62.93 63.24 72.90 62.48 64.44 62.97 60.85
## 252316 77.56 63.06 63.73 80.66 63.07 66.51 63.33 60.34
#差值为全部为0,说明rma方法没有对mm数值进行处理
head(mm(data.raw)-mm(data.rma),2)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501843 0 0 0 0 0 0 0 0
## 252316 0 0 0 0 0 0 0 0
identical(mm(data.raw), mm(data.rma))
## [1] TRUE
归一化处理(normalization)
同一个RNA样品用相同类型的几块芯片进行杂交,获得的结果(信号强度等)都不可能完全相同,甚至差别很大。为了使不同芯片获得的结果具有可比性,必需进行归一化处理。这一步的方法也很多。归一化处理的affy函数为normalize(),以Affybatch对象和处理方法为参数。
线性缩放方法
这是Affy公司在其软件(4.0和5.0版本)中使用的方法。这种方法先选择一个芯片作为参考,将其他芯片和参考芯片逐一做线性回归分析,用回归直线(没有截距)对其他芯片的信号值做缩放。Affy公司的软件做回归分析前去除了2%最强和最弱信号。使用很简单,mas方法做背景消减后做归一化分析的R语句是:
data.mas.ln <- normalize(data.mas, method = "constant")
head(pm(data.mas.ln)/pm(data.mas), 5)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501131 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 251604 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 261891 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 230387 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 217334 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
head(mm(data.mas.ln)/mm(data.mas), 5)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501843 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 252316 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 262603 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 231099 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
## 218046 1 0.8788 1.002 1.141 0.9929 1.029 0.7578 0.8834
可以看出,线性缩放方法以第一块芯片为参考,它的数值没有被处理,而其他芯片都被缩放了。对同一块芯片,不同探针的缩放倍数是一个常数。PM和MM的缩放方法完全一样。
非线性缩放方法
此方法获得的结果比线性方法要好,做非线性拟合时不是取整张芯片而仅取部分(一列)作为基线。
data.mas.nl <- normalize(data.mas, method = "invariantset")
head(pm(data.mas.nl)/pm(data.mas), 5)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501131 1 1.050 1.417 1.359 1.399 2.023 1.0086 1.3124
## 251604 1 1.348 2.279 1.838 1.587 2.431 1.1126 1.3059
## 261891 1 1.564 1.397 1.675 1.497 1.556 1.3167 1.5509
## 230387 1 1.259 1.683 1.390 1.360 1.504 1.0145 1.2239
## 217334 1 1.009 1.127 1.241 1.229 1.263 0.8417 0.9934
可以看到,同一芯片不同探针的信号值的缩放倍数是不一样的。
分位数(quantile)方法
这种方法认为(或假设)每张芯片探针信号的经验分布函数应完全一样,使用任两张芯片的数据做QQ图应该得到一条斜率为1截距为0的直线。
data.mas.qt <- normalize(data.mas, method = "quantiles")
head(pm(data.mas.qt)/pm(data.mas), 5)
## CHIP1 CHIP2 CHIP3 CHIP4 CHIP5 CHIP6 CHIP7 CHIP8
## 501131 0.7176 0.9602 1.140 1.181 1.112 1.187 0.8145 1.0156
## 251604 0.6984 0.9814 1.199 1.209 1.112 1.172 0.8287 1.0143
## 261891 0.6939 0.9956 1.137 1.205 1.111 1.186 0.8389 1.0579
## 230387 0.7653 0.9777 1.171 1.183 1.115 1.185 0.8153 0.9921
## 217334 0.9508 0.9547 1.063 1.162 1.130 1.173 0.7945 0.9266
其他
如循环局部加权回归法(Cyclic loess)和 Contrasts方法。
汇总(Summarization)
最后一步汇总是使用合适的统计方法通过probeset(包含多个探针)的杂交信号计算出计算表达量。affy包的函数为computeExprSet。需要注意的是computeExprSet函数除需要指定统计方法外还需要指定PM校正的方式:computeExprSet(x, pmcorrect.method, summary.method, …)
两个参数可以设定的值可以通过下面函数获得:
pmcorrect.methods()
## [1] "mas" "methods" "pmonly" "subtractmm"
generateExprSet.methods()
## [1] "avgdiff" "liwong" "mas" "medianpolish"
## [5] "playerout"
常用的汇总方法是medianpolish, liwong和mas。liwong方法仅使用PM做背景校正(pmcorrect.method="pmonly")。例如:
eset.rma.liwong <- computeExprSet(data.rma.lo, pmcorrect.method="pmonly",
summary.method="liwong")
最后的结果 ExpressionSet 类型数据
三合一处理和“自动化”函数
mas5函数
eset.mas5 <- mas5(data.raw)
rma函数
也是由affy包提供,其背景处理方法为rma法,归一化处理使用分位数法,而汇总方法使用medianpolish:
eset.rma <- rma(data.raw)
## Background correcting
## Normalizing
## Calculating Expression
gcrma函数
由gcrma包提供。 Affy的软件(比如mas方法)使用MM数据做背景处理,但由于MM出现的问题(上面提到过),这些方法可能高估了背景值。而rma方法在做背景处理时没有使用MM数据,这可能又低估了背景值。MM序列公布后有人对其特异性进行了评估,并使用这些评估结果建立了新方法。gcrma就是这样的方法,也是封装好的三合一方法
library(gcrma)
eset.gcrma <- gcrma(data.raw)
## Adjusting for optical effect........Done.
## Computing affinities.Done.
## Adjusting for non-specific binding........Done.
## Normalizing
## Calculating Expression