DFFITs für die Beta-Regression

Aug 26 2020

Ich versuche, DFFITS für GLM zu berechnen, wobei die Antworten einer Beta-Verteilung folgen. Mit dem betaregR-Paket. Aber ich denke, dieses Paket wird nicht unterstützt, influence.measures()weil durch die Verwendung von dffits() Code

require(betareg)
df<-data("ReadingSkills")
y<-ReadingSkills$accuracy
n<-length(y)

bfit<-betareg(accuracy ~ dyslexia + iq, data = ReadingSkills)
DFFITS<-dffits(bfit, infl=influence(bfit, do.coef = FALSE))
DFFITS

es gibt nach

Fehler in if (Modell $ rank == 0) {: Argument hat die Länge Null

Ich bin ein Neuling in R. Ich weiß nicht, wie ich dieses Problem lösen soll. Bitte helfen Sie, dieses Problem zu lösen, oder geben Sie mir einige Tipps über den R-Code, wie DFFITs manuell berechnet werden. Grüße

Antworten

2 jay.sf Aug 26 2020 at 15:12

dffitssind nicht für "betareg"Objekte implementiert , aber Sie können versuchen, sie manuell zu berechnen.

Entsprechend dieser Frage / Antwort zum Stapelüberlauf könnten wir diese Funktion schreiben:

dffits1 <- function(x1, bres.type="response") {
  stopifnot(class(x1) %in% c("lm", "betareg"))
  sapply(1:length(x1$fitted.values), function(i) { x2 <- update(x1, data=x1$model[-i, ]) # leave one out
    h <- hatvalues(x1)
    nm <- rownames(x1$model[i, ]) num_dffits <- suppressWarnings(predict(x1, x1$model[i, ]) - 
                                     predict(x2, x1$model[i, ])) residx <- if (class(x1) == "betareg") { betareg:::residuals.betareg(x2, type=bres.type) } else { x2$residuals
    }
    denom_dffits <- sqrt(c(crossprod(residx)) / x2$df.residual*h[i])
    return(num_dffits / denom_dffits)
  })
}

Es funktioniert gut für lm:

fit <- lm(mpg ~ hp, mtcars)
dffits1(fit)
stopifnot(all.equal(dffits1(fit), dffits(fit)))

Versuchen wir jetzt betareg:

library(betareg)
data("ReadingSkills")

bfit <- betareg(accuracy ~ dyslexia + iq, data=ReadingSkills)
dffits1(bfit)
#           1           2           3           4           5           6           7 
# -0.07590185 -0.21862047 -0.03620530  0.07349169 -0.11344968 -0.39255172 -0.25739032 
#           8           9          10          11          12          13          14 
#  0.33722706  0.16606198  0.10427684  0.11949807  0.09932991  0.11545263  0.09889406 
#          15          16          17          18          19          20          21 
#  0.21732090  0.11545263 -0.34296030  0.09850239 -0.36810187  0.09824013  0.01513643 
#          22          23          24          25          26          27          28 
#  0.18635669 -0.31192106 -0.39038732  0.09862045 -0.10859676  0.04362528 -0.28811277 
#          29          30          31          32          33          34          35 
#  0.07951977  0.02734462 -0.08419156 -0.38471945 -0.43879762  0.28583882 -0.12650591 
#          36          37          38          39          40          41          42 
# -0.12072976 -0.01701615  0.38653773 -0.06440176  0.15768684  0.05629040  0.12134228 
#          43          44 
#  0.13347935  0.19670715 

Sieht nicht schlecht aus.

Anmerkungen:

  • Auch wenn dies im Code funktioniert, sollten Sie überprüfen, ob es Ihren statistischen Anforderungen entspricht!
  • Ich habe suppressWarningsin Zeilen 5:6von verwendet dffits1. predict(bfit, ReadingSkills)lässt das contrastsirgendwie fallen, wohingegen predict(bfit)nicht (sollte praktisch gleich sein). Die Ergebnisse sind jedoch identisch: all.equal(predict(bfit, ReadingSkills), predict(bfit))Daher ist das Ignorieren der Warnungen sicher.