骨肉瘤组织差异表达基因筛选及其与患者预后的关系

张晋煜

(云南省肿瘤医院，昆明 650000)

骨肉瘤组织差异表达基因筛选及其与患者预后的关系

张晋煜

(云南省肿瘤医院，昆明 650000)

目的采用生物信息学的方法筛选骨肉瘤组织差异表达基因，并探讨其与患者预后的关系。方法收集骨肉瘤组织及其配对癌旁正常组织标本各4例份，采用转录组测序技术筛选两种组织中的差异表达基因，并进行基因本体(GO)功能显著性富集分析。在癌症和肿瘤基因图谱骨肉瘤数据库中选择258例份骨肉瘤组织标本，采用Kaplan-Meier曲线分析富集于生物学过程的差异表达基因对患者生存时间的影响。结果骨肉瘤及其癌旁正常组织存在差异表达基因875个，其中表达上调的基因346个、表达下调的基因529个。GO功能显著性富集分析共识别出14个有意义的类型；在6个生物学过程类型中，差异表达基因主要富集于外源性代谢过程、细胞黏附、细胞外基质组成、类固醇代谢过程等；在8个细胞成分类型中，差异表达基因主要富集于细胞外体、细胞外空间、胶原三聚体、核小体等。258例份骨肉瘤组织中，还原酶1(NQO1)基因高表达233例份、低表达25例份，NQO1基因高表达、低表达患者生存时间分别为(2 763±15)、(1 568±13)d，二者比较P<0.01；醛脱氢酶3族A1(ALDH3A1)基因高表达233例份、低表达25例份，ALDH3A1基因高表达、低表达患者生存时间分别为(2 725±11)、(1 421±14)d，二者比较P<0.01；其余基因表达情况对患者的生存时间均无明显影响(P均>0.05)。结论采用生物信息学的方法筛选出骨肉瘤组织差异表达基因875个，主要富集于生物学过程和细胞成分类型；其中NQO1、ALDH3A1基因低表达提示患者预后不良。

骨肉瘤；基因本体功能显著性富集分析；基因筛选；转录组测序；生物信息学；生存时间

Abstract:To screen out the differentially expressed genes in osteosarcoma tissues by bioinformatics and to investigate its relationship with the prognosis of patients.MethodsWe performed RNA-Seq to screen out the differentially expressed genes (DEGs)in osteosarcoma and matched normal tissues of four osteosarcoma patients followed by the Gene Oncology (GO) enrichment. In addition, we selected 258 cases of osteosarcoma tissue specimens in the Cancer Genome Atlas (TCGA) HCC data to analyze the influence of DEGs on survival time of patients by using based on Kaplan-Meier curve analysis.ResultsIn total, 875 DEGs were identified in four pairs of osteosarcoma specimens, including down-regulated 529 genes and 346 up-regulated genes. According to GO enrichment, 14 significant GO terms in the biological process were identified. In the six biological process types, the DEGs were mainly enriched in exogenous metabolic processes, cell adhesion, extracellular matrix composition, and steroid metabolism process; in the 8 cell composition types, the DEGs were mainly enriched in the extracellular domain, extracellular space, collagen trimer, nucleosome and so on. Two novel genes (NQO1 and ALDH3A1) were significantly associated with poor prognosis in the TCGA osteosarcoma cohort (n=258). The average survival time in 233 cases with over-expression and 25 cases with low-expression of NQO1 was (2763±15) and (1568±13) days, respectively. The average survival time in 233 cases with over-expression and 25 cases with low-expression of ALDH3A1 was (2725±11) and (1421±14) days, respectively. The expression of other genes had no significant effect on the survival time of the patients (allP>0.05).ConclusionWe identified 875 DEGs that mainly enrich in biological process and cellular component and the low-expression of NQO1 and ALDH3A1 is significantly associated with poor prognosis for osteosarcoma.

Keywords: osteosarcoma; Gene Oncology enrichment; gene screening; transcriptome sequencing; bioinformatics; survival time

骨肉瘤好发于青少年，75%的患者发病年龄为15～25岁。目前骨肉瘤的发病机制尚不清楚，但是其遗传易感性研究已有报道。转录组测序技术能够全面快速地获得某一物种特定组织或器官在某一状态下的几乎所有转录本序列信息，可更加清晰地明确等位基因的差别表达。最近有报道对骨肉瘤进行转录组测序，识别了新的基因FOXC2，并发现其参与了阿霉素耐药[1]。本研究采用生物信息学的方法筛选骨肉瘤组织差异表达基因，并探讨其与患者预后的关系。

1 临床资料

