swap: MCMC sampling of multigene models

swapR Documentation

MCMC sampling of multigene models

Description

Given a k-gene model as a starting point, one gene is deleted and another is sampled in its place. This is done using an approximation to the posterior. Then another gene is deleted and another sampled,...

Usage

swap(varcov, invars, rparm, nreps, ana.obj, ...)

Arguments

varcov

The result of make.varcov

invars

Vector of numerical indexes of ana.obj$reg.names telling which variables to start in the model. The first of these is immediately removed, so it is merely a placeholder. The number of genes in the model is therefore k <- length(invars) (except when ana.obj$method=="F2" when it is k <- length(unique(col(ana.obj$reg.names)[invars])))

rparm

Scalar or vector with nrow(varcov$var.x) elements; the 'ridge' parameters for the independent variables - larger values imply more shrinkage or a more concentrated prior for the regresion coefficients.

nreps

How many cycles (of k samples each) to perform.

ana.obj

An analysis.object — see make.analysis.obj

...

Additional arguments override the default choices of candidate loci (locs), prior for locus (locs.prior), or method specified by ana.obj. Also, the default prior for model (combo.prior) when ana.obj$method=="F2" can be overridden. See swapbc1 and swapf2 for details.

Details

An MCMC sampler for loci using the object of make.varcov is executed. This sampler uses the exact posterior probability under the assumed correctness of the regression model using expected genotypes given marker values. This amounts to linearizing the likelihood with respect to the (possibly unknown) locus states. For models where the loci are fully informative markers this is the true posterior.

The chain is implemented as follows: given a set of regressor variables to start, one variable is removed, all regressor variables not in the model are examined to determine the effect of each on the posterior. One variable is sampled. The process is repeated until each variable has been removed and a new one sampled in its place (possibly the same variable that was removed is sampled). And this whole cycle is repeated nreps times.

Value

A list with components:

config

A k by k by nreps array (or, for ana.obj$method=="F2", a 2k by k by nreps array) of the locations (variables) sampled in each iteration.

posteriors

A vector of length k*nreps with the posteriors of the models.

coefs

A k by k matrix of the regression coefficients(or, for ana.obj$method=="F2", a 2k by nreps matrix).

call

The call to swap

cond

The k*nreps posterior probabilities of the k-1 gene models.

marg

The k*nreps marginal posteriors for all k gene models that could be formed using the current k-1 gene model

alt.marginal

A vector with length(locs) (or 2*length(locs)) elements. At each step, the posterior associated with each candidate locus is added to an element of this vector. After all steps are finished, the result is normalized to sum to one. This turns out to be a stable estimate of the marginal posterior.

alt.coef

A vector with length(locs) (or 2*length(locs)) elements. At each step, the product of each posterior times the coefficient(s) associated with a candidate locus is added to an element of this vector. After all steps are finished, the result is normalized by the total marginal posterior. This turns out to be a stable estimate of the marginal (over all models) posterior mean of the regression coefficients.

Author(s)

Charles C. Berry cberry@ucsd.edu

References

Berry C.C. (1998) Computationally Efficient Bayesian QTL Mapping in Experimental Crosses. ASA Proceedings of the Biometrics Section, 164-169.

Examples

data( little.ana.bc )
little.vc <- varcov( bc.phenotype~locus(all), little.ana.bc)
little.4 <- swap( little.vc, c(1,15,55,75), rparm=1, 50, little.ana.bc )
little.4.smry <- summary( little.4 )
print(c("Bayes Factor (3 vs 4)"=little.4.smry$ratio$mean))
par(mfrow=c(3,2))
plot( little.ana.bc, little.4.smry$loc.posterior, type="h",
 ylab="E(genes)" )
rm(little.4,little.vc,little.ana.bc)

bqtl documentation built on Sept. 25, 2024, 1:08 a.m.

Related to swap in bqtl...