这样的免疫浸润套路分析今年难发3+分了
Exploration of the immune-related signature and immune infiltration analysis for breast ductal and lobular carcinoma乳腺导管、小叶癌的免疫相关signature及免疫浸润分析
一、研究背景
乳腺癌（BC）是女性中发生率最高的癌症，而免疫反应在乳腺癌的发生发展中起到了重要作用，但早前的BC免疫相关研究主要集中于差异表达基因和肿瘤免疫环境方面，对肿瘤免疫细胞生物标记的探索仍不够充分。作者希望通过对上述生物标记的探索寻找到可靠性更高的肿瘤免疫情况预测以及免疫治疗的靶点。
二、研究流程
三、结果解析
1.鉴定BC中免疫相关signature
作者首先在TCGA数据库中获得了1178个BC样本，包括1066个肿瘤和112个正常样本，其中临床数据存在缺失的病人被删去，完整基线资料展示于表1。1045个成功匹配到预后以及转录组数据的肿瘤样本被包含在随后的研究中。
表1.TCGA样本队列基线资料表
作者使用R中的limma包来将转录组数据标准化并筛选出了10563个差异表达基因（|log fold change|>1，FDR<0.05），将其绘制成火山图（图1A），整体来看，图上下的基因分布较为对称，没有在图中看出明显的因数据处理不当等原因导致的偏差。
图1A.差异表达基因火山图
获得差异表达基因后，作者又将其与InnateDB数据库获得的4678个免疫相关基因取交集，得到了1399个免疫相关signatures（图1B）。GO及KEGG富集分析结果提示这些基因在几个不同的免疫相关过程存在crosstalk，包括细胞因子-细胞因子受体互作，Th1/Th2细胞分化以及JAK-STAT信号通路（图1C，D）。
图1B.差异表达基因与免疫相关基因的交集
图1C，D.富集分析结果
2.肿瘤相关免疫基因（TAIG）风险评分模型的建立及评估
在通过上述分析得知了肿瘤和正常样本中的免疫相关signatures存在差异后，作者想进一步探索与BC预后相关的基因，并建立一个TAIG风险评分模型以便量化肿瘤相关免疫基因的不同表达情况所对应的风险。
为了达成上述目标，作者首先将基因表达与生存相关数据匹配后进行了单变量Cox回归分析，以p<0.01为阈值获得了136个与预后相关的signature，然后使用LASSO回归进一步筛选。从图1E中我们可以看到，随着LASSO回归中参数λ（可以简单理解为防止过拟合的参数）的增大，被纳入模型中的自变量从136个开始逐渐减少，最终收束到仅一个自变量，而在图1F中我们可以看到随着λ的增大偏似然偏差先减后增。综合自变量收束的速度以及偏似然偏差的大小，作者最终选择了17个自变量，对应17个关键基因（表2）。使用多变量Cox回归（图3A）的结果作为每个基因在模型中的权重，作者构建出了TAIG风险评分模型。
图1E，F.LASSO回归分析
表2.关键基因列表
图3A.17个基因的多变量Cox回归分析
通过该模型，作者计算了队列中每位病人对应的TAIG风险评分值。通过使用评分值的中位数作为分组依据，作者将病人分为了低TAIG组（n=523）和高TAIG组（n=522）。
为了判断高低TAIG评分值的病人之间是否存在差异，作者进行了生存分析，发现TAIG风险评分高的患者的预后风险更高，总生存期往往更短（图2A，B）。此外，两组中的基因表达也有较为显著的差异（图2C）。KM法进行生存分析的结果也进一步证明了高TAIG评分值的病人预后更差（图3C）。基于上述模型对预后进行评估的ROC曲线AUC值为0.793（图3B），表明TIAG风险评分模型对TCGA数据库中BC的预后评估方面表现优秀。
图2A，B.患者的TAIG分数与预后的关系图2C.高TAIG组和低TAIG组的基因表达差异
图3B.模型用于预后评估的ROC曲线图3C.对高/低TAIG评分患者的生存分析
接着，作者将目光放到临床层面，对TAIG评分值与病人的临床特征进行相关分析，也发现高TAIG分数与更差的AJCC-TNM分期和病理学分级相关（图4A，B，C，D，AJCC-TNM是美国癌症联合委员会提出的一种肿瘤分期标准）。
以上结果是基于与模型训练集相同的数据得出的，所以更大的意义可能在于对建模效果的评估，为了验证模型的普适性，作者在另一个独立的METABRIC数据库的数据集对自己的模型进行了检验，发现使用模型对患者3，5，7年的生存情况预测的ROC曲线AUC分别为0.732，0.784和0.829（图4E）。与TCGA中的数据一样，将模型输出的参数与病人预后数据匹配后进行生存分析的结果也提示高TAIG风险分数对应更差的预后（图4F）。以上结果提示作者的模型在预后评估中具有一定的普适性。
图4A，B，C，D.患者的肿瘤分期、分级与TAIG评分值图4E.使用模型对患者预后进行预测的ROC曲线图4F.不同TAIG评分值的预后情况
3.TAIG与肿瘤免疫浸润
因为免疫细胞在肿瘤微环境中起到的重要作用，故有必要对其进行深入研究。因此作者在本节中整合了对免疫signatures的综合分析以及免疫浸润的研究结果。
基于TIMER数据库，作者探索了上文中筛选出的17个主要免疫signatures与免疫细胞之间的关系。作者研究发现，17个signature中大多数都与免疫细胞相关，尤其是PCDHGA2，SPIB，ADRB1，FLT3和NFKBIE这几个基因，从图5的偏相关分析结果中可以看出上述几个基因的表达与肿瘤纯度以及几种肿瘤相关免疫细胞浸润之间存在有一些相关性。
图5.免疫相关基因表达与免疫浸润的偏相关分析
另外，作者使用了CIBERSORT算法计算了TCGA的BC样本中22种免疫细胞的分数，p<0.05的样本被剔除以保证结果的准确性，图6中描绘了各个样本中不同免疫细胞的相对百分比，从图中可以看出，各个样本间免疫细胞的浸润情况是存在有一定差异的。
图6.CIBERSORT算法得到的肿瘤免疫浸润情况
4.两个TAIG组中肿瘤浸润免疫细胞丰度的差异
依据上文中得到的结果，作者推测——因为高低TAIG组中基因表达存在差异，而基因表达又与免疫细胞浸润存在一定的联系，那么高低TAIG组中免疫细胞的分布可能也会不同。通过热图我们可以看到，一些免疫细胞，如M0，M2巨噬细胞在两个TAIG组中的浸润情况有较为明显的差异（图7）。
图7.两个TAIG组中免疫细胞的浸润情况
此外，作者用Wilcoxon秩和检验发现某些免疫细胞在高TAIG组中浸润密度更低，包括CD8+T细胞（p=0.031），静息CD4+T细胞 (P=0.026) ，M0巨噬细胞(P=0.023)，M2巨噬细胞 (P=0.048)等。通过查阅资料，作者发现几个较早的研究中已经证明了更少的免疫细胞浸润对应着BC中更差的预后。
到这里，我们可以梳理一下作者的论证思路——作者先通过差异表达基因（DEG）分析获得了乳腺癌相比于正常组织的DEG，然后将这些DEG与InnateDB数据库获得的免疫相关基因取交集，得到了免疫相关的DEG。其后，作者使用GO和KEGG富集分析，发现这些基因在几个不同的免疫相关过程的crosstalk。再然后，作者进一步精选基因建立了一个TAIG风险评分模型并对模型进行了评估。最后，基于模型输出的参数值，作者探讨了其高低与临床病理、免疫浸润、预后之间的联系。
小结
文章的数据都来自公共数据库TCGA、InnateDB、METABRIC和TIMER，差异基因→风险模型→临床病理、免疫浸润、预后这样的论证流程中规中矩，使用的DEG分析，GO、KEGG富集分析，CIBERSORT算法等也较为常见，可以进行参考学习。