Package {grasps}

Type:	Package
Title:	Groupwise Regularized Adaptive Sparse Precision Solution
Version:	0.1.3
Maintainer:	Shiying Xiao <shiying.xiao@outlook.com>
Description:	Provides a unified framework for sparse-group precision matrix estimation in Gaussian graphical models using combined element-wise and group-wise regularization. It supports multiple sparse-group penalties, including lasso, adaptive lasso, and several non-convex penalties such as SCAD and MCP. Non-convex estimates can be computed using one-step or iterative LLA, with weighted subproblems solved efficiently using ADMM-based optimization. The package is designed for high-dimensional network inference where both sparsity and group structure are present.
License:	GPL (≥ 3)
URL:	https://github.com/Carol-seven/grasps, https://shiying-xiao.com/grasps/
BugReports:	https://github.com/Carol-seven/grasps/issues
Encoding:	UTF-8
Imports:	igraph, ggforce, ggplot2, grDevices, Rcpp, Rdpack, scales
LinkingTo:	Rcpp, RcppArmadillo
RdMacros:	Rdpack
Suggests:	knitr, MASS, quarto, rmarkdown
VignetteBuilder:	knitr, quarto
Config/roxygen2/version:	8.1.1
NeedsCompilation:	yes
Packaged:	2026-11-13 20:05:55 UTC; csxiao
Author:	Shiying Xiao [aut, cre], Jun Yan [aut], Panpan Zhang [aut]
Repository:	CRAN
Date/Publication:	2026-11-14 00:10:13 UTC

Penalty Derivative Computation

Description

Compute one or more derivative values for a given omega, allowing vectorized specifications of penalty, lambda, and gamma.

Usage

compute_derivative(omega, penalty, lambda, gamma = NA)

Arguments

Value

A data frame with S3 class "derivative" containing:

omega

The input omega values.

penalty

The penalty type for each row.

lambda

The regularization parameter used.

gamma

The additional penalty parameter used.

value

The computed derivative value.

The number of rows equals length(omega) multiplied by max(length(penalty), length(lambda), length(gamma)). Each of penalty, lambda, and gamma must have length 1 or the same common length. Arguments of length 1 are recycled to this common length.

References

Candès EJ, Wakin MB, Boyd SP (2008). “Enhancing Sparsity by Reweighted \ell_1 Minimization.” Journal of Fourier Analysis and Applications, 14(5), 877–905. doi:10.1007/s00041-008-9045-x.

Fan J, Feng Y, Wu Y (2009). “Network Exploration via the Adaptive LASSO and SCAD Penalties.” The Annals of Applied Statistics, 3(2), 521–541. doi:10.1214/08-aoas215.

Fan J, Li R (2001). “Variable Selection via Nonconcave Penalized Likelihood and its Oracle Properties.” Journal of the American Statistical Association, 96(456), 1348–1360. doi:10.1198/016214501753382273.

Frank LE, Friedman JH (1993). “A Statistical View of Some Chemometrics Regression Tools.” Technometrics, 35(2), 109–135. doi:10.1080/00401706.1993.10485033.

Friedman J, Hastie T, Tibshirani R (2008). “Sparse Inverse Covariance Estimation with the Graphical Lasso.” Biostatistics, 9(3), 432–441. doi:10.1093/biostatistics/kxm045.

Fu WJ (1998). “Penalized Regressions: The Bridge versus the Lasso.” Journal of Computational and Graphical Statistics, 7(3), 397–416. doi:10.1080/10618600.1998.10474784.

Tibshirani R (1996). “Regression Shrinkage and Selection via the Lasso.” Journal of the Royal Statistical Society: Series B (Methodological), 58(1), 267–288. doi:10.1111/j.2517-6161.1996.tb02080.x.

Wang Y, Fan Q, Zhu L (2018). “Variable Selection and Estimation using a Continuous Approximation to the L_0 Penalty.” Annals of the Institute of Statistical Mathematics, 70(1), 191–214. doi:10.1007/s10463-016-0588-3.

