相關(guān)軟件，比如gapit，rMVP，都會自動出圖，而GEMMA，GCTA則是需要后期自己作圖。
無論是軟件自動出圖，還是需要自己作圖，學(xué)習(xí)根據(jù)GWAS結(jié)果手動作圖都是必須的。
我們一般使用qqman作圖和cmplot兩個(gè)包畫GWAS的QQ圖和曼哈頓圖，后者顏色更漂亮。
這篇博客，介紹一下這兩個(gè)包如何畫GWAS的結(jié)果可視化圖。
第一個(gè)是qqman， 因?yàn)檫@個(gè)軟件函數(shù)很方便。
「安裝qqman軟件」
# Install the stable release from CRAn
install.packages("qqman")
# 直接從Github上安裝
devtools::install_github("stephenturner/qqman")
# 從Github上安裝最新的開發(fā)版
devtools::install_github("stephenturner/qqman",ref="dev")
如果你只知道qqman，而不知道cmplot的話，這說明你還不太專業(yè)。
「安裝cmplot代碼」
install.packages("CMplot")
看一下cmplot包的出圖：
另外，cmplot也可以出這樣的圖：
示例數(shù)據(jù)
使用qqman中的gwasResults數(shù)據(jù)集：
library(qqman)
data("gwasResults")
dat = gwasResults
head(dat)
str(dat)
table(dat$CHR)
數(shù)據(jù)如下：
如果是自己的數(shù)據(jù)，整理成上面的格式，名稱也變?yōu)樯厦娴拿Q。
第一列是SNP的ID，SNP
第二列是染色體，CHR
第三列是物理位置，BP
第四列是P值，P
qqman作圖
「QQ圖繪制」
這里，只需要一列P值即可。
qq(dat$P)
「曼哈頓圖」如果數(shù)據(jù)結(jié)構(gòu)如上所示，直接調(diào)用數(shù)據(jù)即可：
manhattan(dat)
當(dāng)然，更通用的是指定染色體、物理位置、P值：
manhattan(dat,
chr = "CHR", # 染色體
bp = "BP", # 物理位置
p = "P", # P值
)
這里，還可以設(shè)置閾值線，默認(rèn)是5和8：
?suggestiveline
Where to draw a "suggestive" line. Default -log10(1e-5). Set to FALSE to disable.
?
?genomewideline
Where to draw a "genome-wide sigificant" line. Default -log10(5e-8). Set to FALSE to disable.
?
cmplot作圖
我們用同樣的數(shù)據(jù)，使用cmplot作圖。
「qq圖繪制」
CMplot(dat,plot.type = "q",threshold = 0.05)
對比一下cmplot和qqman的QQ圖：可以看到，cmplot的QQ圖更好看，而且還有置信區(qū)間！666
「曼哈頓圖：」代碼：
CMplot(dat,plot.type = "m",threshold = c(0.01,0.05)/nrow(dat),threshold.col=c('grey','black'),
threshold.lty = c(1,2),threshold.lwd = c(1,1), amplify = T,
signal.cex = c(1,1), signal.pch = c(20,20),signal.col = c("red","orange"))
比較一下兩個(gè)軟件的曼哈頓圖：可以看到，cmplot的作圖更好看，不同染色體不同的顏色，不同閾值的結(jié)果不同的顏色。
高級作圖
「SNP密度圖」
CMplot(dat,plot.type = "d",bin.size = 1e6, col = c("darkgreen","yellow","red"))
「環(huán)形曼哈頓圖：」
CMplot(dat,plot.type="c",r=0.5,threshold=c(0.01,0.05)/nrow(pig60K),cex = 0.5,
threshold.col = c("red","orange"), threshold.lty = c(1,2),amplify = T, cir.chr.h = 2,
signal.cex = c(2,2), signal.pch = c(19,20), signal.col=c("red","green"),outward=TRUE)
注意，這個(gè)圖，很多參數(shù)都是默認(rèn)的，具體還可以再設(shè)置美觀一下，但是輪廓圖是有了！
合并密度圖和圓形曼哈頓圖：
CMplot(dat,plot.type="c",r=0.4,col=c("grey30","grey60"),chr.labels=paste("Chr",c(1:22),sep=""),
threshold=c(1e-6,1e-4),cir.chr.h=1.5,amplify=TRUE,threshold.lty=c(1,2),threshold.col=c("red","blue"),signal.line=1,signal.col=c("red","green"),chr.den.col=c("darkgreen","yellow","red"),
bin.size=1e6,outward=FALSE,file="jpg",memo="",dpi=300,file.output=TRUE,verbose=TRUE)
