多算法免疫浸润

library(immunedeconv)library(job)library(tcltk)library(tidyverse)library(ggsci)library(scales)library(ggtext)library(reshape2)library(limma)library(parallel)library(ggplot2)library(ggpubr)library(future.apply)library(car)library(ridge)library(e1071)library(preprocessCore)library(foreach)  library(doParallel)library(data.table)

#生成随机颜色randomColor <- function() {  paste0("#",paste0(sample(c(0:9, letters[1:6]), 6, replace = TRUE),collapse = ""))}
# 生成200个随机颜色randomColors <- replicate(200,randomColor())

# 读入表达数据exprMatrix <- read.table("TCGA.txt",header=TRUE,sep = "\t",as.is = T,row.names = 1)#仅保留肿瘤数据exprMatrix= exprMatrix%>% dplyr::select(str_which(colnames(.), ".01"))%>%as.data.frame()

cancer="BRCA"#TCGA癌症缩写brca_vector <- rep(cancer, ncol(exprMatrix))

exp=exprMatrixres_estimate <- deconvolute(exp, method="estimate",tumor = T,indications = brca_vector)write.table(res_estimate,"estimate.txt",sep="\t",row.names = F)res_estimate$cell_type=paste0(res_estimate$cell_type,"_ESTIMATE")

exp=exprMatrixres_timer <- deconvolute(exp, method="timer",tumor = T,indications = brca_vector)write.table(res_timer,"timer.txt",sep="\t",row.names = F)res_timer$cell_type=paste0(res_timer$cell_type,"_TIMER")

exp=exprMatrixres_abis <- deconvolute(exp, method="abis",tumor = T,indications = brca_vector)write.table(res_abis,"abis.txt",sep="\t",row.names = F)res_abis$cell_type=paste0(res_abis$cell_type,"_ABIS")

exp=exprMatrixres_consensus_tme <- deconvolute(exp, method="consensus_tme",tumor = T,indications = brca_vector)write.table(res_consensus_tme,"consensus_tme.txt",sep="\t",row.names = F)res_consensus_tme$cell_type=paste0(res_consensus_tme$cell_type,"_ConsensusTME")

exp=exprMatrixres_xcell <- deconvolute(exp, method="xcell",tumor = T,indications = brca_vector)write.table(res_xcell,"xcell.txt",sep="\t",row.names = F)res_xcell$cell_type=paste0(res_xcell$cell_type,"_xCell")

exp=exprMatrixres_epic <- deconvolute(exp, method="epic",tumor = T,indications = brca_vector)write.table(res_epic,"epic.txt",sep="\t",row.names = F)res_epic$cell_type=paste0(res_epic$cell_type,"_EPIC")

exp=exprMatrixres_quantiseq <- deconvolute(exp, method="quantiseq",tumor = T,indications = brca_vector) write.table(res_quantiseq,"quantiseq.txt",sep="\t",row.names = F)res_quantiseq$cell_type=paste0(res_quantiseq$cell_type,"_quanTIseq")

data=exprMatrixv=voom(data, plot=F, save.plot=F)out=v$Eout=rbind(ID=colnames(out), out)write.table(out,file="uniq.symbol.txt",sep="\t",quote=F,col.names=F)source("CIBERSORT.R",encoding = "utf-8")results=CIBERSORT("LM22.txt", "uniq.symbol.txt", perm=1000, QN=TRUE)res_CIBERSORT=read.table("CIBERSORT-Results.txt",sep="\t",header=T,row.names=1,check.names=F)res_CIBERSORT=t(res_CIBERSORT)rownames(res_CIBERSORT)=paste0(rownames(res_CIBERSORT),"_CIBERSORT")res_CIBERSORT=data.frame(cell_type=rownames(res_CIBERSORT),res_CIBERSORT)

res_consensus_tme=res_consensus_tme[1:(nrow(res_consensus_tme)-1),]res_xcell=res_xcell[1:(nrow(res_xcell)-3),]res_CIBERSORT=res_CIBERSORT[1:(nrow(res_CIBERSORT)-3),]

res=rbind(res_abis,res_consensus_tme,res_epic,res_quantiseq,res_timer,res_xcell,res_CIBERSORT)res=as.matrix(res)rownames(res)=res[,1]exp=res[,2:ncol(res)]dimnames=list(rownames(exp),colnames(exp))data=matrix(as.numeric(as.matrix(exp)),nrow=nrow(exp),dimnames=dimnames)

#输入要研究的基因gene="TP53"gene_exp=exprMatrix[gene,]gene_exp=t(gene_exp)
samesample=intersect(rownames(gene_exp),colnames(data))gene_exp=gene_exp[samesample,]gene_exp=as.data.frame(gene_exp)colnames(gene_exp)=genedata=data[,samesample]
x=as.numeric(gene_exp[,1])outTab=data.frame()for(i in rownames(data)){  y=as.numeric(data[i,])  if(sd(y)<0.001){next}  corT=cor.test(x, y, method="spearman")  cor=corT$estimate  pvalue=corT$p.value  if(pvalue<0.05){    outTab=rbind(outTab,cbind(immune=i, cor, pvalue))  }}#输出相关性结果write.table(file="corResult.txt", outTab, sep="\t", quote=F, row.names=F)

#绘制气泡图corResult=read.table("corResult.txt", head=T, sep="\t")corResult$Software=sapply(strsplit(corResult[,1],"_"), '[', 2)corResult$Software=factor(corResult$Software,level=as.character(unique(corResult$Software[rev(order(as.character(corResult$Software)))])))b=corResult[order(corResult$Software),]b$immune=factor(b$immune,levels=rev(as.character(b$immune)))colslabels=rep(hue_pal()(length(levels(b$Software))),table(b$Software))     #定义图形颜色#输出图形pdf(file="correlation.pdf", width=9, height=10)ggplot(data=b, aes(x=cor, y=immune, color=Software))+  labs(x="Correlation coefficient",y="Immune cell")+  geom_point(size=4.1)+  theme(panel.background=element_rect(fill="white",size=1,color="black"),        panel.grid=element_line(color="grey75",size=0.5),        axis.ticks = element_line(size=0.5),        axis.text.y = ggtext::element_markdown(colour=rev(colslabels)))dev.off()