Movatterモバイル変換

Type:

Package

Title:

Groupwise Regularized Adaptive Sparse Precision Solution

Version:

0.1.0

Maintainer:

Shiying Xiao <shiying.xiao@outlook.com>

Description:

Provides a unified framework for sparse-group regularization and precision matrix estimation in Gaussian graphical models. It implements multiple sparse-group penalties, including sparse-group lasso, sparse-group adaptive lasso, sparse-group SCAD, and sparse-group MCP, and solves them 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

RoxygenNote:

7.3.3

RdMacros:

Rdpack

Suggests:

knitr, MASS, quarto, rmarkdown

VignetteBuilder:

knitr, quarto

NeedsCompilation:

yes

Packaged:

2025-11-24 04:09:47 UTC; csxiao

Author:

Shiying Xiao

[aut, cre]

Repository:

CRAN

Date/Publication:

2025-11-27 19:10:14 UTC

Penalty Derivative Computation

Description

Compute one or more derivative values for a givenomega, allowingvectorized specifications ofpenalty,lambda, andgamma.

Usage

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

Arguments

omega

A numeric value or vector at which the penalty is evaluated.

penalty

A character string or vector specifying one or more penaltytypes. Available options include:

"lasso": Least absolute shrinkage and selection operator(Tibshirani 1996; Friedman et al. 2008).
"atan": Arctangent type penalty(Wang and Zhu 2016).
"exp": Exponential type penalty(Wang et al. 2018).
"lq": Lq penalty(Frank and Friedman 1993; Fu 1998; Fan and Li 2001).
"lsp": Log-sum penalty(Candès et al. 2008).
"mcp": Minimax concave penalty(Zhang 2010).
"scad": Smoothly clipped absolute deviation(Fan and Li 2001; Fan et al. 2009).

Ifpenalty has length 1, it is recycled to the common lengthdetermined bypenalty,lambda, andgamma.

lambda

A non-negative numeric value or vector specifyingthe regularization parameter.Iflambda has length 1, it is recycled to the common lengthdetermined bypenalty,lambda, andgamma.

gamma

A numeric value or vector specifying the additional parameterfor the penalty function.Iflambda has length 1, it is recycled to the common lengthdetermined bypenalty,lambda, andgamma.The penalty-specific defaults are:

"atan": 0.005
"exp": 0.01
"lq": 0.5
"lsp": 0.1
"mcp": 3
"scad": 3.7

For"lasso",gamma is ignored.

Value

A data frame with S3 class"derivative" containing:

omega: The inputomega 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 equalsmax(length(penalty), length(lambda), length(gamma)).Any ofpenalty,lambda, orgamma with length 1is 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 theL_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 givenomega, allowingvectorized specifications ofpenalty,lambda, andgamma.

Usage

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

Arguments

omega

A numeric value or vector at which the penalty is evaluated.

penalty

A character string or vector specifying one or more penaltytypes. Available options include:

"lasso": Least absolute shrinkage and selection operator(Tibshirani 1996; Friedman et al. 2008).
"atan": Arctangent type penalty(Wang and Zhu 2016).
"exp": Exponential type penalty(Wang et al. 2018).
"lq": Lq penalty(Frank and Friedman 1993; Fu 1998; Fan and Li 2001).
"lsp": Log-sum penalty(Candès et al. 2008).
"mcp": Minimax concave penalty(Zhang 2010).
"scad": Smoothly clipped absolute deviation(Fan and Li 2001; Fan et al. 2009).

Ifpenalty has length 1, it is recycled to the common lengthdetermined bypenalty,lambda, andgamma.

lambda

A non-negative numeric value or vector specifyingthe regularization parameter.Iflambda has length 1, it is recycled to the common lengthdetermined bypenalty,lambda, andgamma.

gamma

"atan": 0.005
"exp": 0.01
"lq": 0.5
"lsp": 0.1
"mcp": 3
"scad": 3.7

For"lasso",gamma is ignored.

Value

A data frame with S3 class"penalty" containing:

omega: The inputomega 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 equalsmax(length(penalty), length(lambda), length(gamma)).Any ofpenalty,lambda, orgamma with length 1is recycled to this common length.

