1. 背景
肾透明细胞癌：肾癌是泌尿系统最常见的恶性肿瘤之一，发病率男性高于女性，且呈现逐年上升的态势[1]。而肾透明细胞癌是来源于肾小管上皮细胞的腺癌，是肾癌中最常见的病理亚型，约占肾癌的60%~85% [2]。流行病学调查发现，大多数的肾细胞癌病例是通过影像学检查偶然发现的，患者的生存率高度依赖于诊断时的分期，然而，当肾癌进展至远处转移期时，5年生存率仅为12% [3]。而早期肾癌往往缺乏临床表现，当出现肾癌的三大典型症状“血尿、腹部包块、腰痛”时，约60%患者已达晚期。因此，临床上亟需找到一种可早期预测肾癌预后的新的诊断方法。
TCGA数据库：于2006年启动，由美国国家癌症和肿瘤研究所和国家人类基因组研究所联合创立的一个癌症数据库，目前已经有20多种组织类型的30多种癌症11,000多个病人的临床与基因表达信息。TCGA数据库中保存有肿瘤患者全部的肿瘤测序数据以及规范的临床数据，并且这些信息处于动态更新中，数据库中的所有信息全部开放供研究者下载使用。作为目前最齐全的基因数据库之一，只要涉及肿瘤的大数据研究，都能看到它的影子。随着TCGA数据库的不断发展，越来越多的研究者选择使用TCGA数据库，在研究肿瘤的同时，也不断丰富着数据库的发展[4]。
肿瘤免疫微环境：是指肿瘤细胞与周围细胞、基质和分子环境的复杂相互作用。它包括免疫细胞、成纤维细胞、血管内皮细胞、基质成分及细胞因子等。肿瘤微环境在肿瘤发生、发展和转移中起关键作用，影响肿瘤细胞的增殖、免疫逃逸和药物反应[5]。而肿瘤免疫微环境则是肿瘤微环境的一个组成部分，它指的是肿瘤微环境中免疫细胞的变化及免疫调节分子与肿瘤细胞相互作用的过程[6]。研究肿瘤的免疫微环境，进而指导肿瘤在临床中的诊断、治疗、预防，正成为炙手可热的研究新向标。
2. 研究方法
2.1. 研究数据的下载及整理
从TCGA官网(https://portal.gdc.cancer.gov/)中下载肾透明细胞癌的蛋白质表达谱、RNA表达谱以及所有病例的临床数据。本研究所需的计算机工作：蛋白质表达数据和临床数据的提取与整理，是使用Strawberry Perl工作站进行的，Perl的版本号为5.30.0.1，分析过程及图表可视化是在RStudio中进行的，R版本号为R-4.3.1。研究中使用的所有数据均来自于公共数据库。
2.2. 建立肾透明细胞癌的预后模型
将525个肾癌患者按照1:1随机分为训练集和测试集。对训练组使用单因素COX回归曲线，筛选出与生存相关的蛋白质。再使用LASSO分析对它们进行筛选，避免变量过拟合，最后进行多因素COX分析，得到与肾癌预后显著相关的蛋白质。以上方法同样作用于测试组，当多因素COX所得结果与训练集相同时，输出最终的预后显著相关蛋白质。多因素COX分析时，获得每种蛋白质的coefficient评分，为每一位患者赋予Risk评分，评分方法为
。其中，exp为蛋白表达量，β是多因素COX得到的风险系数。Risk评分的中位数为截断点，将所有患者分为高低风险组，进而进行后续分析。
2.3. 验证预后模型
使用生存曲线、ROC曲线和COX分析对高低风险组进行组间分析，从多维度验证预后模型的优越性。使用“survival”[7]包对高低风险组勾画总体生存期以及无进展生存期的K-M曲线。根据高低临床分期患者分组，绘制高低风险组的生存曲线，从不同维度判断模型蛋白对患者生存的影响。接着，根据时间以及临床特征分别绘制ROC曲线，从时间纵向维度以及临床特征横向维度判断模型的优越性。最后，使用风险曲线评估高低风险评分对患者生产的影响，并判断模型基因在Risk评分中承担的作用。
2.4. 模型的应用
根据每一种临床特征的表达情况，为患者进行赋分，根据分数的大小，使用Nomo图可以预测患者的生存期[8]。接下来，分别在模型蛋白间以及模型蛋白与所有蛋白间进行pearson相关性分析，所得的相关性结果使用“circle”[9]包和“ggalluvial”[10]进行可视化。为了进一步探究模型蛋白间的关系，我们从STRING官网(https://cn.string-db.org/)下载了模型蛋白间的互作网络图。
2.5. 评估肾癌患者的免疫浸润情况
使用“cibersort”[11]算法根据患者的RNA表达情况赋予打分，用以评估不同免疫细胞的表达情况，22种免疫细胞在高低风险组中的表达情况使用雷达图进行可视化。对具有高低风险组表达差异的免疫细胞绘制箱式图，评估在高低风险组中的表达情况。为了增强研究的说服力，我们进一步使用了“GSVA”算法[12]，使用28种免疫细胞作为评估指标，进行高低风险组免疫细胞表达差异分析。最后，从TIDE (Tumor Immune Dysfunction and Exclusion)官网(http://tide.dfci.harvard.edu/)下载肾透明细胞癌患者的免疫应答评分，并使用“ggplot2”[13]进行可视化。
3. 研究结果
3.1. 单因素COX、LASSO-COX以及多因素COX分析
患者的生存数据与蛋白表达数据合并，使用单因素COX分析，初步筛选出148个与生存相关的蛋白质，使用LASSO对148个变量进行分析，避免数据过拟合，通过缩减变量选取最佳lambda，得到15个与生存显著相关的蛋白质，LASSO筛选过程中的系数路径图1(a)和交叉验证曲线图1(b)如下所示。将单因素COX所得的p值最小的10个蛋白质的结果进行可视化，森林图1(c)如下所示。
Figure 1. (a) LASSO coefficient path diagram; (b) LASSO regularization curve graph; (c) Single factor COX forest plot
图1. (a) LASSO系数路径图；(b) LASSO正则化曲线图；(c) 单因素COX森林图
3.2. 生存分析
根据Risk评分获取高低风险组肾透明细胞癌患者，针对总体生存期和无进展生存期绘制K-M曲线，所得结果如图2(a)和图2(b)所示，高风险组患者的总体生存期和无进展生存期显著低于低风险组，这验证了我们模型的有效性。接下来，我们根据临床分期分别进行生存分析，低临床分期(图2(c))和高临床分期(图2(d))高风险组的生存率均小于低风险组，这表明模型不受临床分期干扰，在肿瘤的早期与晚期均展现了强大的预测能力。
3.3. ROC曲线
根据不同的临床性状：Risk评分、年龄、性别、病理分期和临床分期，绘制ROC曲线，所得结果如图3(A)所示，Risk评分系统的AUC为0.75，相较于其他临床性状，仅次于临床分期，这表明我们的模型在预测准确性上具有良好的表现。然后，根据不同的生存时间，1年、3年、5年生存期绘制时间依赖
Figure 2. (a) KM curve of overall survival for high and low-risk groups; (b) KM curve of progression free survival for high-risk groups; (c) KM curve of high-risk group in low clinical stage T1~2; (d) KM curve of high-risk group in high clinical stage T3~4
图2. (a) 高低风险组总生存期KM曲线；(b) 高低风险组无进展生存期KM曲线；(c) 低临床分期T1~2中高低风险组KM曲线；(d) 高临床分期T3~4中高低风险组KM曲线
Figure 3. (a) ROC curve classified by clinical features; (b) Time-dependent ROC curve
图3. (a) 按临床特征分类的ROC曲线；(b) 时间依赖性ROC曲线
性ROC曲线，所得AUC如图3(b)所示，1年生存期的AUC为0.755、3年生存期的AUC为0.721、5年生存期的AUC为0.755，均高于0.7，这表明随着时间的进展，我们的模型均具有良好的预测潜力。
3.4. 风险曲线
患者的Risk分数从低至高按照从左至右依次排列，如图4(a)所示，左侧为低风险人群，右侧为高风险人群，纵坐标为Risk评分，横坐标为人群数量。将每位患者化成点，横坐标依照图4(a)的横坐标表示，纵坐标为生存时间，红色代表死亡，蓝色代表存活。如图4(b)所示，可以看出高风险组患者的死亡比例较大，且生存时间较短，这进一步证明我们模型的可靠性。最后，将8种模型基因按照高低风险组拆分，如图4(c)所示，横坐标同样按照原种方法，随着Risk评分的升高，SCD1、GAB2、PDK1_pS241、PRAS40_pT246、SHC_pY317的表达逐渐降低，RPA32_pS4_S8、DNMT1、TFRC的表达逐渐升高。这表明SCD1、GAB2、PDK1_pS241、PRAS40_pT246、SHC_pY317是肿瘤进展的保护因素，RPA32_pS4_S8、DNMT1、TFRC是肿瘤进展的危险因素，这与前文生存曲线得出的结果一致，进一步说明了在蛋白质层面与肿瘤预后相关的分子情况。
Figure 4. (a) Risk curves for all patients; (b) Scatter plot of all patients and survival time; (c) Expression heatmap of model genes in patients with high and low-risk groups
图4. (a) 所有患者的风险曲线；(b) 所有患者和生存时间的散点图；(c) 模型基因在高低风险组患者的表达热图
3.5. 列线图
建立模型后，使用列线图预测肾透明细胞癌患者的1，3，5年生存率，每种临床形状赋予打分，所有临床特征所得总分之和用来预测生存期的概率。如图5(a)所示，得到448分的患者的1年生存概率为81.2%，3年生存概率为54.1%，5年生存概率为34.2%。使用此图可以预测肾透明细胞癌患者的生存时间和生存状态。列线图的校准曲线如图5(b)所示，横坐标是预测值，纵坐标是观察值，通过校准曲线可看出我们的列线图具有良好的预测能力。
Figure 5. (a) Column chart predicting patient 1, 3, and 5-year survival based on all clinical features; (b) Calibration curve of column chart
图5. (a) 根据所有临床特征预测患者1，3，5年生存期的列线图；(b) 列线图的校准曲线
3.6. 共表达分析和互作网络
采用pearson相关性分析评估模型蛋白间和模型蛋白与其他蛋白间内在的联系。模型蛋白间的共表达结果如图6(a)所示，SCD1与PRAS40_pT246，DNMT1与TFRC具有内在强正相关性，这表明他们在表达时可能具有协同作用。模型蛋白与其他蛋白间的关系如图6(b)所示，将相关系数绝对值大于0.4以及p值小于0.001的关联关系可视化，模型蛋白与多种蛋白的表达间存在内在的联系。最后，为了探究模型蛋白间的作用网络，我们从STRING官网下载了模型蛋白间的作用结构。如图6(c)所示，图中可以看出，
Figure 6. (a) Circle diagram of co expression analysis between model proteins; (b) Waterfall plot of co expression analysis between model protein and all proteins; (c) Protein interaction network diagram
图6. (a) 模型蛋白间的共表达分析圈图；(b) 模型蛋白与所有蛋白间的共表达分析瀑布图；(c) 蛋白互作网络图
除了RPA2蛋白外，其余模型蛋白间存在强相互作用，这为后续的研究打开了新的方向。
3.7. 评估模型对免疫浸润的影响
22种免疫细胞通过cibersort算法获取表达结果，依据高低风险组进行组间差异分析，所得雷达图如图7(a)所示。*表示p值 < 0.05，**表示p值 < 0.01，***表示p值 < 0.001，静息状态的肥大细胞、M0巨噬细胞具有组间差异。静息状态的肥大细胞在高风险组表达较低，M0巨噬细胞在高风险组表达较高。接下来，使用GSVA算法从另一层面分析高低风险组的免疫浸润差异依据高低风险组不同的免疫肿瘤微环境，所得结果如图7(b)所示，肥大细胞、巨噬细胞同样存在组间差异。箱式图更直观地显示了免疫细胞在高低风险组的表达情况。然后，我们接着进行了免疫应答分析。根据肾透明细胞癌的基因表达情况，从TIDE官网下载了免疫应答分数，使用“ggplot2”[13]进行结果可视化。从图7(c)可看出，肾透明细胞癌的总体免疫应答较好，表明具有良好的免疫治疗潜力。根据高低风险组分别绘制箱式图，在图7(d)中可看出，高风险组的免疫不应答比例高于低风险组，这表明高风险患者对免疫治疗的效果较差。
(a)
(b)
(c)
(d)
Figure 7. (a) Calculation of differential expression of immune cells between high-risk and low-risk groups based on cibersort algorithm; (b) Calculate the differential expression map of immune cells between high-risk and low-risk groups based on the ssGSEA algorithm; (c) Bar chart of immune response for all patients; (d) Bar chart of immune response in high-risk and low-risk groups
图7. (a) 基于cibersort算法计算高低风险组免疫细胞表达差异图；(b) 基于ssGSEA算法计算高低风险组免疫细胞表达差异图；(c) 所有患者免疫应答柱状图；(d) 高低风险组免疫应答条形图
4. 研究结果与讨论
本研究通过TCGA数据库使用蛋白组学建立了一个包含8个蛋白质的肾透明细胞癌的预后模型，并通过生存分析、ROC曲线、和风险曲线的生物信息学方法进行验证，评估了模型的有效性。通过分析模型蛋白的表达情况，发现SCD1、GAB2、PDK1_pS241、PRAS40_pT246、SHC_pY317的高表达是肿瘤进展的保护因素，RPA32_pS4_S8、DNMT1、TFRC的高表达是肿瘤进展的危险因素。接下来，我们进一步开发模型，建立了列线图用来预测患者生存期，表明模型在临床中实践的可能性。
为了阐明模型蛋白的作用机制，使用共表达分析评估了分子层面模型蛋白间的相互作用。发现SCD1与PRAS40_pT246，DNMT1与TFRC具有内在强正相关性，在模型中可能具有正协同作用。同时，共表达分析结果显示了多种与模型蛋白具有协同作用的蛋白质，它们在肿瘤的发生与进展中可能起到了关键性的作用，蛋白相互作用网络分析同样证实了这一点，这有待于进一步实验室结论的证实。
最后，依照预后模型，我们探究了肾透明细胞癌的免疫微环境，发现静息状态的肥大细胞、M0巨噬细胞在肿瘤的发生发展中起到了至关重要的作用。静息状态的肥大细胞高表达是一种保护因素，M0巨噬细胞高表达是肿瘤的危险因素。肾透明细胞癌的免疫微景观揭示其强大的免疫相关性，我们后续的TIDE分析也证实了这一点，高风险组患者具有更差的预后和更差的免疫应答。高风险组患者的免疫治疗效果更差，这表明其可能具有更复杂的免疫景观，进而影响患者的生存期，这有待于研究者后续研究。
尽管如此，研究仍存在一些不足。一、研究主要依赖TCGA数据库，缺乏其他独立数据库的验证，模型的泛化不足，可能存在数据偏倚，未来仍需要大量外部验证集验证。二、免疫浸润分析主要依赖计算工具，可能无法全面反映肿瘤微环境的复杂性。三、研究还未深入探讨免疫浸润与预后的具体机制，限制了结果的生物学解释。最后，尽管指出免疫微环境的重要性，但未进一步探讨如何基于免疫特征优化免疫治疗策略。
综上所述，我们的研究提供了一个新的肾透明细胞癌的预后模型，揭示了肾透明细胞癌的预后相关蛋白标志物，并初步探究了肿瘤的免疫微环境，为肾癌的临床诊断和治疗提供了重要的参考依据。
NOTES
*通讯作者。