Rozumiem, że numpy.linalg.lstsq opiera się na procedurze LAPACK dgelsd .
Problem polega na rozwiązaniu:
minimalizuj ( ponadx )∥ A x - b ∥2)
Oczywiście nie ma to unikalnego rozwiązania dla macierzy A, której ranga jest mniejsza niż długość wektora . W przypadku nieokreślonego systemu udostępnia rozwiązanie , które:bzdgelsd
z
- A z = b
- ∥ z ∥2)≤ ∥ x ∥2) dla wszystkich które spełniają . (tzn. jest minimalnym rozwiązaniem normalnym dla nieokreślonego systemu.xA x = bz
Przykład: jeśli system to , numpy.linalg.lstsq zwróci .x + y= 1x = 0,5 , y= .5
Jak działa dgelsd?
Procedura dgelsd
oblicza rozkład wartości pojedynczej (SVD) A.
Po prostu naszkicuję pomysł użycia SVD do rozwiązania układu liniowego. Rozkład wartości w liczbie pojedynczej jest faktoryzacją gdzie i są macierzami ortogonalnymi, a jest macierzą diagonalną, w której wpisy diagonalne są znane jako wartości osobliwe.UΣ V.′= AUVΣ
Efektywna ranga macierzy będzie liczbą pojedynczych wartości, które są faktycznie niezerowe (tj. Wystarczająco różne od zera w stosunku do precyzji maszyny itp.). Niech będzie macierzą diagonalną niezerowych wartości pojedynczych. SVD jest zatem:AS
A=U[S000]V′
Pseudo odwrotności z jest dana przez:A
A†=V[S−1000]U′
Rozważ rozwiązanie . Następnie:x=A†b
Ax−b=U[S000]V′V[S−1000]U′b−b=U[I000]U′b−b
Są tu w zasadzie dwa przypadki:
- Liczba niezerowych wartości pojedynczych (tj. Rozmiar macierzy ) jest mniejsza niż długość . Rozwiązanie tutaj nie będzie dokładne; rozwiążemy układ liniowy w sensie najmniejszych kwadratów.Ib
- Ax−b=0
Ta ostatnia część jest nieco trudna ... trzeba śledzić wymiary macierzy i używać tego, że jest macierzą ortogonalną.U
Równoważność pseudo-odwrotności
Kiedy ma liniowo niezależne rzędy (np. Mamy macierz tłuszczu), to:
AA†=A′(AA′)−1
W przypadku nieokreślonego systemu można pokazać, że pseudo-odwrotność daje minimalne rozwiązanie normy.
Kiedy ma liniowo niezależne kolumny (np. Mamy cienką matrycę), to:
AA†=(A′A)−1A′