Получение доверительных интервалов для робастного коэффициента регрессии (MASS::rlm)

Question

Есть ли какой-нибудь возможный способ получить 95% ДИ для коэффициентов регрессии из надежной регрессии, как это реализовано в MASS::rlm?

# libraries needed
library(MASS)
library(stats)
library(datasets)

# running robust regression
(x <-
  MASS::rlm(formula = scale(Sepal.Length) ~ scale(Sepal.Width),
            data = iris))
#> Call:
#> rlm(formula = scale(Sepal.Length) ~ scale(Sepal.Width), data = iris)
#> Converged in 5 iterations
#> 
#> Coefficients:
#>        (Intercept) scale(Sepal.Width) 
#>        -0.03728607        -0.14343268 
#> 
#> Degrees of freedom: 150 total; 148 residual
#> Scale estimate: 1.06

# getting confidence interval for the regression coefficient
stats::confint(object = x,
               parm = "scale(Sepal.Width)",
               level = 0.95)
#>                    2.5 % 97.5 %
#> scale(Sepal.Width)    NA     NA

Answer 1

Явный вызов confint.default кажется, дает хорошие результаты, как это:

confint.default(object = x, parm = "scale(Sepal.Width)", level = 0.95)

#                        2.5 %     97.5 %
#scale(Sepal.Width) -0.3058138 0.01894847

---

редактировать

confint использует метод confint.lm когда это пройдено x так как x имеет класс lm (так же как rlm). призвание confint.default явно избегает этого. Эти две функции отличаются только одной строкой кода, как показано ниже:

confint.lm

fac <- qt(a, object$df.residual)

confint.default

fac <- qnorm(a)

Проблема в том, что x$df.residual является NA и следовательно, qt(a, object$df.residual) производит NA в то время как qnorm(a) не имеет этой проблемы

Answer 2

Опоздали на вечеринку, но учтите, что ДИ из нормального распределения будут иметь более низкое, чем ожидалось, покрытие для небольших размеров выборки.

Остаточные степени свободы для объекта rlm можно получить из

      summary(x)$df[2]  # see code for MASS:::summary.rlm

Чтобы написать собственный метод confint для результатов rlm, назначьте df слоту df.residual, а затем вызовите confint.lm:

      confint.rlm <- function(object, ...){
  object$df.residual <- MASS:::summary.rlm(object)$df[2]
  confint.lm(object, ...)
}

Теперь confint ведет себя так, как ожидалось, и также основан на t Стьюдента:

      confint(x)
                        2.5 %     97.5 %
(Intercept)        -0.2004593 0.12588715
scale(Sepal.Width) -0.3071526 0.02028719