Bagaimana seseorang dapat secara empiris menunjukkan dalam R yang mana metode validasi silang yang setara dengan AIC dan BIC?

26

Dalam sebuah pertanyaan di tempat lain di situs ini, beberapa jawaban menyebutkan bahwa AIC setara dengan validasi silang leave-one-out (LOO) dan bahwa BIC setara dengan validasi silang K-fold. Apakah ada cara untuk secara empiris menunjukkan ini dalam R sehingga teknik yang terlibat dalam LOO dan K-fold diperjelas dan didemonstrasikan setara dengan nilai-nilai AIC dan BIC? Kode yang dikomentari dengan baik akan membantu dalam hal ini. Selain itu, dalam menunjukkan BIC silakan gunakan paket lme4. Lihat di bawah untuk dataset sampel ...

library(lme4) #for the BIC function

generate.data <- function(seed)
{
    set.seed(seed) #Set a seed so the results are consistent (I hope)
    a <- rnorm(60) #predictor
    b <- rnorm(60) #predictor
    c <- rnorm(60) #predictor
    y <- rnorm(60)*3.5+a+b #the outcome is really a function of predictor a and b but not predictor c
    data <- data.frame(y,a,b,c) 
    return(data)    
}

data <- generate.data(76)
good.model <- lm(y ~ a+b,data=data)
bad.model <- lm(y ~ a+b+c,data=data)
AIC(good.model)
BIC(logLik(good.model))
AIC(bad.model)
BIC(logLik(bad.model))

Per komentar sebelumnya, di bawah ini saya telah memberikan daftar benih dari 1 hingga 10.000 di mana AIC dan BIC tidak setuju. Ini dilakukan dengan pencarian sederhana melalui benih yang tersedia, tetapi jika seseorang dapat memberikan cara untuk menghasilkan data yang cenderung menghasilkan jawaban yang berbeda dari dua kriteria informasi ini, itu mungkin sangat informatif.

notable.seeds <- read.csv("http://student.ucr.edu/~rpier001/res.csv")$seed

Sebagai tambahan, saya berpikir tentang memesan benih-benih ini sejauh mana AIC dan BIC tidak setuju yang telah saya coba kuantifikasi sebagai jumlah dari perbedaan absolut antara AIC dan BIC. Sebagai contoh,

AICDiff <- AIC(bad.model) - AIC(good.model) 
BICDiff <- BIC(logLik(bad.model)) - BIC(logLik(good.model))
disagreement <- sum(abs(c(AICDiff,BICDiff)))

di mana metrik ketidaksepakatan saya hanya secara wajar berlaku ketika pengamatan penting. Sebagai contoh,

are.diff <- sum(sign(c(AICDiff,BICDiff)))
notable <- ifelse(are.diff == 0 & AICDiff != 0,TRUE,FALSE)

Namun dalam kasus di mana AIC dan BIC tidak setuju, nilai ketidaksepakatan yang dihitung selalu sama (dan merupakan fungsi dari ukuran sampel). Melihat kembali bagaimana AIC dan BIC dihitung, saya dapat melihat mengapa ini mungkin terjadi secara komputasi, tetapi saya tidak yakin mengapa itu akan menjadi kasus secara konseptual. Jika seseorang dapat menjelaskan masalah itu juga, saya akan menghargainya.

russellpierce
sumber
+1 Kode akan mudah ditulis, masih saya sangat tertarik melihat dataset yang jelas dan ilustratif.
Saya tidak yakin apa yang harus ada dalam dataset yang jelas dan ilustratif, tetapi saya telah berupaya memasukkan sampel dataset.
russellpierce
Jadi lihat: apa yang Anda berikan adalah contoh dari set yang tidak berguna, karena BIC dan AIC memberikan hasil yang sama: 340 v. 342 untuk AIC dan 349 v. 353 untuk BIC - sangat bagus. Model menang dalam kedua kasus. Seluruh gagasan dengan konvergensi itu adalah bahwa validasi silang tertentu akan memilih model yang sama dengan IC yang sesuai.
Saya telah membuat pemindaian sederhana dan misalnya untuk seed 76, IC tidak setuju.
1
Wow, ini adalah sesuatu yang akan semakin sulit saya dapatkan; Poin umum saya dalam seluruh diskusi adalah bahwa konvergensi teorema-teorema itu terlalu lemah sehingga perbedaan dapat muncul dari fluktuasi acak. (Dan itu tidak bekerja untuk pembelajaran mesin, tapi saya harap ini jelas.)

Jawaban:

5

Dalam upaya untuk menjawab sebagian pertanyaan saya sendiri, saya membaca deskripsi Wikipedia tentang validasi silang meninggalkan-satu-keluar

melibatkan menggunakan pengamatan tunggal dari sampel asli sebagai data validasi, dan pengamatan yang tersisa sebagai data pelatihan. Ini diulangi sehingga setiap observasi dalam sampel digunakan satu kali sebagai data validasi.

Dalam kode R, saya menduga bahwa itu akan berarti sesuatu seperti ini ...

resid <- rep(NA, Nobs) 
for (lcv in 1:Nobs)
    {
        data.loo <- data[-lcv,] #drop the data point that will be used for validation
        loo.model <- lm(y ~ a+b,data=data.loo) #construct a model without that data point
            resid[lcv] <- data[lcv,"y"] - (coef(loo.model)[1] + coef(loo.model)[2]*data[lcv,"a"]+coef(loo.model)[3]*data[lcv,"b"]) #compare the observed value to the value predicted by the loo model for each possible observation, and store that value
    }

... seharusnya menghasilkan nilai dalam residu yang terkait dengan AIC. Dalam prakteknya jumlah residu kuadrat dari setiap iterasi loop LOO yang dirinci di atas adalah prediktor yang baik dari AIC untuk yang terkenal. Benih, r ^ 2 = 0,9776. Tetapi, di tempat lain seorang kontributor menyarankan bahwa LOO harus setara asymptotically dengan AIC (setidaknya untuk model linier), jadi saya sedikit kecewa bahwa r ^ 2 tidak lebih dekat dengan 1. Jelas ini bukan jawaban - lebih seperti kode tambahan untuk mencoba mendorong seseorang untuk mencoba memberikan jawaban yang lebih baik.

Tambahan: Karena AIC dan BIC untuk model ukuran sampel tetap hanya bervariasi dengan konstanta, korelasi BIC dengan residu kuadrat adalah sama dengan korelasi AIC dengan residu kuadrat, sehingga pendekatan yang saya ambil di atas tampaknya tidak membuahkan hasil.

russellpierce
sumber
perhatikan bahwa ini akan menjadi jawaban yang diterima untuk hadiah (jika Anda tidak memilih jawaban, hadiah akan secara otomatis memilih jawaban dengan poin terbanyak)
robin girard
1
memberikan hadiah itu kepada saya sendiri tampaknya konyol - tetapi tidak ada orang lain yang mengirimkan jawaban.
Russelpierce