尔云间 一个专门做科研的团队

原创 小果 生信果 关注我们

小果今天的项目需要计算不同基因表达量间的相关性，是计算关注的ALKBH5基因与其他基因的相关性，撸起袖子开干

rm(list=ls());gc();
library(data.table)
data <-read.csv("counts478.CSV",header = T,row.names = 1)
#创建空向量
gene_name1<-c()
gene_name2<-c()
cor_r<-c()
pvalue<-c()

#for循环，12076是关注的ALKBH5在数据中的行数
for (r in 1:nrow(data)){
  g1=rownames(data)[12076]
  g2=rownames(data)[r]
  c_r=cor(as.numeric(data[12076,]),as.numeric(data[r,]),method="pearson")
  p=cor.test(as.numeric(data[12076,]),as.numeric(data[r,]),method ="pearson")[[3]]
  ##保存每一步的数据，而不可直接以空向量作为每一步运行的结果
  gene_name1=c(gene_name1,g1)
  gene_name2=c(gene_name2,g2)
  cor_r=c(cor_r,c_r)
  pvalue=c(pvalue,p)
}
#输导出结果文件
data_cor<-data.frame(gene_name1,gene_name2,cor_r,pvalue)
head(data_cor)
dim(data_cor)
write.csv(data_cor,"data_cor.CSV",row.names = T)

输入文件格式如下，列名为样本名，行名为基因名的矩阵

输出文件如下

如果我们需要计算基因表达量两两之间的相关性，把上面的循环换成一个嵌套循环即可，代码如下

#嵌套的for循环
for (i in 1:nrow(data)){
for (r in i:nrow(data)){
g1=rownames(data)[i]
g2=rownames(data)[r]
c_r=cor(as.numeric(data[i,]),as.numeric(data[r,]),method="pearson")
p=cor.test(as.numeric(data[i,]),as.numeric(data[r,]),method ="pearson")[[3]]
gene_name1=c(gene_name1,g1)
gene_name2=c(gene_name2,g2)
cor_r=c(cor_r,c_r)
pvalue=c(pvalue,p)
}
}

至此结束

推荐阅读

• 多种花样都能懂，富集图的看法
  
• 相关性图谱之相关性热图学习
  
• 小果教你快速优雅使用Genecards数据库
  
• 从UCSC下载TCGA数据
  
• 简单高效利用Batch Entrez批量获取基因别名
  

关注小果，小果将会持续为你带来更多生信干货哦。

“生信果”，生信入门、R语言、生信图解读与绘制、软件操作、代码复现、生信硬核知识技能、服务器、生物信息学的教程，以及基于R的分析和可视化等原创内容，一起见证小白和大佬的成长。