我们GWAS常用的软件GAPIT很方便,但是如果表型数据多,将最后的数据进行汇总和处理将也是一件很麻烦的事情。这里,我尝试编写了一个小脚本,可以在GAPIT成功运行之后,将所有的result的表格进行汇总,并自动计算-logP值。可以简化最后的分析过程。
以下是代码,复制到R中即可运行。
myfiles <- Sys.glob("GAPIT*.GWAS.Results.csv")
print(paste(myfiles,"has been found")) #抓取所有名字是results的结果文档
output<-read.csv("origin_mini.csv",sep = ",",header = TRUE)
for(i in myfiles) #循环汇总
{new_col<-read.csv(file=i,header = TRUE,sep = ",")
tmp1<-new_col[,c(1,4)]
tmp1[,2]<-(-log10(tmp1[,2])) #计算-logP
output<-merge(output,tmp1,by="SNP")
output<-rename(output,c(P.value=i))
}
write.csv(output,"output_-log(P).csv") #输出csv格式汇总表格
注意:
运行之前需要进入结果的目录,CD进入目标文件夹即可。
运行还需origin_mini.csv这个初始文件,已经上传到群文件夹。
运行需要reshape这个R包。可以运行以下代码安装并加载这个包。
install.packages("reshape")
library("reshape")
我们下回再见!
网友评论