选择2013年1月～2015年12月云南省肿瘤医院收治的骨肉瘤患者4例，男2例、女2例，年龄8～26(14.6±4.2)岁；均经病理活检确诊。病理类型：软骨母细胞型2例，成纤维型2例；发生部位：上肢2例，下肢2例；转移情况：转移3例，无转移1例。排除标准：①肝、肾、肺功能异常者；②精神异常者；③合并其他肿瘤者；④术前进行放化疗或免疫治疗者。本研究通过医院伦理委员会审核，患者及其家属均知情同意。

2 方法与结果

2.1 骨肉瘤组织差异表达基因筛选 采用转录组测序技术。患者均行手术治疗，术中收集其切除的骨肉瘤组织及配对癌旁正常组织(距离肿瘤组织>2 cm)，保存于液氮中。取部分骨肉瘤组织行常规HE染色，确定肿瘤细胞比例>70%。转录组测序主要步骤：首先将组织从液氮中取出并进行研磨，采用TRIzol试剂提取总RNA。使用带有Oligo(dT)的磁珠富集样品中的mRNA，超声波将mRNA打成短片段，并以其为模板合成双链cDNA。纯化处理，添加A尾，并对接头进行连接测序。利用PCR富集得到cDNA文库，对库检合格的cDNA在HiSeq2000仪器上进行测序。对测序获得的原始数据进行质量评估，包括测序错误率检查、A/T/C/G含量分布检查、测序数据过滤(去掉带接头的以及低质量的原始序列)。使用Tophat软件将不同组织转录组测序与人类参考序列(Human Genome v19)进行比对[2]，记录骨肉瘤组织标本1～4号(S1T～S4T)和癌旁正常组织标本1～4号(S1N～S4N)的总读段数及比对率。结果显示，骨肉瘤组织及其癌旁正常组织每例份样本产生约9 G的数据量，骨肉瘤组织S1T～S4T的比对率分别为91.45%、91.76%、92.56%、92.34%，癌旁正常组织S1N～S4N的比对率分别为91.36%、92.45%、90.29%、88.25%。将骨肉瘤组织及其癌旁正常组织各基因表达量进行比较，筛选出差异表达基因。经过再校正获取数据，参照Trapnell等[3]的方法进行生物学变异排除，最后利用CuffDiff软件进行差异分析，筛选出差异表达基因。结果显示，骨肉瘤组织及其癌旁正常组织的差异表达基因共有875个，其中表达上调的基因346个、表达下调的基因529个。

