gRim (version 0.3.4)

fit_ggm_grips: Fit Gaussian graphical models

Description

Fit Gaussian graphical models using various algorithms.

Usage

fit_ggm_grips(
  S,
  formula = NULL,
  nobs,
  K = NULL,
  maxit = 10000L,
  eps = 0.01,
  convcrit = 1,
  aux = list(),
  method = "ncd",
  print = 0
)

Arguments

S

Sample covariance matrix.

formula

Generators of model; a list of integer vectors or a 2 x p matrix of integers.

nobs

Number of observations

K

Initial value of concentration matrix.

maxit

Maximum number of iterations.

eps

Convergence criterion.

convcrit

Convergence criterions. See section details.

aux

A list of form name=value.

method

Either "ncd" (default), "covips" or "conips".

print

Should output from fitting be printed?

Author

Søren Højsgaard,

Details

Convergence criterion:

• 1: max absolute difference between S and Sigmahat on edges.

• 2: difference in log likelihood divided by number of parameters in the model (number of edges + number of nodes) between successive iterations.

• 3: computed duality gap may turn negative due to rounding error, so its absolute value is returned. This still provides upper bound on error of likelihood function.

Methods:

• "ncd": Neighbour coordinate descent.

• "covips": IPS based on working with the covariance matrix.

• "conips": IPS based on working with the concentration matrix.

ncd is very fast but may fail to converge in rare cases. Both covips and conips are guaranteed to converge provided the maximum likelihood estimate exists, and covips are considerably faster than conips.

Examples

options("digits"=3)
data(math, package="gRbase")

S <- cov(math)
nobs <- nrow(math)
gl <- list(1:3, 3:5)
em <- matrix(c(1,2, 2,3, 1,3, 3,4, 3,5, 4,5), nrow=2)

EPS = 1e-2

fit_cov = fit_ggm_grips(S, gl, nobs=nobs, eps=EPS, method="cov")
fit_con = fit_ggm_grips(S, gl, nobs=nobs, eps=EPS, method="con")
fit_ncd = fit_ggm_grips(S, gl, nobs=nobs, eps=EPS, method="ncd")

K <- solve(S)
(fit_con$K - K)  |> abs() |> max()
(fit_cov$K - K)  |> abs() |> max()
(fit_ncd$K - K)  |> abs() |> max()