References

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 bya stochastic block model (SBM).

Usage

gen_prec_sbm(  d,  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

d

An integer specifying the number of variables (dimensions).

block.sizes

An integer vector (default =NULL) specifyingthe size of each group. IfNULL, thed variables are dividedas evenly as possible acrossK groups.

K

An integer (default = 3) specifying the number of groups.Ignored ifblock.sizes is provided; thenK <- length(block.sizes).

prob.mat

AK \times K symmetric matrix (default =NULL)specifying the Bernoulli rates. Element(i, j) gives the probability ofcreating an edge between vertices from groupsi andj.IfNULL, a matrix withwithin.prob on the diagonal andbetween.prob on the off-diagonal is used.

within.prob

A numeric value in [0, 1] (default = 0.25) specifyingthe probability of creating an edge between vertices within the same group.This argument is used only whenprob.mat = NULL.

between.prob

A numeric value in [0, 1] (default = 0.05) specifyingthe probability of creating an edge between vertices from different groups.This argument is used only whenprob.mat = NULL.

weight.mat

Ad \times d symmetric matrix (default =NULL)specifying the edge weights. IfNULL, weights are generated block-wiseaccording toweight.dists andweight.paras.

weight.dists

A list (default =list("gamma", "unif"))specifying the sampling distribution for each block of weights.Its length determines how the distributions are assigned:

length = 1: Same specification for all blocks.
length = 2: First for within-group blocks, second for between-groupblocks.
length =K + K(K-1)/2: Full specification for each block.The firstK elements correspond to within-group blocks with indices1, \dots, K, and the remainingK(K-1)/2 elements correspond tobetween-group blocks ordered as(1,2),(1,3),(1,4), ...,(1,K),(2,3), ...,(K-1,K).

Each element ofweight.dists can be:

A user-supplied sampling function. The function must accept an argumentn specifying the number of samples.
A character string specifying the distribution family.Accepted distributions (base R samplers in parentheses) include:
- "beta": Beta distribution (rbeta)
- "cauchy": Cauchy distribution (rcauchy).
- "chisq": Chi-squared distribution (rchisq).
- "exp": Exponential distribution (rexp).
- "f": F distribution (rf).
- "gamma": Gamma distribution (rgamma).
- "lnorm": Log normal distribution (rlnorm).
- "norm": Normal distribution (rnorm).
- "t": Student's t distribution (rt).
- "unif": Uniform distribution (runif).
- "weibull": Weibull distribution (rweibull).

weight.paras

A list(default =list(c(shape = 100, rate = 10), c(min = 0, max = 5)))specifying the parameters associated withweight.dists. It must followthe same length rules asweight.dists. Each element should be a namedvector or list suitable for the corresponding sampler.

cond.target

A numeric value > 1 (default = 100) specifying the targetcondition number for the precision matrix. When necessary, a diagonal shiftis applied to ensure positive definiteness and numerical stability.

Details

Edge sampling.Within- and between-group edges are sampled independently according toBernoulli distributions specified byprob.mat, or bywithin.probandbetween.prob ifprob.mat is not supplied.

Weight sampling.Conditional on the adjacency structure, edge weights are sampled block-wisefrom samplers specified inweight.dists andweight.paras.The length ofweight.dists (andweight.paras) determines howweight distributions are assigned:

length = 1: Same specification for all blocks.
length = 2: first for within-group blocks, second for between-groupblocks.
length =K + K(K-1)/2: Full specification for each block.

Block indexing.The order for blocks is:

Within-group blocks: Indices1, \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.The weighted adjacency matrix is symmetrized and used as the precision matrix\Omega_0. Since arbitrary block-structured weights may not be positivedefinite, a diagonal adjustment is applied to control the eigenvalue spectrum.Specifically, let\lambda_{\max} and\lambda_{\min} denotethe largest and smallest eigenvalues of a matrix. A non-negative numericvalue\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 conditionnumber does not exceedcond.target, providing numerical stabilityeven in high-dimensional settings.

Value

An object with S3 class"gen_prec_sbm" containing the followingcomponents:

Omega: The precision matrix with SBM block structure.
Sigma: The covariance matrix, i.e., the inverse ofOmega.
sparsity: Proportion of zero entries inOmega.
membership: An integer vector specifying the group membership.

Examples

library(grasps)## reproducibility for everythingset.seed(1234)## block-structured precision matrix based on SBM#### case 1: base R distributionsim1 <- gen_prec_sbm(d = 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)#### visualizationplot(sim1)#### case 2: user-defined samplermy_gamma <- function(n) {  rgamma(n, shape = 1e4, scale = 1e2)}sim2 <- gen_prec_sbm(d = 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)#### visualizationplot(sim2)

Groupwise Regularized Adaptive Sparse Precision Solution

Description

Provide a collection of statistical methods that incorporate bothelement-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,  nlambda = 10,  lambda.min.ratio = 0.01,  growiter.lambda = 30,  tol.lambda = 0.001,  maxiter.lambda = 50,  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

X

Ann \times d data matrix with sample sizen anddimensiond.
Ad \times d sample covariance matrix with dimensiond.

n

An integer (default =nrow(X)) specifying the sample size.This is only required when the input matrixX is ad \times dsample covariance matrix with dimensiond.

membership

An integer vector specifying the group membership.The length ofmembership must be consistent with the dimensiond.

penalty

A character string specifying the penalty for estimatingprecision matrix. Available options include:

"lasso": Least absolute shrinkage and selection operator(Tibshirani 1996; Friedman et al. 2008).
"adapt": Adaptive lasso(Zou 2006; Fan et al. 2009).
"atan": Arctangent type penalty(Wang and Zhu 2016).
"exp": Exponential type penalty(Wang et al. 2018).
"lq": Lq penalty(Frank and Friedman 1993; Fu 1998; Fan and Li 2001).
"lsp": Log-sum penalty(Candès et al. 2008).
"mcp": Minimax concave penalty(Zhang 2010).
"scad": Smoothly clipped absolute deviation(Fan and Li 2001; Fan et al. 2009).

diag.ind

A logical value (default = TRUE) specifying whether topenalize the diagonal elements.

diag.grp

A logical value (default = TRUE) specifying whether topenalize the within-group blocks.

diag.include

A logical value (default = FALSE) specifying whether toinclude the diagonal entries in the penalty for within-group blocks whendiag.grp = TRUE.

lambda

A non-negative numeric vector specifying the grid forthe regularization parameter. The default isNULL, which generatesits ownlambda sequence based onnlambda andlambda.min.ratio.

alpha

A numeric vector in [0, 1] specifying the grid forthe mixing parameter balancing the element-wise individual L1 penalty andthe block-wise group L2 penalty.An alpha of 1 corresponds to the individual penalty only; an alpha of 0corresponds to the group penalty only.The default value is a sequence from 0.1 to 0.9 with increments of 0.1.

gamma

A numeric value specifying the additional parameter fothe chosenpenalty. The default value depends on the penalty:

"adapt": 0.5
"atan": 0.005
"exp": 0.01
"lq": 0.5
"lsp": 0.1
"mcp": 3
"scad": 3.7

nlambda

An integer (default = 10) specifying the number oflambda values to generate whenlambda = NULL.

lambda.min.ratio

A numeric value > 0 (default = 0.01) specifyingthe fraction of the maximumlambda value\lambda_{max} togenerate the minimumlambda\lambda_{min}.Iflambda = NULL, alambda grid of lengthnlambda isautomatically generated on a log scale, ranging from\lambda_{max}down to\lambda_{min}.

growiter.lambda

An integer (default = 30) specifying the maximumnumber of exponential growth steps during the initial search for anadmissible upper bound\lambda_{\max}.

tol.lambda

A numeric value > 0 (default = 1e-03) specifyingthe relative tolerance for the bisection stopping rule on the interval width.

maxiter.lambda

An integer (default = 50) specifying the maximum numberof bisection iterations in the line search for\lambda_{\max}.

rho

A numeric value > 0 (default = 2) specifying the ADMMaugmented-Lagrangian penalty parameter (often called the ADMM step size).Larger values typically put more weight on enforcing the consensusconstraints at each iteration; smaller values yield more conservative updates.

tau.incr

A numeric value > 1 (default = 2) specifyingthe multiplicative factor used to increaserho when the primalresidual dominates the dual residual in ADMM.

tau.decr

A numeric value > 1 (default = 2) specifyingthe multiplicative factor used to decreaserho when the dual residualdominates the primal residual in ADMM.

nu

A numeric value > 1 (default = 10) controlling how aggressivelyrho is rescaled in the adaptive-rho scheme (residual balancing).

tol.abs

A numeric value > 0 (default = 1e-04) specifying the absolutetolerance for ADMM stopping (applied to primal/dual residual norms).

tol.rel

A numeric value > 0 (default = 1e-04) specifying the relativetolerance for ADMM stopping (applied to primal/dual residual norms).

maxiter

An integer (default = 1e+04) specifying the maximum number ofADMM iterations.

crit

A character string (default = "BIC") specifying the parameterselection criterion to use. Available options include:

"AIC": Akaike information criterion(Akaike 1973).
"BIC": Bayesian information criterion(Schwarz 1978).
"EBIC": extended Bayesian information criterion(Chen and Chen 2008; Foygel and Drton 2010).
"HBIC": high dimensional Bayesian information criterion(Wang et al. 2013; Fan et al. 2017).
"CV": k-fold cross validation with negative log-likelihood loss.

kfold

An integer (default = 5) specifying the number of folds used forcrit = "CV".

ebic.tuning

A numeric value in [0, 1] (default = 0.5) specifyingthe tuning parameter to calculate forcrit = "EBIC".

Value

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

hatOmega: The estimated precision matrix.
lambda: The optimal regularization parameter.
alpha: The optimal mixing parameter.
initial: The initial estimate ofhatOmega when a non-convexpenalty is chosen viapenalty.
gamma: The optimal addtional parameter when a non-convex penaltyis chosen viapenalty.
iterations: The number of ADMM iterations.
lambda.grid: The actual lambda grid used in the program.
alpha.grid: The actual alpha grid used in the program.
lambda.safe: The bisection-refined upper bound\lambda_{\max},corresponding toalpha.grid, whenlambda = NULL.
loss: The optimal k-fold loss whencrit = "CV".
CV.loss: Matrix of CV losses, with rows for parameter combinations andcolumns for CV folds, whencrit = "CV".
score: The optimal information criterion score whencrit is setto"AIC","BIC","EBIC", or"HBIC".
IC.score: The information criterion score for each parametercombination whencrit is set to"AIC","BIC","EBIC", or"HBIC".
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\'emiai Kiad\'o, 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 theL_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.

Examples

library(grasps)## reproducibility for everythingset.seed(1234)## block-structured precision matrix based on SBMsim <- gen_prec_sbm(d = 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)## visualizationplot(sim)## n-by-d data matrixlibrary(MASS)X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)## adapt, HBICres <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "HBIC")## visualizationplot(res)## performanceperformance(hatOmega = res$hatOmega, Omega = sim$Omega)

Performance Measures for Precision Matrix Estimation

Description

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

Usage

performance(hatOmega, Omega)

Arguments

hatOmega

A numericd \times d matrix giving the estimatedprecision matrix.

Omega

A numericd \times d matrix giving the reference(typically true) precision matrix.

Details

Let\Omega_{d \times d} and\hat{\Omega}_{d \times d} bethe reference (true) and estimated precision matrices, respectively, with\Sigma = \Omega^{-1} being the corresponding covariance matrix.Edges are defined by nonzero off-diagonal entries in the upper triangle ofthe precision matrices.

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

"sparsity": Sparsity is computed as the proportion of zero entries amongthe 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}) - d.
"quadratic": Quadratic norm loss= \Vert \Sigma \hat{\Omega} - I_d \Vert_F^2.
"spectral": Spectral (operator) norm loss= \Vert \Omega - \hat{\Omega} \Vert_{2,2} = e_1,wheree_1^2 is the largest eigenvalue of(\Omega - \hat{\Omega})^2.

