基于WGCNA鉴定山羊金黄色葡萄球菌型乳腺炎关键应答基因

王丽敏，靳丽娜，刘亚芳，闫晓茹，郑惠玲

(西北农林科技大学 动物科技学院，陕西 杨凌 712100)

乳腺炎是山羊的常见疾病，具有发病率高、治愈率低、高产奶山羊泌乳盛期常见的特点。该病导致母羊泌乳量和乳品质下降，使母羊乳腺机能减退、使用年限缩短，养殖生产经济损失严重[1]。山羊乳腺炎多由病原微生物入侵引起，金黄色葡萄球菌(Staphylococcusaureus)是主要菌种之一[2]，其所致乳房炎常为慢性，有时也可是急性[3]。乳腺组织由乳腺上皮细胞、基底层结缔组织、乳腺腺泡腔、中性粒细胞、巨噬细胞、淋巴细胞等多种细胞组成。当乳腺组织受到金黄色葡萄球菌感染时，机体的免疫系统使炎症组织中白细胞增多[4]，同时损伤脱落的上皮细胞进入乳中[5]，最终使乳中体细胞数增加。

随着转录组测序技术的发展，以基因芯片为代表的基因表达数据已成广大研究者的研究热点之一，而利用加权基因共表达网络分析(weighted gene co-expression network analysis, WGCNA)等生物信息学方法分析基因芯片测序数据，可在短时间内筛选出与疾病相关的关键基因，对解析疾病发生发展的内在机理发挥着重要作用[6-7]。因此，本研究采用WGCNA对乳腺感染金黄色葡萄球菌24 h山羊的乳汁中体细胞的转录组数据集综合分析，并结合基因本体(gene ontology, GO)富集、蛋白互作分析(protein-protein interaction, PPI)等方法，鉴定金黄色葡萄球菌感染山羊乳腺后的关键应答基因，为揭示山羊金黄色葡萄球菌型乳腺炎发生发展的分子机制，为寻找新的山羊临床乳腺炎的防治方法提供新的理论参考。

1 材料与方法

1.1 数据的获取及预处理

