Li Shen Blog

Stay hungry, stay foolish.

ESTIMATE计算肿瘤免疫得分

沈大力

—— 前言 ——

当评估Gene signature based的预后模型时,通常的做法是查看基于该signature分的高低风险组的肿瘤微环境情况。肿瘤微环境中,免疫细胞和基质细胞是两种主要的细胞类型;因此,对这两种细胞的含量进行评估可以一定程度上反映出肿瘤微环境的状况。基于该背景诞生了几个非常实用的分析工具,ESTIMATE就是其中之一。

 

—— 介绍 ——

ESTIMATE (Estimation of STromal and Immune cells in MAlignant Tumor tissues using Expression data) 是由MDACC的Kosuke Yoshihara等人开发的一款分析肿瘤微环境的R包。用户只需输入表达谱数据(如TCGA中下载的counts数据)即可得到每个样本的肿瘤微环境评估。该工具对肿瘤微环境的评估结果会给出三个分数,即免疫分数、基质分数以及ESTIMATE分数。其中免疫分数是用来表示肿瘤组织中免疫细胞的浸润情况,基质分数则是用来表示肿瘤组织中基质细胞的情况。而最后一个ESTIMATE分数则是综合前面两者给肿瘤纯度的一个打分。在原文中作者使用了以下公式对肿瘤纯度进行了评估:


\(Purity = cos(0.6049872018+0.0001467884 * ESTIMATE)\)

 

但该公式存在争议,即有人认为estimate score 计算肿瘤纯度的公式是根据Affymetrix的芯片数据得出的,是专门针对芯片数据使用,因此不可以用于转录组。建议只计算出estimate score,用这个分数来代替肿瘤纯度的绝对数值用于后续分析。不过,此方法早在15年就已经发表,到现在已经过去了6年,似乎并没有人明确发文反对,因此我个人认为如果你一定想使用这个公式来算纯度的话,应该也不会有太大问题。

 

—— 如何使用 ——

首先准备好相关的包:

library(utils)
rforge <- "http://r-forge.r-project.org"
if(!require("estimate"))install.packages("estimate", repos=rforge, dependencies=TRUE)
library(estimate)

然后准备好表达矩阵(比如TCGA的counts矩阵):

##         TCGA-OR-A5JP-01A TCGA-OR-A5JG-01A TCGA-OR-A5K1-01A
## MT-RNR2           810396          1190259          1206077
## MT-CO1            579888          1298037          1400198
## MT-ND4            623896           768059          1050890

紧接着就可以套用一下的代码进行分析(该段代码采自生信技能树的曾老师):

rd<-read.csv('counts_exp.csv', header = T)
rownames(rd)<-rd[,1]
rd<-rd[,-1]

dat=log2(rd+1)
library(estimate)
estimate <- function(dat,pro){
  input.f=paste0(pro,'_estimate_input.txt')
  output.f=paste0(pro,'_estimate_gene.gct')
  output.ds=paste0(pro,'_estimate_score.gct')
  write.table(dat,file = input.f,sep = '\t',quote = F)
  library(estimate)
  filterCommonGenes(input.f=input.f,
                    output.f=output.f ,
                    id="GeneSymbol")
  estimateScore(input.ds = output.f,
                output.ds=output.ds,
                platform="illumina")   ## 注意这个platform参数
  scores=read.table(output.ds,skip = 2,header = T)
  rownames(scores)=scores[,1]
  scores=t(scores[,3:ncol(scores)])
  return(scores)
}
pro='your_project_name'
scores=estimate(dat,pro)

计算完成后你就会得到如下矩阵:

##                  StromalScore ImmuneScore ESTIMATEScore
## TCGA.OR.A5JP.01A    -773.8226  -1143.9749    -1917.7975
## TCGA.OR.A5JG.01A    -878.7773   -685.5286    -1564.3059
## TCGA.OR.A5K1.01A    -663.8511   -360.2218    -1024.0729
## TCGA.OR.A5JR.01A    -931.0601   -344.3306    -1275.3907
## TCGA.OR.A5KU.01A    -925.6045  -1222.4672    -2148.0717
## TCGA.OR.A5L9.01A    -247.1255    404.5509      157.4254

大功告成!

46人读过
生物信息 R语言
Nov. 12, 2021, 3:17 p.m.

评论:

登录后方可评论,点击登录注册

评论列表:

暂无评论。

苏公网安备 32050602011302号

苏ICP备2020062135号-1
Copyright© Li Shen. All rights reserved.