Gleichzeitige Heteroskedastizität und schwere Schwänze in einem Regressionsmodell

8

Ich versuche, mithilfe der Regression ein Vorhersagemodell zu erstellen. Dies ist das Diagnosediagramm für das Modell, das ich durch die Verwendung von lm () in R erhalte: diagnostische Diagramme von R.

Was ich aus dem QQ-Diagramm gelesen habe, ist, dass die Residuen eine starke Verteilung haben, und das Diagramm Residuen gegen Angepasst scheint darauf hinzudeuten, dass die Varianz der Residuen nicht konstant ist. Ich kann die schweren Schwänze der Residuen mit einem robusten Modell zähmen:

fitRobust = rlm(formula, method = "MM", data = myData)

Aber hier kommen die Dinge zum Stillstand. Das robuste Modell wiegt mehrere Punkte 0. Nachdem ich diese Punkte entfernt habe, sehen die Residuen und die angepassten Werte des robusten Modells folgendermaßen aus:Residuen vs Passend für das robuste Modell

Die Heteroskedastizität scheint immer noch da zu sein. Verwenden von

logtrans(model, alpha) 

α

rlm(formula, method = "MM") 

log(Y+α)X1++Xnα

Residuen vs Angepasst für logarithmisch transformierte Antwort

Es sieht für mich so aus, als ob die Residuen immer noch keine konstante Varianz haben. Ich habe andere Transformationen der Reaktion ausprobiert (einschließlich Box-Cox), aber sie scheinen auch keine Verbesserung zu sein. Ich bin mir nicht einmal sicher, ob die zweite Phase meiner Arbeit (dh das Finden einer Transformation der Antwort in einem robusten Modell) von irgendeiner Theorie unterstützt wird. Ich würde mich über Kommentare, Gedanken oder Vorschläge sehr freuen.

user765195
quelle
2
Ich denke, Sie sind etwas wählerisch in Bezug auf die nicht konstante Varianz. Es scheint mir in Ordnung zu sein. Was ist der Zweck der Regression? Erklärung / Hypothesentest oder Vorhersage?
Wahrscheinlichkeitslogik
@probabilityislogic, danke für deinen Kommentar. Ich schätze es sehr. Mein Ziel ist die Vorhersage. Du hast recht. Ich bin wahrscheinlich zu wählerisch. Gibt es ein Maß für Heteroskedastizität, das ich betrachten kann? Ich dachte daran, Varianz gegen angepasste Werte zu zeichnen, aber es gibt nicht viele Punkte für jeden vorhergesagten Wert, um die Varianz zu berechnen. Ich bin auch neugierig zu verstehen, was die Lösung für dieses Problem im Allgemeinen ist. Sind Box-Cox- und Log-Transformationen auch auf robuste Modelle anwendbar?
user765195
Sie können einen paarweisen Test auf Gleichheit der Varianzen durchführen, indem Sie den F-Test für ein Modell mit Gaußschen Fehlertermen verwenden. Wenn diese eine nicht-Gaußsche Verteilung aufweisen, gibt es robuste Dispersionstests wie den Levene-Test.
Michael R. Chernick
Vielen Dank an MichaelChernick. Ich freue mich sehr über Ihren Kommentar. Ich habe schließlich Koenkers Verallgemeinerung des Breusch-Pagan-Tests auf Heteroskedastizität verwendet, wie er im lmtest-Paket in R ( hosho.ees.hokudai.ac.jp/~kubo/Rdoc/library/lmtest/html/… ) implementiert ist .
user765195

Antworten:

3

Heteroskedastizität und Leptokurtose lassen sich bei der Datenanalyse leicht miteinander verbinden. Nehmen Sie ein Datenmodell, das einen Fehlerterm als Cauchy generiert. Dies erfüllt die Kriterien für Homoskedastizität. Die Cauchy-Verteilung weist eine unendliche Varianz auf. Ein Cauchy-Fehler ist die Methode eines Simulators, einen Ausreißer-Stichprobenprozess einzuschließen.

Mit diesen schwerwiegenden Fehlern führt der Ausreißer zu einem großen Residuum, selbst wenn Sie das richtige mittlere Modell anpassen. Ein Test der Heteroskedastizität hat den Typ-I-Fehler unter diesem Modell stark erhöht. Eine Cauchy-Verteilung hat auch einen Skalierungsparameter. Das Erzeugen von Fehlertermen mit einer linearen Vergrößerung des Maßstabs erzeugt heteroskedastische Daten, aber die Fähigkeit, solche Effekte zu erkennen, ist praktisch null, so dass der Typ-II-Fehler ebenfalls aufgeblasen wird.

Lassen Sie mich dann vorschlagen, dass der richtige datenanalytische Ansatz nicht darin besteht, sich in Tests zu verfangen. Statistische Tests sind in erster Linie irreführend. Nirgendwo ist dies offensichtlicher als bei Tests zur Überprüfung sekundärer Modellierungsannahmen. Sie sind kein Ersatz für gesunden Menschenverstand. Für Ihre Daten sehen Sie deutlich zwei große Residuen. Ihre Auswirkung auf den Trend ist minimal, da nur wenige Residuen in einer linearen Abweichung von der 0-Linie in der Darstellung der Residuen gegenüber der Anpassung versetzt sind. Das ist alles was Sie wissen müssen.

Was dann gewünscht wird, ist ein Mittel zum Schätzen eines flexiblen Varianzmodells, mit dem Sie Vorhersageintervalle über einen Bereich angepasster Antworten erstellen können. Interessanterweise ist dieser Ansatz in der Lage, die meisten gesunden Formen sowohl der Heteroskedastizität als auch der Kurtotis zu behandeln. Verwenden Sie dann einen Glättungs-Spline-Ansatz, um den mittleren quadratischen Fehler abzuschätzen.

Nehmen Sie das folgende Beispiel:

set.seed(123)
x <- sort(rexp(100))
y <- rcauchy(100, 10*x)

f <- lm(y ~ x)
abline(f, col='red')
p <- predict(f)
r <- residuals(f)^2

s <- smooth.spline(x=p, y=r)

phi <- p + 1.96*sqrt(s$y)
plo <- p - 1.96*sqrt(s$y)

par(mfrow=c(2,1))
plot(p, r, xlab='Fitted', ylab='Squared-residuals')
lines(s, col='red')
legend('topleft', lty=1, col='red', "predicted variance")

plot(x,y, ylim=range(c(plo, phi), na.rm=T))
abline(f, col='red')
lines(x, plo, col='red', lty=2)
lines(x, phi, col='red', lty=2)

Gibt das folgende Vorhersageintervall an, das sich "erweitert", um dem Ausreißer Rechnung zu tragen. Es ist immer noch ein konsistenter Schätzer der Varianz und sagt den Leuten sinnvollerweise: "Hey, es gibt diese große, wackelige Beobachtung um X = 4 und wir können dort keine sehr nützlichen Werte vorhersagen."

Geben Sie hier die Bildbeschreibung ein

AdamO
quelle
Würde dies für andere Arten von lms wie gls funktionieren?
user2974951