Tag Archives: confidence

What is a Linear Trend, by the way?

01/01/2017 Arthur Charpentier 16 Comments

I had a very strange discussion on twitter (yes, another one), about regression curves. I think it started with a tweet based on some xkcd picture (just for fun, because it was New Year’s Day)

“don’t trust linear regressions” https://t.co/exUCvyRd1G pic.twitter.com/O6rBJfkULa

— Arthur Charpentier (@freakonometrics) 1 janvier 2017

There were comments on that picture, by econometricians, mainly about ‘significant’ trends when datasets are very noisy. And I mentioned a graph that I saw earlier, a couple of days ago

@AndyHarless @mileskimball actually, all that reminds me of a post by @RogerPielkeJr earlier (not a big fan of the regression line) pic.twitter.com/NQgzgVsBcE

— Arthur Charpentier (@freakonometrics) 1 janvier 2017

Let us reproduce that graph (Roger kindly sent me the dataset)

db=data.frame(year=1990:2016, ratio=c(.23,.27,.32,.37,.22,.26,.29,.15,.40,.28,.14,.09,.24,.18,.29,.51,.13,.17,.25,.13,.21,.29,.25,.2,.15,.12,.12)) library(ggplot2)

The graph is here (with the same aesthetic conventions as Roger’s initial graph, i.e. using some sort of barplot)

ggplot(db, aes(year, ratio)) + geom_bar(stat="identity") + stat_smooth(method = "lm", se = FALSE)

My point was that we miss the ‘confidence band’ of the regression

@freakonometrics @AndyHarless @mileskimball Because it is not a sample. Since 1990 weather losses/global GDP have gone down.

— Roger Pielke Jr. (@RogerPielkeJr) 1 janvier 2017

In R, at least, it is quite natural to get (and actually, it is the default version of the graph function)

ggplot(db, aes(year, ratio)) + geom_bar(stat="identity") + stat_smooth(method = "lm", se = TRUE)

It is hard to claim that the ‘regression line’ is significant (in the sense “significantly non horizontal”). To be more specific, if we look at the output of the regression model, we get

summary(lm(ratio~year,data=db))

Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 9.158531 4.549672 2.013 0.055 . year -0.004457 0.002271 -1.962 0.061 . --- Signif. codes: 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(which is exactly what Roger used in his graph to plot his red straight line). The p-value of the estimator of the slope, in a linear regression model is here 6%. But I found Roger’s point puzzling

@freakonometrics @AndyHarless @mileskimball Disagree. U can create one, of course, but doesnt mean much. These data are not balls from urns.

— Roger Pielke Jr. (@RogerPielkeJr) 1 janvier 2017

Confidence vs. Credibility Intervals

26/11/2014 Arthur Charpentier 13 Comments

Tomorrow, for the final lecture of the Mathematical Statistics course, I will try to illustrate – using Monte Carlo simulations – the difference between classical statistics, and the Bayesien approach.

The (simple) way I see it is the following,

for frequentists, a probability is a measure of the the frequency of repeated events, so the interpretation is that parameters are fixed (but unknown), and data are random
for Bayesians, a probability is a measure of the degree of certainty about values, so the interpretation is that parameters are random and data are fixed

Or to quote Frequentism and Bayesianism: A Python-driven Primer, a Bayesian statistician would say “given our observed data, there is a 95% probability that the true value of $\theta$ falls within the credible region” while a Frequentist statistician would say “there is a 95% probability that when I compute a confidence interval from data of this sort, the true value of $\theta$ will fall within it”.

To get more intuition about those quotes, consider a simple problem, with Bernoulli trials, with insurance claims. We want to derive some confidence interval for the probability to claim a loss. There were $https://latex.codecogs.com/gif.latex?n$ = 1047 policies. And 159 claims.

Consider the standard (frequentist) confidence interval. What does that mean that $\overline{x}\pm\sqrt{\frac{\overline{x}(1-\overline{x})}{n}}$ is the (asymptotic) 95% confidence interval? The way I see it is very simple. Let us generate some samples, of size $n$ , with the same probability as the empirical one, i.e. $\widehat{\theta}$ (which is the meaning of “from data of this sort”). For each sample, compute the confidence interval with the relationship above. It is a 95% confidence interval because in 95% of the scenarios, the empirical value lies in the confidence interval. From a computation point of view, it is the following idea,

