spca.cavi.mvn {VBsparsePCA}

Function for the PX-CAVI algorithm using the multivariate normal slab

Description

This function employs the PX-CAVI algorithm proposed in Ning (2020). The g in the slab density of the spike and slab prior is chosen to be the multivariate normal distribution, i.e., N(0, \sigma^2/\lambda_1 I_r). Details of the model and the prior can be found in the Details section in the description of the 'VBsparsePCA()' function.

Usage

spca.cavi.mvn(
  x,
  r,
  lambda = 1,
  max.iter = 100,
  eps = 1e-04,
  jointly.row.sparse = TRUE,
  sig2.true = NA,
  threshold = 0.5,
  theta.int = NA,
  theta.var.int = NA,
  kappa.para1 = NA,
  kappa.para2 = NA,
  sigma.a = NA,
  sigma.b = NA
)

Arguments

x

Data an n*p matrix.

r

Rank.

lambda

Tuning parameter for the density g.

max.iter

The maximum number of iterations for running the algorithm.

eps

The convergence threshold; the default is 10^{-4}.

jointly.row.sparse

The default is true, which means that the jointly row sparsity assumption is used; one could not use this assumptio by changing it to false.

sig2.true

The default is false, \sigma^2 will be estimated; if sig2 is known and its value is given, then \sigma^2 will not be estimated.

threshold

The threshold to determine whether \gamma_j is 0 or 1; the default value is 0.5.

theta.int

The initial value of theta mean; if not provided, the algorithm will estimate it using PCA.

theta.var.int

The initial value of theta.var; if not provided, the algorithm will set it to be 1e-3*diag(r).

kappa.para1

The value of \alpha_1 of \pi(\kappa); default is 1.

kappa.para2

The value of \alpha_2 of \pi(\kappa); default is p+1.

sigma.a

The value of \sigma_a of \pi(\sigma^2); default is 1.

sigma.b

The value of \sigma_b of \pi(\sigma^2); default is 2.

Value

iter

The number of iterations to reach convergence.

selection

A vector (if r = 1 or with the jointly row-sparsity assumption) or a matrix (if otherwise) containing the estimated value for \boldsymbol \gamma.

theta.mean

The loadings matrix.

theta.var

The covariance of each non-zero rows in the loadings matrix.

sig2

Variance of the noise.

obj.fn

A vector contains the value of the objective function of each iteration. It can be used to check whether the algorithm converges

Examples

#In this example, the first 20 rows in the loadings matrix are nonzero, the rank is 1
set.seed(2021)
library(MASS)
library(pracma)
n <- 200
p <- 1000
s <- 20
r <- 1
sig2 <- 0.1
# generate eigenvectors
U.s <- randortho(s, type = c("orthonormal"))
U <- rep(0, p)
U[1:s] <- as.vector(U.s[, 1:r])
s.star <- rep(0, p)
s.star[1:s] <- 1
eigenvalue <- seq(20, 10, length.out = r)
# generate Sigma
theta.true <- U * sqrt(eigenvalue)
Sigma <- tcrossprod(theta.true) + sig2*diag(p)
# generate n*p dataset
X <- t(mvrnorm(n, mu = rep(0, p), Sigma = Sigma))
result <- spca.cavi.mvn(x = X, r = 1)
loadings <- result$theta.mean