| Title: | Bayesian Inference of Boolean Genetic Networks |
|---|---|
| Description: | Implements a fully Bayesian Markov chain Monte Carlo (MCMC) approach for inferring the topology and Boolean logic transition functions of gene regulatory networks from noisy, binary time-series expression data. Network structure and Boolean rules are sampled jointly from their posterior distribution, providing principled uncertainty quantification rather than a single point estimate. Method described in Han et al. (2014) <doi:10.1371/journal.pone.0115806>. |
| Authors: | Anson Li [aut, cre], Shengtong Han [aut] |
| Maintainer: | Anson Li <[email protected]> |
| License: | BSD_3_clause + file LICENSE |
| Version: | 0.2.0 |
| Built: | 2026-07-17 18:47:53 UTC |
| Source: | https://github.com/anson-li8/bbni |
Randomly generates a valid directed acyclic graph (DAG) topology
and assigns a corresponding Boolean transition function to each node.
The algorithm samples parent set configurations, keeping the constraint
that the maximum in-degree for any node is 2, and further ensures
the resulting structure does not contain directed cyclic loops.
GenerateNetwork(num.node)GenerateNetwork(num.node)
num.node |
An integer representing the total number of genes/nodes in the network. |
A square transition function matrix combining the initial DAG topology with randomly assigned Boolean logic functions. Elements with a value of 0 indicate no directed edge, while positive integers indicate the presence of an edge and specify the defining Boolean function type (1-14).
# Generate a true network topology and Boolean rules for 5 nodes set.seed(123) true_network <- GenerateNetwork(num.node = 5) # Graph the network with built-in exported function plot_network plot_network(true_network)# Generate a true network topology and Boolean rules for 5 nodes set.seed(123) true_network <- GenerateNetwork(num.node = 5) # Graph the network with built-in exported function plot_network plot_network(true_network)
Simulates a synthetic binary observation dataset () based on a given
directed acyclic graph topology with Boolean rules. It starts by running
independent Bernoulli trials on root nodes. The remaining non-root nodes
are calculated based on their assigned Boolean logic functions and parent states.
A pre-generated binary noise matrix is applied via a bitwise XOR operation
to occasionally flip the Boolean outputs, injecting natural biological
noise expected by the model.
GenerateSample( trans_matrix, SampleSize = 50, num.node = nrow(trans_matrix), para = rep(0.5, nrow(trans_matrix)), error = matrix(0, nrow = nrow(trans_matrix), ncol = SampleSize), timeseries = TRUE )GenerateSample( trans_matrix, SampleSize = 50, num.node = nrow(trans_matrix), para = rep(0.5, nrow(trans_matrix)), error = matrix(0, nrow = nrow(trans_matrix), ncol = SampleSize), timeseries = TRUE )
trans_matrix |
A square matrix combining the network topology |
SampleSize |
An integer representing the total number of samples or time points to simulate. Defaults to 50 if not specified. |
num.node |
An integer representing the total number of network nodes. Defaults to |
para |
A numeric vector of baseline success probabilities ( |
error |
A pre-generated binary noise matrix applied to occasionally flip Boolean outputs, injecting natural noise. Defaults to a zero matrix of size |
timeseries |
Logical. If TRUE, simulates time-lagged data where child state at time |
A simulated binary gene expression matrix , where rows represent individual genes/nodes and columns represent samples.
# 1. Generate a 5-node network set.seed(123) num_nodes <- 5 sample_size <- 10 true_network <- GenerateNetwork(num.node = num_nodes) # 2. Set baseline probabilities and simulate zero-noise error matrix root_probs <- rep(0.5, num_nodes) error_matrix <- matrix(0, nrow = num_nodes, ncol = sample_size) # 3. Generate cross-sectional (independent mode) synthetic data dummy_data <- GenerateSample( trans_matrix = true_network, SampleSize = sample_size, para = root_probs, error = error_matrix, timeseries = FALSE ) print(dummy_data)# 1. Generate a 5-node network set.seed(123) num_nodes <- 5 sample_size <- 10 true_network <- GenerateNetwork(num.node = num_nodes) # 2. Set baseline probabilities and simulate zero-noise error matrix root_probs <- rep(0.5, num_nodes) error_matrix <- matrix(0, nrow = num_nodes, ncol = sample_size) # 3. Generate cross-sectional (independent mode) synthetic data dummy_data <- GenerateSample( trans_matrix = true_network, SampleSize = sample_size, para = root_probs, error = error_matrix, timeseries = FALSE ) print(dummy_data)
Visualizes the causal network structure inferred by the BBNI MCMC sampler.
This function takes the marginal posterior probability of each directed edge
and plots the network using the igraph package. Edges with posterior
probabilities below the specified threshold are omitted from the plot.
plot_bbni( results, threshold = 0.5, node_names = NULL, true_network = NULL, ... )plot_bbni( results, threshold = 0.5, node_names = NULL, true_network = NULL, ... )
results |
The list returned by |
threshold |
Numeric. The minimum posterior probability required to draw an edge. Defaults to 0.5. |
node_names |
Character vector. Optional names for the nodes. Defaults to "N1", "N2", etc. |
true_network |
Optional square matrix representing the true network topology. If provided, edges will be color-coded to indicate true positives (along with displaying wrong function inferences), false positives, and false negatives. Purely for simulation purposes. |
... |
Additional graphical parameters passed to |
An invisible igraph object.
# 1. Generate synthetic network and time-series data set.seed(123) true_network <- GenerateNetwork(num.node = 5) dummy_data <- GenerateSample(true_network, SampleSize = 15) # 2. Run BBNI sampler prior_para <- matrix(3, nrow = 6, ncol = 2) prior_para[6, 1] <- 2 prior_para[6, 2] <- 100 results <- run_bbni(dummy_data, prior_para = prior_para, num_update = 100) # 3. Plot inferred network plot_bbni(results, true_network = true_network, threshold = 0.5)# 1. Generate synthetic network and time-series data set.seed(123) true_network <- GenerateNetwork(num.node = 5) dummy_data <- GenerateSample(true_network, SampleSize = 15) # 2. Run BBNI sampler prior_para <- matrix(3, nrow = 6, ncol = 2) prior_para[6, 1] <- 2 prior_para[6, 2] <- 100 results <- run_bbni(dummy_data, prior_para = prior_para, num_update = 100) # 3. Plot inferred network plot_bbni(results, true_network = true_network, threshold = 0.5)
Visualizes a single directed network topology. For instance, a synthetic "true" network
generated by GenerateNetwork().
plot_network(trans_matrix, node_names = NULL, ...)plot_network(trans_matrix, node_names = NULL, ...)
trans_matrix |
A square matrix representing the network topology and Boolean rules (e.g., the output of |
node_names |
Character vector. Optional names for the nodes. Defaults to "N1", "N2", etc. |
... |
Additional graphical parameters passed to |
An invisible igraph object.
# 1. Generate synthetic 5-node network topology set.seed(123) true_network <- GenerateNetwork(num.node = 5) # 2. Plot the generated network plot_network(true_network)# 1. Generate synthetic 5-node network topology set.seed(123) true_network <- GenerateNetwork(num.node = 5) # 2. Plot the generated network plot_network(true_network)
Generates a trace plot of the log-posterior values over iterations of the MCMC to visually assess convergence and stability of the executed Markov chain. Burn-in line is graphed to show when the data started to be utilized for edge-probability calculations.
plot_trace(results)plot_trace(results)
results |
The list returned by |
A base R trace plot.
# 1. Generate synthetic network and time-series data set.seed(123) true_network <- GenerateNetwork(num.node = 5) dummy_data <- GenerateSample(true_network, SampleSize = 15) # 2. Run BBNI sampler prior_para <- matrix(3, nrow = 6, ncol = 2) prior_para[6, 1] <- 2 prior_para[6, 2] <- 100 results <- run_bbni(dummy_data, prior_para = prior_para, num_update = 100) # 3. Visualize MCMC results plot_trace(results)# 1. Generate synthetic network and time-series data set.seed(123) true_network <- GenerateNetwork(num.node = 5) dummy_data <- GenerateSample(true_network, SampleSize = 15) # 2. Run BBNI sampler prior_para <- matrix(3, nrow = 6, ncol = 2) prior_para[6, 1] <- 2 prior_para[6, 2] <- 100 results <- run_bbni(dummy_data, prior_para = prior_para, num_update = 100) # 3. Visualize MCMC results plot_trace(results)
Executes a Metropolis-within-Gibbs Markov chain Monte Carlo (MCMC) algorithm to sample
from the joint posterior distribution of Directed Acyclic Graph (DAG) topologies () and
Boolean logic transition functions (). The algorithm iterates through individual
network nodes and proposes parent set mutations (edge additions, removals, or swaps)
paired with transition function reassignments to one of 14 candidate Boolean rules.
Proposed states transitions are strictly verified to follow the DAG constraint and
evaluated with a Metropolis-Hastings acceptance threshold using log-posterior values.
run_bbni( GeneData, num.node = nrow(GeneData), SampleSize = ncol(GeneData), prior_para = NULL, num_update = 4000, penalty = 0.1, prop.ratio = 0.5, verbose = FALSE, timeseries = TRUE, burn_in = 0.7 )run_bbni( GeneData, num.node = nrow(GeneData), SampleSize = ncol(GeneData), prior_para = NULL, num_update = 4000, penalty = 0.1, prop.ratio = 0.5, verbose = FALSE, timeseries = TRUE, burn_in = 0.7 )
GeneData |
A binary empirical observation matrix ( |
num.node |
An integer representing the total number of network nodes. Defaults to |
SampleSize |
An integer representing the total number of time points or independent samples in the dataset. Defaults to |
prior_para |
A matrix (with dimensions |
num_update |
An integer representing the total number of MCMC iterations to perform. Defaults to 4000 if not specified. |
penalty |
A numeric value representing the structural prior probability per edge used to penalize network complexity |
prop.ratio |
A numeric value between 0 and 1 representing the probability of choosing a uniform proposal distribution over an empirical proposal distribution at each iteration. Defaults to 0.5 |
verbose |
Logical. If TRUE, prints verbose MCMC iteration progress to the console. Default is FALSE. |
timeseries |
Logical. If TRUE, the algorithm assumes a time-series dataset. If FALSE, the algorithm assumes independent samples. Default is TRUE. |
burn_in |
A numeric value between 0 and 1 representing the proportion of initial MCMC samples to discard as burn-in. Defaults to 0.7 if not specified. |
A list containing the full trajectory of the MCMC chain. Specifically, networks (a list of sampled transition function matrices) and log_posterior (a numeric vector of log-posterior scores for each iteration). These represent samples drawn from the marginal posterior distribution used for Bayesian model averaging. Additionally, the post_edge_prob (matrix of marginal posterior edge probabilities) and burn_in ratio are returned in the list.
# 1. Define network parameters set.seed(235) num_nodes <- 8 sample_size <- 50 # 2. Generate true network and simulate data true_network <- GenerateNetwork(num.node = num_nodes) # Set up Beta priors for root-node probabilities and the noise rate prior_para <- matrix(3, nrow = num_nodes + 1, ncol = 2) prior_para[num_nodes + 1, 1] <- 2 prior_para[num_nodes + 1, 2] <- 100 # Simulate parameters para <- numeric(num_nodes + 1) for (i in 1:(num_nodes + 1)) { para[i] <- stats::rbeta(1, prior_para[i, 1], prior_para[i, 2]) } para[num_nodes + 1] <- 0.1 # Fixed noise rate for simulation error_matrix <- matrix(stats::rbinom(num_nodes * sample_size, 1, para[num_nodes + 1]), nrow = num_nodes, ncol = sample_size ) dummy_data <- GenerateSample( trans_matrix = true_network, SampleSize = sample_size, para = para, error = error_matrix ) # 3. Run the MCMC sampler (silently) mcmc_results <- run_bbni( GeneData = dummy_data, prior_para = prior_para, num_update = 80, # Scaled down for example speed prop.ratio = 0.1 ) # 4. Visualize results plot_bbni(mcmc_results, true_network = true_network, threshold = 0.5) plot_trace(mcmc_results)# 1. Define network parameters set.seed(235) num_nodes <- 8 sample_size <- 50 # 2. Generate true network and simulate data true_network <- GenerateNetwork(num.node = num_nodes) # Set up Beta priors for root-node probabilities and the noise rate prior_para <- matrix(3, nrow = num_nodes + 1, ncol = 2) prior_para[num_nodes + 1, 1] <- 2 prior_para[num_nodes + 1, 2] <- 100 # Simulate parameters para <- numeric(num_nodes + 1) for (i in 1:(num_nodes + 1)) { para[i] <- stats::rbeta(1, prior_para[i, 1], prior_para[i, 2]) } para[num_nodes + 1] <- 0.1 # Fixed noise rate for simulation error_matrix <- matrix(stats::rbinom(num_nodes * sample_size, 1, para[num_nodes + 1]), nrow = num_nodes, ncol = sample_size ) dummy_data <- GenerateSample( trans_matrix = true_network, SampleSize = sample_size, para = para, error = error_matrix ) # 3. Run the MCMC sampler (silently) mcmc_results <- run_bbni( GeneData = dummy_data, prior_para = prior_para, num_update = 80, # Scaled down for example speed prop.ratio = 0.1 ) # 4. Visualize results plot_bbni(mcmc_results, true_network = true_network, threshold = 0.5) plot_trace(mcmc_results)
A binary empirical dataset with the gene expression states of 14 yeast cell-cycle genes across various experimental conditions. The dataset is used in the original Han et al. (2014) paper to demonstrate the BBNI method for inferring Boolean gene regulatory networks.
data(yeast_data)data(yeast_data)
A numeric matrix with 14 rows (representing individual genes/nodes) and 385 columns (representing independent samples or sequential time points).
Han, S., Wong, R. K. W., Lee, T. C. M., Shen, L., Li, S.-Y. R., & Fan, X. (2014). A Full Bayesian Approach for Boolean Genetic Network Inference. PLOS ONE, 9(12), e115806. https://doi.org/10.1371/journal.pone.0115806