GWAS

GWAS中的曼哈顿图与QQ图(CMplot)

2022-11-13  本文已影响0人  曹草BioInfo

输入文件为exmma的输出ps文件,整理为ps2文件,格式如下


head common_LD_seed_color.ps2

调用命令为

Rscript script/manhattan_cmplot.R data/common_LD_seed_color.ps2 out/prefix
#!/usr/bin/env Rscript

args=commandArgs(TRUE)
if (length(args) != 2) {
          print ("usage: <gwas file(4 column)>   <out.prefix> ")
  q()
}

input <- args[1]
#max_y <- as.numeric(args[3])
prefix <- args[2]

library(CMplot)

gwas<- read.table(input, header = T);

cutoff <- 0.05/nrow(gwas)

sigSNP <- gwas[gwas[,4] < cutoff,]

write.table( sigSNP, file = paste( prefix, "sigSite.out", sep =  "."), row.names = F, quote = F)


png(paste(prefix, "_manhattan_threshold.png", sep = ""), width=960, height=480)
CMplot(gwas,
       plot.type = "m",
       LOG10 = T,
       col = c("blue4", "orange3"),
       cex = 0.3,
       ylab.pos = 2, 
       cex.axis = 1,
       threshold = c(0.01,0.05)/nrow(gwas), ## 显著性阈值
       threshold.col=c('grey','black'),  ## 阈值线颜色
       threshold.lty = c(1,2), ## 阈值线线型
       threshold.lwd = c(1,1), ## 阈值线粗细
       amplify = T,  ## 放大显著SNP
       signal.cex = c(1,1), ## 点大小
       signal.pch = c(20,20), ## 点形状
       signal.col = c("red","blue"), ## 点颜色
       file.output = F )
dev.off()

png(paste(prefix, "_manhattan.png", sep = ""), width=960, height=480)
CMplot(gwas,plot.type = "m",
       LOG10 = T,
       col = c("blue4", "orange3"),
       cex = 0.3,
       ylab.pos = 2, 
       cex.axis = 1,
       threshold=NULL,
       file.output = F )

dev.off()

png(paste(prefix, "_qqplot.png", sep = ""), width=480, height=480)
CMplot(gwas,
       plot.type = "q", ## 绘制QQplot
       box=T, ## 是否加边框
       conf.int=T, ## 是否绘制置信区间
       conf.int.col=NULL, ## 置信区间颜色
       threshold.col="red", ## 对角线颜色
       threshold.lty=2,  ## 线型
       cex = 0.8,
       ylab.pos = 2, 
       cex.axis = 1,
       main = "QQ-plot",
       file.output = F )
dev.off()


## 曼哈顿图1
pdf(paste(prefix, "_manhattan_threshold.pdf", sep = ""), width=10, height=5)
CMplot(gwas,
       plot.type = "m",
       LOG10 = T,
       col = c("blue4", "orange3"),
       cex = 0.3,
       ylab.pos = 2, 
       cex.axis = 1,
       threshold = c(0.01,0.05)/nrow(gwas), ## 显著性阈值
       threshold.col=c('grey','black'),  ## 阈值线颜色
       threshold.lty = c(1,2), ## 阈值线线型
       threshold.lwd = c(1,1), ## 阈值线粗细
       amplify = T,  ## 放大显著SNP
       signal.cex = c(1,1), ## 点大小
       signal.pch = c(20,20), ## 点形状
       signal.col = c("red","blue"), ## 点颜色
       file.output = F )
dev.off()

## 曼哈顿图2

pdf(paste(prefix, "_manhattan.pdf", sep = ""), width=10, height=5)
CMplot(gwas,plot.type = "m",
       LOG10 = T,
       col = c("blue4", "orange3"),
       cex = 0.3,
       ylab.pos = 2, 
       cex.axis = 1,
       threshold=NULL,
       file.output = F )

dev.off()

## QQplot 
pdf(paste(prefix, "_qqplot.pdf", sep = ""), width=10, height=10)
CMplot(gwas,
       plot.type = "q", ## 绘制QQplot
       box = T, ## 是否加边框
       conf.int=T, ## 是否绘制置信区间
       conf.int.col=NULL, ## 置信区间颜色
       threshold.col="red", ## 对角线颜色
       threshold.lty=2,  ## 线型
       cex = 0.8,
       ylab.pos = 2, 
       cex.axis = 1,
       main = "QQ-plot",
       file.output = F )
dev.off()
上一篇下一篇

猜你喜欢

热点阅读