Package 'MVNGmod'

Title: Matrix-Variate Non-Gaussian Linear Regression Models
Description: An implementation of the expectation conditional maximization (ECM) algorithm for matrix-variate variance gamma (MVVG) and normal-inverse Gaussian (MVNIG) linear models. These models are designed for settings of multivariate analysis with clustered non-uniform observations and correlated responses. The package includes fitting and prediction functions for both models, and an example dataset from a periodontal on Gullah-speaking African Americans, with responses in 'gaad_res', and covariates in 'gaad_cov'. For more details on the matrix-variate distributions used, see Gallaugher & McNicholas (2019) <doi:10.1016/j.spl.2018.08.012>.
Authors: Samuel Soon [aut, cre], Dipankar Bandyopadhyay [aut], Qingyang Liu [aut]
Maintainer: Samuel Soon <[email protected]>
License: MIT + file LICENSE
Version: 0.1.0
Built: 2026-05-25 09:19:00 UTC
Source: https://github.com/soonsk-vcu/mvngmod

Help Index


GAAD Data

Description

These data sets describe periodontal measurements performed on members of the Gullah-Speaking African American community.

Usage

gaad_cov

Format

Each is a list of matrices, with rows denoting tooth sites and columns denoting CAL/PPD response.

Details

gaad_res and gaad_cov contain the response and covariate matrices of the GAAD data.

Examples

gaad_cov
gaad_res

GAAD Data

Description

These data sets describe periodontal measurements performed on members of the Gullah-Speaking African American community.

Usage

gaad_res

Format

Each is a list of matrices, with rows denoting tooth sites and columns denoting CAL/PPD response.

Details

gaad_res and gaad_cov contain the response and covariate matrices of the GAAD data.

Examples

gaad_cov
gaad_res

MVVG Parameter Format

Description

This is an example of the format of input parameter list theta.

Usage

gaad_theta_mvvg

Format

List of model parameters.

Examples

gaad_theta_mvvg

AECM Estimation for Matrix-Variate Normal-Inverse Gaussian Models

Description

This function fits MVNIG linear models for matrix-variate skew data with non-uniform data rows between subjects. Exchangeable observation row correlation and skewness structures are imposed to accommodate the varying row counts across matrices. Note that multiple restarts may be needed to account for unstable local maxima.

Usage

MVNIGmod(Y, X, theta_g = NULL, stopping = 0.001, max_iter = 50)

Arguments

Y

List of ni×pn_i \times p response matrices. Matrices must have same number of columns.

X

List of ni×qn_i \times q design matrices. Matrices must have same number of columns.

theta_g

List of parameters to pass as initial values in the AECM algorithm. If NULL, will be randomly generated. See Details for an in-depth explanation.

stopping

Stopping threshold for the L-infinity norm of differences in consecutive parameter space, evaluated at iteration t+1t+1 as θ^t+1θ^t|\hat\theta^{t+1} - \hat\theta^{t}|_\infty. Default is 0.001

max_iter

Maximum number of iterations, default is 50.

Details

Fits the matrix-variate skew regression model

Yi=XiΘ+Ei,Y_i = X_i \Theta + E_i,

where each response YiY_i is a ni×pn_i \times p matrix that indexes nin_i observations and pp response variables. XiX_i corresponds to a ni×qn_i \times q design matrix, and Θ\Theta corresponds to a q×pq \times p coefficient matrix. EiE_i corresponds to a ni×pn_i \times p error matrix, following a matrix-variate variance-gamma distribution.

The model estimates MVVG parameters Θ,a,r,Ψ,γ~\Theta, \underline{a}, r, \Psi, \tilde\gamma using the alternating expectation conditional maximization (AECM) algorithm, using the density

f(YiMi,a,r,Ψ,γ~,ni,p)=2exp[matlib::tr(Σi1(YiMi)Ψ1AiT)+γ~](2π)nip2+1Σip2Ψni2(δ(Yi;Mi,Σi,Ψ)+1ρ(Ai,Σi,Ψ)+γ~2)(1+nip)4×K(1+nip)2([ρ(Ai,Σi,Ψ)+γ~2][δ(Yi;Mi,Σi,Ψ)+1]),f(Y_i|M_i,\underline{a},r, \Psi, \tilde\gamma, n_i, p) = \dfrac{2 \exp[matlib::tr(\Sigma_i^{-1}(Y_i-M_i)\Psi^{-1}A_i^T) + \tilde\gamma]}{(2\pi)^{\frac{n_ip}{2} + 1} |\Sigma_i|^{\frac{p}{2}} |\Psi|^{\frac{n_i}{2}}} \bigg( \dfrac{\delta(Y_i; M_i, \Sigma_i, \Psi) + 1}{\rho (A_i, \Sigma_i,\Psi) + \tilde\gamma^2} \bigg)^{-\frac{(1+n_ip)}{4}} \\ \times K_{-\frac{(1+n_ip)}{2}} \big( \sqrt{[\rho(A_i, \Sigma_i, \Psi) + \tilde\gamma^2][\delta(Y_i;M_i,\Sigma_i,\Psi) + 1]} \big),