> xbar <- 159
> n <- 1047
> ns <- 100
> M=matrix(rbinom(n*ns,size=1,prob=xbar/n),nrow=n)

I generate 100 samples of size $https://latex.codecogs.com/gif.latex?n$ . For each sample, I compute the mean, and the confidence interval, from the previous relationship

> fIC=function(x) mean(x)+c(-1,1)*1.96*sqrt(mean(x)*(1-mean(x)))/sqrt(n)
> IC=t(apply(M,2,fIC))
> MN=apply(M,2,mean)

Then we plot all those confidence intervals. In red when they do not contain the empirical mean

> k=(xbar/n<IC[,1])|(xbar/n>IC[,2])
> plot(MN,1:ns,xlim=range(IC),axes=FALSE,
+ xlab="",ylab="",pch=19,cex=.7,
+ col=c("blue","red")[1+k])
> axis(1)
> segments(IC[,1],1:ns,IC[,2],1:
+ ns,col=c("blue","red")[1+k])
> abline(v=xbar/n)

Now, what about the Bayesian credible interval ? Assume that the prior distribution for the probability to claim a loss has a $https://latex.codecogs.com/gif.latex?\mathcal{B}(\alpha,\beta)$ distribution. We’ve seen in the course that, since the Beta distribution is the conjugate of the Bernoulli one, the posterior distribution will also be Beta. More precisely

$https://latex.codecogs.com/gif.latex?\mathcal{B}\left(\alpha+\sum%20x_i,\beta+n-\sum%20x_i\right)$

Based on that property, the confidence interval is based on quantiles of that (posterior) distribution

> u=seq(.1,.2,length=501)
> v=dbeta(u,1+xbar,1+n-xbar)
> plot(u,v,axes=FALSE,type="l")
> I=u<qbeta(.025,1+xbar,1+n-xbar)
> polygon(c(u[I],rev(u[I])),c(v[I],
+ rep(0,sum(I))),col="red",density=30,border=NA)
> I=u>qbeta(.975,1+xbar,1+n-xbar)
> polygon(c(u[I],rev(u[I])),c(v[I],
+ rep(0,sum(I))),col="red",density=30,border=NA)
> axis(1)

What does that mean, here, that we have a 95% credible interval. Well, this time, we do not draw using the empirical mean, but some possible probability, based on that posterior distribution (given the observations)

> pk <- rbeta(ns,1+xbar,1+n-xbar)

In green, below, we can visualize the histogram of those values

> hist(pk,prob=TRUE,col="light green",
+ border="white",axes=FALSE,
+ main="",xlab="",ylab="",lwd=3,xlim=c(.12,.18))

And here again, let us generate samples, and compute the empirical probabilities,

> M=matrix(rbinom(n*ns,size=1,prob=rep(pk,
+ each=n)),nrow=n)
> MN=apply(M,2,mean)

Here, there is 95% chance that those empirical means lie in the credible interval, defined using quantiles of the posterior distribution. We can actually visualize all those means : in black the mean used to generate the sample, and then, in blue or red, the averages obtained on those simulated samples,

> abline(v=qbeta(c(.025,.975),1+xbar,1+
+ n-xbar),col="red",lty=2)
> points(pk,seq(1,40,length=ns),pch=19,cex=.7)
> k=(MN<qbeta(.025,1+xbar,1+n-xbar))|
+ (MN>qbeta(.975,1+xbar,1+n-xbar))
> points(MN,seq(1,40,length=ns),
+ pch=19,cex=.7,col=c("blue","red")[1+k])
> segments(MN,seq(1,40,length=ns),
+ pk,seq(1,40,length=ns),col="grey")

More details and exemple on Bayesian statistics, seen with the eyes of a (probably) not Bayesian statistician in my slides, from my talk in London, last Summer,

2011, ACT2040, Computer, GLM, Ratemaking, Regression