Wang Y, Zhu L (2016). “Variable Selection and Parameter Estimation with the Atan Regularization Method.” Journal of Probability and Statistics, 2016, 6495417. doi:10.1155/2016/6495417.

Zhang C (2010). “Nearly Unbiased Variable Selection under Minimax Concave Penalty.” The Annals of Statistics, 38(2), 894–942. doi:10.1214/09-AOS729.

Examples

library(grasps)
library(ggplot2)

deriv_df <- compute_derivative(
  omega = seq(0, 4, by = 0.01),
  penalty = c("atan", "exp", "lasso", "lq", "lsp", "mcp", "scad"),
  lambda = 1)

plot(deriv_df) +
  scale_y_continuous(limits = c(0, 1.5)) +
  guides(color = guide_legend(nrow = 2, byrow = TRUE))

Penalty Function Computation

Description

Compute one or more penalty values for a given omega, allowing vectorized specifications of penalty, lambda, and gamma.

Usage

compute_penalty(omega, penalty, lambda, gamma = NA)

Arguments

Value

A data frame with S3 class "penalty" containing:

omega

The input omega values.

penalty

The penalty type for each row.

lambda

The regularization parameter used.

gamma

The additional penalty parameter used.

value

The computed penalty value.

The number of rows equals length(omega) multiplied by max(length(penalty), length(lambda), length(gamma)). Each of penalty, lambda, and gamma must have length 1 or the same common length. Arguments of length 1 are recycled to this common length.

References

Candès EJ, Wakin MB, Boyd SP (2008). “Enhancing Sparsity by Reweighted \ell_1 Minimization.” Journal of Fourier Analysis and Applications, 14(5), 877–905. doi:10.1007/s00041-008-9045-x.

Fan J, Feng Y, Wu Y (2009). “Network Exploration via the Adaptive LASSO and SCAD Penalties.” The Annals of Applied Statistics, 3(2), 521–541. doi:10.1214/08-aoas215.

Fan J, Li R (2001). “Variable Selection via Nonconcave Penalized Likelihood and its Oracle Properties.” Journal of the American Statistical Association, 96(456), 1348–1360. doi:10.1198/016214501753382273.

Frank LE, Friedman JH (1993). “A Statistical View of Some Chemometrics Regression Tools.” Technometrics, 35(2), 109–135. doi:10.1080/00401706.1993.10485033.

Friedman J, Hastie T, Tibshirani R (2008). “Sparse Inverse Covariance Estimation with the Graphical Lasso.” Biostatistics, 9(3), 432–441. doi:10.1093/biostatistics/kxm045.

Fu WJ (1998). “Penalized Regressions: The Bridge versus the Lasso.” Journal of Computational and Graphical Statistics, 7(3), 397–416. doi:10.1080/10618600.1998.10474784.

Tibshirani R (1996). “Regression Shrinkage and Selection via the Lasso.” Journal of the Royal Statistical Society: Series B (Methodological), 58(1), 267–288. doi:10.1111/j.2517-6161.1996.tb02080.x.

Wang Y, Fan Q, Zhu L (2018). “Variable Selection and Estimation using a Continuous Approximation to the L_0 Penalty.” Annals of the Institute of Statistical Mathematics, 70(1), 191–214. doi:10.1007/s10463-016-0588-3.

Wang Y, Zhu L (2016). “Variable Selection and Parameter Estimation with the Atan Regularization Method.” Journal of Probability and Statistics, 2016, 6495417. doi:10.1155/2016/6495417.

Zhang C (2010). “Nearly Unbiased Variable Selection under Minimax Concave Penalty.” The Annals of Statistics, 38(2), 894–942. doi:10.1214/09-AOS729.

Examples

library(grasps)
library(ggplot2)

pen_df <- compute_penalty(
  omega = seq(-4, 4, by = 0.01),
  penalty = c("atan", "exp", "lasso", "lq", "lsp", "mcp", "scad"),
  lambda = 1)

plot(pen_df, xlim = c(-1, 1), ylim = c(0, 1), zoom.size = 1) +
  guides(color = guide_legend(nrow = 2, byrow = TRUE))

