| Title: | G-Wishart Normalising Constants for Gaussian Graphical Models |
|---|---|
| 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>. Also includes approximations of the ratio of normalising constants, see details at C. Wong, G. Moffa and J. Kuipers (2025), <doi:10.48550/arXiv.2503.13046>. |
| Authors: | Ching Wong [aut], Jack Kuipers [aut, cre] |
| Maintainer: | Jack Kuipers <[email protected]> |
| License: | GPL (>= 2) |
| Version: | 1.0.1 |
| Built: | 2026-05-27 11:44:49 UTC |
| Source: | https://github.com/cran/GWnorm |
this helper function adds information to cliques related to the missing edges
annotate_cliques(cliques, clique_flags, missing_edges, beta)annotate_cliques(cliques, clique_flags, missing_edges, beta)
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 |
annotated list of cliques
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
C_GtoI_G(graph, CGval, delta)C_GtoI_G(graph, CGval, delta)
graph |
An igraph object, corresponding to a prime connected graph |
CGval |
A number log(C_G(delta, D)) |
delta |
A constant > 0 |
A number log(I_G(beta, D))
This function just checks if a graph is connected and prime
check_prime_connected(graph, check_prime = TRUE)check_prime_connected(graph, check_prime = TRUE)
graph |
An igraph object |
check_prime |
Boolean flag to check primality too (default TRUE) |
A boolean, TRUE if connected and prime
Chordal Factor # do not export
chordal_factor(graph)chordal_factor(graph)
graph |
An igraph object, a chordal graph |
A list of maximal cliques and separators
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.
Clique_complete( D, missing_edges, cliques, cliques_ann, prec_count = 0, tol = 1e-12 )Clique_complete( D, missing_edges, cliques, cliques_ann, prec_count = 0, tol = 1e-12 )
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 |
A real symmetric positive definite p by p matrix
Newton-Raphson update for the clique-completion
clique_update_D(D, tau, cliques, cliques_ann, prec_flag = TRUE)clique_update_D(D, tau, cliques, cliques_ann, prec_flag = TRUE)
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 |
the update step
This function checks whether there is a triangle in the graph(graph2) contains at least two edges from the missing edges (graph1).
form_triangle(graph1, graph2)form_triangle(graph1, graph2)
graph1 |
An igraph object, corresponding to the missing edges |
graph2 |
An igraph object, corresponding to the super graph |
TRUE if there is such a triangle
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
I_G_BD(graph, beta, D, n_samp = 1000, C_G_flag = FALSE)I_G_BD(graph, beta, D, n_samp = 1000, C_G_flag = FALSE)
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 |
A number log(I_G(beta, D)) [or log(C_G(delta, D))]
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).
I_G_chordal( graph, beta, D = NULL, cliques = NULL, separators = NULL, ratio = FALSE )I_G_chordal( graph, beta, D = NULL, cliques = NULL, separators = NULL, ratio = FALSE )
graph |
An igraph object, corresponding to a chordal graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix (default, NULL for identity matrix) |
cliques |
A list of cliques (if known) |
separators |
A list of separators (if known) |
ratio |
Whether to output just the ratio compared to the value with I (default FALSE) |
A number log(I_G(beta, D))
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)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)
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)
I_G_complete(p, beta, D = NULL, ratio = FALSE)I_G_complete(p, beta, D = NULL, ratio = FALSE)
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 (default, NULL, if for identity matrix) |
ratio |
Whether to output just the ratio compared to the value with I (default FALSE) |
A number log(I_G(beta, D))
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)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)
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.
I_G_MC( graph, beta, D, n_samp = 1000, nu = 10, err_flag = FALSE, int1d = TRUE, robust = FALSE, quench_k = 0, chordal = NULL )I_G_MC( graph, beta, D, n_samp = 1000, nu = 10, err_flag = FALSE, int1d = TRUE, robust = FALSE, quench_k = 0, chordal = NULL )
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 integration otherwise) |
robust |
flag of whether to use robust means (default not to, but may offer numerical stability at the risk of bias) |
quench_k |
scaling factor of von Mises quenching (default not used, value of 0, but value around 1 may offer numerical stability at the risk of bias) |
chordal |
Optional. An igraph object, corresponding to a chordal completion of graph |
A number log(I_G(beta, D))
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)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)
This function returns the approximation of the ratio of log transformed G-Wishart normalising constants I_G(beta, D) / I_G(beta, I) for graphs G. Graphs are decomposed into connected prime components. For each component if there is no explicit formula, the approximation is used. Note that this is the same as the ratio C_G(delta, D) / C_G(delta, I) with delta = 2*beta + 2
I_G_ratio_approx( graph, beta, D, connected = FALSE, prime = FALSE, chordal = NULL, use_comp = NULL )I_G_ratio_approx( graph, beta, D, connected = FALSE, prime = FALSE, chordal = NULL, use_comp = NULL )
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 |
use_comp |
Whether to approximate using the graph (FALSE), the complement (TRUE) or the more efficient (default NULL) |
A number approximating log(I_G(beta, D)) - log(I_G(beta, I))
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_G_ratio_approx(graph, beta, D) # compare with MC estimates I_G_MC(graph, beta, D, n_samp = 1e5) - I_G_MC(graph, beta, diag(p), n_samp = 1e5)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_G_ratio_approx(graph, beta, D) # compare with MC estimates I_G_MC(graph, beta, D, n_samp = 1e5) - I_G_MC(graph, beta, diag(p), n_samp = 1e5)
This function returns the approximation of the ratio of log transformed G-Wishart normalising constants I_G(beta, D) / I_G(beta, I) for connected prime graphs G. If there is no explicit formula, the approximation is used. Note that this is the same as the ratio C_G(delta, D) / C_G(delta, I) with delta = 2*beta + 2
I_G_ratio_approx_prime(graph, beta, D, chordal = NULL, use_comp = NULL)I_G_ratio_approx_prime(graph, beta, D, chordal = NULL, use_comp = NULL)
graph |
An igraph object, corresponding to a prime connected graph |
beta |
A constant > 0 |
D |
A p by p real symmetric positive definite matrix |
chordal |
Optional. An igraph object, corresponding to a chordal completion of graph |
use_comp |
Whether to approximate using the graph (FALSE), the complement (TRUE) or the more efficient (default NULL) |
A number approximating log(I_G(beta, D)) - log(I_G(beta, I))
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_ratio_approx_prime(graph, beta, D, use_comp = TRUE) # compare with MC estimates I_G_MC(graph, beta, D, n_samp = 1e5) - I_G_MC(graph, beta, diag(p), n_samp = 1e5)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_ratio_approx_prime(graph, beta, D, use_comp = TRUE) # compare with MC estimates I_G_MC(graph, beta, D, n_samp = 1e5) - I_G_MC(graph, beta, diag(p), n_samp = 1e5)
This function returns the log of transformed G-Wishart normalising constant I_G(beta, I) = int_[S^p_++(G)] det(K)^beta * exp(-tr(K)) dK for connected graphs G for special cases where there is an explicit formula, including involving a low-dimensional integral.
I_G_special( graph, beta, connected = FALSE, prime = FALSE, test_6_cycle = TRUE, test_k_partite = TRUE, chordal = NULL )I_G_special( graph, beta, connected = FALSE, prime = FALSE, test_6_cycle = TRUE, test_k_partite = TRUE, chordal = NULL )
graph |
An igraph object, corresponding to a prime connected graph |
beta |
A constant > 0 |
connected |
Whether the graph is connected, if known (default FALSE) |
prime |
Whether the graph is prime, if known (default FALSE) |
test_6_cycle |
Flag whether to test whether the graph is a 6-cycle or its complement (default TRUE) |
test_k_partite |
Flag whether to test whether the graph is complete k partite (default TRUE) |
chordal |
Optional. An igraph object, corresponding to a chordal completion of graph |
A number log(I_G(beta, I))
beta <- 1 p <- 5 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_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) p <- 6 adj <- matrix(0, nrow = p, ncol = p) adj[1,6] <- adj[1,2] <- adj[2,3] <- adj[3,4] <- adj[4,5] <- adj[5,6] <- 1 # 6-cycle graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max")) I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) graph <- igraph::complementer(graph) # complement of 6-cycle I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) adj[upper.tri(adj)] <- 1 adj[2,3] <- adj[1,4] <- adj[5,6] <- 0 # complete graph with some edges missing graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max")) I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5)beta <- 1 p <- 5 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_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) p <- 6 adj <- matrix(0, nrow = p, ncol = p) adj[1,6] <- adj[1,2] <- adj[2,3] <- adj[3,4] <- adj[4,5] <- adj[5,6] <- 1 # 6-cycle graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max")) I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) graph <- igraph::complementer(graph) # complement of 6-cycle I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5) adj[upper.tri(adj)] <- 1 adj[2,3] <- adj[1,4] <- adj[5,6] <- 0 # complete graph with some edges missing graph <- igraph::graph_from_adjacency_matrix(adj, mode = c("max")) I_G_special(graph, beta) # compare with BDgraph gnorm() I_G_BD(graph, beta, diag(p), n_samp = 1e5)
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.
I_Gnorm( graph, beta, D, connected = FALSE, prime = FALSE, chordal = NULL, n_samp = 1000, nu = 10, robust = FALSE, quench_k = 0, err_flag = FALSE, useBDgraph = FALSE )I_Gnorm( graph, beta, D, connected = FALSE, prime = FALSE, chordal = NULL, n_samp = 1000, nu = 10, robust = FALSE, quench_k = 0, err_flag = FALSE, useBDgraph = FALSE )
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) |
quench_k |
scaling factor of von Mises quenching (default value of 0 means not used, but value around 1 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 |
useBDgraph |
flag of whether to use BDgraph instead (or not, default), useful wrapper for comparisons |
A number log(I_G(beta, D))
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 I_Gnorm(graph, beta, D, n_samp = 1e5, useBDgraph = TRUE)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 I_Gnorm(graph, beta, D, n_samp = 1e5, useBDgraph = TRUE)
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
I_GtoC_G(graph, IGval, beta)I_GtoC_G(graph, IGval, beta)
graph |
An igraph object, corresponding to a prime connected graph |
IGval |
A number log(I_G(beta, D)) |
beta |
A constant > 0 |
A number log(C_G(delta, D))
Input an igraph object, this function tells you whether graph is the cycle of length 6, or its complement.
is_6_cycle(graph)is_6_cycle(graph)
graph |
An igraph object |
An integer, 1 means the graph is the cycle of length 6, 2 means the the graph is the complement of the cycle of length 6, 0 means otherwise
# C6 example 5 1 3 2 4 6 is_6_cycle(igraph::make_graph(edges = c(1,3, 1,5, 2,3, 2,4, 5,6, 6,4), n = 6, directed = FALSE)) # C6 complement example is_6_cycle(igraph::make_graph(edges = c(1,2, 1,4, 1,6, 2,5, 2,6, 3,4, 3,5, 3,6, 4,5), n = 6, directed = FALSE)) # not C6 examples is_6_cycle(igraph::make_graph(edges = c(1,3, 1,5, 2,3, 2,4, 5,6, 6,4), n = 7, directed = FALSE)) is_6_cycle(igraph::make_graph(edges = c(1,4, 1,5, 2,3, 2,4, 5,6, 6,4), n = 6, directed = FALSE))# C6 example 5 1 3 2 4 6 is_6_cycle(igraph::make_graph(edges = c(1,3, 1,5, 2,3, 2,4, 5,6, 6,4), n = 6, directed = FALSE)) # C6 complement example is_6_cycle(igraph::make_graph(edges = c(1,2, 1,4, 1,6, 2,5, 2,6, 3,4, 3,5, 3,6, 4,5), n = 6, directed = FALSE)) # not C6 examples is_6_cycle(igraph::make_graph(edges = c(1,3, 1,5, 2,3, 2,4, 5,6, 6,4), n = 7, directed = FALSE)) is_6_cycle(igraph::make_graph(edges = c(1,4, 1,5, 2,3, 2,4, 5,6, 6,4), n = 6, directed = FALSE))
Input an igraph object, this function tells you the size of the partition if the graph is a k partite graph.
is_k_partite(graph)is_k_partite(graph)
graph |
An igraph object |
A vector of the size of the partition. The length of the vector is k. Empty vector means that the graph is not a complete k partite graph.
# not complete k-partite example is_k_partite(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)) # complete 4-partite example, (1,5), (2,4), (3,7), (6) is_k_partite(igraph::make_graph(edges = c(1,2, 1,3, 1,4, 1,6, 1,7, 2,3, 2,5, 2,6, 2,7, 3,4, 3,5, 3,6, 4,5, 4,6, 4,7, 5,6, 5,7, 6,7), n = 7, directed = FALSE))# not complete k-partite example is_k_partite(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)) # complete 4-partite example, (1,5), (2,4), (3,7), (6) is_k_partite(igraph::make_graph(edges = c(1,2, 1,3, 1,4, 1,6, 1,7, 2,3, 2,5, 2,6, 2,7, 3,4, 3,5, 3,6, 4,5, 4,6, 4,7, 5,6, 5,7, 6,7), n = 7, directed = FALSE))
This function returns the Isserlis matrix (rows/columns arranged in some order) with respect to a matrix and a graph for the missing edges (complement of the graph).
Iss_cmat(M, graph)Iss_cmat(M, graph)
M |
An p by p matrix |
graph |
An igraph object, corresponding to the graph |
An m by m matrix, where m is the number of missing edges in the graph
# example code adj <- matrix(0, nrow = 5, ncol = 5) 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")) Iss_cmat(matrix(1:25, nrow = 5, byrow = TRUE), graph)# example code adj <- matrix(0, nrow = 5, ncol = 5) 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")) Iss_cmat(matrix(1:25, nrow = 5, byrow = TRUE), graph)
This function returns the Isserlis matrix (rows/columns arranged in some order) with respect to a matrix and a graph
Iss_mat(M, graph)Iss_mat(M, graph)
M |
An p by p matrix |
graph |
An igraph object, corresponding to the graph |
An m by m matrix, where m is the number of vertices and edges in the graph
# example code adj <- matrix(0, nrow = 5, ncol = 5) 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")) Iss_mat(matrix(1:25, nrow = 5, byrow = TRUE), graph)# example code adj <- matrix(0, nrow = 5, ncol = 5) 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")) Iss_mat(matrix(1:25, nrow = 5, byrow = TRUE), graph)
Compute means and gradients for the Newton-Raphson update for the clique-completion
local_mean_grad(D, tau, cliques, cliques_ann, prec_flag = TRUE)local_mean_grad(D, tau, cliques, cliques_ann, prec_flag = TRUE)
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 |
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)
local_precision(D, tau, cliques, cliques_ann)local_precision(D, tau, cliques, cliques_ann)
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 |
the precision matrix
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.
PD_complete(graph, D, missing_edges = NULL, tol = 1e-12)PD_complete(graph, D, missing_edges = NULL, tol = 1e-12)
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 |
A real symmetric positive definite p by p matrix
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 edgesset.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
Predict Row # do not export
predict_row(i, j, D)predict_row(i, j, D)
i |
A number, the row to predict |
j |
A vector, the columns to predicted |
D |
A real symmetric positive definite matrix |
the predicted values for the vector j
This function returns the prime components and the minimal separators of a connected graph.
prime_decomp(graph)prime_decomp(graph)
graph |
An igraph object, a connected graph |
A list of the prime components and the separators
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)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)