Cannot solve linear system if input matrix has less rows than columns - but in python I could

[heltonbiker]

I am trying to find out why I can solve a given linear system in python (scipy.linalg.lstsq) but not in Mathnet.

Here is the python code:

import numpy as np
import scipy.linalg

inputMatrix = np.array([[1.0, 1.2964161374116825, 1.2146279550953718, 12084.083925925881, 1.6806948013414265, 1.4753210692991645, 146025084.32882026, 1.5746632819369923, 15666.00140740743, 14677.666148148206], [1.0, -1.2190083824510871, 1.1724011058248773, 11951.71503703699, 1.4859814364860158, 1.3745243529393951, 142843492.32653612, -1.4291667755954496, -14569.240814814797, 14012.203925925984], [1.0, -1.2137377805023088, -1.1673040054155028, 11187.74392592588, 1.473159399818671, 1.3625986410590762, 125165614.15209144, 1.4168009727044675, -13578.987481481467, -13059.498296296244], [1.0, 1.3042926271247444, -1.2018810741469206, 12190.3911481481, 1.7011792571719675, 1.4445181163925556, 148605636.34484753, -1.5676046236905967, 15899.837296296313, -14651.40040740735]])

resultMatrix = np.array([[ 20, 30,90],
                        [-20, 30,90],
                        [-20,-30,90],
                        [ 20,-30,90]])

solutionMatrix = scipy.linalg.lstsq(inputMatrix, resultMatrix)[0]

print(solutionMatrix)

And here is the Mathnet code:

using MathNet.Numerics.LinearAlgebra;
using Newtonsoft.Json;

namespace ComparandoMinimosQuadradosPythonMathnet
{
    class Program
    {
        static void Main(string[] args)
        {
            var inputstring = "[[1.0, 1.2964161374116825, 1.2146279550953718, 12084.083925925881, 1.6806948013414265, 1.4753210692991645, 146025084.32882026, 1.5746632819369923, 15666.00140740743, 14677.666148148206], [1.0, -1.2190083824510871, 1.1724011058248773, 11951.71503703699, 1.4859814364860158, 1.3745243529393951, 142843492.32653612, -1.4291667755954496, -14569.240814814797, 14012.203925925984], [1.0, -1.2137377805023088, -1.1673040054155028, 11187.74392592588, 1.473159399818671, 1.3625986410590762, 125165614.15209144, 1.4168009727044675, -13578.987481481467, -13059.498296296244], [1.0, 1.3042926271247444, -1.2018810741469206, 12190.3911481481, 1.7011792571719675, 1.4445181163925556, 148605636.34484753, -1.5676046236905967, 15899.837296296313, -14651.40040740735]]";
            var input = JsonConvert.DeserializeObject<double[,]>(inputstring);
            var inputmatrix = Matrix<double>.Build.DenseOfArray(input);

            var referencia = new double[,]
            {
                { 20,  30, 90},
                {-20,  30, 90},
                {-20, -30, 90},
                { 20, -30, 90}
            };

            var refmatrix = Matrix<double>.Build.DenseOfArray(referencia);

            var outputmatrix = inputmatrix.Solve(refmatrix);
        }
    }
}

Which gives me this error:

System.ArgumentException: 'Matrix dimensions must agree: 4x10.'

There's an almost identical question in the forums, and the answer makes sense, that is, using a larger number of equations (rows) than variables (columns) to have a determinate system.

But the fact is: the python code not only runs, but it works, resulting in a calibration matrix that actually make my data (force applied to a load cell array) be correct.

So how could I do that in Mathnet the same way I did in python, and why the current form doesn't work?

Also, the error message does not feel quite appropriate for this situation, IMO.

[sebhofer]

Have you tried solving the underlying optimization problem? From the scipy documentation

Compute least-squares solution to equation Ax = b. Compute a vector x such that the 2-norm |b - A x| is minimized.

I'm sure you can also do this using Math.Net.

[heltonbiker]

Thanks @sebhofer for the input.

I'm a bit confused by your answser. As I understand, the two sentence quoted from scipy.linalg.lstsq mean the same thing, that is, they are like synonims of one another.

In my python example above, the solution is found performing one step (that is, one single call to lstsq).

If I test with larger input sizes, the same results are obtained from python and from MathNet with one call to either lstsq and Solve, since apparently they are the same thing.

But below a minimum of rows, the python code keeps working, but MathNet doesn't.

[sebhofer]

I agree that the two sentences are synonymous, I just copied them from the scipy documentation.

But then I disagree. How do you know that Solve and lstsq do the same thing? There are many different ways to solve a set of linear equations and many of them won't work on underdetermined systems. I'm assuming that MathNet is trying to use one of the algorithms that can't handle these underdetermined problems. (I have to admit that I don't know MathNet well, so this is speculation.)

My suggestion is: Set up the minimization problem as stated in the scipy documentation and use one of MathNet's minimization routines to solve it. It might not be the same solution as given by scipy, but it should give you one of the many solutions. (Don't forget to check if it actually does solve the original problem!) If you require the exact same solution as the scipy one, you'd need to dig into the scipy source to find out which algorithm (and possibly initial conditions they are using).