Rank-One Nonnegative Matrix Factorization

Author

CVXPY Developers and Balasubramanian Narasimhan

Introduction

The DGP atom library has several functions of positive matrices, including the trace, (matrix) product, sum, Perron-Frobenius eigenvalue, and \((I - X)^{-1}\) (eye-minus-inverse). In this example, we use some of these atoms to approximate a partially known elementwise positive matrix as the outer product of two positive vectors.

We would like to approximate \(A\) as the outer product of two positive vectors \(x\) and \(y\), with \(x\) normalized so that the product of its entries equals \(1\). Our criterion is the average relative deviation between the entries of \(A\) and \(xy^T\), that is,

\[ \frac{1}{mn} \sum_{i=1}^{m} \sum_{j=1}^{n} R(A_{ij}, x_i y_j), \]

where \(R\) is the relative deviation of two positive numbers, defined as

\[ R(a, b) = \max\{a/b, b/a\} - 1. \]

The corresponding optimization problem is

\[ \begin{array}{ll} \mbox{minimize} & \frac{1}{mn} \sum_{i=1}^{m} \sum_{j=1}^{n} R(X_{ij}, x_i y_j) \\ \mbox{subject to} & x_1 x_2 \cdots x_m = 1 \\ & X_{ij} = A_{ij}, \quad \text{for } (i, j) \in \Omega, \end{array} \]

with variables \(X \in \mathbf{R}^{m \times n}_{++}\), \(x \in \mathbf{R}^{m}_{++}\), and \(y \in \mathbf{R}^{n}_{++}\). We can cast this problem as an equivalent generalized geometric program by discarding the \(-1\) from the relative deviations.

Problem Data

We consider the partially known matrix

\[ A = \begin{bmatrix} 1.0 & ? & 1.9 \\ ? & 0.8 & ? \\ 3.2 & 5.9& ? \end{bmatrix}, \]

where the question marks denote the missing entries.

Problem Formulation and Solution

m <- 3
n <- 3

X <- Variable(c(m, n), pos = TRUE)
x <- Variable(m, pos = TRUE)
y <- Variable(n, pos = TRUE)

# Construct the outer product matrix: each row i is x[i] * y^T
outer_product <- vstack(t(x[1] * y), t(x[2] * y), t(x[3] * y))

# Relative deviations (without the -1): max(X/outer_product, outer_product/X)
relative_deviations <- max_elemwise(
  multiply(X, power(outer_product, -1)),
  multiply(power(X, -1), outer_product)
)

Warning: `multiply()` is deprecated. Use `x * y` instead.
This warning is displayed once per session.

objective <- sum_entries(relative_deviations)

constraints <- list(
  X[1, 1] == 1.0,
  X[1, 3] == 1.9,
  X[2, 2] == 0.8,
  X[3, 1] == 3.2,
  X[3, 2] == 5.9,
  x[1] * x[2] * x[3] == 1.0
)

problem <- Problem(Minimize(objective), constraints)
result <- psolve(problem, gp = TRUE)
check_solver_status(problem)

cat(sprintf("Optimal value (average relative deviation): %g\n",
            1.0 / (m * n) * (result - m * n)))
cat("Outer product approximation:\n")
print(value(outer_product))
cat("x:", value(x), "\n")
cat("y:", value(y), "\n")

Optimal value (average relative deviation): 0
Outer product approximation:
          [,1]    [,2]      [,3]
[1,] 1.0000000 1.84375 1.9000000
[2,] 0.4338983 0.80000 0.8244068
[3,] 3.2000000 5.90000 6.0800000
x: 0.8963701 0.3889335 2.868384 
y: 1.115611 2.056907 2.11966

The optimal average relative deviation is essentially zero, meaning the outer product \(xy^T\) matches \(A\) exactly at the known entries. The solver fills in the unknown entries of \(X\) to be consistent with the rank-one outer product structure.

Session Info

R version 4.5.3 (2026-03-11)
Platform: aarch64-apple-darwin20
Running under: macOS Tahoe 26.3.1

Matrix products: default
BLAS:   /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRblas.0.dylib 
LAPACK: /Library/Frameworks/R.framework/Versions/4.5-arm64/Resources/lib/libRlapack.dylib;  LAPACK version 3.12.1

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

time zone: America/Los_Angeles
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] CVXR_1.8.1

loaded via a namespace (and not attached):
 [1] slam_0.1-55       cli_3.6.5         knitr_1.51        ECOSolveR_0.6.1  
 [5] rlang_1.1.7       xfun_0.56         clarabel_0.11.2   otel_0.2.0       
 [9] gurobi_13.0-1     Rglpk_0.6-5.1     highs_1.12.0-3    cccp_0.3-3       
[13] scs_3.2.7         S7_0.2.1          jsonlite_2.0.0    backports_1.5.0  
[17] rprojroot_2.1.1   htmltools_0.5.9   Rmosek_11.1.1     gmp_0.7-5.1      
[21] piqp_0.6.2        rmarkdown_2.30    grid_4.5.3        evaluate_1.0.5   
[25] fastmap_1.2.0     yaml_2.3.12       compiler_4.5.3    codetools_0.2-20 
[29] htmlwidgets_1.6.4 Rcpp_1.1.1        here_1.0.2        osqp_1.0.0       
[33] lattice_0.22-9    digest_0.6.39     checkmate_2.3.4   Matrix_1.7-4     
[37] tools_4.5.3

References

Agrawal, A., Diamond, S., Boyd, S. (2019). Disciplined Geometric Programming. Optimization Letters, 13(5), 961–976.