Étant donné deux tableaux x et y, tous deux de longueur n, j'adapte un modèle y = a + b * x et je veux calculer un intervalle de confiance à 95% pour la pente. C'est (b - delta, b + delta) où b se trouve de la manière habituelle et
delta = qt(0.975,df=n-2)*se.slope
et se.slope est l'erreur standard dans la pente. Une façon d'obtenir l'erreur standard de la pente de R est summary(lm(y~x))$coef[2,2]
.
Supposons maintenant que j'écrive la probabilité de la pente donnée x et y, multipliez cela par un "plat" avant et utilisez une technique MCMC pour tirer un échantillon m de la distribution postérieure. Définir
lims = quantile(m,c(0.025,0.975))
Ma question: est-elle (lims[[2]]-lims[[1]])/2
approximativement égale au delta tel que défini ci-dessus?
Addendum Ci-dessous est un modèle JAGS simple où ces deux semblent différents.
model {
for (i in 1:N) {
y[i] ~ dnorm(mu[i], tau)
mu[i] <- a + b * x[i]
}
a ~ dnorm(0, .00001)
b ~ dnorm(0, .00001)
tau <- pow(sigma, -2)
sigma ~ dunif(0, 100)
}
Je lance ce qui suit dans R:
N <- 10
x <- 1:10
y <- c(30.5,40.6,20.5,59.1,52.5,
96.0,121.4,78.9,112.1,128.4)
lin <- lm(y~x)
#Calculate delta for a 95% confidence interval on the slope
delta.lm <- qt(0.975,df=N-2)*summary(lin)$coef[2,2]
library('rjags')
jags <- jags.model('example.bug', data = list('x' = x,'y' = y,'N' = N),
n.chains = 4,n.adapt = 100)
update(jags, 1000)
params <- jags.samples(jags,c('a', 'b', 'sigma'),7500)
lims <- quantile(params$b,c(0.025,0.975))
delta.bayes <- (lims[[2]]-lims[[1]])/2
cat("Classical confidence region: +/-",round(delta.lm, digits=4),"\n")
cat("Bayesian confidence region: +/-",round(delta.bayes,digits=4),"\n")
Et obtenir:
Région de confiance classique: +/- 4.6939
Région de confiance bayésienne: +/- 5.1605
En réexécutant cette opération plusieurs fois, la région de confiance bayésienne est toujours plus large que la région classique. Est-ce dû aux prieurs que j'ai choisis?
sigma <- pow(tau, -1/2)
ousigma <- 1/sqrt(tau)
Si vous échantillonnez du postérieur de b | y et calculer lims (comme vous le définissez), il doit être identique à (b - delta, b + delta). Plus précisément, si vous calculez la distribution postérieure de b | y sous un a priori plat, il est identique à la distribution d'échantillonnage classique de b.
Pour plus de détails, voir: Gelman et al. (2003). Analyse des données bayésiennes. CRC Press. Section 3.6
Éditer:
Ringold, le comportement que vous observez est conforme à l'idée bayésienne. L'intervalle crédible bayésien (IC) est généralement plus large que les intervalles classiques. Et la raison en est, comme vous l'avez bien deviné, que les hyperpriors ont pris en compte la variabilité à cause des paramètres inconnus.
Pour des scénarios simples comme ceux-ci (PAS EN GÉNÉRAL):
CI baysien> CI bayésien empirique> CI classique; > == plus large
la source
Pour les modèles gaussiens linéaires, il est préférable d'utiliser le package bayesm. Il implémente la famille semi-conjuguée de prieurs, et le prieur de Jeffreys est un cas limite de cette famille. Voir mon exemple ci-dessous. Ce sont des simulations classiques, il n'est pas nécessaire d'utiliser MCMC.
Je ne me souviens pas si les intervalles de crédibilité des paramètres de régression sont exactement les mêmes que les intervalles de confiance habituels des moindres carrés, mais en tout cas ils sont très proches.
la source
Étant donné que la régression linéaire simple est analytiquement identique entre l'analyse classique et bayésienne avec le précédent de Jeffrey, qui sont tous deux analytiques, il semble un peu étrange de recourir à une méthode numérique telle que MCMC pour faire l'analyse bayésienne. MCMC est juste un outil d'intégration numérique, qui permet d'utiliser des méthodes bayésiennes dans des problèmes plus complexes qui sont insolubles analytiquement, tout comme Newton-Rhapson ou Fisher Scoring sont des méthodes numériques pour résoudre des problèmes classiques qui sont insolubles.
La distribution postérieure p (b | y) utilisant le p de Jeffrey antérieur (a, b, s) proportionnel à 1 / s (où s est l'écart-type de l'erreur) est une distribution t de Student avec l'emplacement b_ols, l'échelle se_b_ols (" ols "pour" les moindres carrés ordinaires ") et n-2 degrés de liberté. Mais la distribution d'échantillonnage de b_ols est également un étudiant t avec l'emplacement b, l'échelle se_b_ols et n-2 degrés de liberté. Ainsi, ils sont identiques sauf que b et b_ols ont été échangés, donc quand il s'agit de créer l'intervalle, la "borne est + -" de l'intervalle de confiance est inversée en une "borne est + +" dans l'intervalle crédible.
Ainsi, l'intervalle de confiance et l'intervalle crédible sont analytiquement identiques, et peu importe la méthode utilisée (à condition qu'il n'y ait pas d'informations préalables supplémentaires) - alors prenez la méthode qui est moins coûteuse en termes de calcul (par exemple, celle avec moins d'inversions de matrice). Votre résultat avec MCMC montre que l'approximation particulière utilisée avec MCMC donne un intervalle crédible qui est trop large par rapport à l'intervalle crédible analytique exact. C'est probablement une bonne chose (bien que nous voudrions que l'approximation soit meilleure) que la solution bayésienne approximative apparaisse plus conservatrice que la solution bayésienne exacte.
la source