gaurnmf {rnnmf}R Documentation

gaurnmf .

Description

Additive update Non-negative matrix factorization with regularization, general form.

Usage

gaurnmf(
  Y,
  L,
  R,
  W_0R = NULL,
  W_0C = NULL,
  W_1L = 0,
  W_1R = 0,
  W_2RL = 0,
  W_2CL = 0,
  W_2RR = 0,
  W_2CR = 0,
  tau = 0.1,
  annealing_rate = 0.01,
  check_optimal_step = TRUE,
  zero_tolerance = 1e-12,
  max_iterations = 1000L,
  min_xstep = 1e-09,
  on_iteration_end = NULL,
  verbosity = 0
)

Arguments

Y

an r \times c matrix to be decomposed. Should have non-negative elements; an error is thrown otherwise.

L

an r \times d matrix of the initial estimate of L. Should have non-negative elements; an error is thrown otherwise.

R

an d \times c matrix of the initial estimate of R. Should have non-negative elements; an error is thrown otherwise.

W_0R

the row space weighting matrix. This should be a positive definite non-negative symmetric r \times r matrix. If omitted, it defaults to the properly sized identity matrix.

W_0C

the column space weighting matrix. This should be a positive definite non-negative symmetric c \times c matrix. If omitted, it defaults to the properly sized identity matrix.

W_1L

the \ell_1 penalty matrix for the matrix R. If a scalar, corresponds to that scalar times the all-ones matrix. Defaults to all-zeroes matrix, which is no penalty term.

W_1R

the \ell_1 penalty matrix for the matrix L. If a scalar, corresponds to that scalar times the all-ones matrix. Defaults to all-zeroes matrix, which is no penalty term.

W_2RL

the \ell_2 row penalty matrix for the matrix L. If a scalar, corresponds to that scalar times the identity matrix. Can also be a list, in which case W_2CL must be a list of the same length. The list should consist of \ell_2 row penalty matrices. Defaults to all-zeroes matrix, which is no penalty term.

W_2CL

the \ell_2 column penalty matrix for the matrix L. If a scalar, corresponds to that scalar times the identity matrix. Can also be a list, in which case W_2RL must be a list of the same length. The list should consist of \ell_2 column penalty matrices. Defaults to all-zeroes matrix, which is no penalty term.

W_2RR

the \ell_2 row penalty matrix for the matrix R. If a scalar, corresponds to that scalar times the identity matrix. Can also be a list, in which case W_2CR must be a list of the same length. The list should consist of \ell_2 row penalty matrices. Defaults to all-zeroes matrix, which is no penalty term.

W_2CR

the \ell_2 column penalty matrix for the matrix R. If a scalar, corresponds to that scalar times the identity matrix. Can also be a list, in which case W_2RR must be a list of the same length. The list should consist of \ell_2 column penalty matrices. Defaults to all-zeroes matrix, which is no penalty term.

tau

the starting shrinkage factor applied to the step length. Should be a value in (0,1).

annealing_rate

the rate at which we scale the shrinkage factor towards 1. Should be a value in [0,1).

check_optimal_step

if TRUE, we attempt to take the optimal step length in the given direction. If not, we merely take the longest feasible step in the step direction.

zero_tolerance

values of x less than this will be ‘snapped’ to zero. This happens at the end of the iteration and does not affect the measurement of convergence.

max_iterations

the maximum number of iterations to perform.

min_xstep

the minimum L-infinity norm of the step taken. Once the step falls under this value, we terminate.

on_iteration_end

an optional function that is called at the end of each iteration. The function is called as on_iteration_end(iteration=iteration, Y=Y, L=L, R=R, Lstep=Lstep, Rstep=Rstep, ...)

verbosity

controls whether we print information to the console.

Details

Attempts to factor given non-negative matrix Y as the product LR of two non-negative matrices. The objective function is Frobenius norm with \ell_1 and \ell_2 regularization terms. We seek to minimize the objective

