Ihr Ziel ist es zu lösen:Was ist die genaueste Methode in Python für die Berechnung der Minimum-Norm-Lösung oder die Lösung aus der Pseudo-Inverse?
Kc=y
mit der pseudo-inversen (d.h. Minimalnormlösung):
c=K^{+}y
, so daß das Modell (hoffentlich) hohe Grad f(x) = sum_i c_i x^i Polynom-Modell. Ich interessiere mich besonders für den unterbestimmten Fall, in dem wir mehr polynomische Merkmale als Daten haben (wenige Gleichungen zu viele Variablen/Unbekannte) columns = deg+1 > N = rows. Hinweis K ist die Vandermode-Matrix von Polynom-Features.
Ich war zunächst mit der Python-Funktion np.linalg.pinv, aber dann bemerkte ich etwas funky ging weiter, wie ich hier notiert: Why do different methods for solving Xc=y in python give different solution when they should not?. In dieser Frage verwende ich eine quadratische Matrix, um eine Funktion auf dem Intervall [-1.+1] mit einem Polynom hohen Grades zu lernen. Die Antwort dort schlug mir vor, den Grad des Polynoms zu verringern und/oder die Intervallgröße zu erhöhen. Das Hauptproblem ist, dass es mir nicht klar ist, wie man das Intervall oder den maximalen Grad wählt, bevor die Sache unzuverlässig wird. Ich denke, mein Hauptproblem ist, dass die Wahl eines solchen numerisch stabilen Bereichs von der Methode abhängt, die ich verwenden kann. Am Ende, was wirklich ich wichtig ist, dass
- die Methode, die ich verwenden ist genau (oder ganz in der Nähe) an die pseudo-inverse für dieses Polynomanpassung Problem
- , dass sein numerisch stabil
Idealerweise möchte ich ein großes Polynom versuchen, aber das könnte durch meine Maschinengenauigkeit begrenzt sein. Ist es möglich, die numerische Präzision der Maschine zu erhöhen, indem man etwas genauer als Schwimmer verwendet?
Auch ich wirklich kümmern uns, dass was Funktion von Python Ich benutze es die nächste Antwort auf die bietet Pseudo inverse (und hoffentlich, dass seine numerisch stabil, so kann ich es tatsächlich nutzen). die Antwort für die pseudo-inverse Um zu überprüfen, ich das folgende Skript geschrieben:
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
def l2_loss(y,y_):
N = y.shape[0]
return (1/N)*np.linalg.norm(y-y_)
## some parameters
lb,ub = -200,200
N=100
D0=1
degree_mdl = 120
## target function
freq_cos = 2
f_target = lambda x: np.cos(freq_cos*2*np.pi*x)
## evaluate target_f on x_points
X = np.linspace(lb,ub,N) # [N,]
Y = f_target(X) # [N,]
# get pinv solution
poly_feat = PolynomialFeatures(degree=degree_mdl)
Kern = poly_feat.fit_transform(X.reshape(N,D0)) # low degrees first [1,x,x**2,...]
c_pinv = np.dot(np.linalg.pinv(Kern), Y)
## get polyfit solution
c_polyfit = np.polyfit(X,Y,degree_mdl)[::-1] # need to reverse to get low degrees first [1,x,x**2,...]
##
c_lstsq,_,_,_ = np.linalg.lstsq(Kern,Y.reshape(N,1))
##
print('lb,ub = {} '.format((lb,ub)))
print('differences with c_pinv')
print('||c_pinv-c_pinv||^2 = {}'.format(np.linalg.norm(c_pinv-c_pinv)))
print('||c_pinv-c_polyfit||^2 = {}'.format(np.linalg.norm(c_pinv-c_polyfit)))
print('||c_pinv-c_lstsq||^2 = {}'.format(np.linalg.norm(c_pinv-c_lstsq)))
##
print('differences with c_polyfit')
print('||c_polyfit-c_pinv||^2 = {}'.format(np.linalg.norm(c_polyfit-c_pinv)))
print('||c_polyfit-c_polyfit||^2 = {}'.format(np.linalg.norm(c_polyfit-c_polyfit)))
print('||c_polyfit-c_lstsq||^2 = {}'.format(np.linalg.norm(c_polyfit-c_lstsq)))
##
print('differences with c_lstsq')
print('||c_lstsq-c_pinv||^2 = {}'.format(np.linalg.norm(c_lstsq-c_pinv)))
print('||c_lstsq-c_polyfit||^2 = {}'.format(np.linalg.norm(c_lstsq-c_polyfit)))
print('||c_lstsq-c_lstsq||^2 = {}'.format(np.linalg.norm(c_lstsq-c_lstsq)))
##
print('Data set errors')
y_polyfit = np.dot(Kern,c_polyfit)
print('J_data(c_polyfit) = {}'.format(l2_loss(y_polyfit,Y)))
y_pinv = np.dot(Kern,c_pinv)
print('J_data(c_pinv) = {}'.format(l2_loss(y_pinv,Y)))
y_lstsq = np.dot(Kern,c_lstsq)
print('J_data(c_lstsq) = {}'.format(l2_loss(y_lstsq,Y)))
mit, dass ich es geschafft, zu bemerken, dass selten polyfit tut immer passt die Parameter, die pinv Anwendungen. Ich weiß, Pinv gibt die Pseudoinverse definitiv zurück, also denke ich, wenn mein Hauptziel ist, "sicherzustellen, dass ich die Pseudoinverse verwende", dann ist es eine gute Idee, np.pinv zu verwenden. Allerdings weiß ich auch mathematisch, dass die Pseudo-Inverse immer den kleinsten Fehlerquadrat J(c) = || Kc - y ||^2 egal was (Beweis here Theorem 11.1.2 Seite 446) minimiert. Daher sollte mein Ziel vielleicht sein, einfach die Python-Funktion zu verwenden, die den kleinsten Fehler der kleinsten Quadrate J zurückgibt.So lief ich (im unterbestimmt Fall) einen Vergleich der drei Verfahren
- Polygit,
np.polyfit - pseudo-inverse,
np.linalg.pinv - kleinsten Quadrate,
np.linalg.lstsq
und verglichen Welche Fehler Kleinste Quadrate Fehler sie gaben mir auf die Daten:
Dann besichtigte ich die seltsamen taucht die Funktion zu erleben scheint (die übrigens wie ein völliges Rätsel scheint, warum es taucht, wenn die Algorithmen nicht stochastische sind) und die Zahlen in der Regel kleiner war für polyfit, zum Beispiel:
lb,ub = (-100, 100)
Data set errors
J_data(c_polyfit) = 5.329753025633029e-12
J_data(c_pinv) = 0.06670557822873546
J_data(c_lstsq) = 0.7479733306782645
angesichts dieser Ergebnisse und dass Pseudo-Inverse ist ein Minimizer der kleinsten Quadrate, scheint es, dass das beste Ding ist zu ignorieren np.pinv. Ist das das Beste? Oder fehlt mir etwas Offensichtliches?
Als zusätzliche Note, die ich in polyfit code ging, um zu sehen, was genau macht es besser kleinsten Quadrate Fehler haben (was ich verwende als eine Möglichkeit, jetzt sein die beste Näherung für die pseudo-inverse zu sagen) und es scheint es einige seltsame Zustand/numerische Stabilität Code hat:
# scale lhs to improve condition number and solve
scale = NX.sqrt((lhs*lhs).sum(axis=0))
lhs /= scale
c, resids, rank, s = lstsq(lhs, rhs, rcond)
c = (c.T/scale).T # broadcast scale coefficients
was ich davon ausgehen, ist das, was die zusätzliche Stabilität für den polyfit bringt, dass pinv nicht hat?
Ist dies die richtige Entscheidung, polyfit für meine Aufgabe der linearen Systemannäherung des hohen Gradpolynoms zu verwenden?
auch an diesem Punkt bin ich bereit, andere Software wie Matlab zu verwenden, wenn es mir die richtige pseudo-inverse UND mehr numerische Stabilität (für die meisten Grad und alle Grenzen) zur Verfügung stellt.
Eine weitere zufällige Idee, die ich war gerade hatte, dass vielleicht gibt es eine schöne Art und Weise ist die Funktion zur Probe, so dass die Stabilität der pseudo-inversen gut ist. Meine Vermutung ist, dass mit einem Polynom einen Cosinus annähert, irgendeine Art von Anzahl von Proben oder dem Abstand zwischen ihnen erfordert (wie das Nyquist-Shannon-Abtasttheorem sagt, wenn die Basisfunktionen sinusoidals sind ...)
Es ist festzustellen heraus, dass wahrscheinlich invertieren (oder Pseudo ivnerting) und dann Multiplizieren ist eine schlechte Idee. Siehe:
https://www.johndcook.com/blog/2010/01/19/dont-invert-that-matrix/
, dass man nur spricht über inverse aber ich nehme an, es auch zu Umkehrungen Pseudo erstreckt.
jetzt ist meine Verwirrung, die in der Regel wollen wir die pseudo-inverse nicht wirklich berechnen explizit und tun A^+y=x_min_norm die minimale Normlösung zu erhalten.Allerdings hätte ich gedacht, dass np.lstsq die Antwort liefern würde, die ich wollte, aber sein Fehler unterscheidet sich auch stark von den anderen. Ich finde das extrem verwirrend ... lass mich denken, dass ich den falschen Weg benutze, um die minimale Normlösung in Python zu bekommen.
Ich versuche nicht, eine regulierte Lösung zu bekommen. Ich versuche, die minimale Normlösung und nichts anderes, so numerisch genau wie möglich zu bekommen.
SciPy Dokumente empfehlen SciPy Versionen von linearen Algebra-Routinen, wie [Pinv] (https: // docs. scipy.org/doc/scipy-0.16.1/reference/generated/scipy.linalg.pinv.html). Auch ist die Pseudoinvertierung im Prinzip nicht numerisch stabil; Pseudoinverse ist keine kontinuierliche Funktion seines Arguments. – FTP
@Desire hmm interessant, ich denke, es macht Sinn wegen der '1/singular_value ... aber wenn ich in meinem Beispiel bin von einem Kosinus und approximieren mit einem hohen Grad Polynom, was mich verwirrt ist warum es nicht * immer * voller Rang. Es scheint, dass ein (endliches) Polynom niemals in der Lage sein sollte, einen Kosinus zu approximieren. Selbst wenn der Pseudo-Invserse nicht kontinuierlich ist, ist das Problem vielleicht in der Art, wie ich sample? – Pinocchio
danke @Desire obwohl es scheint, dass es nicht wirklich einen Unterschied gemacht hat, scheint es 'polyfit' ist besser. Der "scipy pinv" erzeugt die gleiche orange Linie wie oben immer noch .... – Pinocchio