Wenn Sie es nicht orthogonal machen können, machen Sie es roh (Polynomregression)

11

Bei der Durchführung einer Polynomregression für Y auf X Menschen manchmal rohe Polynome, manchmal orthogonale Polynome. Aber wenn sie verwenden, was völlig willkürlich erscheint.

Hier und hier werden rohe Polynome verwendet. Aber hier und hier scheinen orthogonale Polynome die richtigen Ergebnisse zu liefern. Was, wie, warum?!

Im Gegensatz dazu, wenn man aus einem Lehrbuch etwas über die Polynomregression lernt (z ISLR ) nicht einmal rohe oder orthogonale Polynome erwähnt - es wird nur das anzupassende Modell angegeben.

Wann müssen wir also was verwenden?
Und warum sind die einzelnen p-Werte für , X 2XX2 usw. stark zwischen diesen beiden Werten?

l7ll7
quelle
1
Sie sollten sich überlegen, welche p-Werte unterschiedlich sind, wenn Sie dasselbe Modell mit rohen und orthogonalen Polynomen und deren Interpretation an dieselben Daten anpassen. Was ist mit den Modellvorhersagen?
Scortchi - Monica wieder einsetzen
@ Scortchi Ich habe die relevanten Informationen zu meiner Frage hinzugefügt.
17ll7
4
Ein weiterer guter Grund für die Verwendung orthogonaler Polynome ist die numerische Stabilität. Die zugehörige Entwurfsmatrix für die Anpassung an die Monombasis kann für eine Anpassung mit hohem Grad ziemlich schlecht konditioniert sein, da die Monome höherer Ordnung "sehr nahezu linear abhängig" sind (ein Konzept, das mathematisch präziser gemacht werden könnte), während die Entwurfsmatrix denn orthogonale Polynome verhalten sich etwas besser. Ich besprach die im gleichen Abstand abscissas (Gram) Fall hier , aber der Deal ist ähnlich in dem nicht-äquidistanten Fall.
JM ist kein Statistiker
(Trotzdem sollte man ohne guten Grund nicht zu hochgradigen Polynomen passen.)
JM ist kein Statistiker

Antworten:

7

Die Variablen und X 2 sind nicht linear unabhängig. Selbst wenn es keinen quadratischen Effekt gibt, ändert das Hinzufügen von X 2 zum Modell den geschätzten Effekt von X.XX2X2X .

Lassen Sie uns einen Blick mit einer sehr einfachen Simulation werfen.

> x <- runif(1e3)
> y <- x + rnorm(length(x))
> summary(lm(y~x))

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -0.03486    0.06233  -0.559    0.576    
x            1.05843    0.10755   9.841   <2e-16 ***

Jetzt mit einem quadratischen Term im Modell passen.

> summary(lm(y~x+I(x^2)))

Coefficients:
            Estimate Std. Error t value Pr(>|t|)
(Intercept)  0.03275    0.09528   0.344    0.731
x            0.65742    0.44068   1.492    0.136
I(x^2)       0.39914    0.42537   0.938    0.348

Natürlich ist der Omnibus-Test immer noch von Bedeutung, aber ich denke, das Ergebnis, das wir suchen, ist nicht dieses. Die Lösung besteht darin, orthogonale Polynome zu verwenden.

 > summary(lm(y~poly(x,2)))

 Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.49744    0.03098  16.059   <2e-16 ***
poly(x, 2)1  9.63943    0.97954   9.841   <2e-16 ***
poly(x, 2)2  0.91916    0.97954   0.938    0.348    

Beachten Sie, dass die Koeffizienten von xim ersten Modell und von poly(x,2)1im zweiten Modell nicht gleich sind und sogar die Abschnitte unterschiedlich sind. Dies liegt daran, dass polyorthonormale Vektoren geliefert werden, die auch orthogonal zum Vektor sind rep(1, length(x)). Also poly(x,2)1nicht x, sondern (x -mean(x))/sqrt(sum((x-mean(x))**2))...

