베타 회귀를위한 DFFIT

Aug 26 2020

응답이 베타 분포를 따르는 GLM에 대한 DFFITS를 계산하려고합니다. betaregR 패키지 를 사용 합니다. 하지만이 패키지는 코드influence.measures() 를 사용하여 지원하지 않는다고 생각합니다.dffits()

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

그것은 양보한다

if (model $ rank == 0) {: 인수의 길이가 0 인 오류

저는 R의 초보자입니다.이 문제를 해결하는 방법을 모르겠습니다. 이 문제를 해결하는 데 도움이되거나 DFFIT를 수동으로 계산하는 방법에 대한 R 코드를 통해 몇 가지 팁을 제공합니다. 문안 인사

답변

2 jay.sf Aug 26 2020 at 15:12

dffits"betareg"개체에 대해 구현되지는 않지만 수동으로 계산할 수 있습니다.

이 Stack Overflow Q / A 에 따르면 다음 함수를 작성할 수 있습니다.

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)
  })
}

다음을 위해 잘 작동합니다 lm.

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

이제 시도해 보겠습니다 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 

나쁘지 않은 것 같습니다.

메모:

  • 이것이 코드에서 작동하더라도 통계 요구 사항을 충족하는지 확인해야합니다!
  • 내가 사용했던 suppressWarnings라인 5:6의 dffits1. 어떻게 든 predict(bfit, ReadingSkills)삭제 contrasts하지만 predict(bfit)그렇지 않습니다 (실제로 동일해야 함). 그러나 결과는 동일 all.equal(predict(bfit, ReadingSkills), predict(bfit))하므로 경고를 무시하는 것이 안전합니다.