Single-marker MLM GWAS using GAPIT

Overview

This example illustrates fitting a single marker GWAS model using GAPIT. The detailed documentation of GAPIT can be found at http://zzlab.net/GAPIT/

Load packages

# clear working environment
rm(list=ls()) 

# install and load support packages
source("http://www.zzlab.net/GAPIT/GAPIT.library.R") 

# load GAPIT function
source("http://www.zzlab.net/GAPIT/gapit_functions.txt")

Maize data

The maize data set used here is from a maize association panel including 281 diverse lines genotyped with 3,093 markers. The three phenotypes include ear height (EarHT), days to pollination (dpoll), and ear diameter (EarDia). In this example, we will only use EarHT.

# phenotypes
myY <- read.table(file = "http://zzlab.net/GAPIT/data/mdp_traits.txt", header = TRUE)

# marker matrix
myGD <- read.table(file = "http://zzlab.net/GAPIT/data/mdp_numeric.txt", header = TRUE)

# map information of markers 
myGM <- read.table(file = "http://zzlab.net/GAPIT/data/mdp_SNP_information.txt", header = TRUE)

Mixed linear model (MLM) GWAS

In this example, we will perform single marker GWAS analysis using a mixed linear model including the first three principle components (PCs) and the genomic relationship matrix. We will estimate the variance components only once prior to GWAS by setting SNP.3D = TRUE.

# GWAS with MLM 
gwasMLMfit <- GAPIT(Y = myY[, c(1:2)], GD = myGD, GM = myGM, model = 'MLM', 
                    PCA.total = 3, SNP.P3D = TRUE, SNP.MAF = 0.05)
names(gwasMLMfit)

Here, GWAS: gwas results, Pred: genomic prediction results, mc: a vector of snp effects, bc: a vector of fixed effects, mp: a vector of p values of markers, h2: heritability, PCA: first three principle components, GD: marker matrix after removing maf = 0.05, GM: marker location information (after removing maf =0.05), and KI: genomic relationship matrix.

Results

Use the read.csv() function to read the output file GAPIT.MLM.EarHT.GWAS.Results.csv. The columns 4, 9 and 10 are p-value, FDR adjusted p-value, and marker effect, respectively.

gwasMLM <- read.csv('GAPIT.MLM.EarHT.GWAS.Results.csv', header = TRUE)
head(gwasMLM)

##           SNP Chromosome  Position      P.value       maf nobs
## 1  PZA03188.4          1 280719882 0.0002445158 0.2831541  279
## 2       an1.5          1 175504926 0.0014496410 0.3817204  279
## 3  PZA00277.9          2 141942814 0.0023157133 0.1935484  279
## 4  PZB01892.1          3 161573186 0.0025736885 0.0483871  279
## 5  PZB00811.1          8 160536300 0.0025829162 0.1827957  279
## 6 PZA02945.10         10  33619582 0.0033319970 0.3960573  279
##   Rsquare.of.Model.without.SNP Rsquare.of.Model.with.SNP
## 1                     0.256831                 0.2941658
## 2                     0.256831                 0.2848032
## 3                     0.256831                 0.2823886
## 4                     0.256831                 0.2818473
## 5                     0.256831                 0.2818290
## 6                     0.256831                 0.2805293
##   FDR_Adjusted_P.values    effect
## 1             0.6511457 -5.815539
## 2             0.9612014  4.740143
## 3             0.9612014  5.123751
## 4             0.9612014 -8.511138
## 5             0.9612014 -7.801263
## 6             0.9612014  3.905459

Manhattan plot

Open the GAPIT.MLM.EarHT.Manhattan.Plot.Genomewise.pdf to view the Manhattan plot.