MR.RGM (version 0.0.4)

NetworkMotif: Estimating the uncertainty of a specified network

Description

The NetworkMotif function facilitates uncertainty quantification. Specifically, it determines the proportion of posterior samples that contains the given network structure. To use this function, users may use the GammaPst output obtained from the RGM function.

Usage

NetworkMotif(Gamma, GammaPst)

Value

The NetworkMotif function calculates the uncertainty quantification for the provided network structure. A value close to 1 indicates that the given network structure is frequently observed in the posterior samples, while a value close to 0 suggests that the given network structure is rarely observed in the posterior samples.

Arguments

Gamma

A matrix of dimension p * p that signifies a specific network structure among the response variables, where p represents the number of response variables. This matrix is the focus of uncertainty quantification.

GammaPst

An array of dimension p * p * n_pst, where n_pst is the number of posterior samples and p denotes the number of response variables. It comprises the posterior samples of the causal network among the response variables. This input might be obtained from the RGM function. Initially, execute the RGM function and save the resulting GammaPst. Subsequently, utilize this stored GammaPst as input for this function.

References

Ni, Y., Ji, Y., & Müller, P. (2018). Reciprocal graphical models for integrative gene regulatory network analysis. Bayesian Analysis, 13(4), 1095-1110. tools:::Rd_expr_doi("10.1214/17-BA1087").

Examples

#' # ---------------------------------------------------------

# Example 1:
# Run NetworkMotif to do uncertainty quantification for a given network among the response variable

# Data Generation
set.seed(9154)

# Number of data points
n = 10000

# Number of response variables and number of instrument variables
p = 3
k = 4

# Initialize causal interaction matrix between response variables
A = matrix(sample(c(-0.1, 0.1), p^2, replace = TRUE), p, p)

# Diagonal entries of A matrix will always be 0
diag(A) = 0

# Make the network sparse
A[sample(which(A!=0), length(which(A!=0))/2)] = 0

# Create D matrix (Indicator matrix where each row corresponds to a response variable
# and each column corresponds to an instrument variable)
D = matrix(0, nrow = p, ncol = k)

# Manually assign values to D matrix
D[1, 1:2] = 1  # First response variable is influenced by the first 2 instruments
D[2, 3] = 1    # Second response variable is influenced by the 3rd instrument
D[3, 4] = 1    # Third response variable is influenced by the 4th instrument


# Initialize B matrix
B = matrix(0, p, k)  # Initialize B matrix with zeros

# Calculate B matrix based on D matrix
for (i in 1:p) {
  for (j in 1:k) {
    if (D[i, j] == 1) {
      B[i, j] = 1  # Set B[i, j] to 1 if D[i, j] is 1
    }
  }
}

Sigma = diag(p)

Mult_Mat = solve(diag(p) - A)

Variance = Mult_Mat %*% Sigma %*% t(Mult_Mat)

# Generate instrument data matrix
X = matrix(rnorm(n * k, 0, 1), nrow = n, ncol = k)

# Initialize response data matrix
Y = matrix(0, nrow = n, ncol = p)

# Generate response data matrix based on instrument data matrix
for (i in 1:n) {

    Y[i, ] = MASS::mvrnorm(n = 1, Mult_Mat %*% B %*% X[i, ], Variance)

}


# Apply RGM on individual level data with Spike and Slab Prior
Output = RGM(X = X, Y = Y, D = D, prior = "Spike and Slab")

# Store GammaPst
GammaPst = Output$GammaPst

# Define a function to create smaller arrowheads
smaller_arrowheads = function(graph) {
    igraph::E(graph)$arrow.size = 1  # Adjust the arrow size value as needed
    return(graph)
}


# Plot the true graph structure between response variables
plot(smaller_arrowheads(igraph::graph_from_adjacency_matrix((A != 0) * 1,
  mode = "directed")), layout = igraph::layout_in_circle, main = "True Graph")

# Start with a random subgraph
Gamma = matrix(0, nrow = p, ncol = p)
Gamma[2, 1] = 1

# Plot the subgraph to get an idea about the causal network
plot(smaller_arrowheads(igraph::graph_from_adjacency_matrix(Gamma,
         mode = "directed")), layout = igraph::layout_in_circle,
            main = "Subgraph")


# Do uncertainty quantification for the subgraph
NetworkMotif(Gamma = Gamma, GammaPst = GammaPst)