Least square Methoden: normale Gleichung vs svd

Question

Ich habe versucht, meinen eigenen Code für die lineare Regression zu schreiben, nach der normalen Gleichung, die beta = inv(X'X)X'Y. Der quadratische Fehler ist jedoch viel größer als die lstsq-Funktion in numpy.linalg. Könnte mir jemand erklären, warum die SVD-Methode (die lstsq verwendet) genauer ist als die normale Gleichung? Danke

Answer 1

Ich vermute, die Matrix X'X für Ihre Daten hat eine hohe condition number. Der Versuch, eine numerische Umkehrung einer solchen Matrix zu berechnen, kann zu großen Fehlern führen. Es ist normalerweise eine schlechte Idee, explizit eine inverse Matrix zu berechnen (siehe beispielsweise http://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/ oder http://epubs.siam.org/doi/abs/10.1137/1.9780898718027.ch14).

Sie können die Zustandsnummer mit numpy.linalg.cond überprüfen.

Hier ist ein Beispiel. Erstellen Sie zunächst X und Y:

In [186]: X = np.random.randn(500, 30) 

In [187]: Y = np.linspace(0, 1, len(X)) 

Für diese zufälligen X, die Konditionszahl ist nicht groß:

In [188]: np.linalg.cond(X.T.dot(X)) 
Out[188]: 2.4456380658308148 

Die normale Gleichung und lstsq das gleiche Ergebnis (nach numpy.allclose bei der Verwendung von Standard der diese Funktion Argumente):

In [189]: betan = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y) 

In [190]: betal, res, rnk, s = np.linalg.lstsq(X, Y) 

In [191]: np.allclose(betan, betal) 
Out[191]: True 

Jetzt zwicken X, indem zwei Säulen fast gleich gemacht werden. Dies macht X'X fast singulär, und gibt es eine große Anzahl Bedingung:

In [192]: X[:,0] = X[:,1] + 1e-8*np.random.randn(len(X)) 

In [193]: np.linalg.cond(X.T.dot(X)) 
Out[193]: 3954529794300611.5 

Nun ist die normale Gleichung ergibt ein anderes Ergebnis als lstsq:

In [194]: betan = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y) 

In [195]: betal, res, rnk, s = np.linalg.lstsq(X, Y) 

In [196]: np.allclose(betan, betal) 
Out[196]: False