Learn R Programming

baycn (version 1.0.0)

mhEdge: mhEdge

Description

The main function for the Metropolis-Hastings algorithm. It returns the posterior distribution of the edge directions.

Usage

mhEdge(data, adjMatrix, burnIn = 0.2, iterations = 1000, nGV,
  pmr = FALSE, prior = c(0.05, 0.05, 0.9), progress = TRUE,
  thinTo = 200)

Arguments

data

A matrix with the nodes across the columns and the observations down the rows.

adjMatrix

An adjacency matrix indicating the edges that will be considered by the Metropolis-Hastings algorithm. This can be the output from another algorithm (e.g., PC). An adjacency matrix is a matrix of zeros and ones. The ones represent an edge and its direction between two nodes.

burnIn

A number between 0 and 1 indicating the percentage of the sample that will be discarded.

iterations

An integer for the number of iterations to run the MH algorithm.

nGV

The number of genetic variants in the graph.

pmr

Logical. If true the Metropolis-Hastings algorithm will use the Principle of Mendelian Randomization (PMR). This prevents the direction of an edge pointing from a gene expression node to a genetic variant node.

prior

A vector containing the prior probability of each edge state.

progress

Logical. If TRUE a progress bar will be printed.

thinTo

An integer indicating the number of observations the chain should be thinned to.

Value

An object of class mcmc containing 9 elements:

  • burnIn -- The percentage of MCMC iterations that will be discarded from the beginning of the chain.

  • chain -- A matrix where each row contains the vector of edge states for the accepted graph.

  • decimal -- A vector of decimal numbers. Each element in the vector is the decimal of the accepted graph.

  • iterations -- The number of iterations for which the Metropolis-Hastings algorithm is run.

  • posteriorES -- A matrix of posterior probabilities for all three edge states for each edge in the network.

  • posteriorPM -- A posterior probability adjacency matrix.

  • likelihood -- A vector of log likelihood values. Each element in the vector is the log likelihood of the accepted graph.

  • stepSize -- The number of MCMC iterations discarded between each accepted graph.

  • time -- The runtime of the Metropolis-Hastings algorithm in seconds.

References

Martin, E. A. and Fu, A. Q. (2019). A Bayesian approach to directed acyclic graphs with a candidate graph. arXiv preprint arXiv:1909.10678.

Examples

Run this code
# NOT RUN {
# Generate data under topology m1_gv.
# Use ?simdata for a description and graph of m1_gv.
data_m1 <- simdata(b0 = 0,
                   N = 200,
                   s = 1,
                   graph = 'm1_gv',
                   ss = 1,
                   p = 0.27)

# Create an adjacency matrix with the true edges.
am_m1 <- matrix(c(0, 1, 0,
                  0, 0, 1,
                  0, 0, 0),
                byrow = TRUE,
                nrow = 3)

# Run the Metropolis-Hastings algorithm on the data from m1_gv using the
# Principle of Mendelian Randomization (PMR) and the true edges as the input.
mh_m1_pmr <- mhEdge(data = data_m1,
                    adjMatrix = am_m1,
                    burnIn = 0.2,
                    iterations = 1000,
                    nGV = 1,
                    pmr = TRUE,
                    prior = c(0.05,
                              0.05,
                              0.9),
                    progress = FALSE,
                    thinTo = 200)

summary(mh_m1_pmr)

# Generate data under topology gn4.
# Use ?simdata for a description and graph of gn4.
data_gn4 <- simdata(b0 = 0,
                    N = 200,
                    s = 1,
                    graph = 'gn4',
                    ss = 1)

# Create an adjacency matrix with the true edges.
am_gn4 <- matrix(c(0, 1, 1, 0,
                   0, 0, 0, 1,
                   0, 0, 0, 0,
                   0, 0, 1, 0),
                 byrow = TRUE,
                 nrow = 4)

# Run the Metropolis-Hastings algorithm on the data from gn4 with the true
# edges as the input.
mh_gn4 <- mhEdge(data = data_gn4,
                 adjMatrix = am_gn4,
                 burnIn = 0.2,
                 iterations = 1000,
                 nGV = 0,
                 pmr = FALSE,
                 prior = c(0.05,
                           0.05,
                           0.9),
                 progress = FALSE,
                 thinTo = 200)

summary(mh_gn4)

# }

Run the code above in your browser using DataLab