Confidence interval for predictions with GLMs

04/11/2011 Arthur Charpentier 5 Comments

Consider a (simple) Poisson regression . Given a sample where , the goal is to derive a 95% confidence interval for given , where is the prediction. Hence, we want to derive a confidence interval for the prediction, not the potential observation, i.e. the dot on the graph below

> r=glm(dist~speed,data=cars,family=poisson)
> P=predict(r,type="response",
+ newdata=data.frame(speed=seq(-1,35,by=.2)))
> plot(cars,xlim=c(0,31),ylim=c(0,170))
> abline(v=30,lty=2)
> lines(seq(-1,35,by=.2),P,lwd=2,col="red")
> P0=predict(r,type="response",se.fit=TRUE,
+ newdata=data.frame(speed=30))
> points(30,P1$fit,pch=4,lwd=3)

i.e.

Let denote the maximum likelihood estimator of . Then

where is Fisher information of (from standard maximum likelihood theory). Recall that

where computation of those values is based on the following calculations

In the case of the log-Poisson regression

Let us get back to our initial problem.

confidence interval for the linear combination

A first idea to get a confidence interval for is to get a confidence interval for (by taking exponential values of bounds, since the exponential is a monotone function). Asymptotically, we know that

thus, an approximation for the variance matrix of will be based on , obtained by plugging estimators of the parameters.
Then, since as an asymptotic multivariate distribution, any linear combination of the parameters will also be normal, i.e.
has a normal distribution, centered on , with variance where is the variance of . All those quantities can be easily computed. First, we can get the variance of the estimators

> i1=sum(predict(reg,type="response"))
> i2=sum(cars$speed*predict(reg,type="response"))
> i3=sum(cars$speed^2*predict(reg,type="response"))
> I=matrix(c(i1,i2,i2,i3),2,2)
> V=solve(I)

Hence, if we compare with the output of the regression,

> summary(reg)$cov.unscaled
(Intercept)         speed
(Intercept)  0.0066870446 -3.474479e-04
speed       -0.0003474479  1.940302e-05
> V
[,1]          [,2]
[1,]  0.0066871228 -3.474515e-04
[2,] -0.0003474515  1.940318e-05

Based on those values, it is easy to derive the standard deviation for the linear combination,

> x=30
> P2=predict(r,type="link",se.fit=TRUE,
+ newdata=data.frame(speed=x))
> P2
$fit
1
5.046034

$se.fit
[1] 0.05747075

$residual.scale
[1] 1

> sqrt(V[1,1]+2*x*V[2,1]+x^2*V[2,2])
[1] 0.05747084
> sqrt(t(c(1,x))%*%V%*%c(1,x))
[,1]
[1,] 0.05747084

And once we have the standard deviation, and normality (at least asymptotically), confidence intervals are derived, and then, taking the exponential of the bounds, we get confidence interval

> segments(30,exp(P2$fit-1.96*P2$se.fit),
+ 30,exp(P2$fit+1.96*P2$se.fit),col="blue",lwd=3)

Based on that technique, confidence intervals are no longer centered on the prediction. But who cares ?

delta method

Actually, those who like to use “more or less” expressions for confidence intervals will not like non centered intervals. So, an alternative is to use the delta method. Instead of writing (again) something on the theory, we can use a package which computes that method,

> estmean=t(c(1,x))%*%coef(reg)
> var=t(c(1,x))%*%summary(reg)$cov.unscaled%*%c(1,x)
> library(msm)
> deltamethod (~ exp(x1), estmean, var)
[1] 8.931232
> P1=predict(r,type="response",se.fit=TRUE,
+ newdata=data.frame(speed=30))
> P1
$fit
1
155.4048

$se.fit
1
8.931232

$residual.scale
[1] 1

The delta method gives us (asymptotic) normality, so once we have a standard deviation, we get the confidence interval.

> segments(30,P1$fit-1.96*P1$se.fit,30,
+ P1$fit+1.96*P1$se.fit,col="blue",lwd=3)

Note that those quantities – obtained with two different approaches – are rather close here