Block-Structured Precision Matrix based on SBM

Description

Generate a precision matrix that exhibits block structure induced by a stochastic block model (SBM).

Usage

gen_prec_sbm(
  p,
  block.sizes = NULL,
  K = 3,
  prob.mat = NULL,
  within.prob = 0.25,
  between.prob = 0.05,
  weight.mat = NULL,
  weight.dists = list("gamma", "unif"),
  weight.paras = list(c(shape = 100, rate = 10), c(min = 0, max = 5)),
  cond.target = 100
)

Arguments

Details

Edge sampling. Within- and between-group edges are sampled independently according to Bernoulli distributions specified by prob.mat, or by within.prob and between.prob if prob.mat is not supplied.

Weight sampling. Conditional on the adjacency structure, edge weights are sampled block-wise from samplers specified in weight.dists and weight.paras. The length of weight.dists (and weight.paras) determines how weight distributions are assigned:

• length = 1: Same specification for all blocks.

• length = 2: first for within-group blocks, second for between-group blocks.

• length = K + K(K-1)/2: Full specification for each block.

Block indexing. The order for blocks is:

• Within-group blocks: Indices 1, \dots, K.

• Between-group blocks: K(K-1)/2 blocks in order (1,2), (1,3), (1,4), ..., (1,K), (2,3), ..., (K-1,K).

Positive definiteness. Motivated by condition-number-regularized covariance estimation (Won et al. 2013), the weighted adjacency matrix is symmetrized and used as the precision matrix \Omega_0. Since arbitrary block-structured weights may not be positive definite, a diagonal adjustment is applied to control the eigenvalue spectrum. Specifically, let \lambda_{\max} and \lambda_{\min} denote the largest and smallest eigenvalues of a matrix. A non-negative numeric value \tau is added to the diagonal so that

