Passen Sie eine robuste Regressionslinie mit einem MM-Schätzer in R an

Kontext. Ich möchte eine Regressionslinie anpassen, um die Beziehung zwischen einer Antwortvariablen und einer kontinuierlichen Kovariate . Aufgrund des Vorhandenseins schlechter Hebelpunkte habe ich mich für einen MM-Schätzer anstelle des üblichen LS-Schätzers entschieden. $y$ $x$

Methodik. Grundsätzlich ist die MM-Schätzung eine M-Schätzung, die von einem S-Schätzer initialisiert wird. Daher müssen zwei Verlustfunktionen ausgewählt werden. Ich habe die weit verbreitete Verlustfunktion von Tukey Biweight gewählt

$\rho ( u ) = \left\{ \begin{array}{ll} 1 - \left[ 1 - \left( \tfrac{u}{k} \right)^{2} \right]^{3} & \textrm{if } | u | \leq k \\ 1 & \textrm{if } | u | > k, \end{array} \right.$

mit am vorläufigen S-Schätzer (was einen Durchschlagspunkt von ergibt ) und mit am M- (um eine Gaußsche Effizienz von zu gewährleisten ). $k = 1.548$ $50 \%$ $k = 2.697$ $70\%$

Ich möchte R verwenden, um meine robuste Regressionslinie anzupassen.

Frage.

library(MASS)
rlm(y~x, 
    method="MM",
    k0=1.548, c=2.697,
    maxit=50)

Stimmt mein Code mit dem vorherigen Absatz überein?
Würden Sie andere optionale Argumente verwenden?

BEARBEITEN. Nach meiner Diskussion mit @Jason Morgan stelle ich fest, dass mein vorheriger Code falsch ist. (@ Jason Morgan: Vielen Dank dafür!) Ich bin jedoch immer noch nicht von seinem Vorschlag überzeugt. Stattdessen schlage ich jetzt Folgendes vor:

library(robustbase)
lmrob(y~x, 
      tuning.chi=1.548, tuning.psi=2.697)

Ich denke, es bleibt jetzt bei der Methodik. Sind Sie einverstanden?

Vielen Dank!

r robust ocram
quelle

Antworten:

Standardmäßig wird in der Dokumentation angegeben, dass Gewichte rlmverwendet psi=psi.huberwerden. Wenn Sie also Tukeys Bisquadrat verwenden möchten, müssen Sie angeben psi=psi.bisquare. Die Standardeinstellungen sind psi.bisquare(u, c = 4.685, deriv = 0), die Sie nach Bedarf ändern können. Zum Beispiel möglicherweise so etwas

rlm(x ~ y, method="MM", psi=psi.bisquare, maxit=50)

Möglicherweise möchten Sie auch untersuchen, ob Sie init="lts"zum Initialisieren Ihrer Startwerte am wenigsten abgeschnittene Quadrate ( ) verwenden sollten. Standardmäßig werden die kleinsten Quadrate verwendet.

Jason Morgan
quelle

@Janson Morgan: Bist du sicher, was du vorbringst? Haben Sie Erfahrung mit dieser Funktion? In meiner Dokumentation (R 2.13.1) heißt es tatsächlich: "Der anfängliche Koeffizientensatz und die endgültige Skala werden von einem S-Schätzer mit k0 = 1,548 ausgewählt. Dies ergibt (für n >> p) einen Aufschlüsselungspunkt von 0,5. Der endgültige Schätzer ist ein M-Schätzer mit Tukeys Bigewicht und fester Skala, der diesen Aufschlüsselungspunkt erbt, vorausgesetzt c> k0; dies gilt für den Standardwert von c, der einer relativen Effizienz von 95% im Normalfall entspricht. "

Ocram

Ich habe diese Modelle in der Vergangenheit geschätzt. Wie in der Dokumentation angegeben, erfolgt der erste Schritt der MM-Schätzung mit Huber-Gewichten, der zweite mit Bisquadratgewichten. In meinen Notizen (von vor ein paar Jahren) heißt es, dass Sie im ersten S-Schritt bisquadratische Gewichte anstelle von Huber-Gewichten verwenden können, wenn Sie dies entsprechend angeben psi. Ich würde wahrscheinlich zunächst die cStandardeinstellung beibehalten (ich werde meine Antwort entsprechend ändern).

Jason Morgan

Ich benutze auch rlm und benutze die bisquare psi-Funktion wegen ihrer wieder absteigenden Eigenschaft. Manchmal gibt es jedoch Konvergenzprobleme, insbesondere bei kleineren Stichproben.

Jbowman