基于加权基因共表达网络识别结直肠癌预后相关的signature基因

目录中文摘要 IAbstract III1、文献综述 51.1研究背景和国内外研究现状 51.2研究内容和研究意义 62、材料与方法 82.1数据基础 82.1.1TGCA数据库 82.1.2GEO数据库 82.2结直肠癌基因的差异表达分析 92.2.1TCGA结直肠癌的数据预处理 92.2.2识别差异表达基因 92.2.3差异表达基因的火山图与热图 102.2.4差异表达基因的GO功能注释和KEGG通路富集分析 102.3结直肠癌加权基因共表达网络分析 112.3.1WGCNA理论依据 112.3.2数据格式处理，去除缺失值与离散样本 122.3.3构建加权基因共表达网络 122.4筛选关键模块 122.4.1筛选显著富集到差异表达基因的模块 122.4.2模块基因的功能与通路富集分析 122.5筛选关键模块中的signature基因 132.5.1临床信息及表达谱数据处理 132.5.2单因素Cox回归分析 132.5.3Lasso回归分析 142.6构建风险评分模型并进行Kaplan-Meier生存分析 142.7绘制ROC曲线 142.8外部数据集验证 152.8.1外部数据集的下载与预处理 152.8.2外部数据集的KM生存分析 152.8.3外数据集的ROC曲线 152.9关键signature基因的数据库验证 163、结果 173.1结直肠癌的差异表达分析 173.1.1识别结直肠癌差异表达基因 173.1.2差异表达基因的功能注释及通路富集分析 193.2结直肠癌的加权基因共表达网络 203.2.1筛选基因与样本 203.2.2选择软阈值 203.2.3识别加权基因共表达模块 213.2.4分析加权基因共表达网络 223.3筛选关键模块 253.3.1基于超几何分析筛选模块 253.3.2识别关键模块 253.4识别关键signature基因 263.4.1基于单因素Cox筛选基因 263.4.2基于Lasso回归进一步筛选基因 263.5构建风险评分模型并进行KM生存分析 283.5.1风险评分模型 283.5.2高低风险组的KM生存分析 283.5.2ROC曲线 293.6外部验证集验证 303.6.1外部验证集的获取 303.6.2外部验证集的KM生存分析 303.6.3外部验证集的ROC曲线 313.7关键signature基因的数据库验证 314、讨论 335、结论 346、致谢 357、参考文献 36中文摘要结直肠癌在中国是最为常见的一种恶性肿瘤，是第三大最常诊断的癌症。虽然分子和临床生物标志物已被用于指导临床决策，但是对结直肠癌从早期发展到转移相关的潜在分子生物标志物和治疗靶点的深入了解大多尚未完成。本课题通过TCGA数据库获得结肠癌和直肠癌表达谱和样本的临床信息并进行整理，共获得了383例癌症样本、51例正常样本。对count表达谱数据进行TPM标准化及数据预处理，并将所得的基因再比对回count值进行差异表达分析。对结直肠癌表达谱进行整理后一共获得17348个基因的表达值，进行差异表达分析后得到共5112个差异表达基因。对癌症患者的TPM数据格式的编码蛋白质的基因表达谱，进行了结直肠癌的加权基因共表达网络的构建，共获得了29个模块。将所获得的29个模块和所获得的5112个差异表达基因，进行超几何分析，从而筛选模块，再对筛选后的模块进行通路与功能的富集分析，将其中的4个显著富集到结直肠癌相关通路的模块作为关键模块，筛选出来。通过使用Cox回归和Lasso回归，我们可以从总体中挑选出17个具有特定signature的基因，并利用它们的Coef结果，来估算各种样本的风险程度。通过将样本按照Kaplan-Meier、ROC等指标进行评估，我们发现，较高的风险水平往往会导致较差的预后。这些预后情况也被Kaplan-Meier、AUC的结果所支持。最后，通过GEO数据库下载GSE103479作为验证集，以此来评估此实验获得的关键signature基因对结直肠癌潜在的治疗以及预后价值的可靠程度。验证集GSE103479的Kaplan-Meier生存曲线与ROC曲线显示了关键signature基因对结直肠癌有潜在治疗和预后价值。综上，我们共识别到结直肠癌中17个关键signature基因，并验证了这些基因显著与结直肠癌的治疗及预后相关，有助于为结直肠癌患者的临床治疗和预后分析提供思路。关键词：WGCNA；结直肠癌；单因素Cox回归；Lasso回归；风险回归模型Theweightedgeneco-expressionnetworkisemployedtopinpointthesignaturegenesassociatedwiththeprognosisofcolorectalcancerAbstractColorectalcancer(CRC)isthemostcommonmalignancyandthethirdmostcommonlydiagnosedcancerinChina.Althoughmolecularandclinicalbiomarkershavebeenusedtoguideclinicaldecisionmaking,anin-depthunderstandingofthepotentialmolecularbiomarkersandtherapeutictargetsassociatedwiththeprogressionfromearlytometastasisofcolorectalcancerhasmostlynotbeencompleted.Inthisstudy,weobtainedandsortedtheclinicalinformationofcolonandrectalcancerexpressionprofilesandsamplesthroughTCGAdatabase.Atotalof383cancersamplesand51normalsampleswereobtained.TPMstandardizationanddatapreprocessingwereperformedonthecountexpressionprofiledata,andtheresultinggeneswerealignedbacktothecountvaluefordifferentialexpressionanalysis.Atotalof17348geneexpressionvalueswereobtainedaftertheexpressionprofileofcolorectalcancer,andatotalof5112differentiallyexpressedgeneswereobtainedafterdifferentialexpressionanalysis.Theweightedgeneco-expressionnetworkofcolorectalcancerwasconstructedbasedontheprotein-codinggenesintheTPMdataformatofcancerpatients,andatotalof29moduleswereobtained.Theobtained29modulesand5112differentiallyexpressedgeneswereanalyzedbyhypergeometricanalysistoscreenmodules.ThenKEGGpathwayenrichmentanalysiswasperformedontheselectedmodules,andfourmodulessignificantlyenrichedincolorectalcancerrelatedpathwayswereselectedasthekeymodules.UnivariateCoxregressionanalysiswasperformedontheobtainedkeymodulestoscreensurvival-relatedgenes,andthenLassoregressionanalysiswasusedforfurtherscreening.Atotalof17keyfeaturegeneswerefinallyscreenedout,andthentheriskscoreofeachsamplewascalculatedaccordingtotheCoefvalueobtainedbyLassoregressionofthese17genes.Groupingthesamplesbymedianriskscore,Kaplan-MeiersurvivalcurvesandROCcurveswereplotted.Verificationofthelowriskscore'sassociationwithafavorableprognosisofcolorectalcancerwasprovidedbytheKaplan-MeiersurvivalcurveandAUCscoreresults.Finally,GSE103479wasdownloadedfromtheGEOdatabaseasavalidationsettoevaluatethereliabilityofthepotentialtherapeuticandprognosticvalueofthekeysignaturegenesinthisexperiment.TheKaplan-MeiersurvivalcurveandROCcurveofGSE103479inthevalidationsetshowedthatthekeysignaturegeneshadpotentialtherapeuticandprognosticvalueforcolorectalcancer.Insummary,weidentified17keysignaturegenesincolorectalcancer.Verificationofasignificantcorrelationbetweenthesegenesandthetreatmentandprognosisofcolorectalcancerhasbeenestablished,whichisofgreatassistanceinfurnishingideasfortheclinicalassessmentofcolorectalcancerpatients.Keywords:WGCNA;CRC;UnivariateCoxregression;Lassoregression;RiskregressionmodelINDEX\h"—M—"\o"S"\c"2"\z"2052"INDEX\o"S"\c"2"\z"2052"1、文献综述1.1研究背景和国内外研究现状结直肠癌是一种危及人类健康的严重疾病，它可能导致直肠部位的癌细胞扩散，并可能引起严重的肠道功能紊乱，如出血、排泄困难、疼痛。结直肠癌的发病率和死亡率在全球范围内居高不下，是目前世界上常见的癌症之一[1]。尽管结直肠癌的研究一直受到广泛关注，但由于其发病机制十分复杂，可能受到多种因素的影响，因此，深入探索结直肠癌的病因仍然是当前研究领域的一大挑战。目前，结直肠癌的诊断技术主要包括内镜检查、肿瘤标志物检测等方法。但这些方法有一定的局限性，因此需要进一步研究新的诊断技术。结直肠癌的诊断和治疗通常依赖于外科技术和药物治疗。不同患者对治疗的反应可能会有所不同，因此需要进一步研究治疗方法的个体化。因此个体化治疗是结直肠癌治疗的热点之一，目前，肿瘤免疫治疗在结直肠癌治疗中的应用备受关注。结直肠癌的发病过程涉及多个基因的突变，包括APC[2]、KRAS[3]、TP53[4,5]、BRAF[6]等。近年来，结直肠癌的治疗受到越来越多的关注，其原因不仅仅是由于其特异性的遗传特征，而且还包括其他的生物学功能。PD-1/PD-L1抑制剂[7]的应用，不仅能够改善结直肠癌患者的症状，而且还能够帮助他们获得靶向治疗的结果，从而推动结直肠癌的治疗取得重大突破。同时，CRISPR/Cas9[8]已被广泛应用于结直肠癌的研究中。通过基因编辑技术，可以模拟和研究结直肠癌的基因突变，从而探索可能的治疗方法和药物靶点。随着科技的进步，精准医疗已成为结直肠癌治疗的重要方向。通过基因测序、组学分析等技术，可以实现对结直肠癌患者的个体化治疗方案设计，从而提高治疗效果和生存率。总的来说，结直肠癌的研究涉及多个方面，从基础研究到临床应用都在不断地取得进展。这些研究成果将为结直肠癌的早期诊断、治疗和预后提供更加有效和个体化的方法和手段。高通量测序技术的快速发展，越来越多的结直肠癌相关数据产生出来，包括转录组在内的多种不同层次的数据。这些数据可以用来进行生物信息学分析，为我们深入了解结直肠癌的病理生理学、诊断和治疗提供了新的视角和思路。其中，转录组学研究主要包括差异表达基因、可变剪切、非编码RNA等方面。这些研究可以为结直肠癌的分子机制提供更加全面的认识。对结直肠癌转录组层面的研究为我们深入了解结直肠癌的病理生理学、诊断和治疗提供了新的视角和思路，同时也为结直肠癌的个体化治疗提供了可能性。WGCNA（WeightedGeneCo-ExpressionNetworkAnalysis）是一种有效的基因表达模式分析工具，它能够有效地将复杂的基因组织成若干个独立的模块，这些模块之间存在着密切的联系。WGCNA可以在整个转录组水平上捕获基因间的关联模式，并将这些模式映射到生物学过程和复杂疾病中，从而揭示相关的生物学机制。在结直肠癌研究中，WGCNA可以用于识别共表达的基因模块，通过WGCNA将基因聚类为若干共表达模块。进而结合超几何检验、功能富集分析和通路富集分析的方法，筛选出与结直肠癌发病机制相关的基因模块，发现与结直肠癌相关的潜在的生物学通路和调控机制，为深入探索结直肠癌的致病机制提供线索，并识别出关键模块用于后续识别关键signature基因。基于已有数据与方法，本课题使用WGCNA研究结直肠癌的转录组学数据，旨在寻找与结直肠癌进展相关的关键signature基因，通过对关键signature基因的研究，深入地了解与结直肠癌密切相关的关键signature基因的作用，从而为靶向用药提供潜在的药物靶点，为结直肠癌的治疗提供价值。1.2研究内容和研究意义本课题通过加权基因共表达网络的构建，识别结直肠癌进展相关的关键signature基因。首先对TCGA数据库中结直肠癌的表达谱以及样本的临床信息进行整理，选择蛋白编码基因作为研究内容，共获得17348个基因在433个样本中的表达值。同时，获得count表达谱对应的TPM格式的数据，对TPM格式的表达谱进行数据预处理。提取预处理后所得基因的count表达谱，利用R语言的DESeq2包对count表达谱进行差异表达分析，得到与结直肠癌显著相关(p<0.05)的差异表达基因。接下来，对预处理后的TPM格式的表达谱中的肿瘤样本数据进行加权基因共表达网络的构建，根据基因表达模式的相似程度进行聚类，再除去不能被归类的基因得到加权基因共表达网络模块。然后利用超几何分析筛选出能显著富集到差异表达基因的模块，再对模块进行功能富集分析和通路富集分析，根据富集出的功能与通路对模块进行进一步筛选，最终确定关键模块。通过采用单因素Cox回归方法，我们可以有效地检测基因在癌症进展中的表现，并且结合Lasso回归的结果，我们可以更加准确地确定最重要的signature基因。而后根据Lasso回归得到的Coef值作为系数，使用关键基因乘以系数相加之和,作为每个样本的风险评分。通过Kaplan-Meier生存曲线的研究，我们发现，在给定的风险水平下，不同的组别具有明显的生存优势。因此，我们可以把样本划分成两组，即高风险组和低风险组。此外，我们还从GEO数据库下载结直肠癌的表达谱，按照样本的风险评分分组并进行Kaplan-Meier生存分析、绘制ROC曲线，也进一步验证了识别到的关键signature基因对结直肠癌有潜在治疗及预后价值。对于结直肠癌的疾病进程而言，我们所识别的关键的signature基因在其中发挥着重要作用，对临床治疗提供方向。结直肠癌的研究和探索可以极大地改善结直肠癌的预后，同时还可以更好地理解其发生的原因，从而给结直肠癌的诊断和治疗带来更多的可能性。