Confusion-matrix counts:

"TP": True positive= number of edgesin both\Omega and\hat{\Omega}.
"TN": True negative= number of edgesin neither\Omega nor\hat{\Omega}.
"FP": False positive= number of edgesin\hat{\Omega} but not in\Omega.
"FN": False negative= number of edgesin\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:

	Predicted Positive	Predicted Negative
Real Positive (P)	True positive (TP)	False negative (FN)	True positive rate (TPR), recall, sensitivity = TP / P = 1 - FNR	False negative rate (FNR) = FN / P = 1 - TPR
Real Negative (N)	False positive (FP)	True negative (TN)	False positive rate (FPR) = FP / N = 1 - TNR	True negative rate (TNR), specificity = TN / N = 1 - FPR
	Positive predictive value (PPV), precision = TP / (TP + FP) = 1 - FDR	False omission rate (FOR) = FN / (TN + FN) = 1 - NPV
	False discovery rate (FDR) = FP / (TP + FP) = 1 - PPV	Negative predictive value (NPV) = TN / (TN + FN) = 1 - FOR

Value

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

measure: The name of each performance metric. The reported metricsinclude: 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.

Examples

library(grasps)## reproducibility for everythingset.seed(1234)## block-structured precision matrix based on SBMsim <- gen_prec_sbm(d = 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)## visualizationplot(sim)## n-by-d data matrixlibrary(MASS)X <- mvrnorm(n = 20, mu = rep(0, 30), Sigma = sim$Sigma)## adapt, BICres <- grasps(X = X, membership = sim$membership, penalty = "adapt", crit = "BIC")## visualizationplot(res)## performanceperformance(hatOmega = res$hatOmega, Omega = sim$Omega)

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

