Menghitung interval prediksi untuk regresi logistik

20

Saya ingin memahami cara membuat interval prediksi untuk estimasi regresi logistik.

Saya disarankan untuk mengikuti prosedur dalam Pemodelan Data Biner Collett , 2nd Ed hal.98-99. Setelah menerapkan prosedur ini dan membandingkannya dengan R predict.glm, saya benar-benar berpikir buku ini menunjukkan prosedur untuk menghitung interval kepercayaan , bukan interval prediksi.

Implementasi prosedur dari Collett, dengan perbandingannya predict.glm, ditunjukkan di bawah ini.

Saya ingin tahu: bagaimana cara saya pergi dari sini untuk menghasilkan interval prediksi daripada interval kepercayaan?

#Derived from Collett 'Modelling Binary Data' 2nd Edition p.98-99
#Need reproducible "random" numbers.
seed <- 67

num.students <- 1000
which.student <- 1

#Generate data frame with made-up data from students:
set.seed(seed) #reset seed
v1 <- rbinom(num.students,1,0.7)
v2 <- rnorm(length(v1),0.7,0.3)
v3 <- rpois(length(v1),1)

#Create df representing students
students <- data.frame(
    intercept = rep(1,length(v1)),
    outcome = v1,
    score1 = v2,
    score2 = v3
)
print(head(students))

predict.and.append <- function(input){
    #Create a vanilla logistic model as a function of score1 and score2
    data.model <- glm(outcome ~ score1 + score2, data=input, family=binomial)

    #Calculate predictions and SE.fit with the R package's internal method
    # These are in logits.
    predictions <- as.data.frame(predict(data.model, se.fit=TRUE, type='link'))

    predictions$actual <- input$outcome
    predictions$lower <- plogis(predictions$fit - 1.96 * predictions$se.fit)
    predictions$prediction <- plogis(predictions$fit)
    predictions$upper <- plogis(predictions$fit + 1.96 * predictions$se.fit)


    return (list(data.model, predictions))
}

output <- predict.and.append(students)

data.model <- output[[1]]

#summary(data.model)

#Export vcov matrix 
model.vcov <- vcov(data.model)

# Now our goal is to reproduce 'predictions' and the se.fit manually using the vcov matrix
this.student.predictors <- as.matrix(students[which.student,c(1,3,4)])

#Prediction:
this.student.prediction <- sum(this.student.predictors * coef(data.model))
square.student <- t(this.student.predictors) %*% this.student.predictors
se.student <- sqrt(sum(model.vcov * square.student))

manual.prediction <- data.frame(lower = plogis(this.student.prediction - 1.96*se.student), 
    prediction = plogis(this.student.prediction), 
    upper = plogis(this.student.prediction + 1.96*se.student))

print("Data preview:")
print(head(students))
print(paste("Point estimate of the outcome probability for student", which.student,"(2.5%, point prediction, 97.5%) by Collett's procedure:"))
manual.prediction
print(paste("Point estimate of the outcome probability for student", which.student,"(2.5%, point prediction, 97.5%) by R's predict.glm:"))    
print(output[[2]][which.student,c('lower','prediction','upper')])
karbokation
sumber
Pertanyaan dasar, mengapa sqrt (jumlah (model.vcov * square.student)) dianggap sebagai kesalahan standar? Bukankah itu standar deviasi dan perlu dibagi dengan sqrt (n)? Jika demikian, n mana yang harus digunakan, n digunakan agar sesuai dengan model atau n dari bingkai data baru yang digunakan untuk memprediksi?
Rafael

Jawaban:

6

0<=y<=1

Greg Snow
sumber
6
Saya mencari interval prediksi 95% dari prediksi yang ada di ruang log-odds. Kemudian saya mengubahnya menjadi ruang probabilitas. Interval prediksi 100% tidak akan pernah menarik untuk prosedur apa pun, bukan? Misalnya, interval prediksi 100% untuk regresi linier akan mencakup -Inf ke Inf ... Bagaimanapun, seperti yang Anda lihat dalam kode saya, interval prediksi dihitung dalam ruang peluang log, yang kemudian diubah menjadi ruang probabilitas nanti . Jadi saya pikir pertanyaan saya tidak ada gunanya.
carbocation
2
Log-odds dapat dikonversi ke probabilitas dan Anda dapat menghitung interval kepercayaan pada probabilitas (atau log-odds). Tetapi interval prediksi ada pada variabel respons yaitu 0 atau 1. Jika hasil Anda bertahan dengan 0 = mati dan 1 = hidup, maka Anda dapat memprediksi probabilitas hidup untuk sekumpulan kovariat tertentu dan menghitung interval kepercayaan pada probabilitas itu. Tetapi hasilnya adalah 0/1, Anda tidak dapat memiliki pasien yang 62% hidup itu harus 0 atau 1, jadi satu-satunya interval prediksi yang mungkin adalah 0-0, 0-1, dan 1-1 (yang merupakan mengapa kebanyakan orang menempel pada interval kepercayaan diri).
Greg Snow
8
Jika Anda memiliki situasi di mana responsnya binomial (yang bisa berupa agregat 0-1s dalam kondisi yang sama), maka interval prediksi mungkin masuk akal.
Glen_b -Reinstate Monica
7
Regresi logistik adalah regresi probabilitas, mencoba memodelkan probabilitas beberapa peristiwa sebagai fungsi dari variabel regressor. Interval prediksi dalam pengaturan ini diambil sebagai interval pada skala probabilitas, atau skala log-odds, sehingga menghasilkan nada yang sempurna.
kjetil b halvorsen
2
@ Cesar, rumus interval prediksi diperoleh dengan mengasumsikan bahwa Y terdistribusi normal tentang garis, tetapi dalam regresi logistik kami tidak memiliki distribusi normal, kami memiliki Bernoulli atau Binomial. Menerapkan rumus pada halaman itu akan mengarah ke interval kepercayaan (sudah dapat melakukan ini) atau interval kepercayaan yang diperluas secara artifisial yang tidak memenuhi definisi interval prediksi (memprediksi hasil aktual pada skala hasil awal). Seperti yang disebutkan Glen_b, interval prediksi mungkin masuk akal jika hasilnya benar-benar binomial.
Greg Snow