\left\{ \begin{array}{l} \dfrac{\lambda_{\max}(\Omega_0 + \tau I)}{\lambda_{\min}(\Omega_0 + \tau I)} \leq \texttt{cond.target} \\[1em] \lambda_{\min}(\Omega_0 + \tau I) > 0 \\[.5em] \tau \geq 0 \end{array} \right.

which ensures both positive definiteness and guarantees that the condition number does not exceed cond.target, providing numerical stability even in high-dimensional settings.

Value

An object with S3 class "gen_prec_sbm" containing the following components:

Omega

The precision matrix with SBM block structure.

Sigma

The covariance matrix, i.e., the inverse of Omega.

sparsity

Proportion of absent edges among all possible undirected edges in Omega.

membership

An integer vector specifying the group membership.

References

Won J, Lim J, Kim S, Rajaratnam B (2013). “Condition-Number-Regularized Covariance Estimation.” Journal of the Royal Statistical Society Series B: Statistical Methodology, 75(3), 427–450. doi:10.1111/j.1467-9868.2012.01049.x.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
#### case 1: base R distribution
sim1 <- gen_prec_sbm(p = 100, K = 5,
                     within.prob = 0.25, between.prob = 0.1,
                     weight.dists = list("gamma", "unif"),
                     weight.paras = list(c(shape = 100, scale = 1e2),
                                         c(min = 0, max = 10)),
                     cond.target = 100)
#### visualization
plot(sim1)

#### case 2: user-defined sampler
my_gamma <- function(n) {
  rgamma(n, shape = 1e4, scale = 1e2)
}
sim2 <- gen_prec_sbm(p = 100, K = 5,
                     within.prob = 0.2, between.prob = 0.05,
                     weight.dists = list(my_gamma, "unif"),
                     weight.paras = list(NULL,
                                         c(min = 0, max = 1)),
                     cond.target = 100)
#### visualization
plot(sim2)

Groupwise Regularized Adaptive Sparse Precision Solution

Description

Provide a collection of statistical methods that incorporate both element-wise and group-wise penalties to estimate a precision matrix.

Usage

grasps(
  X,
  n = nrow(X),
  membership,
  penalty,
  diag.ind = TRUE,
  diag.grp = TRUE,
  diag.include = FALSE,
  lambda = NULL,
  alpha = NULL,
  gamma = NULL,
  initial = "lasso",
  nlambda = 10,
  lambda.min.ratio = 0.01,
  growiter.lambda = 30,
  tol.lambda = 0.001,
  maxiter.lambda = 50,
  lla = FALSE,
  tol.lla = 1e-04,
  maxiter.lla = 10,
  rho = 2,
  tau.incr = 2,
  tau.decr = 2,
  nu = 10,
  tol.abs = 1e-04,
  tol.rel = 1e-04,
  maxiter = 10000,
  crit = "BIC",
  kfold = 5,
  ebic.tuning = 0.5
)

Arguments

Details

For non-convex penalties, the default lla = FALSE computes a one-step LLA estimator. Penalty weights are constructed from the estimator specified by initial and held fixed while the resulting weighted subproblem is solved.

When lla = TRUE, the penalty weights are updated after each weighted subproblem until the relative change between successive precision matrix estimates is no greater than tol.lla, or maxiter.lla is reached. This iterative procedure is intended to approach a stationary point of the original non-convex objective but does not guarantee a global optimum.

For penalty = "adapt", weights are constructed from the estimator specified by initial and remain fixed; therefore, lla is ignored. The converged component concerns convergence of the fixed-weight ADMM subproblems and, when initial = "lasso", the initial lasso fit. It does not indicate stationarity of the original non-convex objective. The lla.converged component concerns the outer iterative LLA procedure.

Value

An object with S3 class "grasps" containing the following components:

hatOmega

The estimated precision matrix.

lambda

The selected regularization parameter.

alpha

The selected mixing parameter.

initial

The initial estimate used to construct penalty weights when penalty != "lasso".

gamma

The additional parameter used for the selected penalty when penalty != "lasso".

iterations

For penalty = "lasso", the number of ADMM iterations. For other penalties, the total number of ADMM iterations across all LLA subproblems, excluding any ADMM iterations used to construct the initial estimate.

converged

A logical value indicating whether all ADMM subproblems and, when initial = "lasso", the initial lasso fit satisfied their convergence criteria before reaching maxiter.

initial.converged

For penalty != "lasso", a logical value indicating whether the initial lasso fit converged. It is NA when initial = "invS" or a matrix.

lla.iterations

The number of LLA steps performed for non-convex penalties. It is 1 for one-step LLA and NA for "lasso" and "adapt".

lla.converged

For iterative LLA, a logical value indicating whether the stopping criterion was satisfied before reaching maxiter.lla. It is NA for one-step LLA, penalty = "lasso", and penalty = "adapt".

lambda.safe

The conservative upper bound used to initialize the line search for \lambda_{\max}, corresponding to each distinct value in alpha.grid, when lambda = NULL.

loss

The minimum mean k-fold loss when crit = "CV".

CV.loss

Matrix of CV losses, with rows for parameter combinations and columns for CV folds, when crit = "CV".

score

The optimal information criterion score when crit is set to "AIC", "BIC", "EBIC", or "HBIC".

IC.score

The information criterion score for each parameter combination when crit is set to "AIC", "BIC", "EBIC", or "HBIC".

converged.grid

The convergence status of the candidate fits. For information criterion selection, this is a logical vector corresponding to lambda.grid and alpha.grid. For cross-validation, this is a logical matrix with parameter combinations in rows and folds in columns. All candidates are retained for model selection regardless of convergence status. Non-finite scores and estimates that are not positive definite are assigned infinite scores. Selection stops with an error if no candidate has a finite score. A warning is issued if the selected fit or its selected CV fits did not satisfy the ADMM convergence criteria.

lla.converged.grid

The iterative LLA convergence status of the candidate fits. For information criterion selection, this is a logical vector corresponding to lambda.grid and alpha.grid. For cross-validation, this is a logical matrix with parameter combinations in rows and folds in columns. Entries are NA when iterative LLA is not used. A warning is issued if the selected fit or its selected CV fits did not satisfy the iterative LLA convergence criterion.

lambda.grid

The actual lambda grid used in the program.

alpha.grid

The actual alpha grid used in the program.

membership

The group membership.

References

Akaike H (1973). “Information Theory and an Extension of the Maximum Likelihood Principle.” In Petrov BN, Csáki F (eds.), Second International Symposium on Information Theory, 267–281. Akadémiai Kiadó, Budapest, Hungary.

Candès EJ, Wakin MB, Boyd SP (2008). “Enhancing Sparsity by Reweighted \ell_1 Minimization.” Journal of Fourier Analysis and Applications, 14(5), 877–905. doi:10.1007/s00041-008-9045-x.

Chen J, Chen Z (2008). “Extended Bayesian Information Criteria for Model Selection with Large Model Spaces.” Biometrika, 95(3), 759–771. doi:10.1093/biomet/asn034.

Fan J, Feng Y, Wu Y (2009). “Network Exploration via the Adaptive LASSO and SCAD Penalties.” The Annals of Applied Statistics, 3(2), 521–541. doi:10.1214/08-aoas215.

Fan J, Li R (2001). “Variable Selection via Nonconcave Penalized Likelihood and its Oracle Properties.” Journal of the American Statistical Association, 96(456), 1348–1360. doi:10.1198/016214501753382273.

Fan J, Liu H, Ning Y, Zou H (2017). “High Dimensional Semiparametric Latent Graphical Model for Mixed Data.” Journal of the Royal Statistical Society Series B: Statistical Methodology, 79(2), 405–421. doi:10.1111/rssb.12168.

Foygel R, Drton M (2010). “Extended Bayesian Information Criteria for Gaussian Graphical Models.” In Lafferty J, Williams C, Shawe-Taylor J, Zemel R, Culotta A (eds.), Advances in Neural Information Processing Systems 23 (NIPS 2010), 604–612. https://dl.acm.org/doi/10.5555/2997189.2997257.

Frank LE, Friedman JH (1993). “A Statistical View of Some Chemometrics Regression Tools.” Technometrics, 35(2), 109–135. doi:10.1080/00401706.1993.10485033.

Friedman J, Hastie T, Tibshirani R (2008). “Sparse Inverse Covariance Estimation with the Graphical Lasso.” Biostatistics, 9(3), 432–441. doi:10.1093/biostatistics/kxm045.

Fu WJ (1998). “Penalized Regressions: The Bridge versus the Lasso.” Journal of Computational and Graphical Statistics, 7(3), 397–416. doi:10.1080/10618600.1998.10474784.

Schwarz G (1978). “Estimating the Dimension of a Model.” The Annals of Statistics, 6(2), 461–464. doi:10.1214/aos/1176344136.

Tibshirani R (1996). “Regression Shrinkage and Selection via the Lasso.” Journal of the Royal Statistical Society: Series B (Methodological), 58(1), 267–288. doi:10.1111/j.2517-6161.1996.tb02080.x.

Wang L, Kim Y, Li R (2013). “Calibrating Nonconvex Penalized Regression in Ultra-High Dimension.” The Annals of Statistics, 41(5), 2505–2536. doi:10.1214/13-AOS1159.

Wang Y, Fan Q, Zhu L (2018). “Variable Selection and Estimation using a Continuous Approximation to the L_0 Penalty.” Annals of the Institute of Statistical Mathematics, 70(1), 191–214. doi:10.1007/s10463-016-0588-3.

Wang Y, Zhu L (2016). “Variable Selection and Parameter Estimation with the Atan Regularization Method.” Journal of Probability and Statistics, 2016, 6495417. doi:10.1155/2016/6495417.

Zhang C (2010). “Nearly Unbiased Variable Selection under Minimax Concave Penalty.” The Annals of Statistics, 38(2), 894–942. doi:10.1214/09-AOS729.

Zou H (2006). “The Adaptive Lasso and Its Oracle Properties.” Journal of the American Statistical Association, 101(476), 1418–1429. doi:10.1198/016214506000000735.

Zou H, Li R (2008). “One-Step Sparse Estimates in Nonconcave Penalized Likelihood Models.” The Annals of statistics, 36(4), 1509–1533. doi:10.1214/009053607000000802.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
sim <- gen_prec_sbm(p = 30, K = 3,
                    within.prob = 0.25, between.prob = 0.05,
                    weight.dists = list("gamma", "unif"),
                    weight.paras = list(c(shape = 20, rate = 10),
                                        c(min = 0, max = 5)),
                    cond.target = 100)
## ground truth visualization
plot(sim)

## n-by-p data matrix
library(MASS)
X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)

