Saya telah meninjau satu set makalah, masing-masing melaporkan rata-rata yang diamati dan SD pengukuran dalam masing-masing sampel dengan ukuran yang diketahui, . Saya ingin membuat tebakan terbaik tentang kemungkinan distribusi ukuran yang sama dalam studi baru yang saya rancang, dan seberapa besar ketidakpastian dalam tebakan itu. Saya senang menganggap ).
Pikiran pertama saya adalah meta-analisis, tetapi model biasanya menggunakan fokus pada estimasi titik dan interval kepercayaan yang sesuai. Namun, saya ingin mengatakan sesuatu tentang distribusi penuh , yang dalam hal ini juga termasuk membuat tebakan tentang varians, .
Saya telah membaca tentang kemungkinan pendekatan Bayeisan untuk memperkirakan set lengkap parameter dari distribusi yang diberikan mengingat pengetahuan sebelumnya. Ini secara umum lebih masuk akal bagi saya, tetapi saya tidak memiliki pengalaman dengan analisis Bayesian. Ini juga tampak seperti masalah sederhana dan relatif mudah untuk memotong gigi saya.
1) Mengingat masalah saya, pendekatan mana yang paling masuk akal dan mengapa? Meta-analisis atau pendekatan Bayesian?
2) Jika menurut Anda pendekatan Bayesian adalah yang terbaik, dapatkah Anda mengarahkan saya ke cara untuk mengimplementasikan ini (lebih disukai dalam R)?
EDIT:
Saya telah mencoba untuk menyelesaikan ini dengan apa yang saya pikir adalah cara Bayesian yang 'sederhana'.
Seperti yang saya nyatakan di atas, saya tidak hanya tertarik pada estimasi rata-rata, , tetapi juga varians, , mengingat informasi sebelumnya, yaitu
Sekali lagi, saya tidak tahu apa-apa tentang Bayeianisme dalam praktiknya, tetapi tidak butuh waktu lama untuk menemukan bahwa posterior dari distribusi normal dengan mean yang tidak diketahui dan varian memiliki solusi bentuk tertutup melalui konjugasi , dengan distribusi gamma normal-invers.
Masalahnya dirumuskan ulang sebagai .
diperkirakan dengan distribusi normal; dengan distribusi gamma terbalik.
Butuh beberapa saat untuk menyiasatinya, tetapi dari tautan ini ( 1 , 2 ) saya dapat, menurut saya, memilah cara melakukan ini di R.
Saya mulai dengan kerangka data yang terdiri dari satu baris untuk masing-masing 33 studi / sampel, dan kolom untuk rata-rata, varians, dan ukuran sampel. Saya menggunakan mean, varians, dan ukuran sampel dari studi pertama, di baris 1, sebagai informasi saya sebelumnya. Saya kemudian memperbarui ini dengan informasi dari studi berikutnya, menghitung parameter yang relevan, dan mengambil sampel dari gamma normal-terbalik untuk mendapatkan distribusi dan . Ini diulang sampai 33 studi telah dimasukkan.
# Loop start values values
i <- 2
k <- 1
# Results go here
muL <- list() # mean of the estimated mean distribution
varL <- list() # variance of the estimated mean distribution
nL <- list() # sample size
eVarL <- list() # mean of the estimated variance distribution
distL <- list() # sampling 10k times from the mean and variance distributions
# Priors, taken from the study in row 1 of the data frame
muPrior <- bayesDf[1, 14] # Starting mean
nPrior <- bayesDf[1, 10] # Starting sample size
varPrior <- bayesDf[1, 16]^2 # Starting variance
for (i in 2:nrow(bayesDf)){
# "New" Data, Sufficient Statistics needed for parameter estimation
muSamp <- bayesDf[i, 14] # mean
nSamp <- bayesDf[i, 10] # sample size
sumSqSamp <- bayesDf[i, 16]^2*(nSamp-1) # sum of squares (variance * (n-1))
# Posteriors
nPost <- nPrior + nSamp
muPost <- (nPrior * muPrior + nSamp * muSamp) / (nPost)
sPost <- (nPrior * varPrior) +
sumSqSamp +
((nPrior * nSamp) / (nPost)) * ((muSamp - muPrior)^2)
varPost <- sPost/nPost
bPost <- (nPrior * varPrior) +
sumSqSamp +
(nPrior * nSamp / (nPost)) * ((muPrior - muSamp)^2)
# Update
muPrior <- muPost
nPrior <- nPost
varPrior <- varPost
# Store
muL[[i]] <- muPost
varL[[i]] <- varPost
nL[[i]] <- nPost
eVarL[[i]] <- (bPost/2) / ((nPost/2) - 1)
# Sample
muDistL <- list()
varDistL <- list()
for (j in 1:10000){
varDistL[[j]] <- 1/rgamma(1, nPost/2, bPost/2)
v <- 1/rgamma(1, nPost/2, bPost/2)
muDistL[[j]] <- rnorm(1, muPost, v/nPost)
}
# Store
varDist <- do.call(rbind, varDistL)
muDist <- do.call(rbind, muDistL)
dist <- as.data.frame(cbind(varDist, muDist))
distL[[k]] <- dist
# Advance
k <- k+1
i <- i+1
}
var <- do.call(rbind, varL)
mu <- do.call(rbind, muL)
n <- do.call(rbind, nL)
eVar <- do.call(rbind, eVarL)
normsDf <- as.data.frame(cbind(mu, var, eVar, n))
colnames(seDf) <- c("mu", "var", "evar", "n")
normsDf$order <- c(1:33)
Berikut adalah diagram jalur yang menunjukkan bagaimana dan berubah saat setiap sampel baru ditambahkan.
Berikut adalah desnitas berdasarkan pengambilan sampel dari estimasi distribusi untuk rata-rata dan varians pada setiap pembaruan.
Saya hanya ingin menambahkan ini jika bermanfaat bagi orang lain, dan agar orang yang tahu bisa memberi tahu saya apakah ini masuk akal, cacat, dll.
sumber
sumber