pmd: Rank-1 Penalized Matrix Decomposition

In schoonees/spca: Lasso-Penalized Singular Value Decomposition and Sparse Principal Components Analysis

Description

Implementation of the rank-1 penalized matrix decomposition of Witten et al. (2009). This function applies lasso penalties to the left and right singular vectors. This function is called by mpmd, which should be used in applications.

Usage

pmd(
  Z,
  c1 = 1,
  c2 = 1,
  maxit = 100,
  eps = sqrt(.Machine$double.eps),
  centre = FALSE
)

Arguments

Z	Matrix to be decomposed
c1	L1-norm bound for u, the left singular vector. Feasible solutions are available when values greater than or equal to 1. For values larger than sqrt(nrow(Z), it has no effect.
c2	L1-norm bound for v, the right singular vector. Feasible solutions are available when values greater than or equal to 1. For values larger than sqrt(ncol(Z), it has no effect.
maxit	Maximum number of iterations
eps	Stopping criterion, an absolute error tolerance on the mean squared reconstruction error
centre	Logical indicating whether to centre the matrix Z using the overall mean before analysis

Value

A list with the penalized singular value decomposition(d, u, v), the vector of errors (error), and the number of iterations (iteration).

References

Witten, D. M., Tibshirani, R., & Hastie, T. (2009). A penalized matrix decomposition, with applications to sparse principal components and canonical correlation analysis. Biostatistics, 10(3), 515-534.

Examples

## Simple random matrix
set.seed(1)
Z <- matrix(rnorm(100), nrow = 20, ncol = 5)

## Default result
pmd(Z)

## Equivalent to rank-1 SVD (almost)
pmd_Z <- pmd(Z, c1 = sqrt(nrow(Z)), c2 = sqrt(ncol(Z)))
all.equal(pmd_Z[-1], svd(Z, nu = 1, nv = 1)[-1], 
          tolerance = (.Machine$double.eps)^0.25)
         
## Increasing c2 incrementally
pts <- seq(from = 1, to = sqrt(ncol(Z)), length.out = 10)
vmat <- matrix(NA, nrow = ncol(Z), ncol = length(pts))
for (i in seq_along(pts)) {
  vmat[, i] <- pmd(Z, c1 = sqrt(nrow(Z)), c2 = pts[i])$v
}
vmat

## Increasing c1 incrementally
pts <- seq(from = 1, to = sqrt(nrow(Z)), length.out = 10)
umat <- matrix(NA, nrow = nrow(Z), ncol = length(pts))
for (i in seq_along(pts)) {
  umat[, i] <- pmd(Z, c1 = pts[i], c2 = sqrt(ncol(Z)))$u
}
umat

## Testing reconstruction error
res <- pmd(Z, c1 = 3, c2 = 1.75)
plot(res$error)

schoonees/spca documentation built on Jan. 31, 2021, 6:21 p.m.