> exp(P2$fit-1.96*P2$se.fit)
1
138.8495
> P1$fit-1.96*P1$se.fit
1
137.8996
> exp(P2$fit+1.96*P2$se.fit)
1
173.9341
> P1$fit+1.96*P1$se.fit
1
172.9101

bootstrap techniques

And a third method (but far from what I expect to teach on that course) is to use bootstrap techniques to about those results based on asymptotic normality (we have only 50 observations). The idea is to sample from out dataset, and to run a log-Poisson regression on those new samples, and to repeat a lot of time,

Time Series

Time horizon in forecasting, and rules of thumb

07/04/2011 Arthur Charpentier 3 Comments

I recently received an email about forecasting and rules of thumb. “Dans la profession […] se transmet une règle empirique qui voudrait que l’on prenne un historique du double de l’horizon de prévision : 20 ans de données pour une prévision à 10 ans, etc… Je souhaite savoir si cette règle n’aurait pas, par hasard, un fondement théorique quitte à ce que le rapport ne soit pas de 2 pour 1, mais de 3 pour 1 ou de 1 pour 1 par exemple.” To summarize briefly, the rule is to consider a 2-1 ratio for the period of observation vs. forecast horizon. And the interesting question is if there are justifications for such a rule…

At first, I remembered a rules of thumb, from the book by Box and Jenkins, which states that it is meaningless to look at autocorrelations when lags exceed the sample size over 6. So with 12 years of data, autocorrelations with a lag higher than two years are useless. But it is not what is mentioned here. So I looked at some dataset, and some standard time series models.

It depends on the series

It might obvious… but if it is the case, it means that it will be difficult to have a general rule of thumb. Consider e.g. the number of airline passengers,

library(forecast)
X = AirPassengers
ETS = ets(X)
plot(forecast(ETS,h=length(X)/2))

or some sales in a big store,

or car casualties in France, or the temperature in Nottingham Castle,

or the water level at Lack Hurron, or the flow of the Nile river,

or see also here for forecasting techniques in demography. Actually, in the case of life insurance, actuaries have to forecast future demography, i.e. try to assess death rates of those who currently purchase retirement contracts, who might be 20 years old. So they have to forecast death rate until 2100, say. One the one hand, it sounds difficult to make forecast over a century (it is already difficult for climate, I guess it is even more complex for human life). On the other hand, a 2-1 ratio means that we have to use data from 1800… Here again, it is difficult to justify that mortality in the 1850 could be interesting to say anything about mortality in 2050. So I guess it will be difficult to justify the use of general rules of thumb….

It depends on the model

Consider the following (simulated) series. Several models can be fitted. And the shape on the forecast (and the forecast error) will depend on the model considered. The benchmark can be the model without any dynamics, i.e. we assume that observations are i.i.d. Or more classically, assume that it is simple a white noise, i.e. an i.i.d centered process. Then the forecast is the following,

With that kind of assumption, we see that the 2-1 ratio is useless since we can get forecasts up to any horizon…. But that does not seem very robust. For instance, if we consider exponential smoothing techniques, we can obtain

Which is rather different. And with the 2-1 ratio, obviously, there is a lot of uncertainty at the end ! It would be even worst if we assume that we look at a random walk. Because actually a dozen models – at least – can be considered, from ARIMA, seasonal ARIMA, Holt Winters, Exponential Smoothing, etc…

http://freakonometrics.blog.free.fr/public/perso2/animationforecast.gif

So I do not see any theoretical justification of that rule of thumb. Obviously, the maximum horizon can not be extremely far away if the series is non-stationary, with a very irregular pattern, and with a lot of noise… So we’re back at the beginning. If anyone is willing to share his or her experience, comments are open.

Statistics

Margin of error, and comparing proportions in the same sample

15/10/2010 Arthur Charpentier Leave a comment