2.2 差异表达基因的基因本体(GO)功能显著性富集分析 选择2.1中筛选出来的875个骨肉瘤组织及其癌旁正常组织差异表达基因，将其上传至DAVID(https://david.ncifcrf.gov/)[4]在线工具，进行GO功能显著性富集分析，识别出有意义的GO功能分类(FDR<0.05为差异有统计学意义)。结果显示，共识别出14个有意义的GO类型，包括6个生物学过程类型和8个细胞成分类型。见表1。

表1 骨肉瘤组织中差异表达基因的GO功能显著性富集分析结果

2.3 骨肉瘤组织差异表达基因与患者预后的关系 选取富集于生物学过程类型(FDR最小的3个)的差异表达基因，即富集于外源性代谢过程、细胞黏附以及细胞外基质组成的基因，在癌症和肿瘤基因图谱(TCGA，https://tcga-data.nci.nih.gov/)中的骨肉瘤数据库进行Kaplan-Meier生存曲线分析。结果显示，258例份骨肉瘤组织中，还原酶1(NQO1)基因高表达233例份、低表达25例份，NQO1基因高表达、低表达患者生存时间分别为(2 763±15)、(1 568±13)d，二者比较P<0.01；醛脱氢酶3族A1(ALDH3A1)基因高表达233例份、低表达25例份，ALDH3A1基因高表达、低表达患者生存时间分别为(2 725±11)、(1 421±14)d，二者比较P<0.01；其余基因表达情况对患者的生存时间均无明显影响(P均>0.05)。

3 讨论

骨肉瘤的发生是由多种基因相互作用的复杂的生物学过程，其发病机制尚不清楚。为了提高骨肉瘤患者的生存率，目前临床上亟需找到一种新的生物学标志物，用于预测患者预后或作为新的药物靶点。随着基因测序技术的快速发展，该技术在疾病的预测和诊断方面均可发挥重要作用，比如在分子层面上利用全基因组关联研究能够将识别的基因与信号通路相结合，以阐释骨肉瘤的发病机理[5,6]。临床上运用全外显子组以及全基因组测序技术发现了一系列与骨肉瘤发生有关的重要基因，如TP53、PTEN及PRB2等[7,8]。此外，CDKN1A/B基因片段的缺失以及MDM2扩增均与骨肉瘤的发生有关[7]。近年来，转录组测序技术在肿瘤的研究中已得到广泛应用。转录组测序技术与之前的芯片技术比较具有明显的优势，其结果更精确、可重复性更好、应用领域更广泛[9]。由于个体化差异，不同患者的基因表达可能存在一些差异，转录组测序技术能够更加直接地对同一基因的表达情况进行比较。本研究结果显示，骨肉瘤组织及其癌旁正常组织的差异表达基因共有875个，其中表达上调的基因346个、表达下调的基因529个。

GO功能显著性富集分析是一种利用高通量测序数据对基因施行功能注释的方法，检测结果分为细胞成分、生物学过程及分子功能[10]。本研究共识别出14个有意义的GO类型；在6个生物学过程类型中，差异表达基因主要富集于外源性代谢过程、细胞黏附、细胞外基质组成、类固醇代谢过程等；在8个细胞成分类型中，差异表达基因主要富集于细胞外体、细胞外空间、胶原三聚体、核小体等。生物学过程阐释了生物体活动的重要过程，人体内生物学过程包含了许多代谢、生化反应等重要途径。此次研究选择了富集于生物学过程类型(FDR最小的3个)的差异表达基因，即富集于外源性代谢过程、细胞黏附以及细胞外基质组成的基因进行生存分析。结果显示，骨肉瘤组织NQO1、ALDH3A1基因高表达者生存时间均长于低表达者，说明NQO1、ALDH3A1基因低表达提示患者预后不良。NQO1是NAD脱氢酶家族的一员，编码一种细胞质电子还原酶，这种蛋白酶能够抑制醌类化合物的电子还原产生自由基。NQO1基因是引起各种肿瘤发生的易感基因，该基因突变与迟发型运动障碍有关，一些长期暴露于苯酚环境中致血液毒性增加的患者也能检测到该基因突变。ALDH3A1基因参与体内生物胺、神经传递素和脂质过氧化作用的代谢过程，具有抗紫外线以及在角膜中抗氧化损伤的作用。既往关于NQO1、ALDH3A1基因对骨肉瘤患者预后的影响鲜见报道，本研究首次发现NQO1、ALDH3A1基因可能参与骨肉瘤的发生、发展，并影响患者预后，其可能成为骨肉瘤诊断及治疗的新靶点。

综上所述，本研究采用生物信息学的方法筛选出骨肉瘤组织差异表达基因875个，主要富集于生物学过程和细胞成分类型；其中NQO1、ALDH3A1基因低表达提示患者预后不良。

[1] Zhang CL, Zhu KP, Ma XL. Antisense lncRNA FOXC2-AS1 promotes doxorubicin resistance in osteosarcoma by increasing the expression of FOXC2[J]. Cancer Lett, 2016(396):66-75.

[2] Trapnell C, Pachter L, Salzberg SL. TopHat: discovering splice junctions with RNA-Seq[J]. Bioinformatics, 2009,25(9):1105-1111.

[3] Trapnell C, Williams BA, Pertea G, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation[J]. Nat Biotechnol, 2010,28(5):516-520.

[4] Dennis G Jr, Sherman BT, Hosack DA, et al. DAVID: database for annotation, visualization, and integrated discovery[J]. Genome Biol, 2003,4(5):3.

[5] Kuijjer ML, Hogendoorn PC, Cleton-Jansen AM. Genome-wide analyses on high-grade osteosarcoma: making sense of a genomically most unstable tumor[J]. Int J Cancer, 2013,133(11):2512-2521.

[6] Savage SA, Mirabello L, Wang Z, et al. Genome-wide association study identifies two susceptibility loci for osteosarcoma[J]. Nat Genet, 2013,45(7):799-803.

[7] Ottaviano L, Schaefer KL, Gajewski M, et al. Molecular characterization of commonly used cell lines for bone tumor research: a trans-European eurobonet effort[J]. Genes Chromosomes Cancer, 2010,49(1):40-51.

[8] Namløs HM, Meza-Zepeda LA, Barøy T, et al. Modulation of the osteosarcoma expression phenotype by microRNAs[J]. PLoS One, 2012,7(10):48086.

[9] Raz T, Kapranov P, Lipson D, et al. Protocol dependence of sequencing-based gene expression measurements[J]. PLoS One, 2011,6(5):e19287.

[10] Ashburner M, Ball CA, Blake JA, et al. Gene ontology: tool for the unification of biology. The gene ontology consortium[J]. Nat Genet, 2000,25(1):25-29.

Screening out differentially expressed genes in osteosarcoma and its relationship with the prognosis

ZHANGJinyu

(YunanCancerHospital,Kunming650000,China)

张晋煜(1972-)，男，副主任医师，研究方向为四肢与脊柱肿瘤切除重建术。E-mail： 903170564@qq.com

10.3969/j.issn.1002-266X.2017.32.003

R738.1

A

1002-266X(2017)32-0010-03

2016-12-07)