DFFITs cho hồi quy beta

Question

Tôi đang cố gắng tính toán DFFITS cho GLM, trong đó các phản hồi tuân theo phân phối Beta. Bằng cách sử dụng betareggói R. Nhưng tôi nghĩ rằng gói này không hỗ trợ influence.measures()vì bằng cách sử dụng dffits() Mã

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

nó mang lại

Lỗi trong if (model $ rank == 0) {: đối số có độ dài bằng 0

Tôi là một người mới trong R. Tôi không biết làm thế nào để giải quyết vấn đề này. Vui lòng giúp đỡ để giải quyết vấn đề này hoặc cho tôi một số mẹo thông qua mã R về cách tính toán DFFITs theo cách thủ công. Trân trọng

Answer 1

dffitskhông được triển khai cho "betareg"các đối tượng, nhưng bạn có thể thử tính toán chúng theo cách thủ công.

Theo Q / A Stack Overflow này, chúng ta có thể viết hàm này:

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

Nó hoạt động tốt cho lm:

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

Bây giờ chúng ta hãy thử 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 

Trông không tệ.

Ghi chú:

• Ngay cả khi điều này hoạt động trong mã, bạn nên kiểm tra xem nó có đáp ứng các yêu cầu thống kê của bạn hay không!
• Tôi đã sử dụng suppressWarningstrong các dòng 5:6của dffits1. predict(bfit, ReadingSkills)giảm contrastsbằng cách nào đó, trong khi predict(bfit)không (thực tế phải giống nhau). Tuy nhiên, kết quả là giống hệt nhau: all.equal(predict(bfit, ReadingSkills), predict(bfit))do đó việc bỏ qua các cảnh báo vẫn an toàn.