GPA: Fit GPA model
In carter-allen/GPA: GPA (Genetic analysis incorporating Pleiotropy and Annotation)

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/GPA.R

Fit GPA model.

GPA( gwasPval, annMat=NULL, pleiotropyH0=FALSE, empiricalNull=FALSE,
	maxIter=2000, stopping="relative", epsStopLL=1e-10,
	initBetaAlpha=0.1, initPi=0.1, initQ=0.75,
	lbPi=NA, lbBetaAlpha=0.001, lbQ=0.001, lbPval=1e-30,
	vDigit=1000, verbose=1 )

`gwasPval`	p-value matrix from GWAS data, where row and column correspond to SNP and phenotype, respectively.
`annMat`	Binary matrix from annotation data, where row and column correspond to SNP and annotation, respectively.
`pleiotropyH0`	Fit GPA under the null hypothesis of pleiotropy test? Possible values are `TRUE` (under the null hypothesis of pleiotropy test) or `FALSE` (usual assumption for GPA). Default is `FALSE`.
`empiricalNull`	Empirically estimate null distribution for GPA? Possible values are `TRUE` (empirical null distribution) or `FALSE` (theoretical null distribution). Default is `FALSE`.
`maxIter`	Maximum number of EM iteration. Default is 2000.
`stopping`	Stopping rule for EM iteration. Possible values are `"absolute"` (based on absolute difference in log likelihood), `"relative"` (based on relative difference in log likelihood), or `"aitken"` (based on Aitken acceleration-based stopping rule). Default is `"relative"`.
`epsStopLL`	Threshold to stop the EM iteration. Default is 1e-100.
`initBetaAlpha`	Initial value for alpha estimate. Default is 0.1.
`initPi`	Initial value for pi estimate. Default is 0.1.
`initQ`	Initial value for q estimate. Default is 0.75.
`lbPi`	Lower bound for pi estimate. If `lbPi=NA`, lower bound is set to 1 / [number of SNPs]. Default is `NA`.
`lbBetaAlpha`	Lower bound for alpha estimate. Default is 0.001.
`lbQ`	Lower bound for q estimate. Default is 0.001.
`lbPval`	Lower bound for GWAS p-value. Any GWAS p-values smaller than `lbPval` are set to `lbPval`. Default is 1e-30.
`vDigit`	Number of digits for reporting parameter estimates. For example, setting it to 1000 means printing out values up to three digits below zero.
`verbose`	Amount of progress report during the fitting procedure. Possible values are 0 (minimal output), 1, 2, or 3 (maximal output). Default is 1.

GPA fits the GPA model. It requires to provide GWAS p-value to gwasPval, while users can also provide annotation data to annMat. It is assumed that number of rows of matrix provided to gwasPval equals to that provided to annMat.

pTest implements the hypothesis testing for pleiotropy. It requires two GPA model fits, one of interest and one under the null hypothesis, and they can be obtained by setting pleiotropyH0=FALSE and pleiotropyH0=TRUE, respectively.

aTest implements the hypothesis testing for annotation enrichment. It requires two GPA model fits, one fitted with using annotation data and one fitted without using annotation data, and they can be obtained by providing annotation data to annMat and not, respectively.

Construct GPA class object.

Dongjun Chung

Chung D*, Yang C*, Li C, Gelernter J, and Zhao H (2014), "GPA: A statistical approach to prioritizing GWAS results by integrating pleiotropy information and annotation data," PLoS Genetics, 10: e1004787. (* joint first authors)

Kortemeier E, Ramos PS, Hunt KJ, Kim HJ, Hardiman G, and Chung D (2018), "ShinyGPA: An interactive and dynamic visualization toolkit for genetic studies," PLOS One, 13(1): e0190949.

assoc, pTest, aTest, GPA.

# simulator function

simulator <- function( risk.ind, nsnp=20000, alpha=0.6 ) {
  
  m <- length(risk.ind)
  
  p.sig <- rbeta( m, alpha, 1 )
  pvec <- runif(nsnp)
  pvec[ risk.ind ] <- p.sig
  
  return(pvec)
}

# run simulation

set.seed(12345)
nsnp <- 1000
alpha <- 0.3
pmat <- matrix( NA, nsnp, 5 )

pmat[,1] <- simulator( c(1:200), nsnp=nsnp, alpha=alpha )
pmat[,2] <- simulator( c(51:250), nsnp=nsnp, alpha=alpha )
pmat[,3] <- simulator( c(401:600), nsnp=nsnp, alpha=alpha )
pmat[,4] <- simulator( c(451:750), nsnp=nsnp, alpha=alpha )
pmat[,5] <- simulator( c(801:1000), nsnp=nsnp, alpha=alpha )

ann <- rbinom(n = nrow(pmat), size = 1, prob = 0.15)
ann <- as.matrix(ann,ncol = 1)

# GPA without annotation data

fit.GPA.noAnn <- GPA( pmat, NULL, maxIter = 100 )
cov.GPA.noAnn <- cov( fit.GPA.noAnn )
		
# GPA with annotation data

fit.GPA.wAnn <- GPA( pmat, ann, maxIter = 100 )
cov.GPA.wAnn <- cov( fit.GPA.wAnn )

# GPA under the null hypothesis of pleiotropy test

fit.GPA.pleiotropy.H0 <- GPA( pmat, NULL, pleiotropyH0=TRUE, maxIter = 100 )