Andreas M.
Andreas M.

Reputation: 153

How can i vectorize a matrix/ the input so that scipy.optimize.minimize can work with it?

I ran into an issue while converting a unconstrainted problem for scipy.optimize.minimize. I want to run the method L-BFGS.

The baseproblem looked like this:

min|| A - XY||

s.t. X,Y

while A is a given Matrix and X $\in \R^{nxl}$ and Y $\in \R^{lxm} Since scipy only accept Vector inputs, i tried to interpret XY as a bigger variable: Z=(X,Y) where i already put the columns of X and Y under each other.

First i tried to program the function that it will convert my input vector. For a base example it worked fine (maybe because the matrix was dense? idk)

Here is my code:

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

Above is just a (seemingly?) working dummy code. It gave me the/a result:

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])

If i run it with a big one it doesnt work at all.

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

without moving further. R is in real a more or less sparse Matrix. But i really dont know WHERE to begin with. Is it my functioncode? Is it the sparsity of R? Both? What is the work around?

Update: Solver works with VERY small dimension. If i go a bit bigger this error occurs:

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.

As u can see at the User time, the problem is quite small and stops working. Running it handwritten (L-BFGS) results in doing no steps/descents at all.

Upvotes: 1

Views: 430

Answers (1)

overfull hbox
overfull hbox

Reputation: 861

If we are trying to solve

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

the solution is well known to be rank k truncated SVD of A.

Instead, we are trying to solve

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

where X has shape n × k and Y has shape k × n (I'm using k because it is easier to read than l)

Claim: these are equivalent problems. To see this we need to show that:

  1. XY is rank <= k (with the above shapes).
  2. any rank k matrix can be written as a product XY (with the above shapes).

Proof:

  1. The first follows from the fact that Y is rank <= k and that the nulls pace of XY contains in the null space of Y

  2. The second follows from writing the SVD of a rank <= k matrix B = U D V* and observing that (UD) has shape n × k and V has shape k × n where we have dropped all but the first k singular values from the decomposition since they are guaranteed to be zero.

Implementation To solve the problem OP stated we need only compute the SVD of A and truncate it to rank k. You can use np.linalg.svd or sp.sparse.linalg.svds depending on if your matrix is sparse. For the numpy version the rank k svd can be computed as:

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))

The syntax of sp.sparse.linalg.svds is almost the same except you are able to specify the rank you want ahead of time.

Upvotes: 1

Related Questions