Mengevaluasi distribusi prediksi posterior dalam regresi linear Bayesian

10

Saya bingung bagaimana cara mengevaluasi distribusi prediksi posterior untuk regresi linier Bayesian, melewati kasus dasar yang dijelaskan di sini pada halaman 3, dan disalin di bawah ini.

p(y~y)=p(y~β,σ2)p(β,σ2y)

Kasus dasar adalah model regresi linier ini:

y=Xβ+ϵ,yN(Xβ,σ2)

Jika kita menggunakan seragam sebelum , dengan skala-Inv sebelum , ATAU gamma normal-terbalik sebelumnya (lihat di sini ) distribusi prediktif posterior analitik dan adalah siswa t. χ 2 σ 2βχ2σ2

Bagaimana dengan model ini?

y=Xβ+ϵ,yN(Xβ,Σ)

Ketika , tetapi diketahui, distribusi prediktif posterior adalah Gaussian multivarian. Biasanya, Anda tidak tahu , tetapi harus memperkirakannya. Mungkin Anda mengatakan diagonal dan menjadikan diagonal fungsi kovariat. Ini dibahas dalam bab regresi linier Analisis Data Bayesian Gelman .Σ ΣyN(Xβ,Σ)ΣΣ

Apakah ada bentuk analitik untuk distribusi prediksi posterior dalam kasus ini? Bisakah saya tancapkan perkiraan saya ke dalam t multivariate student? Jika Anda memperkirakan lebih dari satu varian, apakah distribusinya masih multivarian t?

Saya bertanya karena katakanlah saya memiliki beberapa sudah di tangan. Saya ingin tahu apakah itu lebih mungkin telah diprediksi oleh misalnya regresi linier A, regresi linier B y~

bill_e
sumber
1
Jika Anda memiliki sampel posterior dari distribusi posterior Anda dapat mengevaluasi distribusi prediktif dengan perkiraan Monte Carlo
niandra82
Ah terima kasih, aku selalu bisa melakukan itu. Apakah tidak ada rumus analitik dalam kasus ini?
bill_e
Omong-omong, tautannya rusak. Akan lebih bagus jika Anda memasukkan referensi dengan cara lain.
Maxim.K

Jawaban:

6

Jika Anda menganggap seragam sebelum , maka posterior untuk adalah dengan Untuk menemukan distribusi prediktif, kami memerlukan informasi lebih lanjut. Jika dan tidak tergantung pada diberikan , maka Tetapi biasanya untuk jenis model ini, dan tidak independen secara kondisional, sebaliknya, biasanya kita miliki β β | y ~ N ( β , V β ) . Β = [ X ' Σ - 1 X ] X ' yββ

β|yN(β^,Vβ).
β^=[XΣ-1X]XydanVβ=[XΣ-1X]-1.
y~N(X~β,Σ~)yβ
y~|yN(X~β^,Σ~+Vβ).
yy~˜y| y~N(~Xβ+Σ21Σ-1(y-Xβ),~Σ-Σ21Σ-1Σ12). Σ,Σ12,˜Σ
(yy~)N([XβX~β],[ΣΣ12Σ21Σ~]).
Jika demikian, maka Ini mengasumsikan dan semuanya diketahui. Seperti yang Anda tunjukkan, biasanya mereka tidak dikenal dan perlu diperkirakan. Untuk model umum yang memiliki struktur ini, misalnya seri waktu dan model spasial, umumnya tidak akan ada bentuk tertutup untuk distribusi prediktif.
y~|yN(X~β^+Σ21Σ-1(y-Xβ^),Σ~-Σ21Σ-1Σ12).
Σ,Σ12,Σ~
jaradniemi
sumber
2

Di bawah prior-Wishart non-informatif atau multivariat, Anda memiliki bentuk analitis sebagai distribusi Siswa multivariat, untuk mutlivariat klasik, regresi berganda. Saya kira perkembangan ini dokumen yang terkait dengan pertanyaan Anda (Anda mungkin seperti Lampiran A :-)). Saya biasanya membandingkan hasilnya dengan distribusi prediksi posterior yang diperoleh dengan menggunakan WinBUGS dan bentuk analitis: mereka persis sama. Masalahnya hanya menjadi sulit ketika Anda memiliki efek acak tambahan dalam model efek campuran, terutama dalam desain yang tidak seimbang.

Secara umum, dengan regresi klasik, y dan ỹ bersyarat independen (residual iid)! Tentu saja jika bukan itu masalahnya, maka solusi yang diusulkan di sini tidak benar.

Di R, (di sini, solusi untuk prior uniform), dengan asumsi Anda membuat model lm (bernama "model") dari salah satu respons dalam model Anda, dan menyebutnya "model", berikut adalah cara mendapatkan distribusi prediksi multivarian

library(mvtnorm)
Y = as.matrix(datas[,c("resp1","resp2","resp3")])
X =  model.matrix(delete.response(terms(model)), 
           data, model$contrasts)
XprimeX  = t(X) %*% X
XprimeXinv = solve(xprimex)
hatB =  xprimexinv %*% t(X) %*% Y
A = t(Y - X%*%hatB)%*% (Y-X%*%hatB)
F = ncol(X)
M = ncol(Y)
N = nrow(Y)
nu= N-(M+F)+1 #nu must be positive
C_1 =  c(1  + x0 %*% xprimexinv %*% t(x0)) #for a prediction of the factor setting x0 (a vector of size F=ncol(X))
varY = A/(nu) 
postmean = x0 %*% hatB
nsim = 2000
ysim = rmvt(n=nsim,delta=postmux0,C_1*varY,df=nu) 

Sekarang, kuantil dari ysim adalah interval toleransi beta-ekspektasi dari distribusi prediktif, Anda tentu saja dapat langsung menggunakan distribusi sampel untuk melakukan apa pun yang Anda inginkan.

Pierre Lebrun
sumber