Description

Visualize a precision matrix as a heatmap with dashed boundary linesseparating group blocks. This function is shared by objects returned fromgrasps,gen_prec_sbm, andsparsify_block_banded, all of which inherit fromthe S3 class"blkmat".

Usage

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

Arguments

x

An object inheriting from S3 class"blkmat", typicallyreturned bygrasps,gen_prec_sbmorsparsify_block_banded.

colors

A vector of colors specifying an n-color gradient scale forthe fill aesthetics.

...

Additional arguments passed toggplot.

Value

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

Examples

library(grasps)## reproducibility for everythingset.seed(1234)## block-structured precision matrix based on SBMsim <- gen_prec_sbm(d = 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)## visualizationplot(sim)

Plot Function for S3 Class "penderiv"

Description

Generate a visualization of penalty functions produced bycompute_penalty, or penalty derivatives produced bycompute_derivative.The plot automatically summarizes multiple configurations of penalty type,\lambda, and\gamma. Optional zooming is supported throughfacet_zoom.

Usage

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

Arguments

x

An object inheriting from S3 class"penderiv", typicallyreturned bycompute_penalty, orcompute_derivative.

...

Optional arguments passed tofacet_zoomto zoom in on a subset of the data, while keeping the view of the fulldataset as a separate panel.

