Ajuster une droite de régression robuste à l'aide d'un estimateur MM dans R

Le contexte. Je voudrais ajuster une droite de régression pour étudier la relation entre une variable de réponse et une covariable continue . En raison de la présence de mauvais points de levier, j'ai opté pour un estimateur MM au lieu de l'estimateur LS habituel.$y$$x$

Méthodologie. Fondamentalement, l'estimation MM est une estimation M initialisée par un estimateur S. Par conséquent, deux fonctions de perte doivent être choisies. J'ai choisi la fonction de perte de Tukey Biweight largement utilisée

$\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.$

avec à l'estimateur S préliminaire (qui donne un point de rupture égal à ), et avec à l'étape M-estimation (pour garantir une efficacité gaussienne de ).$k = 1.548$$50 \%$$k = 2.697$$70\%$

Je voudrais utiliser R pour adapter ma robuste ligne de régression.

Question.

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

• Mon code est-il conforme au paragraphe précédent?
• Souhaitez-vous utiliser d'autres arguments facultatifs?

ÉDITER. Suite à ma discussion avec @Jason Morgan, je me rends compte que mon code précédent est incorrect. (@Jason Morgan: Merci beaucoup pour cela!) Cependant, je ne suis toujours pas convaincu par sa proposition. Au lieu de cela, voici ce que je propose maintenant:

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

Je pense que ça colle à la méthodologie maintenant. Êtes-vous d'accord?

Merci!

Par défaut, la documentation indique que rlmutilise des psi=psi.huberpoids. Ainsi, si vous souhaitez utiliser le bisquare de Tukey, vous devez le spécifier psi=psi.bisquare. Les paramètres par défaut sont ceux psi.bisquare(u, c = 4.685, deriv = 0)que vous pouvez modifier à votre guise. Par exemple, peut-être quelque chose comme

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

Vous pouvez également rechercher si vous devez utiliser les carrés les moins coupés ( init="lts") pour initialiser vos valeurs de départ. La valeur par défaut est d'utiliser les moindres carrés.