王丽萍,唐旭清
(江南大学 理学院,中国江苏 无锡 214122)
急性肺损伤(acute lung injury,ALI)是临床上常见的健康问题,是各种直接和间接致伤因素引起的肺泡上皮细胞及毛细血管内皮细胞损伤,可造成弥漫性肺间质及肺泡水肿,导致急性低氧性呼吸功能不全,发展至严重阶段(氧合指数<200)被称为急性呼吸窘迫综合征(acute respiratory distress syndrome,ARDS)。在过去的20年里,ARDS的死亡率一直保持在40%左右[1]。目前,临床上可使用的ALI/ARDS的标志物很少[1],因此,获得肺损伤患者的临床生物学信息,发现新的靶点生物标志物,对于有效治疗疾病至关重要。
疾病的进展是一个动态的过程,Lesterhuis等[2]通过对比复杂疾病的动态生物标志物和静态生物标志物,发现在系统状态中存在“临界点”。Dahlem等[3]研究表明:在某一触发因素下,“临界点”会迅速进入到疾病状态,所以可将临界点视为复杂疾病的早期预警信号。一般来说,复杂疾病的进展过程可以分为3个阶段:正常状态、疾病前状态(或临界状态)和疾病状态。疾病前状态是从正常状态到疾病状态的临界状态,在此阶段,如果采用适当的治疗仍可以恢复到正常状态,且可以收集相关信息以获取疾病的早期预警信号[4]。Chen等[4]提出了动态网络生物标志物(dynamic network biomarkers,DNB)的概念,发现在疾病前状态基因调控网络的1个子模块代表疾病的信息行为,并推导出了基于网络的3个动态标准。生物标志物是生物生理状态的指标,通常用于检查生物学或医学中的器官功能或疾病状态。大多数传统的生物标志物[5~6]是根据疾病状态和正常状态之间信息的差异表达来识别,其目的是将疾病状态与正常状态区分开,而不是预测疾病状态。确定临界点或疾病前状态的生物标志物是医学和生物学中的一个重要挑战,除了可在网络层面了解复杂疾病的分子机制之外,还可以尽早预防和治疗疾病[7]。
前述DNB方法已被多个研究小组应用于复杂疾病和生理过程的分析[8~11]。尽管该方法可以检测复杂疾病的临界状态,但在同一疾病状态下需要多个样本数据,这限制了其临床应用。Liu等[12]提出的单样本动态网络生物标志理论只需要1个病例样本,将其他正常样本作为参考样本,更适用于临床应用。此外,Liu等[13]建立了单样本“landscape”动态网络生物标志物,其将单个样本的分子网络图转化为拓扑图构建的模型可用于预测疾病的早期预警信号。赵宏倩等[14]对乳腺癌数据的所有基因计算局部得分,并通过得分排序选取关键基因,没有使用聚类或其他启发式算法。本文基于小鼠急性肺损伤的高通量表达数据,采用单样本动态网络生物标志物的方法构建复合指标,检测疾病的早期预警信号,模块化样本特异性网络,最大化DNB得分,找出了疾病的临界状态。
本研究的数据来自GEO(Gene Expression Omnibus)数据库(https://www.ncbi.nlm.nih.gov/geo/)。根据以下3个筛选条件得到数据集GSE2565:1)来源于人体组织或动物;2)具有多个时间序列的基因表达谱数据;3)同一时间节点既有患病样本的数据,也有正常样本的数据作为参考。数据来源实验是将雄性小鼠全身暴露于空气中或32 mg/m3的光气中20 min,肺组织在暴露后0.5 h、1 h、4 h、8 h、12 h、24 h、48 h、72 h 收集,以确定光气暴露后基因表达的变化。每个采样点有6个病例样本和6个对照样本[15]。使用R3.2.5(http://www.R-project.org/)对下载的数据进行预处理(矫正、标准化及表达值计算),共得到13 662个基因的表达数据,利用Limma程序包[16]筛选每个时刻的差异表达基因用于后续分析。
识别单个样本的临界状态需要对照样本组作为参考。通常,正常样本可用作参考样本,其表达谱可用作参考数据集,将参考样本和病例样本在每个时间节点的基因表达谱数据进行比较,可以提取单个样本的信息。
1.2.1 构建样本特异性网络
给定n个参考样本,参考样本数据中基因x和y之间的皮尔逊相关系数(Pearson correlation coefficient,记为P)可以计算为:
其中,xi和yi是第i个样本中基因x和y的表达值,和是所有样本中基因x和y的平均表达值,Pn(x,y)是n个参考样本中两个基因(x,y)之间的相关性。
将新的单个样本s添加到参考样本中后,可以基于总n+1个样本通过等式(1)重新计算两个基因之间的相关性(图1A),记为Pn+1(x,y)。两个特定基因(x,y)的单样本相关性(记为sP)定义如下:
sP(x,y)是由添加到参考样本中的新单个样本引起的,因此它表征了该单个样本与参考样本的特定相关信息。由于P遵循正态分布,因此等式(2)中的sP同样遵循具有n个共同样本的差分正态分布,通过单样本网络理论[17]可以准确评估sP的统计显著性。具体地,可以通过等式(3)为每个sP计算“Z”分数,并且可以基于“Z”分数从标准正态累积分布近似地获得每个sP的P值[18]。
通过等式(2)计算差异基因中所有基因对之间的相关性sP,通过等式(3)评估其显著性,如果它们的sP显著,则认为在两个基因之间有边缘连接。以差异基因为顶点,sP的值作为连接边的权重,构建单个样本特异性网络(图1B)。采用k-means聚类算法将网络分解为基于sP的多个模块(图1C)。
1.2.2 量化临界状态的复合指标
s-DNB组内基因的相关性可以建模为组内成员之间成对sP绝对值的均值:
s-DNB组内基因和组外其他基因的相关性可以建模为:
基于网络分解的模块,我们使用等式(7)来评估每一个模块:
等式(7)就是依据系统接近临界状态时出现3种现象来构造的[4]。根据等式(7)计算出来的指标值是每个时刻每一个模块的分值,选择每个时刻得分最高的模块作为当前时刻的候选s-DNB(图1D),并设置s-DNB的评分为α。基于DNB理论,所有时刻得分最高的点视为该样本的临界状态,该时刻对应的候选s-DNB模块为整个过程的s-DNB模块。
图1 识别单个样本中候选s-DNB的算法流程图(A)单个样本基因相关性计算;(B)单个样本特异性网络构建;(C)聚类;(D)复合指标计算。Fig.1 Algorithm flow chart for identifying candidate s-DNB in a single sample(A)Gene correlation calculation of a single sample;(B)Construction of a single sample-specific network;(C)Clustering;(D)Calculating composite indicators.
在每个采样点,对当前时刻基因做差异分析,设置参数P为0.05,log2(FC)(FC:fold change)为0.5。8 个时间点分别获得了 53、135、721、896、439、475、714、625 个差异表达基因,取其并集共得到2 308个基因,差异基因的火山图见图2。通过对差异基因进行1.2中的操作,最终确定急性肺损伤的临界状态在8 h,由于个体差异,不同样本的s-DNB不完全相同,根据试验的经验值,设置s-DNB的评分为α=4.0,分别得到6个病例样本的s-DNB模块。生物实验(GSE2565)发现,最突出的生理影响发生在接触后8 h内,肺水肿增加,最终存活率下降。在因氯化碳吸入引起急性肺损伤的小鼠中,12 h后观察到50%~60%的死亡率,24 h后观察到60%~70%的死亡率[15]。具体地,图3显示了复合指标F的变化情况,6个病例样本在临界状态的得分分别为6.333 3、5.853 7、5.294 6、5.739 5、4.592 8、4.979 8。从图 3 可以看到,第 4个时间点(8 h)的值最大,且超过给定阈值4.0,表明基于s-DNB的预测与实际疾病发展一致。
图2 差异基因的火山图P<0.05,log2(FC)>0.5,红色表示上调基因,蓝色表示下调基因。Fig.2 Volcano map of differential genesP<0.05,log2(FC)>0.5,red means up-regulated genes,and blue means down-regulated genes.
图3 复合指标变化曲线横坐标表示肺组织暴露时间,纵坐标表示复合指标,6条折线代表6个病例样本复合指标的变化情况。Fig.3 The changing curves of composite indicatorThe abscissa represents the lung tissue exposure time,the ordinate represents the composite index,and the six broken lines represent the changes in the composite indexes of six case samples.
s-DNB是代表疾病从正常状态过渡到疾病状态的重要网络,因此它们与发病机理中涉及的基因相关联。本研究采用生物信息学数据库DAVID(https://david.abcc.ncifcrf.gov/)[19]中的 GO注释和KEGG通路分析来研究s-DNB的生物功能行为。由于s-DNB模块越大,包含的冗余基因可能越多,所以我们对个数为130的s-DNB模块进行了GO分析,发现所识别的s-DNB基因与炎症反应、细胞趋化性、细胞增殖、调亡的负调控、细胞黏附连接、氧化应激反应等有关,具体见表1。文献研究表明:细胞凋亡的失调在急性肺损伤和其他相关疾病的发生中起着至关重要的作用[20~21];急性肺损伤是肺部炎症反应的广泛表现[22]。KEGG通路分析结果显示,s-DNB中的基因与细胞衰老、凋亡、免疫、氧化应激反应等有关,具体见表2。文献研究表明:内质网中的部分蛋白质加工有助于肺纤维化的发生[23];p53通路基因在癌症易感性位点中显著富集,肿瘤易感基因通常在癌症中发生突变[24]。
表1 GO功能分析Table 1 GO function analysis
表2 KEGG通路分析Table 2 KEGG pathway analysis
为了分析s-DNB的动力学分子机制,我们进一步构建了s-DNB模块的蛋白质-蛋白质相互作用(protein-protein interaction,PPI)网络图(图 4)。PPI网络从系统的角度论述了疾病的分子机制,网络中包含104个s-DNB基因节点和300条相互作用关系。使用CytoHubba插件中的MCC(maximal clique centrality)算法计算网络中每个节点的最大团中心性,筛选出最大团中心度排名前10 的关键基因:HSPA5(heat shock protein 5)、HSPA9(heat shock protein 9)、HSPA1A(heat shock protein 1A)、HSPA1B(heat shock protein 1B)、HSPB1(heat shock protein 1)、HSPB8(heat shock protein 8)、HSPH1(heat shock 105 kDa/110 kDa protein 1)、HSP90AA1[heat shock protein 90,alpha(cytosolic),class A member 1]、HSP90AB1[heat shock protein 90 alpha(cytosolic),class B member 1]、DNAJB1[DnaJ(Hsp40)homolog,subfamily B,member 1],这10个关键基因的热图展示见图5,可以看到所选基因在病例样本中的表达值均高于参考样本的表达值,说明这些基因在疾病的发展进程中起着正调控的作用。
图4 PPI网络图颜色表示节点度的变化,颜色越红,度越大。Fig.4 PPI network diagramThe colors represent the change in the degree of nodes.The redder the color,the greater the degree.
图5 10个关键基因在所有样本中的热图Fig.5 Heatmap of 10 key genes in all samples
HSP90AA1和HSP90AB1同属于HSP90家族,在信号转导、蛋白质折叠、蛋白质降解和形态演变中具有关键作用[25]。HSPA5、HSPA9、HSPA1A和HSPA1B都是HSP70家族成员。HSPA5参与内质网中蛋白质的折叠和组装,Shen等[26]发现,HSPA5与ATF6相结合在响应内质网应激时彼此分离,并通过抑制高尔基体定位信号和细胞分裂来保留内质网中的ATF6,在内质网应激期间HSPA5的解离使ATF6可以转运至高尔基体。HSPA9主要位于线粒体,也存在于内质网、质膜和细胞质囊泡中,在细胞增殖、应激反应和线粒体维持中起作用。HSPA1A和HSPA1B又称HSP70-1和HSP70-2,研究表明二者的同时消耗可抗癌细胞增殖[27]。Choi等[28]发现HSPB1可以抵抗细胞压力,并与癌症进展和肺纤维化有关;Li等[29]分析了HSPB1基因多态性与肺癌患者放射性损伤风险之间的关系,发现HSPB1的RS2868371基因型可能与辐射引起的食道损害有关。HSPB8在多种癌症中发挥作用:通过激活ERK-CREB途径促进癌细胞的生长,并可能成为胃癌患者的潜在预后因素[30];通过抑制磷酸肌醇3-激酶(PI3K)/AKT途径减少肝癌细胞的迁移[31];可以调节乳腺癌细胞的增殖和迁移[32]等。Liang等[33]的研究表明,用白介素-1β抗体或HSPH1抑制剂治疗可减轻急性肺损伤大鼠的肺损伤;Lenna等[34]发现DNAJB1与肺动脉高压的严重程度(通过肺动脉压测量)呈正相关。这些已有的研究进一步验证了我们方法的有效性。
与一般生物标志物相比,动态网络生物标志物更适合于表征系统状态的转移。本研究基于差异关联的信息预测疾病状态,构造用于量化单个样本的早期预警信号。我们使用肺组织暴露在光气和空气中的小鼠急性肺损伤生物数据,基于差异基因之间的差分皮尔逊相关系数构建样本特异性网络并模块化,根据系统达到临界状态分子网络的3种变化构建早期预警信号,识别出了s-DNB模块,找到了疾病的临界点。指标变化显示,在8 h系统达到了临界状态,这与原始数据[暴露后4~12 h,光气暴露小鼠的谷胱甘肽S-转移酶(glutathione S-transferase,GST)水平明显高于空气暴露小鼠,且小鼠在接触光气8 h内,肺水肿增加,最终存活率下降]相吻合。针对关键基因的功能分析、PPI网络分析、热图展示以及相关文献报道,都进一步验证了我们方法的有效性。此外,文中用于检测临界状态的s-DNB评分的阈值是本研究中的经验值,但它对结果没有显著影响。由于个体差异性,同一疾病的每个个体都有不同的DNB,要确定每种疾病所有个体的共同临界阈值可能需要整个群体的数据,以系统有效的方式识别s-DNB阈值是我们未来的重要工作。