一阶非负矩阵分解
DGP原子库提供了几个处理正矩阵的函数，包括迹，（矩阵）乘积，求和，Perron-Frobenius特征值和:math:`(I - X)^{-1}`（单位矩阵减逆矩阵）。在这个笔记本中，我们使用其中的一些原子来近似一个仅部分知道的逐元素正矩阵，使其成为两个正向量的外积。
我们想要将 \(A\) 近似为两个正向量 \(x\) 和 \(y\) 的外积，其中 \(x\) 的各个元素相乘得到的积等于 1。我们的准则是 \(A\) 和 \(xy^T\) 的各个元素之间的平均相对偏差，即
\[\frac{1}{mn} \sum_{i=1}^{m} \sum_{j=1}^{n} R(A_{ij}, x_iy_j),\]
其中 \(R\) 是两个正数的相对偏差，定义为
\[R(a, b) = \max\{a/b, b/a\} - 1.\]
The corresponding optimization problem is
\[\begin{split}\begin{equation}
\begin{array}{ll}
\mbox{minimize} & \frac{1}{mn} \sum_{i=1}^{m} \sum_{j=1}^{n} R(X_{ij}, x_iy_j)
\\
\mbox{subject to} & x_1x_2 \cdots x_m = 1 \\
& X_{ij} = A_{ij}, \quad \text{for } (i, j) \in \Omega,
\end{array}
\end{equation}\end{split}\]
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.
The below code constructs and solves this optimization problem, with specific problem data
\[\begin{split}A = \begin{bmatrix}
1.0 & ? & 1.9 \\
? & 0.8 & ? \\
3.2 & 5.9& ?
\end{bmatrix},\end{split}\]
import cvxpy as cp
m = 3
n = 3
X = cp.Variable((m, n), pos=True)
x = cp.Variable((m,), pos=True)
y = cp.Variable((n,), pos=True)
outer_product = cp.vstack([x[i] * y for i in range(m)])
relative_deviations = cp.maximum(
cp.multiply(X, outer_product ** -1),
cp.multiply(X ** -1, outer_product))
objective = cp.sum(relative_deviations)
constraints = [
X[0, 0] == 1.0,
X[0, 2] == 1.9,
X[1, 1] == 0.8,
X[2, 0] == 3.2,
X[2, 1] == 5.9,
x[0] * x[1] * x[2] == 1.0,
]
problem = cp.Problem(cp.Minimize(objective), constraints)
problem.solve(gp=True)
print("最优值:\n", 1.0/(m * n) * (problem.value - m * n), "\n")
print("外积近似:\n", outer_product.value, "\n")
print("x: ", x.value)
print("y: ", y.value)
Optimal value:
1.7763568394002505e-14
Outer product approximation
[[1. 1.84375 1.9 ]
[0.43389831 0.8 0.82440678]
[3.2 5.89999999 6.07999999]]
x: [0.89637009 0.38893346 2.86838428]
y: [1.11561063 2.0569071 2.1196602 ]