Value

An object of classggplot.

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))

Groupwise Block-Banded Sparsifier

Description

Make a precision-like matrix block-banded according to group membership,keeping only entries within specified group neighborhoods.

Usage

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

Arguments

mat

Ad \times d precision-like matrix specifying the basematrix to be masked.

membership

An integer vector specifying the group membership.The length ofmembership must be consistent with the dimensiond.

neighbor.range

An integer (default = 1) specifying the neighbor range,where groups whose labels differ by at mostneighbor.range areconsidered neighbors and kept in the mask.

Value

An object with S3 class"sparsify_block_banded" containingthe following components:

Omega: The masked precision matrix.
Sigma: The covariance matrix, i.e., the inverse ofOmega.
sparsity: Proportion of zero entries inOmega.
membership: An integer vector specifying the group membership.

Examples

library(grasps)## reproducibility for everythingset.seed(1234)## precision matrix estimationX <- matrix(rnorm(200), 10, 20)membership <- c(rep(1,5), rep(2,5), rep(3,4), rep(4,6))est <- grasps(X, membership = membership, penalty = "lasso", crit = "BIC")## default: keep blocks within ±1 of each groupres1 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 1)plot(res1)## wider band: keep blocks within ±2 of each groupres2 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 2)plot(res2)## special case: block-diagonal matrixres3 <- sparsify_block_banded(est$hatOmega, membership, neighbor.range = 0)plot(res3)