## precision matrix: adaptive lasso; BIC
prec <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "BIC")

## precision matrix visualization
plot(prec)

## performance
performance(hatOmega = prec$hatOmega, Omega = sim$Omega)

## adjacency matrix: diagonal = 0; raw partial correlations;
##                   no thresholding; weighted network
adj <- prec_to_adj(prec$hatOmega,
                   diag.zero = TRUE, absolute = FALSE,
                   threshold = NULL, weighted = TRUE)

## adjacency matrix visualization
plot(adj)

Performance Measures for Precision Matrix Estimation

Description

Compute a collection of loss-based and structure-based measures to evaluate the performance of an estimated precision matrix.

Usage

performance(hatOmega, Omega)

Arguments

Details

Let \Omega_{p \times p} and \hat{\Omega}_{p \times p} be the reference (true) and estimated precision matrices, respectively. When \Omega is positive definite, \Sigma = \Omega^{-1} is the corresponding covariance matrix. Edges are defined by nonzero off-diagonal entries in the upper triangle of the precision matrices.

The input matrices are not required to be positive definite. If either matrix is not positive definite, the function issues a warning but continues computing the performance measures. In this case, covariance- and likelihood-based measures may not have their usual statistical interpretation. If \Omega is singular, the KL and quadratic losses are returned as NA, while the remaining measures are still computed. For non-positive-definite but invertible matrices, the reported KL value uses the log absolute determinant and need not be non-negative. A classification-based measure is returned as NA when its denominator is zero.