2、材料与方法2.1数据基础2.1.1TGCA数据库TCGA数据库（TheCancerGenomeAtlas）[9]是一个重要的研究工具，它可以帮助我们更好地了解人类肿瘤的遗传特征，并为临床治疗提供有价值的信息。其中包含了来自多种类型肿瘤的临床信息、基因组学数据、转录组数据等信息。RNA-Seq是TCGA数据库中常见的一种数据类型，用于检测特定样本中所有RNA分子的表达情况。RNA-Seq数据可以提供大量的转录组信息，如基因的表达水平、可变剪切形式等。使用RNA-Seq数据可以帮助研究者了解肿瘤基因表达的变化，并探索不同基因、信号通路的作用和相互作用。通过对TCGA数据库中的RNA-Seq数据进行分析，可以探索肿瘤基因表达和信号通路的变化，为肿瘤研究和治疗提供有力的支持。2.1.2GEO数据库美国国家医学图书馆（NationalLibraryofMedicine）拥有一个庞大的数据库（GeneExpressionOmnibus）[10]，它汇集了来自世界各地的大量基因表达数据、芯片芯片数据、序列标记数据以及基因组数据，以满足科学家们对基因表达谱及表观遗传学的需求，并且能够有效地推动这一领域的发展。GEO数据库中包含了来自各种不同物种的基因表达数据，包括人类、小鼠、大鼠、果蝇等。GEO数据库提供了一个全面的视角，可以收集多种类型的基因表达数据，以深入了解不同物种在不同环境、病理和治疗条件下的表观遗传学变化，进而更好地掌握基因调控和信号传导的规律，从而更好地探索生物学过程及其相关的疾病机制。2.2结直肠癌基因的差异表达分析2.2.1TCGA结直肠癌的数据预处理从TCGA数据库中总共下载了433个结直肠癌样本数据集，并带有相应的临床和生存注释：TCGA-COAD，TCGA-READ。对于TCGA结直肠癌数据，需要进行一系列的数据预处理才能进行后续的分析。根据基因表达矩阵，使用TPM方法进行数据标准化，使不同样本之间的表达值可比较。从TCGA数据库下载结直肠癌的数据并进行整合，得到结直肠癌的表达谱文件、注释信息文件以及样本的临床信息文件。对结直肠癌的TPM格式的表达谱文件进行筛选，去除掉了表达不足30%的基因后，得到17348个蛋白质编码基因在433个样本中的TPM表达值。将这17348个蛋白质编码基因比对到count表达值，用于差异分析识别差异表达基因。此套数据共有382个肿瘤样本，51个正常样本。2.2.2识别差异表达基因结直肠癌的发生与许多基因的异常表达有关，因此差异表达分析是研究结直肠癌基因调控机制的重要手段。通过对差异表达基因的筛选和分析，可以进一步揭示结直肠癌的致病机制，为结直肠癌的诊断和治疗提供重要的理论支持。R语言的DESeq2包的DESeq2是基于负二项分布模型的差异表达分析方法，适用于RNA-Seq数据的差异分析，对小样本的数据也有较好的处理能力。通过应用R语言的DESeq2函数，我们可以比较结直肠癌患者和健康人的表达谱数据，并将其中的p值低于0.05，而log2foldchange值高于1，低于-1的基因作为差异表达的标志。在选出的差异表达基因中，log2foldchange值大于1的为表达上调的基因，小于-1的为表达下调的基因。识别差异表达基因是识别与结直肠癌进展相关的关键signature基因的一个重要步骤，它可以帮助我们发现在癌症样本与正常样本中表达存在差异的基因。对差异表达的基因进行进一步的注释和功能分析，例如，KEGG通路分析，以揭示结直肠癌生物学意义和与相关的生物学通路。火山图和热图是常用的差异表达基因分析结果展示方式。火山图展示了基因表达的差异水平和显著性，通常横坐标表示基因的表达变化水平（例如对数倍数变化），纵坐标表示基因的差异显著性（例如负对数的p值），可以根据显著性水平和方向将基因分为差异表达基因和非差异表达基因，在火山图中使用不同颜色进行展示。热图则通常用来展示基因表达的模式和趋势，可以将样本和基因进行聚类，展示它们的相似性和差异性。热图一般用颜色表示基因在不同样本中的表达量或相对表达量，可以快速识别表达量高低的差异，并确定表达模式或趋势。在热图中，通常也会使用颜色条示意图来表示基因表达量或相对表达量的范围和大小。火山图和热图结合起来可以更全面地展示差异表达基因的分析结果，从更直观的观察基因表达的变化，从而发掘其生物学意义。功能注释和KEGG通路分析是两种常见的富集分析方法。GO（GeneOntology）是一种广泛使用的基因功能分类系统，通过将基因的生物学功能描述为三个不同的方面（分子功能、细胞组分和生物过程），可以帮助我们更好地理解基因在生物学过程中的作用。GO功能注释可以帮助我们更好地理解基因表达的差异性，它可以将基因表达的差异性映射到GO数据库中的不同分类，从而更准确地评估出各个分类的富集程度。通过GO功能注释可以帮助我们了解差异表达基因在不同功能分类中的富集情况，进而推测差异表达基因在生物学过程中的作用和相关机制。利用KEGG（KyotoEncyclopediaofGenesandGenomes）这一全球性的基因表达资讯库，我们能够更加深入地探索差异表达基因的功能、表达水平及其对生物体内其他生物学行为的影响，从而更好地理解差异表达基因的分布特征，并且更有效地预测它们对生物体内其他生理行为的影响。GO功能注释和KEGG通路富集分析是一种有效的研究方法，它可以深入探究基因表达水平的差异，从而更好地理解基因的生物学意义及其相关功能。综合使用GO功能注释和KEGG通路分析可以帮助我们更全面地理解差异表达基因在生物学过程中的作用和机制，为后续的功能研究和生物实验提供重要的指导意义。2.3结直肠癌加权基因共表达网络分析2.3.1WGCNA理论依据WGCNA（WeightedGeneCo-ExpressionNetworkAnalysis）ADDINEN.CITEADDINEN.CITE.DATA[11]是一种基于基因共表达网络分析的方法，用于发现基因在生物学过程中的关联性及其与表型的相关性。其理论依据是基因共表达网络分析，该分析方法是通过基因表达量的相关性来构建一个基因共表达网络。这种方法在全基因组水平上对基因进行聚类和模块化，可以发现生物过程中的基因调控模块，并找到与特定表型相关的模块和关键基因。WGCNA中的权重指的是节点之间的相似度或相关性度量。当我们构建一个基因共表达网络时，我们需要确定每个基因之间的相互关联，并将这种关联转化为一种可以衡量它们相似程度的指标，比如Pearson和Spearman指标。然后，可以使用一种基于相似性的聚类算法（如层次聚类算法）将基因分组成模块。每个模块中的基因具有相似的表达模式，且这些模式在不同的生物条件下保持稳定。这种模块化方法可以有效地降低数据的维数，从而帮助我们发现复杂的基因调控网络。WGCNA方法可以应用于各种类型的生物数据，例如RNA-seq、microarray等。它可以帮助我们发现生物过程中的关键基因、功能模块和调节网络，为理解生物学机制提供了一种有效的工具。2.3.2数据格式处理，去除缺失值与离散样本在进行WGCNA分析之前，需要对数据的格式进行处理，读入表达矩阵数据和样本信息数据。提取TPM数据格式的癌症样本的表达数据转置，使得行代表基因，列代表样本。同时需要去除缺失值和离散样本。首先，对于缺失值的处理，可以使用R语言中的na.omit()函数将缺失值所在的样本或基因从数据集中删除。其次，对样本进行层次聚类，若存在离散样本，可以将其从数据集中删除，可以减少噪音和提高聚类的准确性。2.3.3构建加权基因共表达网络通过将R中的WGCNA包[12]加载到blockwiseModules函数中，可以构建一个加权的基因共表达网络，其中包括基因相关性矩阵、经过加权优化的相关性邻接矩阵以及拓扑重叠矩阵（TOM），以此来实现更有效的基因表达。相关性邻接矩阵只考虑了两个基因之间的线性关系，但两个基因之间还可以通过若干个中间关联的基因而形成联系，这个联系的强弱经过TOM的构建也会一定程度上被考虑进划分模块的依据中。构建网络后共得到29个基因共表达模块。2.4筛选关键模块2.4.1筛选显著富集到差异表达基因的模块在进行WGCNA分析后，可以将基因分成多个模块，并且每个模块都有一组共同表达的基因。接下来，需要将这些模块与差异表达分析得到的基因列表进行比较，以筛选出与差异表达基因显著富集的模块。使用超几何分布ADDINEN.CITEADDINEN.CITE.DATA[13]计算模块内基因为差异表达基因时的p值，选择p值小于0.05的模块作为显著富集到差异表达基因的模块。2.4.2模块基因的功能与通路富集分析通过利用富集分析，我们可以对基因的功能和通路进行深入的研究，以期更加准确地确定这些基因是否与结直肠癌有关联。该课题使用在线网站KOBAS()ADDINEN.CITEADDINEN.CITE.DATA[14]，进行基因列表富集及KEGG通路结果的可视化。该网站的Gene-listEnrichment功能可以将差异表达基因和经WGCNA分析并进行筛选后的基因模块进行KEGG与GO富集分析，以发现与结直肠癌相关的生物功能与KEGG生物学通路。该网站默认的统计学方法为hypergeometrictest/Fisher'sexacttest，使用多重检验校正，该网站默认为Benjamini-Hochberg校正，以此来确保结果的准确性。筛选P值小于0.01的生物功能与KEGG生物学通路，作为差异表达基因与基因模块主要参与的功能与通路。进而可以根据富集分析的结果来筛选关键的共表达基因模块。根据模块的功能与通路富集分析，进一步筛选出与癌症进展相关的关键模块，一共筛选出了black、pink、turquoise、yellow这4个模块，作为最终的关键模块进行下一步的关键基因的筛选。2.5筛选关键模块中的signature基因2.5.1临床信息及表达谱数据处理从TCGA数据库获取患者的临床信息，主要为生存状态、生存时间。整合样本的生存信息，加入表达谱中（行表示样本，列为临床信息与基因）。2.5.2单因素Cox回归分析通过单因素Cox回归分析，我们可以有效地筛选出与结直肠癌有关的重要signature基因，并且深入探索它们在疾病发生和发展过程中的生物学作用。使用单因素Cox回归分析来分析基因表达水平与肿瘤生存率之间的关系，并对筛选出来的差异表达基因进行生物信息学分析，以便更好地理解这些基因在肿瘤发生和发展过程中的作用。单因素Cox回归分析是一种常用的生存分析方法，可以帮助我们筛选和分析与疾病相关的关键基因和生物学机制，从而更好地理解和预测疾病的发生和发展。2.5.3Lasso回归分析Lasso回归用于特征选择和正则化的线性回归方法。我们通常需要从大量基因表达数据中筛选出与生物过程相关的关键基因。使用Lasso回归可以帮助我们识别那些与响应变量相关性最强的基因，进而加深我们对生物过程的理解。在Lasso回归中，参数的选择对筛选出的基因集合有很大的影响。一般情况下，我们可以使用交叉验证的方法来选择最优的参数值，从而得到最佳的筛选结果。2.6风险评分模型的构建与验证根据Lasso回归分析筛选的基因和这些基因所对应的Coef值，建立风险评分模型，将癌症样本根据风险评分的中值划分为高风险组和低风险组。KM（Kaplan-Meier）生存分析ADDINEN.CITE<EndNote><Cite><Author>Goel</Author><Year>2010</Year><RecNum>16</RecNum><DisplayText>[15]</DisplayText><record><rec-number>16</rec-number><foreign-keys><keyapp="EN"db-id="209sespa0d5zdaev0pqpxp5jf22apw2p5fet"timestamp="1684048254">16</key></foreign-keys><ref-typename="JournalArticle">17</ref-type><contributors><authors><author>Goel,M.K.</author><author>Khanna,P.</author><author>Kishore,J.</author></authors></contributors><auth-address>DepartmentofCommunityMedicine,PostGraduateInstituteofMedicalScience,Rohtak,Haryana,India.</auth-address><titles><title>Understandingsurvivalanalysis:Kaplan-Meierestimate</title><secondary-title>IntJAyurvedaRes</secondary-title></titles><periodical><full-title>IntJAyurvedaRes</full-title></periodical><pages>274-8</pages><volume>1</volume><number>4</number><keywords><keyword>Kaplan-Meierestimate</keyword><keyword>Survivalanalysis</keyword></keywords><dates><year>2010</year><pub-dates><date>Oct</date></pub-dates></dates><isbn>0974-925X(Electronic) 0974-7788(Print) 0974-7788(Linking)</isbn><accession-num>21455458</accession-num><urls><related-urls><url>/pubmed/21455458</url></related-urls></urls><custom1>ConflictofInterest:Nonedeclared</custom1><custom2>PMC3059453</custom2><electronic-resource-num>10.4103/0974-7788.76794</electronic-resource-num><remote-database-name>PubMed-not-MEDLINE</remote-database-name><remote-database-provider>NLM</remote-database-provider></record></Cite></EndNote>[15]是用于评估两组或多组患者生存差异的统计方法。它可以用于评估不同治疗方案、疾病预后、基因表达谱分类等方面的生存差异。KM生存分析基于生存函数（SurvivalFunction）和累积风险（CumulativeHazard）的概念。生存函数指的是在任意时间点上，某个患者或某个群体中有多少比例的患者还活着，而累积风险则指在某个时间点上，某个患者或某个群体中有多少比例的患者已经死亡或发生了事件。在KM生存分析中，首先需要根据样本的生存时间和事件情况（如死亡或疾病复发）计算生存函数和累积风险。然后，将样本按照某个变量（如治疗方案、疾病分期、基因表达谱等）进行分组，利用KM曲线绘制不同分组之间的生存函数曲线，比较不同组之间的生存时间和生存率差异。本课题对高、低风险组绘制KM生存曲线，以此来探索高低风险与疾病进展与预后的关系。ROC（ReceiverOperatingCharacteristic）曲线ADDINEN.CITE<EndNote><Cite><Author>Park</Author><Year>2004</Year><RecNum>17</RecNum><DisplayText>[16]</DisplayText><record><rec-number>17</rec-number><foreign-keys><keyapp="EN"db-id="209sespa0d5zdaev0pqpxp5jf22apw2p5fet"timestamp="1684048349">17</key></foreign-keys><ref-typename="JournalArticle">17</ref-type><contributors><authors><author>Park,S.H.</author><author>Goo,J.M.</author><author>Jo,C.H.</author></authors></contributors><auth-address>DepartmentofRadiology,SeoulNationalUniversityCollegeofMedicineandInstituteofRadiationMedicine,SNUMRC,Seoul,Korea.</auth-address><titles><title>Receiveroperatingcharacteristic(ROC)curve:practicalreviewforradiologists</title><secondary-title>KoreanJRadiol</secondary-title></titles><periodical><full-title>KoreanJRadiol</full-title></periodical><pages>11-8</pages><volume>5</volume><number>1</number><keywords><keyword>AreaUnderCurve</keyword><keyword>Humans</keyword><keyword>*ROCCurve</keyword><keyword>Radiography/*statistics&numericaldata</keyword><keyword>Software</keyword><keyword>SolitaryPulmonaryNodule/diagnosticimaging</keyword><keyword>Statistics,Nonparametric</keyword></keywords><dates><year>2004</year><pub-dates><date>Jan-Mar</date></pub-dates></dates><isbn>1229-6929(Print) 2005-8330(Electronic) 1229-6929(Linking)</isbn><accession-num>15064554</accession-num><urls><related-urls><url>/pubmed/15064554</url></related-urls></urls><custom2>PMC2698108</custom2><electronic-resource-num>10.3348/kjr.2001</electronic-resource-num><remote-database-name>Medline</remote-database-name><remote-database-provider>NLM</remote-database-provider></record></Cite></EndNote>[16]是一种常用的用于评估分类模型性能的图形化工具。ROC曲线描述了二分类模型在不同分类阈值下的敏感性和特异性之间的权衡关系，以及不同阈值下的分类准确率和误差率之间的关系。ROC曲线的表现优劣可以通过AUC（AreaUnderCurve）值来评估，AUC值介于0.5到1之间，越接近1表示分类器性能越好，越接近0.5表示分类器性能越差。当AUC等于0.5时，分类器的性能等同于随机猜测。ROC曲线常用于评估分类模型在样本分类、基因表达谱分类等问题中的性能，可以帮助我们确定最佳分类阈值和选择最佳模型。例如，在肿瘤分类问题中，ROC曲线可以帮助我们确定基于基因表达谱数据的分类模型在不同阈值下的诊断准确性，从而提高癌症诊断的准确性和效率。本课题使用R语言中的timeROC包绘制ROC曲线，分别绘制了1年生存、3年生存及5年生存的ROC曲线，并在图中显示对应的AUC值。根据相应的结果来验证预测模型的性能，进而验证用于构建风险模型的基因与疾病的关系。2.7外部数据集验证2.7.1外部数据集的下载与预处理在GEO数据库中下载结直肠癌的GSE103479的表达谱数据及对应的临床生存信息，根据TCGA数据库的结直肠癌数据处理所得到的模型，使用GSE103479数据集进行验证。GSE103479共包含156例结直肠癌样本，使用的处理平台为GPL23985，根据GPL23985对GSE103479表达谱中的探针名改为gene_Symbol。2.7.2外部数据集的KM生存分析通过对GSE103479数据中的关键signature基因表达谱的分析，结合TCGA数据库中的结直肠癌数据，我们可以建立一个风险模型，以评估样本的风险水平，并将其按照中值进行分组，最终实现KM生存率的预测。2.7.3外数据集的ROC曲线使用R语言中的timeROC包绘制ROC曲线，分别绘制了1个月生存、3个月生存率及6个月生存的ROC曲线，并在图中显示对应的AUC值。2.8关键signature基因的数据库验证使用GEPIA(/)ADDINEN.CITEADDINEN.CITE.DATA[17]这个新开发的交互式网站的单基因分析与多基因分析功能对该实验筛选出关键signature基因进行验证与评估。根据模型所使用的关键signature基因，在PubMed()数据库中搜索与之相关的文献，挖掘关键signature基因对结直肠癌的影响。

3、结果3.1结直肠癌的差异表达分析3.1.1识别结直肠癌差异表达基因经R包的DESeq2方法进行差异表达分析，当p值低于0.05，而log2foldchange值高于1时，可以确定结直肠癌患者的基因表达水平有所上调，而当log2foldchange值低于0.05时，则可以确定患者的基因表达水平有所降低，(表3-1和图3-1、图3-2)。在差异表达基因中选取最显著的前20个基因在火山图中标注出。对显著表达的前100上调基因与下调基因的表达情况进行可视(图3-2)，可以直观地观察到差异表达基因在正常样本中与癌症样本中表达的差异。表3-1TCGA结直肠癌差异表达基因数目Table3-1TheNumberofDifferentExpressionGenesinTCGACRCUpDownTotalGeneCount244426685112图3-1结直肠癌的差异表达基因火山图展示。图中蓝点表示下调基因，红点表示上调基因，进行了基因名标注的是Top20差异表达的基因。Figure3-1.VolcanoplotofdifferentexpressiongenesinCRC.Inthefigure,thebluedotsrepresentdown-regulatedgenes,thereddotsrepresentup-regulatedgenes.ThegenenamemarkedistheTop20differentexpressiongenes.图3-2结直肠癌的基因差异表达分析热图。(A)前一百个上调的mRNA绘制的热图；(B)前一百个下调的mRNA绘制的热图。Figure3-2.HeatmapplotofdifferentexpressiongenesinCRC.(A)Heatmapoftop100up-regulatedmRNA;(B)Heatmapoftop100down-regulatedmRNA.3.1.2研究差异表达基因的功能，并探究它们的信号传导途径选取了前200个上调和前200个下调的差异表达基因，使用在线网站KOBAS，识别并注释这些基因。筛选p值小于0.01的KEGG通路，可以发现前200个上调的差异表达信号通路与Wnt通路显著相关(表3-2)，说明这些基因有可能是通过影响结直肠癌患者的Wnt通路从而在结直肠癌进展中发挥作用。Wnt信号通路ADDINEN.CITEADDINEN.CITE.DATA[18]（也称为Wnt/β-连环蛋白信号通路）是公认的结肠癌驱动因素，也是最具代表性的信号通路之一。作为Wnt信号传导的功能效应分子，β连环蛋白的修饰和降解是Wnt信号通路和结肠癌发展进展的关键事件。Wnt信号通路对于研究结直肠癌（CRC）的发病机制具有至关重要的意义，这一点已经得到了科学家们的广泛认可。表3-2前200个上调基因显著富集的KEGG通路Table3-2KEGGpathwayssignificantlyenrichedforthetop200up-regulatedgenesTermsInputnumberBackgroundnumberP-ValueRibosomebiogenesisineukaryotes81050.0000001Alanine,aspartateandglutamatemetabolism4360.0000460Wntsignalingpathway61600.0001860而前200个下调的差异表达基因富集的功能与代谢通路显著相关(表3-3)，说明这些基因可能是通过影响结直肠癌患者的代谢通路，进而在结直肠癌进展中发挥作用，最终影响结直肠癌的发生与发展。表3-3前200个下调基因显著富集的KEGG通路Table3-3KEGGpathwayssignificantlyenrichedforthetop200down-regulatedgenesTermsInputnumberBackgroundnumberP-ValueMetabolicpathways1714330.00000013.2结直肠癌的加权基因共表达网络3.2.1筛选基因与样本结直肠癌的TPM格式的表达谱转置后，进行样本的筛选，用R包stats中的hclust函数对样本做一个简单的均值聚类分析(图3-4)，可以看出样本间的离散程度不算太大，因此我们对这17348个基因在433个样本中的表达矩阵进行加权基因共表达网络构建。图3-3结直肠癌样本的聚类分析。图3-3是对结直肠癌433个样本的简单聚类，展示了样本间的离散程度。Figure3-3.Clusteranalysisofthe434samplesinCRC.Figures3-3showasimpleclusteringof433samplesforCEC,demonstratingthedegreeofdispersionbetweensamples.3.2.2选择软阈值使用R语言WGCNA包中的pickSoftThreshold函数计算相似度经power值加权后网络的无尺度网络适应系数(图3-4A)和网络的平均连通性(图3-4B)。一般来说，无标度拓扑拟合指数这个图是用来选择软阈值的一个根据。例如下图是1到30。一般选择在0.9以上的，且为第一个达到0.9以上数值。下图的6是第一个达到0.9的数值，可以考虑6作为软阈值。从下图的B图可以看出，当数值为6的时候，已经开始持平，则软阈值为6时，网络的连通性好。AB图3-4网络的无尺度网络适应系数及网络中所有节点的平均连通数。图3-4是对结直肠癌433个样本的简单聚类，展示了样本间的离散程度。(A)纵坐标是相似度经power值加权后网络的无尺度网络适应系数；(B)纵坐标是相似度经power值加权后网络中所有节点的平均连通数。Figure3-4.Thescale-freenetworkadaptationcoefficientofthenetworkandtheaverageconnectednumberofallnodesinthenetwork.(A)Theordinate,ascale-freenetworkadaptationcoefficient,isdeterminedbythepowervalueweightedbysimilarity;(B)Thisordinateistheaverageconnectednumberofallnodesinthenetwork.3.2.3识别加权基因共表达模块使用R语言WGCNA包中的blockwiseModules函数一步法构建加权基因共表达网络，共识别到29个基因共表达模块(图3-5)。图3-5加权基因共表达网络模块的层次聚类图。图中左侧灰色部分是未分配到共表达模块中的基因（在后续实验中不必再分析），彩色部分是29个共表达模块的层次聚类展示。Figure3-5.Hierarchicalclusteringmapofweightedgeneco-expressionnetworkmodules.Theleftgraypartofthefigureisthegenesnotassignedtotheco-expressionmodule(whichdonothavetobeanalyzedagaininsubsequentexperiments),andthecoloredpartisthehierarchicalclusteringpresentationofthe29co-expressionmodules.3.2.4分析加权基因共表达网络通过R语言WGCNA包的TOMsimilarityFromExpr函数，为了更加准确地描述拓扑重叠矩阵，我们从400个基因中抽样，并利用这些信息，构建出一张具有较高精度的基因相关性热图（图3-6）与此同时，使用R语言WGCNA包中的moduleEigengenes函数获得每个模块的表达谱数据对应的一维向量，根据每个模块内基因的表达进行聚类可以得到模块间的相关性热图并对应着它们的聚类情况(图3-7)。图3-6随机选取的400个基因的模块内相关性热图。图中左侧灰色部分是未分配到共表达模块中的基因（在后续实验中不必再分析），彩色部分是29个共表达模块的层次聚类展示。本图展示的是随机选取的400个基因绘制的模块内相关性热图，左侧和上方的图例展示的是这29个模块的层次聚类情况。热图的深浅代表这29个模块之间相关性的大小。Figure3-6.Heatmapofwithin-modulecorrelationsfor400randomlyselectedgenes.Theleftgraypartofthefigureisthegenesnotassignedtotheco-expressionmodule(whichdonothavetobeanalyzedagaininsubsequentexperiments),andthecoloredpartisthehierarchicalclusteringpresentationofthe29co-expressionmodules.Thisfigureshowstheheatmapofthecorrelationwithinthemoduledrawnbythe400genesrandomlyselected,andthelegendontheleftandtopshowsthehierarchicalclusteringofthese29modules.Thedepthoftheheatmaprepresentsthemagnitudeofthecorrelationbetweenthese29modules.图3-7模块间相关性热图以及聚类情况。该图上方展示的是29个模块的聚类情况，下方的热图展示的是两两模块之间的相关性。Figure3-7.Heatmapofcorrelationbetweenmodulesandclustering.Thetopofthefigureshowstheclusteringofthe29modules,andthebottomheatmapshowsthecorrelationbetweenthetwomodules.3.3筛选关键模块3.3.1基于超几何分析筛选模块使用R语言自带的phyper函数对从WGCNA获得的模块与差异表达基因进行超几何分析，筛选p值小于0.01的基因模块，p值小于0.01说明模块显著富集到差异表达基因，共从29个基因共表达模块中筛选出了11个显著富集到差异表达基因的共表达模块(表3-4)。表3-4超几何分析结果Table3-4ResultsofHypergeometricAnalysismoduleDEGRatiocountintersect_countblack0.5105382195magenta0.3817338129midnightblue0.9658146141pink0.4108353145purple0.5802324188red0.4890409200saddlebrown0.55004022tan0.5787235136turquoise0.456934561579white0.88374338yellow0.3781804304注：表3-4是超几何分析计算完后选取的p<0.05的模块。Module是模块名，DEGratio是模块中差异基因的占比。3.3.2识别关键模块对11个模块中的差异表达基因分别进行功能与通路的富集分析。根据模块基因富集到的通路，再结合差异表达基因在模块基因中的占比来对模块做进一步筛选。选取4个显著富集的模块，分别为black、pink、turquoise、yellow模块。black模块可以显著地检测到细胞周期和p53信号通路，同时还可以检测出具有促进细胞分裂的特定基因。pink模块显著富集的通路有病毒性致癌、癌症的转录失调。turquoise模块显著富集的通路有PI3K-AktADDINEN.CITEADDINEN.CITE.DATA[19]信号通路、抗原加工与呈现等。从相关研究可知，PI3K-Akt信号通路由受体酪氨酸激酶激活，在调节血管生成方面很重要。yellow模块可以有效地检测出Wnt信号传导途径的异常，这种传导途径可能会导致多种不同类型的疾病，其中尤以恶性肿瘤为例。因此我们选择这4个模块作为关键模块，对这4个模块内的基因做进一步分析。3.4识别关键signature基因3.4.1基于单因素Cox筛选基因对筛选出的4个模块的基因与差异表达基因取交集，得到2223个基因作为待分析的基因。单因素Cox回归，选择p值小于0.01的基因进行后续分析。将基因从2223个基因筛选到99个基因。3.4.2基于Lasso回归进一步筛选基因将经过单因素Cox回归分析得到的99个基因进行Lasso回归分析（图3-8），共得到了17个基因，将该17个基因作为关键signature基因（表3-5）。表3-5Lasso回归筛选的关键signature基因Table3-5KeysignaturegenesscreenedbyLassoregressionGeneCoef1CDC25C-0.042606812HMMR-0.0323328363TERT0.2913231064GRPEL2-0.4039481415ELFN10.0730851386ATP6V1B10.0277187747SIX40.0630044478SHISAL2A0.2081825929PANX20.30501077110ZNF6761.0733473811GABRD0OXC60.10664178113SLC35F30.1624459614TIMP10.03378586615EEPD10LAC10.165000321续表3-5Lasso回归筛选的关键signature基因GeneCoef17ACOX1-0.262655831图3-8Lasso回归结果图。通过Lasso回归分析有效地筛选出有用的基因。(A)系数分布图，就像数值分叉一样，一条线就代表一个基因，这些基因的末尾会指向一个纵坐标，这个纵坐标就是系数，LASSO会为每一个基因算一个系数；(B)图中有两条虚线，一般会选择左面那条虚线对应的参数，这个参数被称为lambda.min，从图上可以看到，上图对应的上面的数值是17，说明有17个基因的系数可以留下来。Figure3-8.PlotofresultsfromLassoregression.Figure3-8showsthatLassoregressionanalysiswasperformedon99genestofurtherscreengenes.(A)Coefficientdistribution,likeanumericalbifurcation,wherealinerepresentsagene,andtheendofeachgenepointstoaordinate,whichisthecoefficient,andLASSOcalculatesacoefficientforeachgene;(B)Therearetwodashedlinesinthefigure,andtheparametercorrespondingtothedashedlineontheleftisgenerallyselected,whichiscalledlambda.min.Ascanbeseenfromthefigure,thecorrespondingvalueonthetopofthefigureaboveis17,indicatingthatthecoefficientsof17genescanbeleft.3.5构建风险评分模型并进行KM生存分析3.5.1风险评分模型通过Lasso回归模型的筛选，最终将17个关键signature基因(CDC25C,HMMR,TERT,GRPEL2,ELFN1,ATP6V1B1,SIX4,SHISAL2A,PANX2,ZNF676,GABRD,HOXC6,SLC35F3,TIMP1,EEPD1,PLAC1,ACOX1)纳入到生存预后模型。根据该模型中的入选的17个基因来构建模型，为每个样本进行风险打分，并根据中位风险打分将样本分为高低风险。根据各个样本的风险评分并依次划分为高低两组。3.5.2高低风险组的KM生存分析使用高低风险组的生存信息，进行KM生存分析,结果表明低风险评分与良好的预后相关（图3-9）。由此可证明，筛选到的17个基因可能为潜在的影响结直肠癌进展的关键signature基因。图3-9高低风险组的KM生存曲线。Figure3-9.KMsurvivalcurvesforthehighandlowriskgroups.3.5.2ROC曲线通过使用R语言的timeROC包，我们可以绘制出样本在该风险模型下的ROC曲线（图3-10），从而可以清楚地识别出不同时间段患者的生存状况，其中AUC值分别为0.777、0.753和0.811。图3-10ROC曲线。曲线反应模型的准确性，AUC值均大于0.75，说明该模型具有较高的准确性。Figure3-10.ROCcurve.TheAUCvaluesofthecurveresponsemodelwereallgreaterthan0.75,indicatingthatthemodelhadhighaccuracy.3.6外部验证集验证3.6.1外部验证集的获取使用处理好的结直肠癌的GSE103479的表达谱数据及对应的临床生存信息，提取GSE103479数据所包含的关键signature基因的表达谱作为外部验证集。3.6.2外部验证集的KM生存分析使用上述风险模型，评估样本的风险评分，并根据中值分组，分为高低风险两组，并进行KM生存分析（图3-11）。图3-11验证集高低风险组的KM生存曲线。Figure3-9.TheKMsurvivalcurvesofthehighandlowriskgroupsinthevalidationset.3.6.3外部验证集的ROC曲线同样使用R语言timeROC包绘制验证集在该风险模型下的ROC曲线（图3-12）。由于生存数据时间不够长，区分为1个月、3个月、5个月患者的生存状态，其AUC值分别为0.688、0.726、0.779，说明该模型具有较好的区分性能。图3-12ROC曲线。曲线反应模型的准确性，AUC值均大于0.68，说明该模型具有较高的准确性。Figure3-12.ROCcurve.TheAUCvaluesofthecurveresponsemodelwereallgreaterthan0.68,indicatingthatthemodelhadhighaccuracy.3.7关键signature基因的数据库验证经GEPIA数据库检索，我们发现CDC25C、TERT和TIMP1三个基因都可以参与到肿瘤的治疗之中。其中，CDC25C基因的表达水平极其稳定，它可以有效地控制细胞的生长和凋亡，从而参与到肿瘤的治疗之中。酪氨酸磷酸酶编码的蛋白质被称为Cdc25，这一酶类可以促使B结合CDC2，从而实现去磷酸化，从而激活细胞的生长和繁殖过程。它也被认为可以抑制p53诱导的生长停滞。TERT，端粒酶逆转录酶,是一种核糖核蛋白聚合酶体细胞中端粒酶表达的失调可能与肿瘤发生有关TIMP1属于一个多样性的蛋白质编码群，它们不仅可以抑制大多数MMP，而且也可以激活多种不同的细胞生长，甚至可以抵御衰老。MMP的编码蛋白质不仅仅被认为是一种酶，它们也被认为是一种调节机体免疫力的重要物质。这个基因的表达可以显著地影响许多生物体内的生理过程，包括调节肽水平。在PubMed数据库中进一步搜索，发现已有对这三个基因研究的文章，且与结直肠癌相关。