Ein wichtiger Punkt ist, dass die Wald-Tests in diesem letzten Modell unabhängig sind. Sie können orthogonale Polynome verwenden, um zu entscheiden, bis zu welchem ​​Grad Sie gehen möchten, indem Sie sich den Wald-Test ansehen: Hier entscheiden Sie sich, aber nicht X 2 beizubehaltenXX2 . Natürlich würden Sie das gleiche Modell finden, wenn Sie die ersten beiden Modelle vergleichen, aber auf diese Weise ist es einfacher - wenn Sie überlegen, in höhere Grade aufzusteigen, ist es wirklich viel einfacher.

Sobald Sie entschieden haben, welche Begriffe beibehalten werden sollen, möchten Sie möglicherweise zur Interpretierbarkeit oder zur Vorhersage zu den Rohpolynomen und X 2 zurückkehren .XX2

Elvis
quelle
+1 Endlich eine klare Antwort! Vielen Dank! Können Sie mir bitte sagen, bevor ich akzeptierte, gibt es noch andere Statistiken wie R ^ 2 oder die F-Statistik, die ich besser als die rohe Zusammenfassung der orthogonalen Darstellung lesen sollte? Ist die Anpassung unter Verwendung von Rohpolynomen neben dem Zeichnen der Variablen für alles andere in diesem Szenario gut?
17ll7
Und wenn ich mehrere Prädiktoren habe, gilt das auch?
17ll7
Wie würden Sie "orthogonale Polynome verwenden, um zu entscheiden, ob Sie einen quadratischen Term einfügen möchten oder nicht"?
Scortchi - Monica wieder einsetzen
1
Der Punkt ist, dass der Test des Effekts höchster Ordnung, in diesem Fall des Quadrats, der gleiche ist, unabhängig davon, ob Sie rohe oder orthogonale Polynome verwenden. Warum also mit orthogonalen Polynomen arbeiten?
Scortchi - Monica wieder einsetzen
4
Natürlich sollten Sie diese Randtests in diesem Modell einfach nicht durchführen. Sie sollten nach dem Verwerfen des Effekts höchster Ordnung erneut anpassen. Orthogonale Polynome ersparen Ihnen die Mühe und ermöglichen ein einfaches Herabsetzen - vielleicht könnten Sie dies mit einem kubischen Begriff veranschaulichen.
Scortchi - Monica wieder einsetzen
3

Um eine naive Einschätzung der Situation zu geben:

allgemein: Angenommen, Sie haben zwei verschiedene Basisfunktionssysteme sowie { ˜ p } n = 1 für einen Funktions- (Hilbert-) Raum, übliches L 2 ( [ a , b ] ) dh der Raum aller quadratisch integrierbaren Funktionen.{pn}n=1{p~}n=1L2([a,b])

Dies bedeutet, dass jede der beiden Basen verwendet werden kann, um jedes Element von zu erklären , dh für y L 2 ( [ a , b ] ) haben Sie für einige Koeffizienten θ n und ˜ θ nR , n = 1 , 2 , ... (in der L 2 -SENSE): & Sigma; n = 1 ~ θ nL2([a,b])yL2([a,b])θnθ~nRn=1,2,L2

n=1θ~np~n=y=n=1θnpn.

k<

{pn}n=1k
{p~}n=1k,
L2([a,b])

{p~}n=1, is just an orthogonalization of the other basis, {pn}n=1, the overall prediction of y will be the same for each truncated model ({p}n=1k and their orthogonalized counterpart will describe the same k-dimensional subspace of L2([a,b])).

But each individual basis function from the two "different" bases will yield a different contribution to this predcition (obviously as the functions/predictors are different!) resulting in different p-values and coefficients.

Hence, in terms of prediction there is (in this case) no difference.

From a computational point of view a model matrix consisting of orthogonal basis functions have nice numerical/computational properties for the least squares estimator. While at the same time from the statistical point of view, the orthogonalization results in uncorrelated estimates, since var(θ~^)=Iσ² under the standard assumptions.


The natural question arises if there is a best truncated basis system. However the answer to the question is neither simple nor unique and depends for example on the definition of the word "best", i.e. what you are trying to archive.

chRrr
quelle
1
(+1) No difference in terms of prediction; & it might be said no difference in terms of any meaningful inference.
Scortchi - Reinstate Monica