Sparsity is treated as a structural summary, while the remaining measures are grouped into loss-based measures, raw confusion-matrix counts, and classification-based (structure-recovery) measures.

"sparsity": Sparsity is computed as the proportion of zero entries among the off-diagonal elements in the upper triangle of \hat{\Omega}.

Loss-based measures:

• "Frobenius": Frobenius (Hilbert-Schmidt) norm loss = \Vert \Omega - \hat{\Omega} \Vert_F.

• "KL": Kullback-Leibler divergence = \mathrm{tr}(\Sigma \hat{\Omega}) - \log\det(\Sigma \hat{\Omega}) - p.

• "quadratic": Quadratic norm loss = \Vert \Sigma \hat{\Omega} - I_p \Vert_F^2.

• "spectral": Spectral (operator) norm loss = \Vert \Omega - \hat{\Omega} \Vert_{2,2} = e_1, where e_1^2 is the largest eigenvalue of (\Omega - \hat{\Omega})^2.

Confusion-matrix counts:

• "TP": True positive = number of edges in both \Omega and \hat{\Omega}.

• "TN": True negative = number of edges in neither \Omega nor \hat{\Omega}.

• "FP": False positive = number of edges in \hat{\Omega} but not in \Omega.

• "FN": False negative = number of edges in \Omega but not in \hat{\Omega}.

Classification-based (structure-recovery) measures:

• "TPR": True positive rate (TPR), recall, sensitivity = \mathrm{TP} / (\mathrm{TP} + \mathrm{FN}).

• "FPR": False positive rate (FPR) = \mathrm{FP} / (\mathrm{FP} + \mathrm{TN}).

• "F1": F_1 score = 2\,\mathrm{TP} / (2\,\mathrm{TP} + \mathrm{FN} + \mathrm{FP})