Irecently tried to answer a simple question, asked by @adelaigue. Actually, I thought that the answer would be obvious… but it is a little bit more compexe than what I thought. In a recent survey about elections in Brazil, it was mentionned in a French newspapper that “Mme Rousseff, 62 ans, de 46,8% des intentions de vote et José Serra, 68 ans, de 42,7%” (i.e. proportions obtained from the survey). It is also mentioned that “la marge d’erreur du sondage est de 2,2% ” i.e. the margin of error is 2.2%, which means (for the journalist) that there is a “grande probabilité que les 2 candidats soient à égalité” (there is a “large probability” to have equal proportions).
Usually, in sampling theory, we look at the margin of error of a single proportion. The idea is that the variance of $https://latex.codecogs.com/gif.latex?%20\widehat{p}$ , obtained from a sample of size $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial15.png$ is

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m201.png$

thus, the standard error is

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m202.png$

The standard 95% confidence interval, derived from a Gaussian approximation of the binomial distribution is

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m203.png$

The largest value is obtained when p is 1/2, and then we have a worst case confidence interval (an upper bound) which is

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m204.png$

So with a margin of error $https://perso.univ-rennes1.fr/arthur.charpentier/latex/m205.png$ means that $https://perso.univ-rennes1.fr/arthur.charpentier/latex/m206.png$ . Hence, with a 5% margin of error, it means that n=400. While 2.2% means that n=2000:
> 1/.022^2
[1] 2066.116
Classically, we compare proportions between two samples: surveys at two different dates, surveys in different regions, surveys paid by two different newpapers, etc. But here, we wish to compare proportions within the same sample. This has been consider in an “old” paper published in 1993 in the American Statistician,

It contains nice figures to illustrate the difference between the standard approach,

and the one we would like to study here.

This point is mentioned in the book by Kish, survey sampling (thanks Benoit for the reference),

Let $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial05.png$ and $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial06.png$ denote empirical frequencies we have obtained from the sample, based on $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial15.png$ observations. Then since

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial07.png$

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial08.png$

and

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial09.png$

we have

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial11.png$

Thus, a natural margin of error on the difference between the two proportion is here

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m207.png$

which is here 4 points
> n=2000
> p1=46.8/100
> p2=42.7/100
> 1.96*sqrt((p1+p2)-(p1-p2)^2)/sqrt(n)
[1] 0.04142327
Which is exactly the difference we have here ! Hence, the probability of reaching such a value is quite small (2%)
> s=sqrt(p1*(1-p1)/n+p2*(1-p2)/n+2*p1*p2/n)
> (p1-p2)/s
[1] 1.939972
> 1-pnorm(p1-p2,mean=0,sd=sqrt((p1+p2)-(p1-p2)^2)/sqrt(n))
[1] 0.02619152

Actually, we can compare the three margin of errors we have so far,

the upper bound

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m208.png$

the “average one”

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m209.png$

where

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m212.png$

the more accurate one we just obtained,

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/m213.png$

where $https://perso.univ-rennes1.fr/arthur.charpentier/latex/m214.png$ .
> p=seq(0,.5,by=.01)
> ic1=rep(1.96/sqrt(4*n),length(p))
> ic2=1.96*sqrt(p*(1-p))/sqrt(n)
> delta=.01
> ic31=1.96*sqrt(2*p-delta^2)/sqrt(n)
> delta=.2
> ic32=1.96*sqrt(2*p-delta^2)/sqrt(n)
> plot(p,ic32,type=”l”,col=”blue”)
> lines(p,ic31,col=”red”)
> lines(p,ic2)
> lines(p,ic1,lty=2)
So on the graph below, the dotted line is the standard upper bound, the plain line in black being a more accurate one when the probability is $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial19.png$ (the x-axis). The red line is the true margin of error with a large difference between candidates (20 points) and the blue line with a small difference (1 point).

Remark: an alternative is to consider a chi-square test, comparering two multinomial distributions, with probabilities $https://perso.univ-rennes1.fr/arthur.charpentier/latex/m215.png$ and $https://perso.univ-rennes1.fr/arthur.charpentier/latex/m216.png$ where $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial19.png$ is the average proportion, i.e. 44.75%. Then

$https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial21.png$