where Ai=1ni×aTA_i = \underline{1}_{n_i} \times \underline{a}^T, Σi=Ini+r(1ni1niTIni)\Sigma_i = I_{n_i} + r(\underline{1}_{n_i}\underline{1}_{n_i}^T - I_{n_i}), δ(X;M,Σ,Ψ)=matlib::tr(Σ1(XM)Ψ1(XM)T)\delta(X;M, \Sigma, \Psi) = matlib::tr(\Sigma^{-1}(X-M)\Psi^{-1}(X-M)^T), ρ(A,Σ,Ψ)=matlib::tr(Σ1AΨ1AT)\rho(A, \Sigma, \Psi) = matlib::tr(\Sigma^{-1}A\Psi^{-1}A^T), and Kν(x)K_{\nu}(x) is the modified Bessel function of the second kind.

The structure of theta_g and parameter estimates returned by the function must be in the form of a list with the following named elements:

Theta:

q×pq \times p coefficient matrix

a:

p×1p \times 1 skewness vector

rho:

Compound symmetry parameter for row correlation matrix

Psi:

p×pp \times p column covariance matrix

tgamma:

Univariate mixing parameter

Value

MVNIGmod returns a list with the following elements:

Iteration:

Number of iterations taken to convergence. Inf if convergence not reached.

⁠Starting Value⁠:

List of initial parameter values.

⁠Final Value⁠:

List of final parameter estimates.

⁠Stopping Criteria⁠:

Vector of θ^t+1θ^t|\hat\theta^{t+1} - \hat\theta^{t}|_\infty at each iteration.

AIC:

Model AIC

BIC:

Model BIC

Author(s)

Samuel Soon

Dipankar Bandyopadhyay

Qingyang Liu

Examples

MVNIGmod(Y,X,theta_mvnig)

set.seed(1234)
# num response variables
p <- ncol(gaad_res[[1]])
# num covariates
q <- ncol(gaad_cov[[1]])
# generate initial value to input, then run AECM with MVVG distribution
initial_mvnig_theta <- list(Theta = matrix(stats::rnorm(p*q), nrow = q, ncol = p),
                      A = rep(1,p),
                     rho = 0.3,
                     Psi = diag(p),
                     tgamma = 3)
MVNIGmod(gaad_res[1:10], gaad_cov[1:10], initial_mvnig_theta)

AECM Estimation for Matrix-Variate Variance Gamma (MVVG) Models

Description

This function fits MVVG linear models for matrix-variate skew data with non-uniform data rows between subjects. Exchangeable observation row correlation and skewness structures are imposed to accommodate the varying row counts across matrices. Note that multiple restarts may be needed to account for unstable local maxima.

Usage

MVVGmod(Y, X, theta_g = NULL, stopping = 0.001, max_iter = 50)

Arguments

Y

List of ni×pn_i \times p response matrices. Matrices must have same number of columns.

X

List of ni×qn_i \times q design matrices. Matrices must have same number of columns.

theta_g

List of parameters to pass as initial values in the AECM algorithm. If NULL, will be randomly generated. See Details for an in-depth explanation.

stopping

Stopping threshold for the L-infinity norm of differences in consecutive parameter space, evaluated at iteration t+1t+1 as θ^t+1θ^t|\hat\theta^{t+1} - \hat\theta^{t}|_\infty. Default is 0.001

max_iter

Maximum number of iterations, default is 50.

Details

Fits the matrix-variate skew regression model

Yi=XiΘ+Ei,Y_i = X_i \Theta + E_i,

where each response YiY_i is a ni×pn_i \times p matrix that indexes nin_i observations and pp response variables. XiX_i corresponds to a ni×qn_i \times q design matrix, and Θ\Theta corresponds to a q×pq \times p coefficient matrix. EiE_i corresponds to a ni×pn_i \times p error matrix, following a matrix-variate variance-gamma distribution.

The model estimates MVVG parameters Θ,a,r,Ψ,γ\Theta, \underline{a}, r, \Psi, \gamma using the alternating expectation conditional maximization (AECM) algorithm, using the density

