Calculer l'intervalle de confiance pour la moyenne d'une distribution bêta

12

Considérons une distribution bêta pour un ensemble donné de notations dans [0,1]. Après avoir calculé la moyenne:

μ=αα+β

Existe-t-il un moyen de fournir un intervalle de confiance autour de cette moyenne?

dominic
la source
1
Dominique - vous avez défini la moyenne de la population . Un intervalle de confiance serait basé sur une estimation de cette moyenne. Quel exemple de statistiques utilisez-vous?
Glen_b -Reinstate Monica
Glen_b - Salut, j'utilise un ensemble de notes normalisées (d'un produit) dans l'intervalle [0,1]. Ce que je recherche, c'est une estimation d'un intervalle autour de la moyenne (pour un niveau de confiance donné), par exemple: moyenne + - 0,02
dominic
2
dominic: Laisse-moi réessayer. Vous ne connaissez pas la moyenne de la population . Si vous voulez qu'une estimation se situe au milieu de votre intervalle ( estimation demi-largeur , comme dans votre commentaire), vous auriez besoin d'un estimateur pour cette quantité dans l'ordre moyen pour placer un intervalle autour d'elle. Qu'utilisez-vous pour cela? Plausibilité maximum? Méthode des moments? autre chose? ±
Glen_b -Reinstate Monica
Glen_b - merci pour votre patience. Je vais utiliser MLE
dominic
2
dominic; dans ce cas, pour un grand on utiliserait les propriétés asymptotiques des estimateurs du maximum de vraisemblance; l'estimation ML de sera distribuée normalement de manière asymptotique avec la moyenne et l'erreur standard qui peut être calculée à partir des informations de Fisher . Dans de petits échantillons, on peut parfois calculer la distribution du MLE (bien que dans le cas de la bêta, je semble me rappeler que c'est difficile); une alternative consiste à simuler la distribution à la taille de votre échantillon pour comprendre son comportement. nμμ
Glen_b -Reinstate Monica

Réponses:

22

Bien qu'il existe des méthodes spécifiques pour calculer les intervalles de confiance pour les paramètres d'une distribution bêta, je décrirai quelques méthodes générales, qui peuvent être utilisées pour (presque) toutes sortes de distributions , y compris la distribution bêta, et sont facilement implémentées dans R .

Intervalles de confiance de vraisemblance du profil

Commençons par l'estimation du maximum de vraisemblance avec les intervalles de confiance de vraisemblance de profil correspondants. Nous avons d'abord besoin de quelques exemples de données:

# Sample size
n = 10

# Parameters of the beta distribution
alpha = 10
beta = 1.4

# Simulate some data
set.seed(1)
x = rbeta(n, alpha, beta)

# Note that the distribution is not symmetrical
curve(dbeta(x,alpha,beta))

Fonction de densité de probabilité pour la distribution bêta.

La moyenne réelle / théorique est

> alpha/(alpha+beta)
0.877193

Nous devons maintenant créer une fonction pour calculer la fonction de vraisemblance logarithmique négative pour un échantillon à partir de la distribution bêta, avec la moyenne comme l'un des paramètres. Nous pouvons utiliser la dbeta()fonction, mais comme cela n'utilise pas de paramétrage impliquant la moyenne, nous devons exprimer ses paramètres ( α et β ) en fonction de la moyenne et d'un autre paramètre (comme l'écart type):

# Negative log likelihood for the beta distribution
nloglikbeta = function(mu, sig) {
  alpha = mu^2*(1-mu)/sig^2-mu
  beta = alpha*(1/mu-1)
  -sum(dbeta(x, alpha, beta, log=TRUE))
}

Pour trouver l'estimation du maximum de vraisemblance, nous pouvons utiliser la mle()fonction dans la stats4bibliothèque:

library(stats4)
est = mle(nloglikbeta, start=list(mu=mean(x), sig=sd(x)))

Ignorez simplement les avertissements pour l'instant. Ils sont causés par les algorithmes d'optimisation essayant des valeurs invalides pour les paramètres, donnant des valeurs négatives pour α et / ou β . (Pour éviter l'avertissement, vous pouvez ajouter un lowerargument et modifier l'optimisation methodutilisée.)

Nous avons maintenant à la fois des estimations et des intervalles de confiance pour nos deux paramètres:

> est
Call:
mle(minuslogl = nloglikbeta, start = list(mu = mean(x), sig = sd(x)))

Coefficients:
        mu        sig 
0.87304148 0.07129112

> confint(est)
Profiling...
         2.5 %    97.5 %
mu  0.81336555 0.9120350
sig 0.04679421 0.1276783

Notez que, comme prévu, les intervalles de confiance ne sont pas symétriques:

par(mfrow=c(1,2))
plot(profile(est)) # Profile likelihood plot

Tracé de vraisemblance du profil pour la distribution bêta.

(Les deuxièmes lignes magenta externes montrent l'intervalle de confiance à 95%.)

Notez également que même avec seulement 10 observations, nous obtenons de très bonnes estimations (un intervalle de confiance étroit).

Comme alternative à mle(), vous pouvez utiliser la fitdistr()fonction du MASSpackage. Cela calcule également l'estimateur du maximum de vraisemblance et présente l'avantage que vous n'avez besoin que de fournir la densité, pas la probabilité logarithmique négative, mais ne vous donne pas d'intervalles de confiance de vraisemblance de profil, seulement des intervalles de confiance asymptotiques (symétriques).

Une meilleure option est mle2()(et les fonctions associées) du bbmlepackage, qui est un peu plus flexible et puissant que mle(), et donne des tracés légèrement plus agréables.

Intervalles de confiance Bootstrap

Une autre option consiste à utiliser le bootstrap. Il est extrêmement facile à utiliser dans R, et vous n'avez même pas à fournir de fonction de densité:

> library(simpleboot)
> x.boot = one.boot(x, mean, R=10^4)
> hist(x.boot)                # Looks good
> boot.ci(x.boot, type="bca") # Confidence interval
BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 10000 bootstrap replicates

CALL : 
boot.ci(boot.out = x.boot, type = "bca")

Intervals : 
Level       BCa          
95%   ( 0.8246,  0.9132 )  
Calculations and Intervals on Original Scale

Le bootstrap a l'avantage supplémentaire de fonctionner même si vos données ne proviennent pas d'une distribution bêta.

Intervalles de confiance asymptotiques

Pour les intervalles de confiance sur la moyenne, n'oublions pas les bons vieux intervalles de confiance asymptotiques basés sur le théorème central limite (et la distribution t ). Tant que nous avons soit une grande taille d'échantillon (donc le CLT s'applique et la distribution de la moyenne de l'échantillon est approximativement normale), soit de grandes valeurs de α et β (de sorte que la distribution bêta elle-même soit approximativement normale), cela fonctionne bien. Ici, nous n'avons ni l'un ni l'autre, mais l'intervalle de confiance n'est toujours pas trop mauvais:

> t.test(x)$conf.int
[1] 0.8190565 0.9268349

Pour des valeurs de n légèrement plus larges (et des valeurs pas trop extrêmes des deux paramètres), l'intervalle de confiance asymptotique fonctionne extrêmement bien.

Karl Ove Hufthammer
la source
Merci Karl. Question rapide: comment avez-vous déterminé votre alpha et bêta? J'ai utilisé la variance et la moyenne de l'échantillon pour obtenir l'alpha et la bêta, mais je pense que j'ai peut-être confondu la moyenne de l'échantillon avec la moyenne de la population, donc je ne suis pas sûr de l'avoir fait correctement ... voir le commentaire de Glen_b ci-dessus .
dominic
Pour déterminer α et β en tant que fonctions de la moyenne et de l'écart-type, je viens d'inverser les fonctions de la moyenne et des écarts-types en tant que fonctions de α et β (mais je suis sûr que vous pouvez également le consulter sur le net).
Karl Ove Hufthammer
α,β
0

Découvrez la régression bêta. Une bonne introduction à la façon de le faire en utilisant R peut être trouvée ici:

http://cran.r-project.org/web/packages/betareg/vignettes/betareg.pdf

Une autre façon (très facile) de construire un intervalle de confiance serait d'utiliser une approche boostrap non paramétrique. Wikipédia a de bonnes informations:

http://en.wikipedia.org/wiki/Bootstrapping_%28statistics%29

Vidéo aussi sympa ici:

http://www.youtube.com/watch?v=ZCXg64l9R_4

Rasmus Bååth
la source