ssa (version 1.3.0)

nebula.chisq.bin.train: Nonparametric empirical Bayes classifier using latent annotations: chi-square test statistics and binary indicators; training

Description

Assumes that chi-square test statistics for each SNP are available from another study, and binary indicators for each SNP are available as well. Treats the true control and case minor allele frequencies and the chi-square non-centrality parameters as random triples from a bivariate prior distribution G0 for SNPs with indicators 0, and from G1 for SNPs with indicators equal to 1. Estimates the optimal Bayesian classifier given G. Nonparametric maximum likelihood is used as a plug-in estimator for G.

Usage

nebula.chisq.bin.train(pi0, pi1, n0, n1, T, I, d = 25, maxit = 200,
  tol = 1e-04, verbose = FALSE)

Arguments

pi0, pi1

p x 1 vectors of control and case minor allele frequencies, respectively; IMPORTANT: must be relative to the same allele in both cases and controls

n0, n1

number of controls and number of cases, respectively

T

p x 1 vector of chi-square test statistics

I

p x 1 vector of binary indicators

d

if a single number, G0 and G1 are estimated on d x d x d grids; if a three-component vector (d0,d1,dt), G0 and G1 are estimated on d0 x d1 x dt grids

maxit

maximum number of EM iterations

tol

error tolerance

verbose

TRUE to print the error attained by each EM iteration

Value

nebula0

output of nebula.chisq.train using only SNPs with I==0

nebula1

output of nebula.chisq.train using only SNPs with I==1

I

binary indicator

Examples

# NOT RUN {
p <- 1000; ## number of snps
I <- rep(0,p); I[1:10] <- 1; ## which snps are causal
set.seed(1); pi0 <- runif(p,0.1,0.5); ## control minor allele frequencies
set.seed(1); ors <- runif(sum(I),-1,1); ## odds ratios
pi1 <- pi0;
pi1[I==1] <- expit(ors+logit(pi0[I==1]));
set.seed(1); lam <- rep(0,p); lam[I==1] <- rchisq(sum(I==1),1,25); ## ncps
## training data
n0 <- 100; ## number of controls
X0 <- t(replicate(n0,rbinom(p,2,pi0))); ## controls
n1 <- 50; ## number of cases
X1 <- t(replicate(n1,rbinom(p,2,pi1))); ## cases
T <- rchisq(p,1,lam); ## chi-square statistics
nebula <- nebula.chisq.bin.train(colMeans(X0)/2,colMeans(X1)/2,n0,n1,T,I,d=c(10,15,20));
par(mfrow=c(2,3));
contour(nebula$nebula0$Pi0,nebula$nebula0$Pi1,apply(nebula$nebula0$g,c(1,2),sum));
points(pi0[I==0],pi1[I==0]);
contour(nebula$nebula0$Pi0,nebula$nebula0$Lam,apply(nebula$nebula0$g,c(1,3),sum));
points(pi0[I==0],lam[I==0]);
contour(nebula$nebula0$Pi1,nebula$nebula0$Lam,apply(nebula$nebula0$g,c(2,3),sum));
points(pi1[I==0],lam[I==0]);
contour(nebula$nebula1$Pi0,nebula$nebula1$Pi1,apply(nebula$nebula1$g,c(1,2),sum));
points(pi0[I==1],pi1[I==1]);
contour(nebula$nebula1$Pi0,nebula$nebula1$Lam,apply(nebula$nebula1$g,c(1,3),sum));
points(pi0[I==1],lam[I==1]);
contour(nebula$nebula1$Pi1,nebula$nebula1$Lam,apply(nebula$nebula1$g,c(2,3),sum));
points(pi1[I==1],lam[I==1]);

# }