Three-Layer Non-negative Matrix Factorization (NMF-AE)

nmfae fits a three-layer nonnegative matrix factorization model \(Y_1 \approx X_1 \Theta X_2 Y_2\), where \(X_1\) is a decoder basis (column sum 1), \(\Theta\) is a bottleneck parameter matrix, \(X_2\) is an encoder basis (row sum 1), and \(Y_2\) is the input matrix.

When Y2 = Y1, the model acts as a non-negative autoencoder. When Y1 != Y2, it acts as a heteroencoder.

Initialization uses a three-step NMF procedure via nmfkc: (1) nmfkc(Y1, rank=Q) to obtain \(X_1\), (2) nmfkc(Y1, A=Y2, rank=Q) with fixed \(X_1\) to obtain \(C = \Theta X_2\), (3) nmfkc(Y2, rank=R) to factor \(C\) into \(\Theta\) and \(X_2\).

Usage

nmfae(
  Y1,
  Y2 = Y1,
  rank = 2,
  rank.encoder = rank,
  epsilon = 1e-04,
  maxit = 5000,
  verbose = FALSE,
  ...
)

Source

Satoh, K. (2025). Applying Non-negative Matrix Factorization with Covariates to Multivariate Time Series. Japanese Journal of Statistics and Data Science.

Arguments

Y1

Output matrix \(Y_1\) (P1 x N). Non-negative. May contain NAs (handled via Y1.weights).

Y2

Input matrix \(Y_2\) (P2 x N). Non-negative. Default is Y1 (autoencoder).

rank

Integer. Rank of the decoder basis \(X_1\) (P1 x Q). Default is 2. For backward compatibility, Q is accepted via ....

rank.encoder

Integer. Rank of the encoder basis \(X_2\) (R x P2). Default is rank. For backward compatibility, R is accepted via ....

epsilon

Positive convergence tolerance. Default is 1e-4.

maxit

Maximum number of multiplicative update iterations. Default is 5000.

verbose

Logical. If TRUE, prints progress messages during fitting. Default is FALSE.

...

Additional arguments:

Y1.weights: Weight matrix (P1 x N) or vector for \(Y_1\). 0 indicates missing/ignored elements. Default: auto-detect NAs.
C.L1: L1 regularization parameter for \(C\). Default is 0.
X1.L2.ortho: L2 orthogonality regularization for \(X_1\) columns. Default is 0.
X2.L2.ortho: L2 orthogonality regularization for \(X_2\) rows. Default is 0.
seed: Integer seed for reproducibility. Default is 123.
print.trace: Logical. If TRUE, prints progress. Default is FALSE.

Value

An object of class "nmfae", a list with components:

X1: Decoder basis matrix (P1 x Q), column sum 1.
C: Parameter matrix (Q x R).
X2: Encoder basis matrix (R x P2), row sum 1.
Y1hat: Fitted values \(X_1 \Theta X_2 Y_2\) (P1 x N).
rank: Named integer vector c(Q, R).
objfunc: Final objective value.
objfunc.iter: Objective values by iteration.
r.squared: Coefficient of determination \(R^2\).
niter: Number of iterations performed.
runtime: Elapsed time as a difftime object.
n.missing: Number of missing (or zero-weighted) elements in \(Y_1\).
n.total: Total number of elements in \(Y_1\) (P1 x N).

Lifecycle

This function is experimental. The interface may change in future versions.

References

Lee, D. D. and Seung, H. S. (2001). Algorithms for Non-negative Matrix Factorization. Advances in Neural Information Processing Systems, 13.

Saha, S. et al. (2022). Hierarchical Deep Learning Neural Network (HiDeNN): An Artificial Intelligence (AI) Framework for Computational Science and Engineering. Computer Methods in Applied Mechanics and Engineering, 399.

Examples

# Autoencoder example
Y <- matrix(c(1,0,1,0, 0,1,0,1, 1,1,0,0), nrow=3, byrow=TRUE)
res <- nmfae(Y, rank=2, rank.encoder=2)
res$r.squared
#> [1] 0.7480689

# Heteroencoder example
Y1 <- matrix(c(1,0,0,1), nrow=2)
Y2 <- matrix(runif(8), nrow=4)
res2 <- nmfae(Y1, Y2, rank=2, rank.encoder=2)