Single Consolidated Output for Multiple Traits in GAPIT GWAS

I like GAPIT for some occasional GWAS I need to run.

Why?
Because it’s one of the more user-friendly packages for GWAS work.
Because the degrees of customization it offers within a single package are great.
And also because it outputs every single thing. Some unnecessary things too, sure, but the convenience of having every single output is better than not.

Though, if you have several traits for a set of genotypes, you know that you will get a different file for each trait. Consolidating the gwas results for all traits into one single file could be painful if you have say, 3 or more traits. I usually have 10+ traits (different traits from multiple locations and years/seasons) so opening each output file and carrying out copy-pasta is uncool and therefore not for me. Hence, this post! 😉

Not sure if you source the gapit functions script from zzlab: http://zzlab.net/GAPIT/gapit_functions.txt
I don’t. I save the script as ‘gapit_source_code.txt’ locally and call this instead.

Once you save the text file, go to line 5205.
Alternatively, search for this line:

Memory=GAPIT.Memory(Memory=Memory,Infor="Extract GWAS start")

as line numbers may get changed in future versions of the text file.

Then right before the line where it says

if(!byPass){

insert this snippet:

new.PWI.Filtered=cbind(GI,ps,maf,nobs,rsquare_base,rsquare)[index,]
colnames(new.PWI.Filtered)=c("SNP","Chromosome","Position ",
                              paste(name.of.trait,"_P.value",sep=''),
                              paste(name.of.trait,"_maf",sep=''),
                              paste(name.of.trait,"_nobs",sep=''),
                              paste(name.of.trait,"_Rsquare.of.Model.without.SNP",sep=''),
                              paste(name.of.trait,"_Rsquare.of.Model.with.SNP",sep=''))
assign(paste('trait_',name.of.trait,'_gwas_results',sep=''),
              new.PWI.Filtered, envir=globalenv())

Save and close the text file, or save it under a different name but be sure to call the altered source file for your work!

Alternatively, download this file and copy the code, paste into text editor and save as a txt file. In this file, I have only added the few lines of code shown above; none of the functions are altered.

Now, in your R script — the one you are using to import genotype, phenotype, etc — insert these two lines AFTER you run the GAPIT function.

traits_list <- mget(ls(pattern = "trait_.*"))
all_traits_GWAS <- Reduce(function(...) merge(..., all = T), traits_list)

For example, this is how my script looks:

library(multtest)
library(gplots)
library(LDheatmap)
library(genetics)
library(compiler) #this library is already installed in R

source("C:\\Users\\jujumaan\\Desktop\\gapit\\gapit_emma.txt")
source("C:\\Users\\jujumaan\\Desktop\\gapit\\biocLite.R")
source("C:\\Users\\jujumaan\\Desktop\\gapit\\gapit_source_code.txt") 
#Remember, this is the altered source code

myY <- read.table("pheno.txt", head = T)
myG <- read.table("geno.txt", head = F)

myGAPIT <- GAPIT(Y=myY, G=myG, PCA.total=0, Model.selection = F)
#myGAPIT <- GAPIT(Y=myY, G=myG, PCA.total=10, Model.selection = T)

traits_list <- mget(ls(pattern = "trait_.*")) #make a list of all traits
all_traits_GWAS <- Reduce(function(...) merge(..., all = T), traits_list) #merge all lists into a single data frame
write.table(all_traits_GWAS, "all_traits_GWAS.txt", sep = "\t", quote = F, row.names = F)

The “all_traits_GWAS.txt” file is the consolidated output with gwas results from each individual trait, which is much more convenient and an efficient way of doing things.

If this isn’t working for you, feel free to reach out.
Cheers!

Advertisements

याँ मुन्तिर कमेन्ट ठोक्नुस्

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google+ photo

You are commenting using your Google+ account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

w

Connecting to %s