qtlnet (version 1.5.4)

mcmc.qtlnet: Sample genetic architecture and QTL network

Description

Use MCMC to alternatively sample genetic architecture and QTL network as directed acyclic graphs (DAGs).

Usage

mcmc.qtlnet(cross, pheno.col, threshold, addcov = NULL, intcov = NULL,
  nSamples = 1000, thinning = 1, max.parents = 3, M0 = NULL,
  burnin = 0.1, method = "hk", random.seed = NULL, init.edges = 0,
  saved.scores = NULL, rev.method = c("nbhd", "node.edge", "single"),
  verbose = FALSE, …)
init.qtlnet(pheno.col, max.parents, init.edges)

Arguments

cross

Object of class cross. See read.cross.

pheno.col

Phenotype identifiers from cross object. May be numeric, logical or character.

threshold

Scalar or list of thresholds, one per each node.

addcov

Additive covariates for each phenotype (NULL if not used). If entered as scalar or vector (same format as pheno.col), then the same addcov is used for all phenotypes. Altenatively, may be a list of additive covariate identifiers.

intcov

Interactive covariates, entered in the same manner as addcov.

nSamples

Number of samples to record.

thinning

Thinning rate. Number of MCMC samples is nSamples*thinning.

max.parents

Maximum number of parents to a node. This reduces the complexity of graphs and shortens run time. Probably best to consider values of 3-5.

M0

Matrix of 0s and 1s with initial directed graph of row->col if (row,col) entry is 1. Cycles are forbidden (e.g. 1s on diagonal or symmetric 1s across diagonal). Default (if NULL) is sampled by a call to init.qtlnet; all 0s if init.edges = 0 (default).

burnin

Proportion of MCMC samples to use as burnin. Default is 0.1 if burnin is TRUE. Must be between 0 and 1.

method

Model fitting method for scanone.

random.seed

Initialization seed for random number generator. Must be NULL (no reset) or positive numeric. Used in Random.

init.edges

Initial number of edges for M0, to be sampled using {init.qtlnet}. Chosen uniformly from 0 to the number of possible edges if set to NULL.

saved.scores

Updated scores, typically pre-computed by bic.qtlnet.

rev.method

Method to use for reversing edges. See details.

verbose

Print iteration and number of models fit.

…

Additional arguments. Advanced users may want to supply pre-computed saved.scores to speed up calculations.

Value

List of class qtlnet

post.model

Model code (see details).

post.bic

Posterior BIC

Mav

Model average of M across MCMC samples.

freq.accept

Frequency of acceptance M-H proposals.

saved.scores

Saved LOD score for each phenotype and all possible sets of the other phenotypes as parent nodes.

all.bic

cross

The cross object with calculated genotype probabilities.

In addition, a number of attributes are recorded:

M0

Initial network matrix.

threshold

threshold list

nSamples

Number of samples saved

thinning

Thinning rate

pheno.col

Phenotype columns.

pheno.names

Phenotype names

addcov

Additive covariate columns.

intcov

Interactive covariate columns.

burnin

Burnin proportion

method

Method used for scanone.

random.seed

Initial random number generator seed.

random.kind

Random number generator kind from Random.

Details

Models are coded compactly as (1)(2|1)(3|1,2,4,5)(4|2)(5|2). Each parenthetical entry is a of form (node|parents); these each require a model fit, for now with scanone.

The scanone routine is run on multiple phenotypes in the network that could all have the same parents. For instance, for 5 phenotypes, if (1|2,4) is sampled, then do scanone of this model as well as (3|2,4) and (5|2,4). Setting the hidden parameter scan.parents to a value smaller than length(pheno.col) - 1 (default) disallows multiple trait scanning with more than that number of parents.

The saved.scores parameter can greatly reduce MCMC run time, by supplying pre-computed BIC scores. See bic.qtlnet. Another option is to capture saved.scores from a previous mcmc.qtlnet run with the same phenotypes (and covariates). Caution is advised as only a modest amount of checking can be done.

The init.qtlnet routine can be used to randomly find an initial causal network M0 with up to init.edges edges.

MCMC updates include delete, add or reverse edge direction. The early version of this method only considered the edge on its own (rev.method = "single"), while the neighborhood method (rev.method = "nbhd") uses the update

References

Chiabub Neto E, Keller MP, Attie AD, Yandell BS (2010) Causal graphical models in systems genetics: a unified framework for joint inference of causal network and genetic archicecture for correlated phenotypes. Ann Appl Statist 4: 320-339. http://dx.doi.org/10.1214/09-AOAS288

Grzegorczyk and Husmeier (2008) Improving the structure MCMC sampler for Bayesian networks by introducing a new edge reversal move. Mach Learn 71: 265-305. http://dx.doi.org/10.1007/s10994-008-5057-7

See Also

read.cross, scanone, Random, bic.qtlnet.

Examples

# NOT RUN {
data(Pscdbp)
# }
# NOT RUN {
  ## Run of subset of traits. Still takes some time.
  Pscdbp.qtlnet <- mcmc.qtlnet(Pscdbp, pheno.col = c(1,2,4,5,6),
                               threshold = 3.83,
                               nSamples = 1000, thinning = 20, 
                               random.seed = 92387475, verbose = TRUE)
  save(Pscdbp.qtlnet, file = "Pscdbp.qtlnet.RData", compress = TRUE)
# }
# NOT RUN {
data(Pscdbp.qtlnet)

# }
# NOT RUN {
  out.qtlnet <- mcmc.qtlnet(Pscdbp, pheno.col = 1:13,
                            threshold = 3.83,
                            nSamples = 1000, thinning = 20, 
                            random.seed = 92387475, verbose = TRUE,
                            saved.scores = Pscdbp.bic)
# }