4、讨论经过精心分析，我们从TCGA数据库中提取出结肠癌和直肠癌的基因，并以TPM格式构建基因共表达网络，以便更好地探索它们之间的差异。此外，我们还采用超几何分析，从结直肠癌基因集中筛选出4个关键模块，以深入探索它们的功能和通路。使用单因素Cox回归分析，对关键模块中的差异表达基因进行筛选，识别到了99个和结直肠癌进展相关的基因。最后，使用Lasso回归分析，将基因进一步筛选，最终获得了CDC25C、TERT、TIMP1等17个与疾病进展相关的关键signature基因。在这些基因中，我们发现CDC25C参与调节G2/M进展和介导DNA损伤修复，CDC25CADDINEN.CITEADDINEN.CITE.DATA[20]表达的变化与肿瘤发生和肿瘤发展密切相关，可以作为癌症治疗的潜在靶点。研究表明，TERT基因的端粒酶复合物[21]具有重要的生理功能，它们不仅能够影响肿瘤的生长，而且能够提供准确的预测信息，并且提供了一系列全面的抗肿瘤药物。实验得到的结果在文献中得到印证，说明了我们识别的关键signature基因在结直肠癌疾病进程中发挥重要作用。另外，我们发现GRPELADDINEN.CITEADDINEN.CITE.DATA[22]在调节线粒体蛋白输入和氧化还原稳态中起着至关重要的作用，在胶质母细胞瘤细胞中的验证表明，GRPEL2敲低有效抑制细胞增殖，GRPEL2抑制诱导亚G1期的细胞周期停滞。GRPEL2是一种新型的线粒体生物能量调节剂，也是未来治疗胶质母细胞瘤的潜在靶标。那么，我们猜测GRPEL2可能也会是未来治疗结直肠癌的潜在靶标。根据这17个关键signature基因所对应的Geof值计算了383个肿瘤样本的风险评分。根据风险评分将肿瘤样本分为高低风险两组，通过KM生存曲线和ROC曲线，我们发现高低风险组之间的生存存在显著差异，且低风险评分组的结直肠癌患者有良好的预后，证明了所识别的关键signature基因具有潜在的预