本研究所使用的原始数据来源于NCBI的GEO数据库中GSE33894[8]数据集(https://www.ncbi.nlm.nih.gov/geo)。该数据集是以金黄色葡萄球菌侵染10只健康法国高山山羊的乳房，分别于感染前及感染后24 h采集乳汁体细胞样本的微阵列转录组测序数据。基于平台文件GPL14856对原始数据集进行注释，将探针名转换为基因名，对于2个探针对应同1个基因名的情况取平均值，对于缺失值采用K临近(k-nearest neighbor，K-NN)算法进行补充。

1.2 共表达模块的构建

使用R语言WGCNA包[9]的hclust函数对样本进行聚类，剔除离群样本。使用pick-Soft Threshold函数计算，以无标度拓扑拟合指数(R2)大于0.8且网络具有较好的平均连通性为条件选择适当的软阈值。使用adjacency函数将原始矩阵转换为拓扑重叠矩阵(TOM)，依据TOM矩阵相异度(dissTOM)，使用hclust函数对基因进行聚类，得到多个基因共表达模块。采用主成分分析法计算每个模块的特征值(module eigengene，ME)，根据ME对模块聚类，对相似度大于0.75的共表达模块进行合并。

1.3 目标模块的选择与分析

将合并后的模块与表型信息矩阵进行关联分析，选择与感染阶段相关的模块为目标模块(依据相关系数r及P值进行选择，P<0.05为差异显著)。计算目标模块内每个基因的模块身份(moudle membership，MM)和基因显著性(gene significance，GS)的关系，验证所选模块的可靠性。选择目标模块内基因连接度排名前30的基因为枢纽基因(hubgene)。通过DAVID[10]网站对目标模块基因进行GO功能富集分析。

1.4 关键应答基因的选择

通过String在线网站[11]对目标模块基因构建PPI网络，利用Cytoscape[12]软件的CytoHubba插件计算PPI网络中节点的度值(degree)，选择排名前30的基因为核心基因(core gene)，并与模块内的枢纽基因取交集得到关键应答基因。

2 结果与分析

2.1 共表达模块的构建

数据预处理后得到8 032个基因用于后续分析。对20个样本进行聚类分析，未发现离群样本。以R2>0.8且具有一定连通性为标准选择软阈值为5，根据ME进行模块聚类(图1)，图中红色线的高度是0.25，在红色线条以下的模块是相似度较大需要合并的模块。以相似度大于0.75合并共表达模块，得到11个模块(图2)，图中每个颜色代表1个模块，第1行的颜色表示首次聚类的结果，第2行的颜色是模块经过合并后的结果。其中灰色模块的基因属于共表达趋势不明显的基因，单独形成1个模块。

图1 模块特征基因的聚类

图2 基因动态剪切聚类树

2.2 目标模块的选择与分析

依据每个模块的特征值(ME值)，将模块与表型信息矩阵进行关联分析，得到每个模块与不同感染阶段的相关性及P值(图3)，选择与感染24 h相关性最强且P<0.05的blue模块(r=0.69，P=7×10-4)和green模块(r=0.68，P=0.001)为目标模块。

图3 模块与感染阶段关联分析

GO富集是对基因在不同维度和不同层次的描述，包括3部分：生物过程(biological process，BP)、细胞组分(cellular component，CC)和分子功能(molecular function，MF)。分别表示基因参与的生物学过程、基因存在的位置及基因在分子层面的功能。通过对所选目标模块的基因进行GO(gene ontology)富集(每部分富集结果取基因count值前5展示，不足5个的全部展示)。

结果显示，金黄色葡萄球菌感染24 h的blue模块基因主要具有蛋白同源二聚化、相同蛋白识别、ATP结合等分子功能；green模块基因主要参与细胞内蛋白水解、泛素依赖的蛋白酶体降解过程，主要具有泛素结合、SMAD蛋白结合、分子伴侣HSP90蛋白结合等分子功能。

图4 Blue模块GO富集

图5 green模块GO富集

2.3 关键应答基因的筛选

利用Cytoscape软件构建核心基因的PPI网络见图6、图7，图中红色越深表示基因连接度越高。依据目标模块基因连接度筛选出的枢纽基因表示基因在转录水平发挥关键作用，依据PPI网络筛选出的核心基因代表基因在蛋白水平发挥关键作用[13]。而二者的交集基因(图中菱形表示的基因)即为本研究所选择的关键应答基因。因此，山羊乳腺感染金黄色葡萄球菌24 h的4个关键应答基因为NFKB1、MRPL16、TMED10和ERLEC1。

图6 Blue模块PPI网络核心基因互作

图7 Green模块PPI网络核心基因互作

3 讨 论

金黄色葡萄球菌感染山羊乳腺组织的乳汁体细胞转录组数据来源于GSE33894数据集，Cremonesi等[8]于2012年通过利用Limma函数对数据进行差异表达分析，发现金黄色葡萄球菌感染山羊乳腺组织24 h后乳汁体细胞中有300个基因呈现差异表达，其中251个基因表达上调，49个基因表达下调。本研究以样本数据基因表达量为基础，利用WGCNA技术对数据构建基因共表达网络，使表达趋势相同的基因聚为同一模块，表达趋势不同的基因聚为不同模块，将模块与表型(感染时间)相关联后，结合PPI方法，依次获得模块核心基因和模块关键基因，最终交集获得山羊乳腺感染金黄色葡萄球菌24 h的4个关键应答基因(NFKB1、MRPL16、TMED10和ERLEC1)。

NF-κB是真核生物细胞调节炎症和对外来刺激免疫应答的重要响应蛋白[14]，对维持免疫系统的稳态十分重要。NF-κB1和RelA是NF-κB家族的部分成员[15]。IKBKE是NF-κB激活因子，其可通过解除IκB对NF-κB的抑制，使NF-κB入核并激活下游靶基因的转录[16]。王晋鹏等[17]研究表明乳腺炎发生后，细胞通过激活NF-κB信号通路启动机体的先天免疫防御。NF-κB1是NF-κB信号通路的关键节点，是炎症反应和免疫反应的热门靶点[18]。李蕊等[19]研究证实NF-κB1是奶牛金黄色葡萄球菌型乳腺炎的关键抗性候选蛋白。杨星哲等[20]研究发现NF-κB1是丙酸杆菌诱导的兔耳痤疮炎症模型的关键应答因子。他们的研究结果与本研究结果相符。提示虽然感染菌种属和受菌物种存在差异，但NF-κB1在多种炎症免疫模型中具有重要作用。

蛋白质在细胞生长发育、增殖分化、物质代谢等方面有极其复杂的调控过程。本研究中green和blue 2个目标模块是与金黄色葡萄球菌感染24 h最为相关的2个模块，GO富集分子功能中，大多数模块基因聚焦于蛋白同源二聚化活性、同质蛋白质结合、ATP结合、泛素结合、Hsp90蛋白结合、SMAD结合等分子功能。这正说明金黄色葡萄球菌侵染乳腺后，细胞调动并激活了多种蛋白信号传导途径以启动免疫防御。蛋白质泛素化通过将泛素分子递续连接到目标蛋白上，对蛋白进行降解，调节细胞增殖、凋亡和免疫信号传导等过程。Ebstein等[21]通过炎症诱导物刺激树突状细胞后，许多泛素相关基因表达上升。而高慧杰[22]的研究中，以生物信息学挖掘到泛素化相关的基因(UBE2L6、UBA7、USP18、ISG15)，是奶牛金黄色葡萄球菌乳房炎的关键应答基因。本研究中一些泛素相关基因(USP5、UBE2M、UBE2D1、UBE2S)，在金黄色葡萄球菌感染山羊乳腺24 h的green模块中也被鉴定到。同时green模块GO富集到的生物过程集中于细胞蛋白水解过程、泛素依赖的蛋白酶解过程。表明细胞中蛋白质降解系统的激活与重塑在炎症反应中发挥了重要作用。

钙结合蛋白S100A8在炎症相关研究中已广见报道，通常被认为是启动炎症通路、促进炎症发展进程的关键炎症介质[23]。在Cremonesi等[8]的研究中，S100A8被定义为金黄色葡萄球菌感染山羊乳腺后的抵抗因子和促炎介质。同时其也是本研究中感染金黄色葡萄球菌24 h最相关的blue模块的枢纽基因，但该基因并未出现在PPI网络中，可能是因为在炎症条件下S100A8mRNA表达和蛋白表达水平存在差异。关键基因MRP16编码线粒体核糖体蛋白，在恶性肿瘤中通常高表达，具有促进癌细胞增殖、正向调控恶性肿瘤的功能[24]。但MRP16在乳腺炎相关研究中鲜有报道。关键基因TMED10编码跨膜转运蛋白，通过内质网-高尔基体系统参与蛋白质的囊泡运输，在早期胚胎发育及免疫反应中发挥重要作用[25]。关键基因ERLEC1编码常驻内质网蛋白，具有识别n-多聚糖的功能，在乳腺炎相关研究尚未见报道。

乳腺炎的发生发展涉及多通路的调控网络，其内在机理繁多复杂。目前，关于乳腺炎机制的研究报道较多，但本次鉴定出的一些关键基因在乳腺炎领域鲜有报道，其在炎症发展进程的分子作用机制尚待探究。

4 结 论

本研究通过利用WGCNA、GO富集、PPI方法对样本数据进行综合分析发现，金黄色葡萄球菌感染山羊乳腺24 h，多种蛋白参与细胞信号转导(包括蛋白泛素化、分子伴侣蛋白HSP90结合、受体信号转导蛋白Smad结合、蛋白二聚化等)；4个关键基因NFKB1、MRPL16、TMED10、ERLEC1在炎症发展进程中发挥核心调控作用。