I have this line of MATLAB code:
a/b
I am using these inputs:
a = [1,2,3,4,5,6,7,8,9,1,2,3,4,5,6,7,8,9]
b = ones(25, 18)
This is the result (a 1x25 matrix):
[5,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]
What is MATLAB doing? I am trying to duplicate this behavior in Python, and the mrdivide
documentation in MATLAB was unhelpful. Where does the 5 come from, and why are the rest of the values 0?
I have tried this with other inputs and receive similar results, usually just a different first element and zeros filling the remainder of the matrix. In Python when I use linalg.lstsq(b.T,a.T)
, all of the values in the first matrix returned (i.e. not the singular one) are 0.2. I have already tried right division in Python and it gives something completely off with the wrong dimensions.
I understand what a least square approximation is, I just need to know what mrdivide
is doing.
Per this handy "cheat sheet" of numpy for matlab users,
linalg.lstsq(b,a)
--linalg
is numpy.linalg.linalg, a light-weight version of the fullscipy.linalg
.a/b finds the least square solution to the system of linear equations bx = a
if b is invertible, this is a*inv(b), but if it isn't, the it is the x which minimises norm(bx-a)
You can read more about least squares on wikipedia.
according to matlab documentation, mrdivide will return at most k non-zero values, where k is the computed rank of b. my guess is that matlab in your case solves the least squares problem given by replacing b by b(:1) (which has the same rank). In this case the moore-penrose inverse
b2 = b(1,:); inv(b2*b2')*b2*a'
is defined and gives the same answerMRDIVIDE or the
/
operator actually solves thexb = a
linear system, as opposed to MLDIVIDE or the\
operator which will solve the systembx = a
.To solve a system
xb = a
with a non-symmetric, non-invertible matrixb
, you can either rely onmridivide()
, which is done via factorization ofb
with Gauss elimination, orpinv()
, which is done via Singular Value Decomposition, and zero-ing of the singular values below a (default) tolerance level.Here is the difference (for the case of
mldivide
): What is the difference between PINV and MLDIVIDE when I solve A*x=b?In your example:
the system is underdetermined, and the two different solutions will be:
In both cases the approximation error of
xb-a
is non-negligible (non-exact solution) and the same, i.e.norm(x1*b-a)
andnorm(x2*b-a)
will return the same result.What is MATLAB doing?
A great break-down of the algorithms (and checks on properties) invoked by the '\' operator, depending upon the structure of matrix
b
is given in this post in scicomp.stackexchange.com. I am assuming similar options apply for the/
operator.For your example, MATLAB is most probably doing a Gaussian elimination, giving the sparsest solution amongst a infinitude (that's where the 5 comes from).
What is Python doing?
Python, in
linalg.lstsq
uses pseudo-inverse/SVD, as demonstrated above (that's why you get a vector of 0.2's). In effect, the following will both give you the same result as MATLAB'spinv()
: