問題如下:
下載TCGA腫瘤樣本揖庄,但該數(shù)據(jù)中無(wú)正常組織栗菜,合并GETx數(shù)據(jù)進(jìn)行免疫lncRNA提取,代碼部分如下:
rm(list = ls())
library(limma)
corFilter=0.4????
pvalueFilter=0.001??
setwd("C:\\Users\\CXH\\Desktop\\120irlncRNA\\09.irRNA")?
#讀取lncRNA蹄梢,并對(duì)數(shù)據(jù)進(jìn)行處理
rt = read.table("lncRNA.txt", header=T, sep="\t", check.names=F)
rt=as.matrix(rt)
rownames(rt)=rt[,1]
exp=rt[,2:ncol(rt)]
dimnames=list(rownames(exp),colnames(exp))
data=matrix(as.numeric(as.matrix(exp)),nrow=nrow(exp),dimnames=dimnames)
data=avereps(data)
data=data[rowMeans(data)>0.5,]
#刪掉正常樣品
group=sapply(strsplit(colnames(data),"\\-"),"[",4)
class(group)
list(group)
group=sapply(strsplit(group,""), "[", 1)
group=gsub("2","1",group)
lncRNA=data[,group==0]
#讀取免疫基因表達(dá)文件疙筹,并對(duì)數(shù)據(jù)進(jìn)行處理
rt = read.table("immGeneExp.txt", header=T, sep="\t", check.names=F)
rt=as.matrix(rt)
rownames(rt)=rt[,1]
exp=rt[,2:ncol(rt)]
dimnames=list(rownames(exp),colnames(exp))
immuneGene=matrix(as.numeric(as.matrix(exp)),nrow=nrow(exp),dimnames=dimnames)
immuneGene=avereps(immuneGene)
immuneGene=immuneGene[rowMeans(immuneGene)>0.5,]
#刪掉正常樣品
group=sapply(strsplit(colnames(immuneGene),"\\-"),"[",4)
group=sapply(strsplit(group,""),"[",1)
group=gsub("2","1",group)
immuneGene=immuneGene[,group==0]
#相關(guān)性檢驗(yàn)
outTab=data.frame()
for(i in row.names(lncRNA)){
?if(sd(lncRNA[i,])>0.5){
??for(j in row.names(immuneGene)){
???x=as.numeric(lncRNA[i,])
???y=as.numeric(immuneGene[j,])
???corT=cor.test(x,y)
???cor=corT$estimate
???pvalue=corT$p.value
???if((cor>corFilter) & (pvalue<pvalueFilter)){
????outTab=rbind(outTab,cbind(immuneGene=j,lncRNA=i,cor,pvalue,Regulation="postive"))
???}
???if((cor< -corFilter) & (pvalue<pvalueFilter)){
????outTab=rbind(outTab,cbind(immuneGene=j,lncRNA=i,cor,pvalue,Regulation="negative"))
???}
??}
?}
}
但是運(yùn)行相關(guān)性檢驗(yàn)后,結(jié)果:Error in if ((cor > corFilter) & (pvalue < pvalueFilter)) { :?
?需要TRUE/FALSE值的地方不可以用缺少值
此外: Warning message:
In cor(x, y) : 標(biāo)準(zhǔn)差為零
反復(fù)排查禁炒,初步考慮刪掉正常樣品部分的代碼需要修改而咆,但怎么嘗試都不行,請(qǐng)大神指點(diǎn):
#刪掉正常樣品
group=sapply(strsplit(colnames(immuneGene),"\\-"),"[",4)
group=sapply(strsplit(group,""),"[",1)
group=gsub("2","1",group)
immuneGene=immuneGene[,group==0]
因?yàn)門CGA命名規(guī)律:TCGA-20-0987-01A
而GETx命名規(guī)律:GTEX-PWCY-1326-SM-48TCU
還請(qǐng)大神們不惜賜教幕袱,謝謝暴备!謝謝!