相比于網(wǎng)頁工具,使用編程語言處理科研數(shù)據(jù)的一大優(yōu)勢诅迷,在于高度的定制化教翩,以及批量處理數(shù)據(jù)的快捷性和高效性
目錄
批量處理——for循環(huán)批量計算組間差異
批量處理——apply批量計算組間差異
批量處理——for循環(huán)畫圖
批量處理——for循環(huán)遷移文件
對于批量處理數(shù)據(jù)的方法逮刨,之前使用for循環(huán)和apply語句進(jìn)行處理過秽晚,但是不夠系統(tǒng)瓦糟,學(xué)習(xí)果子生信課程后有一個清晰的認(rèn)識,寫下來赴蝇,一是可以調(diào)用方便菩浙,二是自己寫過之后,才能算是完全掌握句伶。當(dāng)然一切以解決問題為主劲蜻,不陷于技術(shù)深究。
批量計算基因和基因之間的相關(guān)性考余,也是一項(xiàng)很好的應(yīng)用先嬉。
場景
計算RNA-seq得到幾個基因與免疫細(xì)胞相關(guān)性
數(shù)據(jù)準(zhǔn)備
所需要的兩項(xiàng)數(shù)據(jù),一個是基因表達(dá)數(shù)據(jù)楚堤,另一個是免疫細(xì)胞浸潤的矩陣疫蔓,行名要一致,都是樣本名稱
1. 單個基因和免疫細(xì)胞相關(guān)性分析
這個和之前提到的單個基因和近2萬個基因求相關(guān)性是一樣的
1.1 for循環(huán)計算
gene <- "KLF5"
y <- as.numeric(expr_data[,gene])
### 批量操作的具體實(shí)現(xiàn)過程:
### 1.設(shè)定容器,最終生成的數(shù)據(jù)放在什么地方身冬?
correlation <- data.frame()
### 2.批量把數(shù)據(jù)導(dǎo)出到容器
for(i in 1:length(colnames(immu_data))){
## 1.指示
print(i)
## 2.計算
dd = cor.test(as.numeric(immu_data[,i]),y,method="spearman")
## 3.填充
correlation[i,1] = colnames(immu_data)[i]
correlation[i,2] = dd$estimate
correlation[i,3] = dd$p.value
}
### 修改名稱
colnames(correlation) <- c("cell","cor","p.value")
計算結(jié)果
1.2 lapply 計算
從批量處理——基因之間的相關(guān)性復(fù)制相關(guān)代碼衅胀,修改個別名稱
yourgene <- 'DDR1'
genelist <- colnames(immu_data)
# 寫一個函數(shù)
mycor <- function(x){
dd = cor.test(expr_data[,yourgene], immu_data[,x], method = 'spearman')
data.frame(gene = yourgene, cell = x, R = dd$estimate, p.value = dd$p.value)
}
# 測試一下
mycor(genelist[1])
# 批量
lapplylist <- lapply(genelist, mycor)
# do.call
cor_data <- do.call(rbind, lapplylist)
# 整理成一個函數(shù)
cor_data <- do.call(rbind, lapply(genelist, function(x){
dd = cor.test(COAD_data[,yourgene], COAD_data[,x], method = 'spearman')
data.frame(gene1 = yourgene, gene2 = x, R = dd$estimate, p.value = dd$p.value)
}))
得到的結(jié)果是這樣
結(jié)果一致
多個基因怎么計算
當(dāng)然可以接著寫一個循環(huán)函數(shù),其實(shí)已經(jīng)有人幫我們做好輪子酥筝,拿來用就好拗小。使用psych這個包(這個包應(yīng)該也是將計算內(nèi)容向量化,提高運(yùn)算速度)
sig_gene <- c("APOE","CASR","CTLA4", "CXCL8", "F2","IL6","MMP9")
library(psych)
x <- expr_data[,sig_gene]
y <- immu_data
library(psych)
d <- corr.test(x,y,use="complete",method = 'spearman')
r <- d$r
p <- d$p
library(ggcorrplot)
ggcorrplot(t(d$r), show.legend = T,
p.mat = t(d$p.adj), digits = 2, sig.level = 0.05,insig = 'blank',lab = T)
得到的結(jié)果如下
結(jié)果還是不錯的樱哼,速度快哀九,有可視化。這算是包裝過后的批量分析搅幅。
后記:研究問題阅束,一個是量的問題,一個是分組的問題茄唐。涉及到量的問題息裸,就會有差異性分析,相關(guān)性分析沪编。這是總體的相關(guān)性分析呼盆,接著將兩兩相關(guān)性圖形展現(xiàn)。