2

I need to solve an equation AX = B using Python where A, X, B are matrices and all values of X must be non-negative.

The best solution I've found is

X = np.linalg.lstsq(A, B, rcond=None)

but as a result X contains negative values. Is it possible to get a solution without negative values? Thanks in advance!

Андрій Немченко
  • 417
  • 6
  • 19
  • Maybe https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.lsq_linear.html#scipy.optimize.lsq_linear – user2699 May 20 '19 at 20:45
  • @user2699 `B` has to be one-dimensional with `scipy.optimize.lsq_linear` and `scipy.optimize.nnls` unfortunately. – user2653663 May 20 '19 at 20:54
  • https://stackoverflow.com/questions/10697995/linear-programming-in-python – Reiner Czerwinski May 20 '19 at 21:04
  • @user2653663 Why is that an issue? A wrapper can easily solve for each column in B, which reference solution in the question using `np.linalg.lstsq` already does. – user2699 May 21 '19 at 12:14

2 Answers2

6

In general, this is not mathematically possible. Given the basic requirements of A and B being invertible, X is a unique matrix. If you don't like the elements that X has, you can't simply ask for another solution: there isn't one. You'll have to change A or B to get a different result.

Prune
  • 76,765
  • 14
  • 60
  • 81
  • A and B do not have to be invertible - that is not a constraint given in the question or even implied. The least squares solver given in the question finds the solution that gives the minimum mean squared error - but there are other definitions of error (or constraints on B) that result in unique solutions for non-invertible A. – user2699 May 21 '19 at 12:20
2

You could solve it with cvxpy:

import cvxpy

def solve(A, B):
    """
    Minimizes |AX - B|**2, assuming A and B are 
    square matrices for simplicity. If this optimized 
    error is zero, this corresponds to solving AX = B.
    """
    n = A.shape[0]
    X = cvxpy.Variable((n,n))
    # Set objective
    obj_fun = cvxpy.sum_squares(A*X - B)
    objective = cvxpy.Minimize(obj_fun)
    # Set constraints
    constraints = [X >= 0]
    prob = cvxpy.Problem(objective, constraints)
    result = prob.solve(solver = "ECOS")
    return X.value

EDIT: The answer by Prune is correct I believe. You can check if the error in the numerical solver is non-zero by inspecting results.

user2653663
  • 2,818
  • 1
  • 18
  • 22
  • Thank you for the answer, but I got "'shape' elements cannot be negative" error with A.shape=(128,1025) and B.shape=(128,146). When I tried the solution on arrays with less dimensios, it works fine. – Андрій Немченко May 20 '19 at 21:56