i.e. $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial22.png$ =3.71
> p=(p1+p2)/2
> (x2=n*((p1-p)^2/p+(p2-p)^2/p))
[1] 3.756425
> 1-pchisq(x2,df=1)
[1] 0.05260495
Under the null hypothesis, $https://perso.univ-rennes1.fr/arthur.charpentier/latex/multinomial22.png$ should have a chi-square distribution, with one degree of freedom (since the average is fixed here). Here the probability to reach that level is around 5% (which can be compared with the 2% we add before).

So finally, I would think that here, stating that there is a “large probability” is not correct…

Bayesian, Inference, Statistics, Tests

Intervalle de confiance pour une proportion

02/02/2009 Arthur Charpentier 1 Comment

Je voulais poster un petit billet puisque c’est la deuxième fois que l’on me pose une question similaire. Pour déterminer un intervalle de confiance pour une proportion, on connaît tous la formule

La question qui m’avait été posée est simple: que se passe-t-il si n est petit, et si – empiriquement – on a trouvé une proportion de 100% ? Dit autrement: “on fait 5 lancers de “pile” ou “face” et que on obtient 5 fois “pile”: que peut-on dire sur p, la probabilité de tomber sur “pile” ?”

On utilisant la formule précédente, on est un peu bloqué…

La première réponse pourrait être fréquentiste, et consiste à proposer une lecture duale de l’intervalle de confiance, et des abaques.

Le graphique de droite montre, en trait violet, l’intervalle de confiance approché, c’est à dire avec l’intervalle de confiance indiqué ci-dessus. La courbe bleue, encadrant la zone verte correspond à l’intervalle de confiance calculé à partir des quantiles d’une loi binomiale de paramètre la fréquence empirique observée. Ici, on aurait envie de dire que p a 95% de chance d’être supérieure à 47.8%.

En fait, cette “lecture duale” des abaques revient à construire – me semble-t-il – l’intervalle de confiance de Clopper-Pearson. L’intervalle de confiance a ici la forme suivante

Notons que c’est cette technique qui est utilisée sous R, via la fonction binom.test du package de base. Ici, on obtient

> binom.test(x, n, p = 0.5,conf.level = 0.95)
Exact binomial test
data: x and n
number of successes = 5, number of trials = 5, p-value = 0.0625
alternative hypothesis: true probability of success is not equal to 0.5
95 percent confidence interval:
0.4781762 1.0000000
sample estimates:
probability of success
1

…. ce qui correspond exactement à l’intervalle de confiance calculé précédemment. On retrouve d’ailleurs cette fonction dans binconf i.e.

> binconf(x = 5, n = 5, alpha = .05, method = "exact")
PointEst Lower Upper
1 0.4781762 1

qui donne (c’est rassurant) toujours la même chose. On peut aussi programme l’intervalle de confiance proposé par Heldge Blacker, qui donne ici un intervalle de confiance de la forme [50%; 100%]. Pour aller plus loin, je pourrais renvoyer au papier de Lawrence Brown, Tony Cai et Anirban DasGupta, dans Statistical Science sur ce sujet.

Une seconde réponse pourrait être une approche bayésienne, qui d’ailleurs correspond exactement à l’expérience de Bayes (1763) ou Laplace (1786).

On suppose que p est une variable aléatoire. On se donne une information a priori, et on peut se demander ce que devient la loi a posteriori de p, sachant que 5 “piles” ont été observés. On suppose que les tirages sont indépendants, et donc la probabilité d’avoir k fois “piles” sur n tirages suit une loi binomiale, conditionnellement à p, i.e.

Le plus simple est de considérer comme loi a priori la loi conjuguée de la loi binomiale, à savoir une loi beta. Dans ce cas, si p suit a priori une loi B(a,b), alors la loi a posteriori de p, sachant que x “pile” ont été observés est une loi B(a+x,b+n-x).

expérience décrite

La courbe en rouge est la densité de la loi a priori (uniforme ou beta de paramètres 2 et 2), et la courbe bleue, la loi a posteriori, avec un intervalle de confiance à 95%. Sur l’exemple de droite, par exemple, on a 95% de chance que p soit compris entre 47,3% et 96,8%.

Freakonometrics