Dlaczego nie skorzystać z „równań normalnych”, aby znaleźć proste współczynniki najmniejszych kwadratów?

Widziałem tę listę tutaj i nie mogłem uwierzyć, że istnieje tak wiele sposobów rozwiązania najmniejszych kwadratów. „Normalne równania” na Wikipedii wydawał się być dość prosty sposób do

Dlaczego więc ich nie użyć? Zakładam, że musi istnieć problem obliczeniowy lub precyzyjny, biorąc pod uwagę, że w pierwszym linku powyżej Mark L. Stone wspomina, że ​​SVD lub QR są popularnymi metodami w oprogramowaniu statystycznym i że równania normalne są „OKAZISTE z punktu widzenia niezawodności i dokładności numerycznej”. Jednak w poniższym kodzie równania normalne dają mi dokładność do ~ 12 miejsc po przecinku w porównaniu z trzema popularnymi funkcjami pytona: polyfit numpy ; regres linowy Scipy'ego ; i scikit-learn'S regresja liniowa .

Co ciekawsze, metoda równania normalnego jest najszybsza, gdy n = 100000000. Czasy obliczeniowe dla mnie wynoszą: 2,5 s dla regresji liniowej; 12,9s dla polyfit; 4.2s dla regresji liniowej; i 1,8 s dla równania normalnego.

Kod:

import numpy as np
from sklearn.linear_model import LinearRegression
from scipy.stats import linregress
import timeit

b0 = 0
b1 = 1
n = 100000000
x = np.linspace(-5, 5, n)
np.random.seed(42)
e = np.random.randn(n)
y = b0 + b1*x + e

# scipy                                                                                                                                     
start = timeit.default_timer()
print(str.format('{0:.30f}', linregress(x, y)[0]))
stop = timeit.default_timer()
print(stop - start)

# numpy                                                                                                                                      
start = timeit.default_timer()
print(str.format('{0:.30f}', np.polyfit(x, y, 1)[0]))
stop = timeit.default_timer()
print(stop - start)

# sklearn                                                                                                                                    
clf = LinearRegression()
start = timeit.default_timer()
clf.fit(x.reshape(-1, 1), y.reshape(-1, 1))
stop = timeit.default_timer()
print(str.format('{0:.30f}', clf.coef_[0, 0]))
print(stop - start)

# normal equation                                                                                                                            
start = timeit.default_timer()
slope = np.sum((x-x.mean())*(y-y.mean()))/np.sum((x-x.mean())**2)
stop = timeit.default_timer()
print(str.format('{0:.30f}', slope))
print(stop - start) 

$Ax \approx b$$A$$A^TA$$log_{10}(cond)$$A^TA$$A^TAx = A^Tb$cyfr dokładności. Tzn., Utworzenie równania normalnego podwoiło liczbę utraconych cyfr dokładności, od samego początku.$log_{10}(cond(A^TA)) = 2*log_{10}(cond(A))$

$10^8$$10^{16}$

Czasami unika się normalnych równań, a czasem nie.

Jeśli musisz rozwiązać tylko ten jeden zmienny problem, skorzystaj z formuły. Nie ma w tym nic złego. Widziałem na przykład, jak piszesz kilka wierszy kodu w ASM dla urządzenia osadzonego. W rzeczywistości użyłem tego rodzaju rozwiązania w niektórych sytuacjach. Oczywiście nie trzeba przeciągać dużych bibliotek statystycznych, aby rozwiązać ten jeden mały problem.

Niestabilność numeryczna i wydajność są problemami większych problemów i ogólnych ustawień. Jeśli rozwiążesz wielowymiarowe najmniejsze kwadraty itp. W przypadku ogólnego problemu oczywiście tego nie użyłbyś.

Żaden nowoczesny pakiet statystyczny nie rozwiązałby regresji liniowej za pomocą równań normalnych. Równania normalne istnieją tylko w księgach statystycznych.

Równań normalnych nie należy używać, ponieważ obliczanie odwrotności macierzy jest bardzo problematyczne.

Po co używać opadania gradientu do regresji liniowej, gdy dostępne jest rozwiązanie matematyczne w formie zamkniętej?

... chociaż dostępne jest bezpośrednie równanie normalne. Zauważ, że w równaniu normalnym należy odwrócić macierz. Teraz odwracanie macierzy kosztuje O (N3) do obliczeń, gdzie N jest liczbą wierszy w macierzy X, tj. Obserwacjami. Co więcej, jeśli X jest źle warunkowany, spowoduje to błędy obliczeniowe w oszacowaniu ...