基因共表达网络分析及蛋白质相互作用筛选核分裂周期蛋白80与肺腺癌分期及预后的相关性

杨芳，李济民，袁伟奇，罗以勤

肺癌作为全球发病率和死亡率最高的恶性肿瘤之一[1]，严重威胁着人类的生命与健康。肺癌分为小细胞肺癌和非小细胞肺癌，其中肺腺癌是非小细胞肺癌最常见的组织学类型之一，约占40%[2-3]。尽管在肺癌的诊治方面取得了不错的进步，但是肺腺癌的预后仍然较差，5年生存率较低[4-5]。所以，寻找更有效的生物学标志物对改善肺腺癌病人的生存和预后意义深远。

随着基因芯片和RNA测序技术的发展，生物学信息分析方法已经被广泛用于识别与癌症相关的生物标志物[6]。具有类似表达模式的基因在功能上可能是相关的，然而这些方法仅仅是筛选出差异基因，却忽视了基因间可能的相关性。加权基因共表达网络分析（weighted gene co-expression network analaysis，WGCNA）可以通过构建基因共表达网络来探索不同基因模块与临床特征之间的相关性[7]。

本研究采用WGCNA构建肺腺癌基因共表达网络，识别共表达模块，并根据相应临床信息探讨其与基因模块间的相关性，筛选出与肺腺癌临床特征高关联的模块，最后对这个模块基因编码蛋白进行蛋白质相互作用（protein-protein interaction，PPI）分析，筛选出枢纽基因核分裂周期蛋白80（nuclear division cycle 80，NDC80）。

1 资料与方法

1.1 差异表达基因的筛选收集TCGA数据库（https：//cancergenome.nih.gov/）中497例肺腺癌癌组织和54例正常肺组织的RNA-seq表达矩阵及相应的临床信息。利用R语言的DESeq2软件，对整理好的RNA-seq数据进行基因表达差异分析，以对数化表达变化倍数（|log2FoldChange|）>1，及校正后的P<0.05作为筛选条件，得到差异表达基因（differentially expressed genes，DEGs），并制作火山图。

1.2 WGCNA

1.2.1 构建基因共表达网络 在R语言下运行WGCNA包建立基因共表达网络，剔除离群值，保证此网络的稳定性。依据无尺度网络的标准挑选出合适的软阈值β[8]，并应用β将相关矩阵转为邻接矩阵，然后将邻接矩阵转换成拓扑矩阵（topological overlap matrix，TOM），并计算1-TOM值。接着采用动态剪切树法进行模块的识别，最小模块尺寸设置为20。最后采用层次聚类法对特征向量基因（module eigengene，ME）进行聚类分析，同时制作出模块树状图。

1.2.2 鉴别与临床信息相关的模块 计算各模块的ME与外部每个临床特征信息的Pearson相关系数及P值，绘制模块-特征的相关性热图，选取Pearson相关系数最高的基因模块作为目的基因模块。

1.3 枢纽基因的筛选将目的基因模块里的所有基因编码蛋白上传到STRING数据库（https：//stringdb.org/），建立PPI网络。综合得分>0.4作为截止标准。并采用Cytoscape软件对PPI网络进行可视化。最后使用其中的MCODE插件鉴定出枢纽基因（seed基因）。

1.4 枢纽基因的表达水平与肺腺癌分期、预后的关系使用基于TCGA的GEIPA数据库（http：//gepia.cancer-pku.cn/）分析枢纽基因在肺腺癌组织中的表达水平及其与肿瘤分期、预后的关系。

1.5 统计学方法统计分析均在R（3.5.3）语言上进行，其中使用的R语言程序包包括DESeq2 1.26.0、WGCNA 1.69。枢纽基因的表达差异及其表达量与分期的关系采用单因素方差分析法；在肺腺癌中表达量与预后的关系采用Log-rank检验，以P<0.05表示差异有统计学意义。

2 结果

2.1 差异表达基因的筛选使用R软件读取芯片数据得到13 289个肺腺癌相关基因，以（|log2FoldChange|）>1，校正后的P<0.05 为筛选标准，得到11 904个DEGs，包括1 112个上调基因和792个下调基因。DEGs火山图见图1。

图1 肺腺癌差异表达基因火山图

