Is mldivide always the same as OLS in MATLAB?

Question

I am doing a comparison of some alternate linear regression techniques.

Clearly these will be bench-marked relative to OLS (Ordinary Least Squares).

But I just want a pure OLS method, no preconditioning of the data to uncover ill-conditioning in the data as you find when you use regress().

I had hoped to simply use the classic (XX)^-1XY expression? However this would necessitate using the inv() function, but in the MATLAB guide page for inv() it recommends that you use mldivide when doing least squares estimation as it is superior in terms of execution time and numerical accuracy.

However, I'm concerned as to whether it's okay to use mldivide to find the OLS estimates? As an operator it seems I can't see what the function is doing by "stepping-in" in the debugger.

Can I be assume that mldivide will produce the same answers as theoretical OLS under all conditions (including in the presence of) singular/i-ll conditioned matrices)?

If not what is the best way to compute pure OLS answers in MATLAB without any preconditioning of the data?

Answer 1

The short answer is:

When the system A*x = b is overdetermined, both algorithms provide the same answer. When the system is underdetermined, PINV will return the solution x, that has the minimum norm (min NORM(x)). MLDIVIDE will pick the solution with least number of non-zero elements.

As for how mldivide works, MathWorks also posted a description of how the function operates.

However, you might also want to have a look at this answer for the first part of the discussion about mldivide vs. other methods when the matrix A is square.

Depending on the shape and composition of the matrix you would use either Cholesky decomposition for symmetric positive definite, LU decomposition for other square matrix or QR otherwise. Then you can can hold onto the factorization and use linsolve to essentially just do back-substitution for you.

As to whether mldivide is preferable to pinv when A is either not square (overspecified) or is square but singular, the two options will give you two of the infinitely many solutions. According to those docs, both solutions will give you exact solutions:

Both of these are exact solutions in the sense that norm(A*x-b) and norm(A*y-b)are on the order of roundoff error.

Answer 2

According to the help page pinv gives a least squares solution to a system of equations, and so to solve the system Ax=b, just do x=pinv(A)*b.