• "MCC": Matthews correlation coefficient (MCC) = (\mathrm{TP}\times\mathrm{TN} - \mathrm{FP}\times\mathrm{FN}) / \sqrt{(\mathrm{TP}+\mathrm{FP})(\mathrm{TP}+\mathrm{FN}) (\mathrm{TN}+\mathrm{FP})(\mathrm{TN}+\mathrm{FN})}

The following table summarizes the confusion matrix and related rates:

Value

A data frame of S3 class "performance", with one row per performance metric and two columns:

measure

The name of each performance metric. The reported metrics include: sparsity, Frobenius norm loss, Kullback-Leibler divergence, quadratic norm loss, spectral norm loss, true positive, true negative, false positive, false negative, true positive rate, false positive rate, F1 score, and Matthews correlation coefficient.

value

The corresponding numeric value. Undefined measures are returned as NA.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
sim <- gen_prec_sbm(p = 30, K = 3,
                    within.prob = 0.25, between.prob = 0.05,
                    weight.dists = list("gamma", "unif"),
                    weight.paras = list(c(shape = 20, rate = 10),
                                        c(min = 0, max = 5)),
                    cond.target = 100)
## ground truth visualization
plot(sim)

## n-by-p data matrix
library(MASS)
X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)

## precision matrix: adaptive lasso; BIC
prec <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "BIC")

## precision matrix visualization
plot(prec)

## performance
performance(hatOmega = prec$hatOmega, Omega = sim$Omega)

## adjacency matrix: diagonal = 0; raw partial correlations;
##                   no thresholding; weighted network
adj <- prec_to_adj(prec$hatOmega,
                   diag.zero = TRUE, absolute = FALSE,
                   threshold = NULL, weighted = TRUE)

## adjacency matrix visualization
plot(adj)

Plot Function for S3 Class "adjmat"

Description

Visualize an adjacency matrix as a heatmap. This function is shared by objects returned from prec_to_adj.

Usage

## S3 method for class 'adjmat'
plot(x, ...)

Arguments

Value

A heatmap of class ggplot showing the matrix entries. The plot title also reports matrix dimension and sparsity.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
sim <- gen_prec_sbm(p = 30, K = 3,
                    within.prob = 0.25, between.prob = 0.05,
                    weight.dists = list("gamma", "unif"),
                    weight.paras = list(c(shape = 20, rate = 10),
                                        c(min = 0, max = 5)),
                    cond.target = 100)
## ground truth visualization
plot(sim)

## n-by-p data matrix
library(MASS)
X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)

## precision matrix: adaptive lasso; BIC
prec <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "BIC")

## precision matrix visualization
plot(prec)

## performance
performance(hatOmega = prec$hatOmega, Omega = sim$Omega)

## adjacency matrix: diagonal = 0; raw partial correlations;
##                   no thresholding; weighted network
adj <- prec_to_adj(prec$hatOmega,
                   diag.zero = TRUE, absolute = FALSE,
                   threshold = NULL, weighted = TRUE)

## adjacency matrix visualization
plot(adj)

Plot Function for Block-Structured Precision Matrices (Visualize a Matrix with Group Boundaries)

Description

Visualize a precision matrix as a heatmap with dashed boundary lines separating group blocks. This function is shared by objects returned from grasps, gen_prec_sbm, and sparsify_block_banded, all of which inherit from the S3 class "blkmat".

Usage

## S3 method for class 'blkmat'
plot(x, colors = NULL, ...)

Arguments

Value

A heatmap of class ggplot showing the matrix entries. Dashed lines indicate group boundaries. The plot title also reports matrix dimension and sparsity.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
sim <- gen_prec_sbm(p = 100, K = 10,
                    within.prob = 0.2, between.prob = 0.05,
                    weight.dists = list("gamma", "unif"),
                    weight.paras = list(c(shape = 100, scale = 10),
                                        c(min = 0, max = 5)),
                    cond.target = 100)

## visualization
plot(sim)

Plot Function for S3 Class "penderiv"

Description

Generate a visualization of penalty functions produced by compute_penalty, or penalty derivatives produced by compute_derivative. The plot automatically summarizes multiple configurations of penalty type, \lambda, and \gamma. Optional zooming is supported through facet_zoom.