f(YiXiΘ,a,r,Ψ,γ,ni,p)=2γγexp[matlib::tr(Σi1(YiXiΘ)Ψ1AiT)](2π)nip/2Σip/2Ψni/2Γ(γ)(δ(Yi;XiΘ,Σi,Ψ)ρ(Ai,Σi,Ψ)+2γ)(γnip/2)/2×K(γnip/2)([ρ(Ai,Σi,Ψ)+2γ][δ(Yi;XiΘ,Σi,Ψ)]),f(Y_i| X_i\Theta,\underline{a},r, \Psi, \gamma, n_i, p) = \dfrac{2\gamma^\gamma \exp[matlib::tr(\Sigma_i^{-1}(Y_i- X_i\Theta)\Psi^{-1}A_i^T)]}{(2\pi)^{n_ip/2} |\Sigma_i|^{p/2} |\Psi|^{n_i/2} \Gamma(\gamma)} \bigg( \dfrac{\delta(Y_i; X_i\Theta, \Sigma_i, \Psi)}{\rho (A_i, \Sigma_i,\Psi) + 2\gamma} \bigg)^{(\gamma - n_ip/2)/2} \\ \times K_{(\gamma - n_ip/2)} \big( \sqrt{[\rho(A_i, \Sigma_i, \Psi) + 2\gamma][\delta(Y_i; X_i\Theta,\Sigma_i,\Psi)]} \big),

where Ai=1ni×aTA_i = \underline{1}_{n_i} \times \underline{a}^T, Σi=Ini+r(1ni1niTIni)\Sigma_i = I_{n_i} + r(\underline{1}_{n_i}\underline{1}_{n_i}^T - I_{n_i}), δ(X;M,Σ,Ψ)=matlib::tr(Σ1(XM)Ψ1(XM)T)\delta(X;M, \Sigma, \Psi) = matlib::tr(\Sigma^{-1}(X-M)\Psi^{-1}(X-M)^T), ρ(A,Σ,Ψ)=matlib::tr(Σ1AΨ1AT)\rho(A, \Sigma, \Psi) = matlib::tr(\Sigma^{-1}A\Psi^{-1}A^T), and Kν(x)K_{\nu}(x) is the modified Bessel function of the second kind.

The structure of theta_g and parameter estimates returned by the function must be in the form of a list with the following named elements:

Theta:

q×pq \times p coefficient matrix

a:

p×1p \times 1 skewness vector

rho:

Compound symmetry parameter for row correlation matrix

Psi:

p×pp \times p column covariance matrix

gamma:

Univariate mixing parameter

Value

MVVGmod returns a list with the following elements:

Iteration:

Number of iterations taken to convergence. Inf if convergence not reached.

⁠Starting Value⁠:

List of initial parameter values.

⁠Final Value⁠:

List of final parameter estimates.

⁠Stopping Criteria⁠:

Vector of θ^t+1θ^t|\hat\theta^{t+1} - \hat\theta^{t}|_\infty at each iteration.

AIC:

Model AIC

BIC:

Model BIC

Author(s)

Samuel Soon

Dipankar Bandyopadhyay

Qingyang Liu

Examples

MVVGmod(Y,X,theta_mvvg)

set.seed(1234)
# num response variables
p <- ncol(gaad_res[[1]])
# num covariates
q <- ncol(gaad_cov[[1]])
# generate initial value to input, then run AECM with MVVG distribution
initial_gaad_theta_mvvg <- list(Theta = matrix(stats::rnorm(p*q), nrow = q, ncol = p),
                      A = rep(1,p),
                     rho = 0.3,
                     Psi = diag(p),
                     gamma = 4)
MVVGmod(gaad_res[1:10], gaad_cov[1:10], initial_gaad_theta_mvvg)

MVNG Model Prediction

Description

Predicts response values given a list of covariate matrices and a model output from either MVVGmod or MVNIGmod.

Usage

predict(mod, X)

Arguments

mod

object outputted by either MVVGmod or MVNIGmod

X

Inputted covariate matrix

Value

Returns a list of predicted response matrices

Author(s)

Samuel Soon

Dipankar Bandyopadhyay

Qingyang Liu

Examples

set.seed(1234)
# num response variables
p <- ncol(gaad_res[[1]])
# num covariates
q <- ncol(gaad_cov[[1]])
# generate initial value to input, then run AECM with MVVG distribution
initial_mvnig_theta <- list(Theta = matrix(stats::rnorm(p*q), nrow = q, ncol = p),
                      A = rep(1,p),
                     rho = 0.3,
                     Psi = diag(p),
                     tgamma = 4)
mvnig_mod <- MVNIGmod(gaad_res[1:10], gaad_cov[1:10], initial_mvnig_theta)

predict(mvnig_mod, gaad_cov[1:10])

Toy Response Initial Parameter (MVNIG)

Description

Part of toy dataset for examples.

Usage

theta_mvnig

Format

List of parameters for input to MVNIGmod function

Examples

theta_mvvg

Toy Response Initial Parameter (MVVG)

Description

Part of toy dataset for examples.

Usage

theta_mvvg

Format

List of parameters for input to MVVGmod function

Examples

theta_mvvg

Toy Covariate Matrices

Description

Part of toy dataset for examples.

Usage

X

Format

List of covariate matrices for individual subjects

Examples

X

Toy Response Matrices

Description

Part of toy dataset for examples.

Usage

Y

Format

List of response matrices for individual subjects

Examples

Y