时间序列研究的是基因表达的动态行为,测量的是一系列和时间点之间有强烈相关性的过程。和针对某一时间点的基因表达进行差异分析不同,时间序列更加关注是发现基因表达的趋势,以有助于理解生物学动态变化过程(比如对刺激的反应、发育过程、周期行为等)。也就是说,时间序列关注的是整体变化趋势而不是某特异表达。
基因表达有时间依赖性,蛋白质会根据不同的功能需要进行合成,即使稳定的状态下,mRNA不断的被转录,蛋白不断的合成与降解,这个过程被高度调控。当细胞面对新的状况,比如,饥饿、感染、应激等,一些调控因子会通过调控自身或其它基因的表达来启动或抑制转录,甚至激活新的表达模式。很多情况下,这种表达模式通过激活一些转录因子开始,这些转录因子又会反过来调控其它的基因,而这些基因几乎都是对新情况的反应。通过时间序列分析,可以鉴定只在一些特定或新的状况下特异表达的部分基因。然而,为了确定这些状况下表达的完整的基因集,进而确定它们之间的相互关系,时间序列的数据分析就尤为重要。也就是说,来确定的不是新情况下稳定状态的那些通路或基因,而是为了达到这种状态(比如肝脏重建)被激活的那些通路或基因,甚至网络。对于静态实验的差异基因的选取已经有很多方法的报道,但因为基因表达是一个动态的过程,尤其对大鼠肝脏从开始切除到最后重建,还有肝癌的生成等过程,所以鉴定并找出那些表达随时间的变化而变化的基因非常重要。也就是说更关注的是整个过程中的总体趋势而不是某特异的表达水平,筛选的是表达模式类似的差异共表达基因。不同样本组中的差异共表达基因更可能是调控因素,也就更能解释表型之间的差异。
这样就有几个挑战,一是要分析的数据量会很大,二是实验条件变多,三是要发掘实验动态变化本质,传统的统计学方法比如t-tests就无能为力了,需要运用新的统计学方法,四是样本间的时间间隔并不总是相等。
主要针对配对比较分析的SAM,LIMMA等方法在对变化趋势的分析上无能为力。而对时间序列的数据处理,有不少报道,比如等级聚类、基于主成份分析的聚类等,虽然这些聚类可以鉴定并可视化共调节的基因,但基因数目多的时候难以解释,还有一个不足就是,不能得到随时间变化有统计学意义的基因。
而对时间序列的分析,需要:首先,可以使用统计学程序来鉴定显著表达变化的基因;第二,把随时间变化发生显著表达变化的基因进行聚类并且可视化。这个可以通过回归来解决,其中,时间被视为数量变量,实验条件视为分类变量,这样就可以分析趋势变化。
maSigPro的全称是Microarray Significant Profiles,采用2步回归策略。第一步选择基因,第二步选择变量。并且,可以调整模型参数更拟合数据,使用虚拟变量代表实验条件。数据需要经过预处理才可以由maSigPro分析,包括背景矫正,log2 ratios计算,lowess标准化,一般的芯片数据处理方法都可以,比如RMA, MAS5等。下面简要概述maSigPro分析的步骤及原理。
在用maSigPro进行分析时,。一般情况都是两个或两个以上的感兴趣的变量,其中一个典型的就是时间变量,另外一个通常都是分类变量,代表实验组别(比如不同的处理,细胞株,组织等)。模型如下:
其中,i=实验组别
J=时间点
r=重复
εijr=随机变量
D=虚拟二进制变量(实验条件)
T=时间
yijr=标准化后的表达值
β,δ,γ,λ=回归系数
β0,δ0,γ0,...,λ0是对应于对照组的回归系数。βi,δi,γi,...,λi解释第(i+1)组和对照组之间的特定差异(线性,二次方,三次方等)的回归系数。注意其中的虚拟二进制变量,一共I-1个,以此区分每组和对照组,见表1。这个模型隐形定义了和实验组一样多的模型。例如,因为第一组中的虚拟变量是0,所以第一组的模型是y1jr=β0+δ0T1jr+γ0T21jr+…+λ0TJ-1 1jr+ε1jr,而对于第二组来说其虚拟变量是1,其模型为y2jr=(β0+β1)+(δ0+δ1)T2jr+(γ0+γ1)T22jr+…+(λ0+λ1)TJ-12jr+ε2jr。
maSigPro分析的第一步是应用最小二乘法来估算每个基因的上面所描述的回归模型的参数,选出有统计学意义的回归模型。第二步是选择变量。根据第一步挑选出的基因,由逐步回归产生新模型,采取向后剔除方法(backward)相对较多。
maSigPro包得到的时间序列数据所有差异表达基因表达模式的动态变化聚类图
maSigPro包得到的时间序列数据差异表达基因表达模式变化
Cluster1-9中的代表基因随时间变化的RMA表达值