2.2 基因共表达网络的建立在R语言下运行WGCNA包建立基因共表达网络，发现无需剔除离群值，故将1 904个DEGs全部用于WGCNA。考虑到无尺度网络的构建及平均连接度的适度保留，本研究选择β=5构建共表达网络（相关系数等于0.98作为标准），见图2。最终得到6个基因模块，见图3。

图2 WGCNA软阈值（β）的确定：A为不同β下计算的无尺度拟合指数分析；B为不同β下计算的平均连接度；C为当β=5时连接度分布直方图；D为当β=5时无尺度拓扑检测

图3 基因聚类树状图

2.3 鉴别与临床信息相关的模块根据各个模块的特征向量，利用R语言下的WGCNA包分别计算这些模块与肺腺癌样本临床信息（肿瘤病人的生存时间、生存状态、年龄、性别及肿瘤的TNM分期及远处转移）之间的相关性。结果表明，棕色模块与肿瘤病人生存时间的相关性较高（P=9e-04），蓝色模块与肿瘤病人性别的相关性较高（P=7e-04），见图4。

图4 肺腺癌不同临床特征与模块特征值相关性热图

2.4 枢纽基因的筛选为了筛选出枢纽基因，本研究将棕色模块里的128个基因编码蛋白和蓝色模块里的517个基因编码蛋白分别上传到STRING数据库构建PPI网络。综合得分>0.4时，棕色模块里的基因编码蛋白构建了包含46个节点和48个边的网络，蓝色模块里的基因编码蛋白构建了包含290个节点和2347个边的网络。考虑到棕色模块的PPI网格节点的度数普遍≤5，所以蓝色模块被确定为枢纽模块。最后使用Cytoscape软件中的MCODE筛选出蓝色模块中的核心模块进行后续分析，并发现NDC80是其枢纽基因（Seed基因）。

2.5 NDC80基因的表达与肺腺癌预后的关系。本研究使用GEIPA数据库验证结果显示，与正常肺组织样本相比，NDC80基因的表达水平在肺腺癌组织中明显升高（|log2FC|>1，P<0.05），且与肺腺癌的分期相关（F=3.58，P<0.05）；预后结果显示NDC80基因的高表达与肺腺癌的总生存期（Overall Survival，OS）较差明显相关（HR=1.6，P<0.05），见图5。

图5 GEPIA数据库中NDC80基因表达水平与肺腺癌分期及预后的关系：A为NDC80基因表达水平（P<0.05）；B为NDC80表达水平与肿瘤分期的关系（P<0.05）；C为NDC80基因表达水平与肿瘤总生存期的关系（P<0.05）

3 讨论

本研究利用基因共表达网络分析及蛋白质相互作用方法，筛选出肺腺癌差异表达基因NDC80，并通过GEIPA数据库验证NDC80与肺腺癌分期及预后的关系。

NDC80复合体是一种异四聚体蛋白复合物，对细胞有丝分裂至关重要[9]。复合体异常表达会引起染色体异常分离，导致染色体不稳定，而染色体不稳定是肿瘤细胞的一个共同特征，是肿瘤形成的重要机制[10]。NDC80作为复合体的重要组成部分之一，其表达量在S期至M期显著升高，且在G2/M期发挥着重要的调节作用[11]。研究表明NDC80基因在肝癌、乳腺癌、胃癌和胰腺癌组织中表达量均明显升高，且NDC80基因的高表达与胰腺癌的预后有关[12-15]。Yuan 等[16]筛选了 13 个 GSE数据集，发现NDC80等7个关键基因与ZW10基因相互作用，参与了肺癌的形成过程。本研究发现肺腺癌组织的NDC80基因的表达量明显高于正常肺组织，且与肿瘤的分期相关。提示NDC80基因可能在肺腺癌的发生发展中发挥着重要的作用，这与前面的文献报道是一致的。其次NDC80基因表达含量高的肺腺癌病人的总生存期明显较短，这暗示了NDC80基因有望成为肺腺癌预后的生物学标志物。

综上所述，本研究通过WGCNA构建共表达网络及PPI网络分析，筛选出与肺腺癌的分期及预后有关的基因NDC80，这些结论揭示了NDC80基因可能在肺腺癌病人的发生发展过程中发挥着重要的作用。