| Type: | Package |
| Title: | G-Wishart Normalising Constants for Gaussian Graphical Models |
| Version: | 1.0 |
| Date: | 2026-03-26 |
| Author: | Ching Wong [aut], Jack Kuipers [aut, cre] |
| Maintainer: | Jack Kuipers <jack.kuipers@bsse.ethz.ch> |
| Description: | Computes G-Wishart normalising constants through a Fourier approach. Either exact analytical results, numerical integration or Monte Carlo estimation are employed. Details at C. Wong, G. Moffa and J. Kuipers (2024), <doi:10.48550/arXiv.2404.06803>. |
| License: | GPL-2 | GPL-3 [expanded from: GPL (≥ 2)] |
| Depends: | R (≥ 4.0.0) |
| Imports: | igraph, BDgraph, CholWishart, MASS, mvtnorm, Rcpp (≥ 1.1.1) |
| LinkingTo: | Rcpp, RcppEigen |
| RoxygenNote: | 7.3.3 |
| Encoding: | UTF-8 |
| NeedsCompilation: | yes |
| Packaged: | 2026-04-10 09:22:10 UTC; jkuipers |
| Repository: | CRAN |
| Date/Publication: | 2026-04-16 09:42:03 UTC |
This function returns the log of the G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK from the transformed version log(C_G(delta, D)) = int_[S^p_++(G)] det(K)^((delta-2)/2) * exp(-tr(KD)/2) dK
Description
This function returns the log of the G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK from the transformed version log(C_G(delta, D)) = int_[S^p_++(G)] det(K)^((delta-2)/2) * exp(-tr(KD)/2) dK
Usage
C_GtoI_G(graph, CGval, delta)
Arguments
graph |
An igraph object, corresponding to a prime connected graph |
CGval |
A number log(C_G(delta, D)) |
delta |
A constant > 0 |
Value
A number log(I_G(beta, D))
Clique-completion
Description
This function returns the completion of a matrix with respect to a graph. The entries of the given matrix are such that in the clique decomposition of the chordal completion there is no linear term in the determinant expansion.
Usage
Clique_complete(
D,
missing_edges,
cliques,
cliques_ann,
prec_count = 0,
tol = 1e-12
)
Arguments
D |
A real symmetric positive definite p by p matrix |
missing_edges |
argument of edges to complete |
cliques |
list of cliques |
cliques_ann |
list of clique annotations |
prec_count |
Use Hessian after a this number of iterations (default is to always use) |
tol |
A tolerance to stop the numerical iterations |
Value
A real symmetric positive definite p by p matrix
This function is a wrapper for BDgraph to compare and to avoid the NOTE in the package checks since we only had BDgraph in the examples
Description
This function is a wrapper for BDgraph to compare and to avoid the NOTE in the package checks since we only had BDgraph in the examples
Usage
I_G_BD(graph, beta, D, n_samp = 1000, C_G_flag = FALSE)
Arguments
graph |
An igraph object, corresponding to a prime connected graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
n_samp |
number of sample points used in MC integration |
C_G_flag |
boolean whether to return C_G instead of I_G |
Value
A number log(I_G(beta, D)) [or log(C_G(delta, D))]
G Wishart normalising constant through MC integration
Description
This function returns the log of transformed G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK for connected graphs G. If there is no explicit formula, or the integral cannot be reduced into lower dimensions, MC integration is used.
Usage
I_G_MC(
graph,
beta,
D,
n_samp = 1000,
nu = 10,
err_flag = FALSE,
int1d = TRUE,
robust = FALSE,
chordal = NULL
)
Arguments
graph |
An igraph object, corresponding to a prime connected graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
n_samp |
number of sample points used in MC integration |
nu |
degree of freedom of the MC importance distribution (nu=0 is Gaussian) |
err_flag |
flag of whether to return the log relative error estimate (or not, default) |
int1d |
flag of whether to integrate numerically in 1d (default, MC integation otherwise) |
robust |
flag of whether to use robust means (default not to, but may offer numerical stability at the risk of bias) |
chordal |
Optional. An igraph object, corresponding to a chordal completion of graph |
Value
A number log(I_G(beta, D))
Examples
set.seed(42)
p <- 5
data <- matrix(runif(10*p), ncol = p)
D <- t(data) %*% data
beta <- 1
adj <- matrix(0, nrow = p, ncol = p)
adj[1,5] <- adj[1,2] <- adj[2,3] <- adj[3,4] <- adj[4,5] <- 1 # 5-cycle
graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max"))
I_G_MC(graph, beta, D)
# compare with BDgraph gnorm()
I_G_BD(graph, beta, D, n_samp = 1e5)
G Wishart normalising constant for chordal graphs
Description
This function returns the log of transformed G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK for chordal graphs G, using explicit formula (cliques / separators).
Usage
I_G_chordal(graph, beta, D, cliques = NULL, separators = NULL)
Arguments
graph |
An igraph object, corresponding to a chordal graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
cliques |
A list of cliques (if known) |
separators |
A list of separators (if known) |
Value
A number log(I_G(beta, D))
Examples
set.seed(42)
p <- 5
data <- matrix(runif(10*p), ncol = p)
D <- t(data) %*% data
beta <- 1
adj <- matrix(0, nrow = p, ncol = p)
adj[1,c(2:5)] <- adj[2,3] <- adj[3,4] <- adj[4,5] <- 1 # chordal completion of 5-cycle
graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max"))
I_G_chordal(graph, beta, D)
# compare with BDgraph gnorm()
I_G_BD(graph, beta, D, n_samp = 1e4)
G Wishart normalising constant for complete graphs
Description
This function returns the log of the transformed G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK for complete graphs G, using explicit formula det(D)^(-beta-(p+1)/2) Gamma_n(beta+(p+1)/2)
Usage
I_G_complete(p, beta, D)
Arguments
p |
A positive integer, indicating the number of vertices of the graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
Value
A number log(I_G(beta, D))
Examples
set.seed(42)
p <- 6
data <- matrix(runif(10*p), ncol = p)
D <- t(data) %*% data
beta <- 1
I_G_complete(6, beta, D)
# compare with BDgraph gnorm()
adj <- matrix(1, nrow = 6, ncol = 6)
graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max"))
I_G_BD(graph, beta, D)
G Wishart normalising constant
Description
This function returns the log of transformed G-Wishart normalising constant I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK for graphs G. Graphs are decomposed into connected prime components. For each component if there is no explicit formula, or the integral cannot be reduced into one dimension, MC integration is used.
Usage
I_Gnorm(
graph,
beta,
D,
connected = FALSE,
prime = FALSE,
chordal = NULL,
n_samp = 1000,
nu = 10,
robust = FALSE,
err_flag = FALSE
)
Arguments
graph |
An igraph object, corresponding to a prime connected graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
connected |
Whether the graph is connected, if known (default FALSE) |
prime |
Whether the graph is prime, if known (default FALSE) |
chordal |
Optional. An igraph object, corresponding to a chordal completion of graph |
n_samp |
number of sample points used in MC integration |
nu |
degree of freedom of the MC importance distribution (nu=0 is Gaussian) |
robust |
flag of whether to use robust means (default not to, but may offer numerical stability at the risk of bias) |
err_flag |
flag of whether to return the log relative error estimate (or not, default), only considered for prime connected graphs |
Value
A number log(I_G(beta, D))
Examples
set.seed(42)
p <- 11
data <- matrix(runif(10*p), ncol = p)
D <- t(data) %*% data
beta <- 1
graph <- igraph::make_graph(edges = c(1,2, 2,4, 2,5, 1,3, 3,5, # disconnected
5,6, 6,7, 5,7, 4,7, # with prime decomposition too
8,9, 8,11, 9,10, 10,11), n = 11, directed = FALSE)
I_Gnorm(graph, beta, D)
# compare with BDgraph gnorm()
# BDgraph gives an error when the graph is disconnected, so we do the components separately
C_G1 <- I_G_BD(igraph::subgraph(graph, 1:7), beta, D[1:7,1:7], n_samp = 1e5, C_G_flag = TRUE)
C_G2 <- I_G_BD(igraph::subgraph(graph, 8:11), beta, D[8:11,8:11], n_samp = 1e5, C_G_flag = TRUE)
C_GtoI_G(graph, C_G1 + C_G2, 2*beta+2)
This function returns the log of the G-Wishart normalising constant log(C_G(delta, D)) = int_[S^p_++(G)] det(K)^((delta-2)/2) * exp(-tr(KD)/2) dK from the transformed version I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK
Description
This function returns the log of the G-Wishart normalising constant log(C_G(delta, D)) = int_[S^p_++(G)] det(K)^((delta-2)/2) * exp(-tr(KD)/2) dK from the transformed version I_G(beta, D) = int_[S^p_++(G)] det(K)^beta * exp(-tr(KD)) dK
Usage
I_GtoC_G(graph, IGval, beta)
Arguments
graph |
An igraph object, corresponding to a prime connected graph |
IGval |
A number log(I_G(beta, D)) |
beta |
A constant > 0 |
Value
A number log(C_G(delta, D))
PD-completion
Description
This function returns the PD-completion of a matrix with respect to a graph. The entries of the given matrix corresponding to the missing edges are modified such that the inverse has zeros in the places of missing edges.
Usage
PD_complete(graph, D, missing_edges = NULL, tol = 1e-12)
Arguments
graph |
An igraph object, corresponding to the graph |
D |
A real symmetric positive definite p by p matrix |
missing_edges |
optional argument of edges to complete |
tol |
A tolerance to stop the numerical iterations |
Value
A real symmetric positive definite p by p matrix
Examples
set.seed(42)
p <- 5
data <- matrix(runif(10*p), ncol = p)
D <- t(data) %*% data
adj <- matrix(0, nrow = p, ncol = p)
adj[1,5] <- adj[1,2] <- adj[2,3] <- adj[3,4] <- adj[4,5] <- 1 # 5-cycle
graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max"))
(D_tilde <- PD_complete(graph, D))
(round(solve(D_tilde), 6)) # has zeros at the missing edges
this helper function adds information to cliques related to the missing edges
Description
this helper function adds information to cliques related to the missing edges
Usage
annotate_cliques(cliques, clique_flags, missing_edges, beta)
Arguments
cliques |
list of cliques |
clique_flags |
indicator of whether the elements are cliques (TRUE) or separators (FALSE) |
missing_edges |
array of missing edges |
beta |
A constant > 0 |
Value
annotated list of cliques
This function just checks if a graph is connected and prime
Description
This function just checks if a graph is connected and prime
Usage
check_prime_connected(graph, check_prime = TRUE)
Arguments
graph |
An igraph object |
check_prime |
Boolean flag to check primality too (default TRUE) |
Value
A boolean, TRUE if connected and prime
Chordal Factor # do not export
Description
Chordal Factor # do not export
Usage
chordal_factor(graph)
Arguments
graph |
An igraph object, a chordal graph |
Value
A list of maximal cliques and separators
Newton-Raphson update for the clique-completion
Description
Newton-Raphson update for the clique-completion
Usage
clique_update_D(D, tau, cliques, cliques_ann, prec_flag = TRUE)
Arguments
D |
A real symmetric positive definite p by p matrix |
tau |
number of missing edges |
cliques |
list of cliques |
cliques_ann |
list of clique annotations |
prec_flag |
whether to use the full Hessian |
Value
the update step
Compute means and gradients for the Newton-Raphson update for the clique-completion
Description
Compute means and gradients for the Newton-Raphson update for the clique-completion
Usage
local_mean_grad(D, tau, cliques, cliques_ann, prec_flag = TRUE)
Arguments
D |
A real symmetric positive definite p by p matrix |
tau |
number of missing edges |
cliques |
list of cliques |
cliques_ann |
list of clique annotations |
prec_flag |
compute the full Hessian |
Value
the mean and gradient for the update step note this does not take into account the covariance/precision for numerical speed
Compute second moment of the log expansion of the determinant (assuming the mean is 0 from Clique_completion)
Description
Compute second moment of the log expansion of the determinant (assuming the mean is 0 from Clique_completion)
Usage
local_precision(D, tau, cliques, cliques_ann)
Arguments
D |
A real symmetric positive definite p by p matrix |
tau |
number of missing edges |
cliques |
list of cliques |
cliques_ann |
list of clique annotations |
Value
the precision matrix
Predict Row # do not export
Description
Predict Row # do not export
Usage
predict_row(i, j, D)
Arguments
i |
A number, the row to predict |
j |
A vector, the columns to predicted |
D |
A real symmetric positive definite matrix |
Value
the predicted values for the vector j
Prime Decomposition
Description
This function returns the prime components and the minimal separators of a connected graph.
Usage
prime_decomp(graph)
Arguments
graph |
An igraph object, a connected graph |
Value
A list of the prime components and the separators
Examples
graph <- igraph::make_graph(edges = c(1,2, 2,4, 2,5, 1,3, 3,5,
5,6, 6,7, 5,7, 4,7), n = 7, directed = FALSE)
prime_decomp(graph)