我如何向量化矩阵/输入以便 scipy.optimize.minimize 可以使用它？

Question

我在为 scipy.optimize.minimize 转换无约束问题时遇到了问题。我想运行 L-BFGS 方法。

基本问题如下所示：

分钟|| A - XY||

s.t. X,Y

而 A 是一个给定的矩阵和 X $\in \R^{nxl}$ 和 Y $\in \R^{lxm} 由于 scipy 只接受向量输入，我试图将 XY 解释为一个更大的变量：Z=(X,Y) 其中我已经将 X 和 Y 的列放在一起。

首先我尝试编写函数来转换我的输入向量。对于一个基本示例，它工作得很好（也许是因为矩阵很密集？idk）

这是我的代码：

import numpy as np
from scipy.optimize import minimize
R=np.array(np.arange(12)).reshape(3, 4)
Z0 = np.array(np.random.random(14))
#X=3x2 = 6
#Y=2x4 = 8
def whatineed (Z):
    return np.linalg.norm(R - np.dot(Z[:6].reshape(3,2),Z[6:].reshape(2,4)))

A = minimize(fun=whatineed, x0=Z0, method='L-BFGS-B', options={'disp': 1})

#print A

上面只是一个（看似？）工作的虚拟代码。它给了我/一个结果：

x: array([ 1.55308851, -0.50000733,  1.89812395,  1.44382572,  2.24315938, 3.38765876,  0.62668062,  1.23575295,  1.8448253 ,  2.45389762, 1.94655245,  1.83844053,  1.73032859,  1.62221667])

如果我用一个大的运行它，它根本不起作用。

RUNNING THE L-BFGS-B CODE

           * * *

Machine precision = 2.220D-16
 N =       377400     M =           10
 This problem is unconstrained.

At X0         0 variables are exactly at the bounds

没有进一步移动。 R 实际上是一个或多或少的稀疏矩阵。但我真的不知道从哪里开始。是我的功能码吗？是R的稀疏性吗？两个都？解决方法是什么？

更新：求解器适用于非常小的尺寸。如果我再大一点，就会出现这个错误：

ABNORMAL_TERMINATION_IN_LNSRCH                              

 Line search cannot locate an adequate point after 20 function
  and gradient evaluations.  Previous x, f and g restored.
 Possible causes: 1 error in function or gradient evaluation;
                  2 rounding error dominate computation.

 Cauchy                time 0.000E+00 seconds.
 Subspace minimization time 0.000E+00 seconds.
 Line search           time 0.000E+00 seconds.

 Total User time 0.000E+00 seconds.

正如您在用户时间看到的，问题非常小并且停止工作。 以手写方式运行 (L-BFGS) 导致根本不执行任何步骤/下降。

Answer 1

如果我们试图解决

min_{B is rank <= k} ‖ A - B ‖_2

众所周知，solution 是 A 的秩 k 截断 SVD。

相反，我们正在尝试解决

min_{X,Y} ‖ A - XY ‖_2

其中 X 的形状为 n × k 而 Y 的形状为 k × n （我正在使用k 因为它比 l 更容易阅读）

声明：这些是等价的问题。要看到这一点，我们需要证明：

1. XY 是等级 k（具有上述形状）。
2. 任何秩 k 矩阵都可以写成乘积 XY（具有上述形状）。

证明：

1. 第一个原因是 Y 是排名 k 并且 XY 的空值速度包含在Y

的零空间

2. 第二个是写一个秩 k 矩阵 B = U D V* 的 SVD 并观察 (UD) 的形状 n × k 和 V 的形状为 k × n，我们从分解中删除了除第一个 k 个奇异值之外的所有奇异值，因为它们保证为零。

实施 为了解决 OP 所说的问题，我们只需要计算 A 的 SVD 并将其截断以排名 k。 您可以使用np.linalg.svd 或sp.sparse.linalg.svds，具体取决于您的矩阵是否稀疏。对于 numpy 版本，排名 k svd 可以计算为：

m,n = 10,20
A = np.random.randn(m,n)

k = 6
u,s,vt = np.linalg.svd(A)

X = u[:,:k]*s[:k]
Y = vt[:k]

print(X.shape, Y.shape)
print(np.linalg.norm(A-X@Y,2))

sp.sparse.linalg.svds 的语法几乎相同，只是您可以提前指定所需的排名。