Usage

## S3 method for class 'penderiv'
plot(x, ...)

Arguments

Value

An object of class ggplot.

Examples

library(grasps)
library(ggplot2)

pen_df <- compute_penalty(
  omega = seq(-4, 4, by = 0.01),
  penalty = c("atan", "exp", "lasso", "lq", "lsp", "mcp", "scad"),
  lambda = 1)

plot(pen_df, xlim = c(-1, 1), ylim = c(0, 1), zoom.size = 1) +
  guides(color = guide_legend(nrow = 2, byrow = TRUE))

deriv_df <- compute_derivative(
  omega = seq(0, 4, by = 0.01),
  penalty = c("atan", "exp", "lasso", "lq", "lsp", "mcp", "scad"),
  lambda = 1)

plot(deriv_df) +
  scale_y_continuous(limits = c(0, 1.5)) +
  guides(color = guide_legend(nrow = 2, byrow = TRUE))

Adjacency Matrix from Precision Matrix

Description

Convert a precision matrix to a partial-correlation-based adjacency matrix.

Usage

prec_to_adj(
  prec.mat,
  diag.zero = TRUE,
  absolute = FALSE,
  threshold = NULL,
  weighted = TRUE
)

Arguments

Details

For a precision matrix \Omega, the partial correlation between nodes i and j is computed as

\rho_{ij} = - \frac{\Omega_{ij}}{\sqrt{\Omega_{ii}\Omega_{jj}}}.

The input matrix is not required to be positive definite. If prec.mat is not positive definite, the function issues a warning but continues the conversion. In that case, the result may not represent valid partial correlations.

Value

A numeric adjacency matrix with S3 class "adjmat".

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## block-structured precision matrix based on SBM
sim <- gen_prec_sbm(p = 30, K = 3,
                    within.prob = 0.25, between.prob = 0.05,
                    weight.dists = list("gamma", "unif"),
                    weight.paras = list(c(shape = 20, rate = 10),
                                        c(min = 0, max = 5)),
                    cond.target = 100)
## ground truth visualization
plot(sim)

## n-by-p data matrix
library(MASS)
X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)

## precision matrix: adaptive lasso; BIC
prec <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "BIC")

## precision matrix visualization
plot(prec)

## performance
performance(hatOmega = prec$hatOmega, Omega = sim$Omega)

## adjacency matrix: diagonal = 0; raw partial correlations;
##                   no thresholding; weighted network
adj <- prec_to_adj(prec$hatOmega,
                   diag.zero = TRUE, absolute = FALSE,
                   threshold = NULL, weighted = TRUE)

## adjacency matrix visualization
plot(adj)

Groupwise Block-Banded Sparsifier

Description

Apply a groupwise block-banded mask to a symmetric matrix, keeping only entries within specified group neighborhoods.

Usage

sparsify_block_banded(mat, membership, neighbor.range = 1)

Arguments

Details

Entry (i, j) is retained when the absolute difference between the group labels of variables i and j does not exceed neighbor.range.

Applying the block-banded mask does not necessarily preserve positive definiteness. If the masked matrix is not positive definite, the function issues a warning and returns Sigma = NULL.

Value

An object with S3 class "sparsify_block_banded" containing the following components:

Omega

The block-banded matrix.

Sigma

The inverse of Omega if Omega is positive definite; otherwise, NULL.

sparsity

Proportion of absent edges among all possible undirected edges in Omega.

membership

An integer vector specifying the group membership.

Examples

library(grasps)

## reproducibility for everything
set.seed(1234)

## precision matrix estimation
X <- matrix(rnorm(200), 10, 20)
membership <- rep(1:5, each = 4)
est <- grasps(X, membership = membership, penalty = "lasso",
              lambda = 0.1, alpha = 0.5)

## default: keep blocks within ±1 of each group
res1 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 1)
plot(res1)

## wider band: keep blocks within ±2 of each group
res2 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 2)
plot(res2)

## special case: block-diagonal matrix
res3 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 0)
plot(res3)