CODEWITHSUNDEEP
pythonnumpyscipylinear-regression
M Waz
M Waz

Reputation: 775

What do extra results of numpy.polyfit mean?

When creating a line of best fit with numpy's polyfit, you can specify the parameter full to be True. This returns 4 extra values, apart from the coefficents. What do these values mean and what do they tell me about how well the function fits my data?

https://docs.scipy.org/doc/numpy/reference/generated/numpy.polyfit.html

What i'm doing is:

bestFit = np.polyfit(x_data, y_data, deg=1, full=True)

and I get the result:

(array([ 0.00062008,  0.00328837]), array([ 0.00323329]), 2, array([ 
1.30236506,  0.55122159]), 1.1102230246251565e-15)

The documentation says that the four extra pieces of information are: residuals, rank, singular_values, and rcond.

Edit: I am looking for a further explanation of how rcond and singular_values describes goodness of fit.

Thank you!

Upvotes: 2

Views: 3235

Answers (1)

user6655984
user6655984

Reputation:

how rcond and singular_values describes goodness of fit.

Short answer: they don't.

They do not describe how well the polynomial fits the data; this is what residuals are for. They describe how numerically robust was the computation of that polynomial.

rcond

The value of rcond is not really about quality of fit, it describes the process by which the fit was obtained, namely a least-squares solution of a linear system. Most of the time the user of polyfit does not provide this parameter, so a suitable value is picked by polyfit itself. This value is then returned to the user for their information.

rcond is used for truncation in ill-conditioned matrices. Least squares solver does two things:

  1. Finds x that minimizes the norm of residuals Ax-b
  2. If multiple x achieve this minimum, returns x with the smallest norm among those.

The second clause occurs when some changes of x do not affect the right-hand side at all. But since floating point computations are imperfect, usually what happens is that some changes of x affect the right hand side very little. And this is where rcond is used to decide when "very little" should be considered as "zero up to noise".

For example, consider the system 

x1                     = 1
x1 + 0.0000000001 * x2 = 2

This one can be solved exactly: x1 = 1 and x2 = 10000000000. But... that tiny coefficient (that in reality, came after some matrix manipulations) has some numeric error in it; for all we know it could be negative, or zero. Should we let it have such huge influence on the solution?

So, in such a situation the matrix (specifically its singular values) gets truncated at level rcond. This leaves 

x1 = 1
x1 = 2

for which the least-squares solution is x1 = 1.5, x2 = 0. Note that this solution is robust: no huge numbers from tiny fluctuations of coefficients.

Singular values

When one solves a linear system Ax = b in the least squares sense, the singular values of A determine how numerically tricky this is. Specifically, large disparity between largest and smallest singular values is problematic: such systems are ill-conditioned. An example is 

0.835*x1 + 0.667*x2 = 0.168
0.333*x1 + 0.266*x2 = 0.0067

The exact solution is (1, -1). But if the right hand side is changed from 0.067 to 0.066, the solution is (-666, 834) -- totally different. The problem is that the singular values of A are (roughly) 1 and 1e-6; this magnifies any changes on the right by the factor of 1e6.

Unfortunately, polynomial fit often results in ill-conditioned matrices. For example, fitting a polynomial of degree 24 to 25 equally spaced data points is unadvisable. 

import numpy as np
x = np.arange(25)
np.polyfit(x, x, 24, full=True)

The singular values are 

 array([4.68696731e+00, 1.55044718e+00, 7.17264545e-01, 3.14298605e-01,
        1.16528492e-01, 3.84141241e-02, 1.15530672e-02, 3.20120674e-03,
        8.20608411e-04, 1.94870760e-04, 4.28461687e-05, 8.70404409e-06,
        1.62785983e-06, 2.78844775e-07, 4.34463936e-08, 6.10212689e-09,
        7.63709211e-10, 8.39231664e-11, 7.94539407e-12, 6.32326226e-13,
        4.09332903e-14, 2.05501534e-15, 7.55397827e-17, 4.81104905e-18,
        8.98275758e-20]),

which, with the default value of rcond (5.55e-15 here), gets four of them truncated to 0.

The difference in magnitude between smallest and largest singular values indicates that perturbing the y-values by numbers of size 1e-15 can result in changes of about 1 to the coefficients. (Not every perturbation will do that, just some that happen to align with a singular vector for a small singular value).

Rank

The effective rank is just the number of singular values above the rcond threshold. In the above example it's 21. This means that even though the fit is for 25 points, and we get a polynomial with 25 coefficients, there are only 21 degrees of freedom in the solution.

Upvotes: 7

Related Questions