\frac{1}{2}tr((Y-LR)' W_{0R} (Y-LR) W_{0C}) + tr(W_{1L}'L) + tr(W_{1R}'R) + \frac{1}{2} \sum_j tr(L'W_{2RLj}LW_{2CLj}) + tr(R'W_{2RRj}RW_{2CRj}),

subject to L \ge 0 and R \ge 0 elementwise, where tr(A) is the trace of A.

The code starts from initial estimates and iteratively improves them, maintaining non-negativity. This implementation uses the Lee and Seung step direction, with a correction to avoid divide-by-zero. The iterative step is optionally re-scaled to take the steepest descent in the step direction.

Value

a list with the elements

L

The final estimate of L.

R

The final estimate of R.

Lstep

The infinity norm of the final step in L

.

Rstep

The infinity norm of the final step in R

.

iterations

The number of iterations taken.

converged

Whether convergence was detected.

Note

This package provides proof of concept code which is unlikely to be fast or robust, and may not solve the optimization problem at hand. User assumes all risk.

Author(s)

Steven E. Pav shabbychef@gmail.com

References

Merritt, Michael, and Zhang, Yin. "Interior-point Gradient Method for Large-Scale Totally Nonnegative Least Squares Problems." Journal of Optimization Theory and Applications 126, no 1 (2005): 191–202. https://scholarship.rice.edu/bitstream/handle/1911/102020/TR04-08.pdf

Pav, S. E. "An Iterative Algorithm for Regularized Non-negative Matrix Factorizations." Forthcoming. (2024)

Lee, Daniel D. and Seung, H. Sebastian. "Algorithms for Non-negative Matrix Factorization." Advances in Neural Information Processing Systems 13 (2001): 556–562. http://papers.nips.cc/paper/1861-algorithms-for-non-negative-matrix-factorization.pdf

See Also

aurnmf

Examples


 nr <- 20
 nc <- 5
 dm <- 2
 
 randmat <- function(nr,nc,...) { matrix(pmax(0,runif(nr*nc,...)),nrow=nr) }
 set.seed(1234)
 real_L <- randmat(nr,dm+2)
 real_R <- randmat(ncol(real_L),nc)
 Y <- real_L %*% real_R
 gram_it <- function(G) { t(G) %*% G }
 W_0R <- gram_it(randmat(nr+5,nr))
 W_0C <- gram_it(randmat(nc+5,nc))
 
 wt_objective <- function(Y, L, R, W_0R, W_0C) { 
   err <- Y - L %*% R
   0.5 * sum((err %*% W_0C) * (t(W_0R) %*% err))
 }
 matrix_trace <- function(G) {
   sum(diag(G))
 }
 wt_objective(Y,real_L,real_R,W_0R,W_0C)
 
 L_0 <- randmat(nr,dm)
 R_0 <- randmat(dm,nc)
 wt_objective(Y,L_0,R_0,W_0R,W_0C)
 out1 <- gaurnmf(Y, L_0, R_0, W_0R=W_0R, W_0C=W_0C, 
         max_iterations=1e4L,check_optimal_step=FALSE)
 wt_objective(Y,out1$L,out1$R,W_0R,W_0C)
 
 W_1L <- randmat(nr,dm)
 out2 <- gaurnmf(Y, out1$L, out1$R, W_0R=W_0R, W_0C=W_0C, W_1L=W_1L, 
         max_iterations=1e4L,check_optimal_step=FALSE)
 wt_objective(Y,out2$L,out2$R,W_0R,W_0C)
 
 W_1R <- randmat(dm,nc)
 out3 <- gaurnmf(Y, out2$L, out2$R, W_0R=W_0R, W_0C=W_0C, W_1R=W_1R, 
         max_iterations=1e4L,check_optimal_step=FALSE)
 wt_objective(Y,out3$L,out3$R,W_0R,W_0C)


# example showing how to use the on_iteration_end callback to save iterates.
 max_iterations <- 1e3L
 it_history <<- rep(NA_real_, max_iterations)
 on_iteration_end <- function(iteration, Y, L, R, ...) {
   it_history[iteration] <<- wt_objective(Y,L,R,W_0R,W_0C)
 }
 out1b <- gaurnmf(Y, L_0, R_0, W_0R=W_0R, W_0C=W_0C, 
   max_iterations=max_iterations, on_iteration_end=on_iteration_end, check_optimal_step=FALSE)


# should work on sparse matrices too.
if (require(Matrix)) { 
 real_L <- randmat(nr,dm,min=-1)
 real_R <- randmat(dm,nc,min=-1)
 Y <- as(real_L %*% real_R, "sparseMatrix")
 L_0 <- as(randmat(nr,dm,min=-0.5), "sparseMatrix")
 R_0 <- as(randmat(dm,nc,min=-0.5), "sparseMatrix")
 out1 <- gaurnmf(Y, L_0, R_0, max_iterations=1e2L,check_optimal_step=TRUE)
}


[Package rnnmf version 0.3.0 Index]