R语言制作GWAS-曼哈顿图-自编教程

自己制作的教程,简单易懂易操作!

R语言制作曼哈顿图-自编教程

一、TXT文件格式

（1）在excel格式（带表）头如下：

（2）“另存为”“文本文件（制表符分割）（*.txt）”

二、R语言制图

1、改变所在R语言工作目录（所在文件目录）

2、更改所编程语言的检验水平线（P=0.05和P=0.01）

自己制作的教程,简单易懂易操作!

在“1”和“2”处的取值分别代表着显著水平（P<0.05）和极显著水平（P<0.01），分别取值为：0.05÷SNP数目。

3、更改导入的“TXT”文件和保存的“pdf”文件名：

程序语言如下：（只输入双虚线中间的语言部分）

==========================================================

manhattan<- function(dataframe, pchsize=19, colors=c("gray10", "gray100"), ymax="max", limitchromosomes=1:19, suggestiveline=-log10(1.78916E-08), genomewideline=-log10(8.94582E-08), annotate=NULL, ...) {

d=dataframe

if (!("CHR" %in% names(d) & "BP" %in% names(d) & "P" %in% names(d))) stop("Make sure your data frame contains columns CHR, BP, and P")

if (any(limitchromosomes)) d=d[d$CHR %in% limitchromosomes, ]

d=subset(na.omit(d[order(d$CHR, d$BP), ]), (P>0 & P<=1)) # remove na's, sort, and keep only 0<P<=1

d$logp = -log10(d$P)

d$pos=NA

ticks=NULL

lastbase=0

colors<- rep(colors,max(d$CHR))[1:max(d$CHR)]

if (ymax=="max") ymax<-ceiling(max(d$logp))

if (ymax<8) ymax<-8

numchroms=length(unique(d$CHR))

if (numchroms==1) {

d$pos=d$BP

ticks=floor(length(d$pos))/2+1

自己制作的教程,简单易懂易操作!

} else {

for (i in unique(d$CHR)) {

if (i==1) {

d[d$CHR==i, ]$pos=d[d$CHR==i, ]$BP

} else {

lastbase=lastbase+tail(subset(d,CHR==i-1)$BP, 1)

d[d$CHR==i, ]$pos=d[d$CHR==i, ]$BP+lastbase

}

ticks=c(ticks, d[d$CHR==i, ]$pos[floor(length(d[d$CHR==i, ]$pos)/2)+1])

}

}

if (numchroms==1) {

with(d, plot(pos, logp, ylim=c(0,ymax), ylab=expression(-log[10](italic(p))), xlab=paste("Chromosome",unique(d$CHR),"position"), ...))

} else {

with(d, plot(pos, logp, ylim=c(0,ymax), ylab=expression(-log[10](italic(p))), xlab="Chromosome", xaxt="n", type="n", ...))

axis(1, at=ticks, lab=unique(d$CHR), ...)

icol=1

for (i in unique(d$CHR)) {

with(d[d$CHR==i, ],points(pos,pch=pchsize, logp, col=colors[icol], ...))

icol=icol+1

}

}

if (!is.null(annotate)) {

d.annotate=d[which(d$SNP %in% annotate), ]

with(d.annotate, points(pos, logp, col="green3", ...))

}

if (suggestiveline) abline(h=suggestiveline, col="blue")

if (genomewideline) abline(h=genomewideline, col="red")

}

d=read.table(file="GLM_眼肌面积_胴体重协变量.txt", header=TRUE)

pdf("GLM_眼肌面积_胴体重协变量.pdf");

manhattan(d, limitchromosomes=1:30,colors=c("#234589", "red","#431287", "green"), pchsize=19);

dev.off();

==========================================================