1. 引言
骨肉瘤(Osteosarcoma, OS)是影响儿童和青少年的最常见的恶性骨肿瘤类型[1]。OS既可累及中轴骨骼,也可累及四肢骨骼,通常发生于负重长骨,最常见的部位是股骨远端和胫骨近端。尽管手术提高了OS患者的存活率,但转移后的5年生存率不到20% [2]。长链非编码RNAs (long non-coding RNA, lncRNAs)是长度超过200个核苷酸的RNA转录物,缺乏蛋白质编码能力[3]。有报道称lncRNAs调节许多生物学过程,包括细胞增殖、凋亡、迁移、侵袭和血管生成[4] [5]。此外,lncRNAs也可以作为致癌基因或抑癌基因,从而在骨肉瘤的发展中发挥关键作用[6] [7]。最近的一项研究表明,二硫化物诱导的细胞死亡(二硫化物中毒),在癌症治疗中发挥着关键作用[8]。冯等人进行了一项研究,强调了二硫化物变性作为一种新的分类系统的重要性,该系统能够准确预测甲状腺癌患者的临床预后和药物敏感性[9]。这种方法涉及鉴定与二硫蛋白变性相关的特定基因,有助于癌症患者预后预测和潜在治疗靶点的鉴定[10]。几种与双硫死亡相关lncRNAs (disulfidptosis-related lncRNAs, DRLs)被证实与预后和免疫微环境有关,可能影响人类癌症的免疫治疗[11]。因此,急需探索DRLs与OS的潜在关系,以确定OS有前景的预后标志物。
在本研究中,基于TCGA数据库中获得的数据集对DRLs进行了探索,并基于筛选出的DRLs建立了预后模型并进行了验证。此外,使用qRT-PCR分析验证关键DRLs在OS细胞及组织中的表达,并分析其表达在患者预后中的价值,以期研究一种基于DRLs的有价值的OS预后模型以及相关预后标志物。
2. 材料与方法
2.1. 数据来源
从UCSC Xena [12]数据库中下载TCGA中骨肉瘤RNA-seq表达水平数据(标准化后的log (FPKM + 1, 2)表达水平)和临床信息,样本共计85例。并将Wang等人[10]获得的双硫死亡相关基因集(包括SLC7A11、INF2、ACTN4、IQGAP1、TLN1、DSTN、NADPH、CD2AP、MYH9、FLNA、MYL6、CAPZB、PDLIM1、MYH10、FLNB和ACTB)纳入当前分析。
2.2. 双硫死亡相关lncRNAs鉴定和功能富集分析
提取Disulfidptosis基因的表达量矩阵,使用Pearson计算Disulfidptosis与lncRNAs之间的关系,根据相关系数|r| > 0.4和P < 0.001,筛选出双硫死亡基因相关DRLs (Disulfidptosis Relative lncRNAs),并通过lncSEA数据库[13]对DRLs进行功能富集分析。
2.3. 预后显著相关DRLs基因识别及预后模型的独立分析
为了评估上述DRLs是否和患者的生存相关,随机将有临床表达数据的85例患者样本突变数据按7:3的比例分成训练集和验证集,利用训练集进行单因素Cox比例风险回归分析,筛选出单因素P < 0.05的基因,进一步通过LASSO回归[14]对特征lncRNAs进行筛选,最后进行多因素Cox回归分析,最终确定的基因为预后生物标志物,用于构建风险模型。随后,通过模型基因的表达量,计算85名OS患者样本的中位风险评分,根据风险评分的中位值将病人分为高低风险两组。使用R包survival [15]对高低风险组进行生存差异分析,绘制Kaplan-Meier曲线。最后,利用上述验证集进行风险模型验证,进一步用所有样本验证该模型的效果。预后模型的价值由1年、2年和3年生存率的曲线下面积(AUC)值决定。
为研究临床病理特征与风险模型的预后情况,通过R包rms构建基于临床信息(Gender, Age at Diagnosis, Year of Diagnosis, Disease at diagnosis, Time to first relapse)和风险评分(RiskScore)的列线图,将临床病理因素纳入Cox独立预后分析评估每个变量对患者生存的贡献。
2.4. 材料及试剂
本研究纳入我院接受初次手术的86例骨肉瘤患者。所有患者术前未行靶向治疗、放化疗或免疫治疗。通过手术切除从每位骨肉瘤患者获得肿瘤组织和癌旁正常组织,并立即在液氮中冷冻。然后,将样品保存在−80℃下进行后续实验。从病历中收集临床病理特征。所有患者出院后通过电话或门诊复查方式进行为期3年的随访,患者出现死亡时则随访结束,记录患者总生存期(overall survival, OS),总生存期定义为从手术日至随访截止日期或死亡时间。如果个体的死亡与骨肉瘤无关,则将其排除在本研究之外。所有患者均签署了参与本研究的书面知情同意书。样本经我院医学伦理委员会批准使用。
三种骨肉瘤细胞系(MG-63,Saos-2和U-2 OS)和人正常成骨细胞(hFOB1.19)从美国ATCC公司获得,DMEM培养基、胎牛血清购自美国Gibco公司,LipofectamineTM2000和TRIzol试剂购于美国Invitrogen公司,CCK-8试剂盒购于中国上海碧云天生物技术有限公司,PCR引物、si-ASB16.AS1、si-NC购自广州锐博生物,实时定量PCR试剂盒购自日本TaKaRa公司。
2.5. 细胞培养及转染
所有细胞均培养在含10%胎牛血清及青链霉素双抗含量为1%的DMEM培养基中。取对数生长期细胞接种于6孔板,当细胞密度达到80%融合时,进行传代。选用MG-63细胞,使用Lipofectamine 2000 (Invitrogen;赛默飞世尔科学公司)进行细胞转染。si-ASB16-AS1为实验组、si-NC组为对照组,转染48 h后收集细胞,qRT-PCR检测ASB16-AS1表达情况。
2.6. qRT-PCR
TRIzol试剂(Invitrogen)提取组织及细胞内总RNA;采用NanoDrop™2000分光光度计(Thermo Fisher Scientific, Inc.)对提取的总RNA进行质量测定。将RNA逆转录成cDNA,使用ABI 7500系统(Applied Biosystem, USA)进行PCR扩增。GAPDH用作检测lncRNAs的内部对照。引物序列如表1所示。所有基因表达均采用量2-ΔΔCt法方法进行分析。
Table 1. The detail information for all primers
表1. 所有引物的详细信息
Primer |
Sequence (5’-3’) |
RP11.304F15.6-F |
GGAGCATTTCTCGCGCTACA |
RP11.304F15.6-R |
ACAGGCAATTCTTGTCGCAAA |
RP11.750H9.5-F |
CCCCACCCTTTGTTCAG |
RP11.750H9.5-R |
TCAGCCTTGTCCCTCCT |
RP11.313F23.4-F |
ATGAGTGCTGCAGTGACTGC |
RP11.313F23.4-R |
TTAGGAGGCCTTGAGACGGT |
RP11.46C24.7-F |
GCCAAGCTAACCAAAGCTC |
RP11.46C24.7-R |
CATAAAGAGGCTACCATAA |
ASB16.AS1-F |
ATTGCAGGTGATGGCTCAGT |
ASB16.AS1-R |
CTGCTCAACAGCATCTATCG |
RP11.452C13.1-F |
GCGTTCTCTTCGTCCTCATC |
RP11.452C13.1-R |
GTGTGAGTCTCCCAGGATGC |
GAPDH-F |
CAAGGTCATCCATGACAACTTCG |
GAPDH-R |
GTCCACCACCCTGTTGCTGTAG |
2.7. CCK-8检测细胞增殖
培养24 h后收集转染上述分子产物的细胞,以2 × 104个/mL的密度制备悬液。将100 μL细胞悬液接种于96孔板的每孔中。每组设5个重复。细胞在37℃、5% CO2的培养箱中培养0、24、48和72 h,此时加入10 µL CCK-8 (Dojindo Laboratories Co. Ltd, Kumamoto, Japan)。孵育4 h后,使用酶标仪(Bio-Rad, Hercules, CA, USA)在450 nm记录光密度(OD)。
2.8. 统计分析
采用SPSS 22.0软件,R包及其相应软件包和GraphPad Prism 8.0软件进行数据处理和统计学分析。采用单因素Cox回归和LASSO分析构建预后模型。采用单因素和多因素回归分析确定独立的预后危险因素。采用ROC曲线分析评估临床风险模型。数据值以均数 ± 标准差(SD)表示。并通过t检验估计骨肉瘤组织与癌旁正常组织之间的差异。用卡方检验评估ASB16-AS1表达与临床病理参数之间的关系。采用Kaplan-Meier法进行生存分析。以P < 0.05为差异有统计学意义。
3. 结果
3.1. 双硫死亡相关lncRNAs鉴定和功能富集分析
本研究共发现411个与双硫死亡相关基因显著相关的lncRNAs (图1(A))。此外,功能富集分析发现这些lncRNAs发现与急性骨髓性白血病和急性血液性白血病这两种癌症类型相关性较大(图1(B))。
Figure 1. Disulfidptosis relative lncRNAs (DRLs) and associated enrichment analysis. A: the volcano plot for DRLs revealed in current study (X-axis represented the Pearson correlation coefficient; the Y-axis represented the value of P; the red node represented the DRLs a positive correlation with immune related genes; the green node represented the DRLs a negative correlation with immune related genes); B: the result of enrichment analysis for DRLs
图1. 双硫死亡相关lncRNAs (DRLs)和富集分析。A:DRLs火山图(X轴表示Pearson相关系数;Y轴表示P的值;红色节点表示DRLs与免疫相关基因呈正相关;绿色节点表示DRL与免疫相关基因呈负相关);B:DRLs富集分析
3.2. 预后显著相关DRLs基因识别
基于上述得到的411个DRLs,结合它们在训练集TCGA-OS样本中的表达值和生存信息,通过单因素Cox回归得到27个lncRNAs,进一步通过LASSO进行特征基因的筛选,筛选出12个lncRNAs用于后续分析(图2(A)、图2(B))。最终确定6个lncRNAs (RP11.304F15.6, RP11.750H9.5, RP11.313F23.4, RP11.46C24.7, ASB16.AS1, RP11.452C13.1) (图2(C)),建立DRLs相关预后模型。
注:A和B中的Y轴表示变量的系数,X轴表示log(lambda)的值。
Figure 2. Selected characteristic genes. A: Characteristic genes screened by LASSO; B: Cross validation to reveal optimal number of genes; C: The multivariate Cox regression revealed 6 prognostic DRLs
图2. 筛选的特征基因。A:LASSO筛选的特征基因;B:交叉验证以揭示最佳基因数量;C:多变量Cox回归揭示了6个预后DRLs
在训练集(图3(A))和验证集(图3(B))中,高危组的生存结局均比低危组差。在训练集(图3(C))和验证集(图3(D))中,1年、2年和3年生存率的AUC值均大于0.77。此外,风险曲线分析表明,风险模型在训练集(图3(E))和验证集(图3(F))均能有效区分不同风险评分和状态的患者。总之,这些发现表明,当前模型预测预后效果良好。
Figure 3. The prognostic model established in current study. A-B: The Kaplan Meier survival curve based on Riskscore model in TCGA-OS training dataset and validation dataset; C-D: the ROC curve of 1-year, 2-year and 3-year survival based on samples in TCGA-OS training dataset and validation dataset; E-F: survival distribution of samples from different risk groups in TCGA-OS training dataset and validation dataset
图3. 当前研究中建立的预后模型。A-B:TCGA-OS训练数据集和验证数据集基于Riskscore模型的Kaplan-Meier生存曲线;C-D:基于TCGA-OS训练数据集和验证数据集中的样本的1年、2年和3年生存率的ROC曲线;E-F:TCGA-OS训练数据集和验证数据集中不同风险组样本的生存分布
3.3. 独立预后分析
将临床病理因素纳入Cox独立预后分析评估每个变量对患者生存的影响,探究风险模型及临床病理因素的独立预后情况。Cox回归结果显示,Disease at diagnosis和riskScore与预后相关(图4(A)、图4(B))。我们根据Age、riskScore进行打分,每一种因素对应一个评分,各因素总评分相加对应总评分(Total Point),然后根据总评分预测1、3、5年生存率,分数越高,生存率越低。最后,我们将预测结果绘制成列线图(图4(C))。通过列线图的准确性分析,证明了评分系统在预后概率估计中的可靠性和有效性(图4(D))。
Figure 4. Independent prognostic analysis. A: The forest map of for all prognostic DRLs reveled based on univariate Cox regression analysis; B: the forest map of for all prognostic DRLs reveled based on multivariate Cox regression analysis; C: nomogram model predicting the 1-, 2- and 3-year OS in patients; D: calibration curve analysis of nomogram
图4. 独立预后分析。A:根据单因素Cox回归分析得出所有预后DRLs的森林图;B:基于多因素Cox回归分析的所有预后DRLs的森林图;C:预测患者1、2、3年OS的列线图模型;D:列线图校准曲线分析
3.4. DRLs在骨肉瘤细胞系中的表达
采用qRT-PCR分析6个预后DRLs (RP11.304F15.6、RP11.750H9.5、RP11.313F23.4、RP11.46C2.47、ASB16.AS1和RP11.452C13.1)在三种骨肉瘤细胞系(MG-63,Saos-2和U-2OS)和一种正常成骨细胞(hFOB1.19)中的表达。结果显示,6个预后相关DRLs在骨肉瘤细胞系中的表达均显著高于正常成骨细胞(均P < 0.05) (图5)。
3.5. si-ASB16.AS1抑制骨肉瘤细胞增殖
在三种骨肉瘤细胞系中表现出相对较高ASB16-AS1表达的MG-63细胞系被选择用于后续的细胞实验。为了揭示ASB16-AS1在骨肉瘤细胞中的可能调节作用,通过将si-ASB16-AS1转染进MG-63细胞,使得ASB16-AS的表达被抑制。qRT-PCR分析验证了ASB16-AS1的表达(图6(A))。采用CCK-8法测定si-ASB16-AS1对MG-63细胞增殖的影响。结果表明,与si-NC相比,转染si-ASB16-AS1的MG-63细胞的增殖能力显著降低(图6(B))。
**,P < 0.01;***,P < 0.001 (下同);**,P < 0.01;***,P < 0.001 (the same below)。
Figure 5. Expression of 6 lncRNAs in osteosarcoma cell lines
图5. 6个lncRNAs在骨肉瘤细胞系中的表达
Figure 6. Effect of ASB16-AS1 on osteosarcoma cell proliferation. A: The transfection efficiency of si-ASB16-AS1 was detected by qRT-PCR; B: Effect of si-ASB16-AS1 on the proliferation of MG-63 cells
图6. ASB16-AS1对骨肉瘤细胞增殖的影响。A:qRT-PCR检测si-ASB16-AS1的转染效率;B:si-ASB16-AS1对MG-63细胞增殖的影响
3.6. ASB16.AS1在人骨肉瘤组织中表达上调
qRT-PCR用于检测86对骨肉瘤组织和癌旁正常组织中ASB16.AS1的表达水平。如图7所示,与相应的癌旁组织相比,ASB16.AS1在骨肉瘤组织中的表达显著升高(P < 0.001)。
Figure 7. ASB16. AS1 expression in osteosarcoma tissues
图7. ASB16.AS1在骨肉瘤组织中的表达
3.7. ASB16.AS1表达与骨肉瘤患者临床病理特征的相关性
为了评估ASB16.AS1在骨肉瘤发展中的潜在作用,以ASB16.AS1表达水平的中位数为截断值,将骨肉瘤患者分为低ASB16.AS1表达组(n = 41)和高ASB16.AS1表达组(n = 45)。然后分析ASB16.AS1表达与临床病理特征的相关性。如表2所示,ASB16.AS1表达与Enneking分期(P = 0.003),远处转移(P = 0.012)和复发(P = 0.014)显著相关,但与年龄,性别,肿瘤大小和肿瘤位置无关(all, P > 0.05)。
Table 2. The association between ASB16.AS1 expression and clinical features of OS patients
表2. ASB16.AS1表达与骨肉瘤患者临床特征的关系
特征 |
病例 (n = 86) |
ASB16.AS1 |
P values |
低表达(n = 41) |
高表达(n = 45) |
年龄(岁) |
|
|
|
0.797 |
<20 |
47 |
23 |
24 |
≥20 |
39 |
18 |
21 |
性别 |
|
|
|
0.780 |
男 |
49 |
24 |
25 |
女 |
37 |
17 |
20 |
肿瘤直径(cm) |
|
|
|
0.105 |
≤8 |
51 |
28 |
23 |
>8 |
35 |
13 |
22 |
肿瘤位置 |
|
|
|
0.322 |
胫骨/股骨 |
61 |
27 |
34 |
其他 |
25 |
14 |
11 |
Enneking分期 |
|
|
|
0.003 |
I-II A |
34 |
23 |
11 |
II B-III |
52 |
18 |
34 |
远处转移 |
|
|
|
0.012 |
无 |
65 |
36 |
29 |
有 |
21 |
5 |
16 |
复发 |
|
|
|
0.014 |
无 |
58 |
33 |
25 |
有 |
28 |
8 |
20 |
3.8. ASB16.AS1高表达与骨肉瘤患者预后不良有关
为了研究ASB16.AS1表达在骨肉瘤中的预后价值,进行了long-rank test和Kaplan-Meier分析。结果显示,与ASB16.AS1低表达的患者相比,ASB16.AS 1高表达的骨肉瘤患者与较短的总生存期(P = 0.03,图8)显著相关。此外,单因素和多因素Cox回归分析的结果进一步表明,ASB16.AS1表达水平可能是影响OS患者总生存期(表3)的独立预后因素(P = 0.004)。
Figure 8. Correlation between ASB16.AS1 expression level and overall survival time in osteosarcoma patients
图8. ASB16.AS1表达水平与骨肉瘤患者总生存时间的相关性
Table 3. Univariate and multivariate analyses for OS in the osteosarcoma patients
表3. 骨肉瘤患者OS的单变量和多变量分析
Factors |
Univariate analysis |
Multivariate analysis |
HR (95% CI) |
P values |
HR (95% CI) |
P values |
ASB16.AS1 expression |
2.620 (1.355~5.066) |
0.004 |
2.620 (1.355~5.066) |
0.004 |
Age |
- |
0.786 |
- |
0.855 |
Gender |
- |
0.903 |
- |
0.993 |
Tumor size |
- |
0.942 |
- |
0.683 |
Tumor location |
- |
0.964 |
- |
0.837 |
Enneking stage |
- |
0.103 |
- |
0.558 |
Distant metastasis |
- |
0.507 |
- |
0.685 |
Recurrence |
- |
0.873 |
- |
0.392 |
4. 讨论
骨肉瘤是一种常见的原发性恶性骨肿瘤,常见于儿童和青少年,其特点是治疗困难,复发和转移,预后不良[16]。尽管人们已经认识到硫死亡及其相关基因在各种人类肿瘤预后中的重要性,但OS中DRLs的确切分子机制仍然鲜为人知。本研究希望通过基于DRLs的预后特征分析和临床验证,揭示DRLs在OS进展中的潜在机制,以期为OS治疗和预后的研究工作奠定基础。
二硫化物变性是指以二硫化物积累为特征的特定程序性细胞死亡过程[17]。据报道,二硫化物变性在癌症的免疫反应和炎症中起着至关重要的作用[18] [19]。目前已证明与双硫死亡相关的预后特征参与了人类癌症的进展[20]。这些特征可以为开发相应的靶向药物或进一步了解癌症的病理机制提供重要参考[21]。Zhang等人[22]建立了一个与二硫化物变性相关的lncRNAs模型,该模型可用于预测肺腺癌患者的预后、肿瘤突变负担、免疫细胞浸润情况以及对免疫疗法和靶向治疗的反应。众所周知,lncRNAs在生物过程中起着至关重要的作用。最近的一项研究表明,通过基因差异表达分析和KM生存分析,四种lncRNAs被确定为OS的独立预后标志[23]。本研究最终确定6个lncRNAs (RP11.304F15.6, RP11.750H9.5, RP11.313F23.4, RP11.46C24.7, ASB16.AS1, RP11.452C13.1),并建立DRLs相关预后模型。我们推测包括这6个DRLs可能是OS的新的双硫死亡相关lncRNAs。
进一步验证发现,6个预后DRLs在骨肉瘤细胞系(MG-63,Saos-2和U-2OS)中的表达均显著高于正常成骨细胞(hFOB1.19)。先前研究发现,ASB16-AS1在人类多种癌症中表达上调,并参与结直肠癌、卵巢癌和食管癌细胞的增殖、迁移、侵袭和凋亡等生物学行为[24]-[26]。因此,我们选用ASB16-AS1为研究对象进一步分析。经增殖分析发现,与si-NC相比,si-ASB16-AS1转染的MG-63细胞的增殖能力显著降低。Tan等人[26]研究发现,ASB16-AS1的沉默减弱了NSCLC的增殖能力,并使细胞周期停滞在G0/1期,并加速了细胞凋亡率。
许多研究分析指出,ASB16-AS1在肿瘤组织中高表达,提示其在多种癌症中预测预后不良。ASB16-AS1在大多数癌症类型中与B细胞、T细胞CD4+和T细胞CD8+呈正相关,在一些癌症类型中与巨噬细胞、树突状细胞和中性粒细胞呈负相关。此外,ASB16-AS1与分子特征之间存在不同的相互作用模式,提示ASB16-AS1是一种潜在的全癌预后标志物,它和多种癌症类型的免疫浸润有关[27] [28]。本研究中,qRT-PCR分析表明,ASB16.AS1在骨肉瘤组织中的表达显著高于相应的癌旁组织;进一步分析发现ASB16.AS1表达与Enneking分期,远处转移和复发显著相关。Kaplan-Meier分析结果显示,与ASB16.AS1低表达的患者相比,ASB16.AS1高表达的骨肉瘤患者总生存期更短;Cox回归分析进一步表明,ASB16.AS1表达水平可能是影响OS患者总生存期的独立预后因素。张等人鉴定lncRNA ASB16-AS1在组织样本中具有显著不同的表达模式,并与肿瘤分期和分级密切相关[29]。Elabd等人[30]研究发现,与良性和对照组相比,结直肠癌患者的组织和血浆样本中ASB16-AS1的表达水平显著升高,与结直肠癌患者的低表达相比,ASB16-AS1在组织和血浆中的高表达水平与总生存率和无进展率降低显著相关,这与我们的研究结果相一致。
值得注意的是本研究有几点局限性。首先,我们只检测了ASB16.AS1对OS细胞增殖的影响,并未进行其他相关细胞生物学分析;其次,在临床预后验证中纳入的样本较少,这可能影响结果的准确性;最后,ASB16.AS1在OS中的作用机制尚未探索。因此,下一步我们将扩大样本量,更好地设计实验,进一步探索ASB16.AS1在OS中的分子机制。
5. 结论
综上所述,本研究确定了6个DRLs (RP11.304F15.6、RP11.750H9.5、RP11.313F23.4、RP11.46C2.47、ASB16.AS1和RP11.452C13.1)作为OS预后的标志,构建了一个有价值的OS预后模型;并发现ASB16.AS1可参与OS细胞的增殖,是OS患者一个独立的预后生物标志物。
基金项目
骨肉瘤中异常表达的IncRNAs的生物学功能及临床意义的研究(WFWSJK-2021-107)。
NOTES
*通讯作者。