CAPÍTULO 2: DISTRIBUCIONES CON-

JUGADAS

Para leer

Gelman et al (1995), Capı́tulo 2.

Lee (1997), Capı́tulo 3, Secciones 3.1–3.5.
Bernardo y Smith (1994), Sección 5.2.

Ejemplo 11 Retomando el Ejemplo 6, sea la

distribución a priori una Beta, B(α, β). Entonces
la distribución a posteriori será

f (θ|x) ∝ f (θ)l(θ|x)
1
∝ θ α−1(1 − θ)β−1 ×
B(α, β)


12
θ 9(1 − θ)3
9
∝ θ α+9−1(1 − θ)β+3−1
θ|x ∼ B(α + 9, β + 3).

55
La distribución a posteriori es de la misma fa-
milia como la distribución a priori.

Generalizando, dada la distribución a priori B(α, β),

si se tira la moneda n veces observando x cruces
y n − x caras, entonces la distribución a poste-
riori es B(α + x, β + n − x).

Entonces para cualquier distribución de mues-

treo de tipo Bernoulli, binomial, geométrica o
binomial negativa, si se utiliza una distribución
a priori beta, luego la distribución a posteriori
es de la misma familia.

Se dice que la distribución beta es conjugada

(Raiffa y Schlaifer 1961) con las distribuciones
de muestreo Bernoulli, binomial, etc..

56
La idea de la distribución conjugada es que la
distribución a posteriori tiene la misma forma
algebráica como la distribución a priori.

Para formalizarla, recordamos la Definición 1

de una estadı́stica suficiente t. Si t es suficiente
para θ , entonces

f (x|θ ) = f (t|θ )f (x|t)

l(θ |x) ∝ f (t|θ ).

Luego, esta claro que si se elige una distribu-

ción a priori f (θ ) con la misma estructura de
f (t|θ ), pensando en esto como función de θ ,
entonces, la distribución a posteriori tendrá la
misma forma.

Sigue una definición formal.

57
Definición 3 (Bernardo y Smith 1994) La fa-
milia conjugada de densidades a priori para θ
con respeto a la verosimilitud l(θ |x) con es-
tadı́stico suficiente t = t(x) = {n, s(x)} (de
dimensión fija k, independiente de x) es

f (θ |τ = (τ0, τ1, . . . , τk ) ∈ T )
donde
  
T = τ : f (s = (τ1, . . . , τk )|θ , n = τ0) dθ < ∞
Θ
y
f (s = (τ1, . . . , τk )|θ , n = τ0)
f (θ |τ ) =  .
Θ f (s = (τ1, . . . , τk )|θ , n = τ0) dθ

58
Ejemplo 12 Sea X|θ ∼ P(θ), una densidad
Poisson.

Dada una muestra de datos x, la verosimilitud

es
n
 θ xi e−θ
l(θ|x) = ∝ θ nx̄e−nθ
i=1 xi !
donde x̄ es suficiente.

Entonces, una distribución conjugada serı́a de

forma
f (θ) ∝ θ τ1 e−τ0θ
es decir una distribución gamma, por ejemplo
G(α, β) (donde α = τ1 + 1 y β = τ0).

Luego si θ ∼ G(α, β), dada una muestra de

tamaño n, se tiene θ|x ∼ G(α + nx̄, β + n).

La distribución gamma es conjugada con la dis-

tribución muestral Poisson.
59
Familias Exponenciales y Distribuciones Con-
jugadas

Dada la definición de una distribución conju-

gada, es obvio que está relacionado con el
concepto de familias exponenciales de distribu-
ciones.

Definición 4 La familia de distribuciones f (x|·)

con densidades de forma
 
k
f (x|θ ) = C(θ )h(x) exp  Ri(θ )si (x)
i=1
se llama una familia exponencial (con k parámet-
ros).

φ = R(θ ) = (R1(θ ), . . . , Rk (θ )) se llama el

parámetro natural de la familia.

60
Si el soporte de X no depende de θ , la familia
se llama regular. Si el soporte de X depende
de θ , la familia es irregular.

Aunque la mayorı́a de las distribuciones co-

munes son de familias exponenciales regulares,
existen algunas excepciones:

Ejemplo 13 La distribución uniforme X|θ ∼

U (0, θ) es una familia exponencial irregular.

Ejemplo 14 La distribución de Cauchy

1
f (x|θ) ∝
1 + (x − θ)2
no pertenece a ninguna familia exponencial.

61
Una familia conjugada con una familia ex-
ponencial regular

Si f (x|θ ) es de una familia exponencial, con

densidad de la forma de la Definición 4, está claro
que existe una familia conjugada.

Teorema 2 Sea f (x|θ ) una familia exponen-

cial como en Definición 4. Entonces la distribu-
ción a priori f (θ ) = g(θ |τ ) con
 
k
g(θ |τ ) = K(τ )C(θ )τ0 exp  Ri(θ )τi
i=1
es conjugada y dada una muestra de tamaño
n, la distribución a posteriori es

f (θ |x) = g(θ |τ + t)
n n 
donde t(x) = n, i=1 s1(xi), . . . , i=1 sk (xi)
es el estadı́stico suficiente.

62
Igualmente, si la distribución predictiva a priori
de unos futuros datos y es

f (y) = f (y|θ )f (θ ) dθ

= f (y|θ )g(θ |τ ) dθ
= h(y|τ )
entonces la distribución predictiva a posteriori
es
f (y|x) = h(y|τ + t).

Demostración

Dada una muestra x, la verosimilitud es

 
n k
l(θ |x) ∝ C(θ ) exp 
n
Ri (θ )si (xj )
j=1 i=1
n n 
y t(x) = n, i=1 s1(xi), . . . , i=1 sk (xi ) es su-
ficiente.
63
Entonces, dada la distribución a priori

f (θ ) = g(θ |τ )
 
k
∝ C(θ ) τ0
exp  Ri(θ )τi
i=1
la distribución a posteriori es
 
k n
f (θ |x) ∝ C(θ )τ0 +n exp  Ri (θ )(τi + si(xj ))
i=1 j=1
= g(θ |τ + t).

Igualmente

f (y|x) = f (y|θ )f (θ |x) dθ

= f (y|θ )g(θ |τ + t) dθ
= h(τ + t).

64
El teorema demuestra que la familia conjugada
de distribuciones está cerrada bajo muestreo
(Weatherill 1961).

La información en la distribución a priori es

equivalente a la información contenida en una
muestra con estadı́stico suficiente τ .

Si se quiere ser objetivo, una posibilidad es con-

siderar la distribución a priori derivada cuando
τ → 0). En este caso, la distribución a priori
puede ser impropia.

65
Ejemplo 15 Sea X|θ ∼ E(θ). Entonces

f (x|θ) = θe−xθ
es una famila exponencial y una distribución a
priori conjugada para θ es

f (θ) ∝ θ τ0 exp(−βτ1)
es decir una distribución gamma G(α, β) donde
α = τ0 + 1 y β = τ1.

Entonces, la distribución a priori gamma es

conjugada con la distribución exponencial y da-
da una muestra de tamaño n, la distribución a
n
posteriori es G(α + n, β + i=1 xi).

Dejando τ0, τ1 → 0, este implica que la dis-

tribución a priori es (G(−1, 0)) f (θ) ∝ 1, es
decir una distribución impropia. En este caso,
θ|x ∼ G(n, nx̄) y la media a posteriori coincide
con el EMV.

66
Ejemplo 16 Sea X un ensayo de Bernoulli con
parámetro θ. Luego

P (X = x|θ) = θ x(1 − θ)1−x

 
θ
= (1 − θ) exp x log
1−θ

Una distribución conjugada es

 
θ
f (θ) ∝ (1 − θ)τ0 exp τ1 log
1−θ
∝ θ τ1 (1 − θ)τ0−τ1
es decir una distribución beta, B(α, β) donde
α = τ1 + 1 y β = τ0 − τ1 + 1.

67
Dada una muestra de tamaño n, el estadı́sti-
n
co suficiente es t(x) = {n, i=1 xi} y luego la
distribución a posteriori es B(α∗, β ∗) donde
n
α∗ = τ1 + xi + 1
i=1
= α + # cruces
n
β ∗ = τ0 + n − τ1 − xi + 1
i=1
= β + # caras

Se ha visto en el Ejemplo 11 que la distribución

beta es conjugada con la binomial.

68
Interpretación de la información propor-
cionada por una distribución a priori con-
jugada

Se ha visto que el parámetro τ0 juega el mis-

mo papel en la distribución a priori como el
tamaño de la muestra, y los demás parámet-
ros (τ1, . . . , τk ) equivalen a los estadı́sticos su-
ficientes. Entonces, se puede interpretar la in-
formación proporcionada por la distribución a
priori como la información de un experimento
de tamaño n donde se observan estadı́sticos
suficientes (τ1, . . . , τk ).

Ejemplo 17 Volviendo al Ejemplo 11 se ve

que dada una distribución a priori B(α, β), la
distribución a posteriori es
B(α + # cruces visto, β + # caras visto)

La información contenida en la distribución a

priori es equivalente a una muestra de α + β
tiradas de la moneda, α de ellas siendo cruces.
69
La media a posteriori como media ponder-
ada

Se puede expresar la media a posteriori como

una media ponderada de la media a priori y el
EMV.

Ejemplo 18 En el Ejemplo 12, tenemos

α + nx̄
E[θ|x] =
β+n
α
= w + (1 − w)x̄
β
β
donde 0 ≤ w = β+n ≤ 1.

La media a posteriori es una media ponderada

con pesos proporcionales al número de obser-
vaciones equivalentes en la distribución inicial
(β) y al número de datos (n).
70
 Formalización

Supongamos una famı́lia exponencial con par-

ámetro natural θ (de dimensión k),
 
k
f (x|θ ) ∝ D(θ ) exp  θisi (x) .
i=1

La distribución a priori conjugada es

 
k
f (θ ) ∝ D(θ )τ0 exp  θiτi
i=1
 
k
= exp  θiτi − τ0G(θ ) donde
i=1
G(θ ) = − log D(θ )
y la distribución a posteriori es
 
k n
f (θ |x) ∝ exp  θi(τi + si (xj )) − (τ0 + n)G(θ ) .
i=1 j=1

71
Luego se tiene el siguiente teorema.

Teorema 3

E[∇G(θ )] = ws̄ + (1 − w)τ donde

∂
[∇G(θ )]i = G(θ )
∂θi
n
w = y
τ0 + n
1
τ = (τ1, . . . , τk )
τ0 
1 n n
s̄ = s1(xj ), . . . , sk (xj )
n j=1 j=1

El teorema demuestra que la media a posteriori

de ∇G(θ ) es una media ponderada de la media
a priori y la media muestral de los estadı́sticos
 ) donde θ
suficientes s̄, es decir ∇G(θ  es el EMV

de θ con ponderaciones proporcionales a τ0 y

el tamaño de la muestra.

72
Demostración

Abusando de la notación, sea τ = (τ1, . . . , τk ).

Luego

τ0 (τ − E [∇G(θ )]) = τ0 (τ − ∇G(θ )) f (θ ) dθ

= ∇f (θ ) dθ donde
∂
∇ [f (θ )]i = f (θ ) que implica
∂θi
E [∇G(θ )] = τ

Entonces, a posteriori, mediante el mismo ar-

gumento, se tiene:
n
E[∇G(θ )|x] = τ + s(xj )
j=1
donde s(x) = (s1(x), . . . , sk (x)). Luego
τ0τ + ns
E[∇G(θ )|x] =
τ0 + n
= ws + (1 − w)τ

73
Ejemplo 19 Retomamos el Ejemplo 18. Se tiene
X|θ ∼ P(θ). Entonces
f (x|θ) ∝ exp(x log θ − θ)
y reparameterizando, φ = log θ,
 
φ
f (x|φ) = exp φx − e

Una distribución a priori conjugada para φ es

 
φ
f (φ) ∝ exp φτ1 − τ0e

y entonces, la media a posteriori de ∂φ ∂ eφ =

eφ = θ es
τ τ + nx
E[θ|x] = 0 1 donde τ 1 = ττ1 .
τ0 + n 0

Haciendo la transformación al revés, θ = log φ,

observamos que la distribución a priori impli-
cada para θ es
f (θ) ∝ θ τ1−1 exp(−τ0θ) es decir
θ ∼ G(τ1, τ0)
y sustituyendo α = τ1 y β = τ0, se ha de-
mostrado el resultado del Ejemplo 18.
74
En el caso más general, las ponderaciones no
necesariamente son proporcionales a n y τ0.

Ejemplo 20 Volvemos al Ejemplo 15 con datos

exponenciales y distribución a priori gamma.
En este caso, la distribución a posteriori es
G(α + n, β + nx̄) y la media a posteriori es
α+n 1 α
E[θ|x] = = w + (1 − w)
β + nx̄ x̄ β
x̄ .
donde w = x̄+β

No obstante, se puede utilizar el teorema para

demostrar que la media a posteriori de 1θ es una
β
combinación lineal de la media a priori α−1 y
x̄:
 
1 β
E | x = wx̄ + (1 − w)
θ α−1
donde w = n+α−1n .

75
Inferencia para la media de una distribu-
ción normal

Ejemplo 21 X|µ ∼ N (µ, σ 2) (σ 2 conocido).

 
1
f (x|µ) ∝ exp − 2 (x − µ)2
 2σ   
1 2 µ
∝ exp − 2 µ exp 2
x
2σ σ
es decir,
 una
 familia exponencial con C(µ) =
exp − 2σ 1 µ2 , s (x) = x y R (µ) = µ .
2 1 1 σ2

Entonces, una distribución a priori conjugada

será de forma
    
1 2 τ0 µ
f (µ) ∝ exp − 2 µ exp 2
τ1

2σ
σ
τ τ
∝ exp − 02 µ2 − 2 1 µ
2σ τ0
y completando el cuadrado, se ve que la dis-
tribución a priori conjugada será normal

µ ∼ N (m, δ 2)
donde m = τ1/τ0 y δ 2 = σ 2/τ0.
76
 La distribución a posteriori. Cálculo directo
mediante el teorema de Bayes.

Dados los datos x, la verosimilitud es

 
n
1
l(µ|x) ∝ (σ 2)−n/2 exp − (x i − µ)2
2σ 2 i=1
  
1
∝ (σ 2)−n/2 exp − 2 (n − 1)s2 + n(µ − x̄)2
2σ
1 n
donde s = n−1 i=1(xi − x̄)2 es el estimador
2
clásico de la varianza.

Entonces, dada una distribución a priori nor-

mal, la distribución a posteriori es
  
1 (µ − m)2 n(µ − x̄)2
f (µ|x) ∝ exp − +
2 δ2 σ2
     
1 1 n m nx̄
∝ exp − + µ 2
−2 + 2 µ
2 δ2 σ2 δ2 σ
 2 
  1 n
1 1 n m + σ2 x̄
∝ exp − 2
+ 2
µ − δ2
1 n

2 δ σ δ2
+ σ2

 −1 
1 n
m + σ2 x̄ 1 n
µ|x ∼ N δ2
1 n
, +
δ 2 + σ2
δ2 σ2

77
Cálculo usando las propiedades de la famil-
ia exponencial.

Se ha visto que la distribucı́on a priori conjuga-

da es N (m, δ 2) donde m = τ1/τ0 y δ 2 = σ 2/τ0
y que
n n
s1(xj ) = xj = nx̄.
j=1 j=1

 %

Entonces, la distribución a posteriori es N m%, δ 2
donde
n
τ1 + i=1 xi
m% =
τ0 + n
2 % σ2
δ =
τ0 + n
y entonces, sustituyendo por τ0 y τ1,

τ0 = σ 2/δ 2
τ1 = mσ 2/δ 2.
78
Luego, se tiene

2% σ2
δ =
σ 2/δ 2 + n
 
1 n −1
= 2
+ 2
δ σ
2 2 n
mσ /δ + i=1 xi
m% = y dividiendo por σ 2
σ 2/δ 2 + n
1 m + n x̄
2 σ2
m% = δ 1 n
δ2
+ σ2
es decir el mismo resultado como anteriormente.

79
Propiedades de la distribución a posteriori

Observación 17

E[µ|x] = wE[µ] + (1 − w)µ̂

1
donde µ̂ = x̄ es el EMV de µ y w = δ2
1+ n . Las
δ2 σ2
ponderaciones son proporcionales a las preci-
siones del estimador a priori (m) y el EMV.

Observación 18 Un intervalo de credibili-

dad a posteriori de 95 % para µ es

m + n x̄  −1
δ 2 σ 2 1 n
1 + n
± 1,96 2
+ 2
δ σ
δ2 σ2

80
Observación 19 Si y sólo si δ12 = 0, el inter-
valo será igual al intervalo clásico de confianza.
σ
x̄ ± 1,96 √
n
En esta situación la distribución a priori es im-
propia

µ ∼ lı́m N (m, δ 2) ⇒ f (µ) ∝ c

δ 2→∞

Observación 20 A menudo se expresa la dis-

tribución a priori de otra forma, escribiendo
δ 2 = σ 2/α. Entonces, la distribución a priori es


σ2
µ∼N m,
α
y la distribución a posteriori reduce a


αm + nx̄ σ2
µ|x ∼ N , .
α+n α+n
Es posible derivar este resultado directamente
a través del Teorema 3.
81
Mixturas de distribuciones conjugadas

La familia de mixturas de distribuciones con-

jugadas es también conjugada. Sea la distribu-
ción a priori
r
p(θ ) = wifi(θ )
i=1
donde fi(·) es conjugada con una distribución
muestral f (x|θ ) en el sentido de Definición 3.
Entonces
r
f (θ |x) = wi∗fi∗(θ )
i=1
donde fi∗(·) es también conjugada.

Observación 21 Es posible aproximar cualquier

distribución f (·) con una mixtura suficiente-
mente grande de densidades conjugadas.

82
 Ejemplo 22 Volviendo al Ejemplo 11, supong-
amos que la distribución a priori es una mixtura
de una distribución uniforme y una distribución
más informativa:
θ ∼ 0,25B(1, 1) + 0,75B(5, 5).


1
f (θ|x) ∝ 0,25 + 0,75 θ 5−1(1 − θ)5−1
B(5, 5)
θ 9(1 − θ)3
∝ 0,25θ 10−1(1 − θ)4−1 +
0,75 14−1
θ (1 − θ)8−1
B(5, 5)
1
∝ 0,25B(10, 4) θ 10−1(1 − θ)4−1 +
B(10, 4)
0,75B(14, 8) 1
θ 14−1(1 − θ)8−1
B(5, 5) B(14, 8)
= wB(10, 4) + (1 − w)B(14, 8)

0,25B(10,4)
donde w = .
0,25B(10,4)+ 0,75B(14,8)
B(5,5)

83
Inferencia para la distribución uniforme

Ejemplo 23 Aunque la distribución uniforme

(X ∼ U (0, θ)) es una familia exponencial irreg-
ular, si existe una distribución conjugada.

Dados los datos x, la verosimilitud es

l(θ|x) = θ −n por θ > máx{x1, . . . , xn}.

Supongamos una distribución a priori Pareto;

θ ∼ PA(α, β). Entonces

f (θ) = βαβ θ −β−1

donde θ > α.

La distribución a posteriori es

f (θ|x) ∝ θ −β−n−1
para θ > α∗ = máx{α, x1, . . . , xn} es decir que
θ|x ∼ PA(α∗, β + n).
84
Distribuciones conjugadas no siempre son
fáciles de utilizar

Si se conocen las propiedades de la familia con-

jugada de distribuciones, entonces la inferencia
es sencilla. No obstante, existen algunas excep-
ciones.

Ejemplo 24 Supongamos que X|θ ∼ B(α, θ)

donde α es conocido. Entonces
Γ(α + θ)
f (x|θ) ∝ (1 − x)θ
Γ(θ)
y una distribución a priori conjugada será

a
Γ(α + θ)
f (θ) ∝ bθ
Γ(θ)
cuando

a+n  n θ
Γ(α + θ) 
f (θ|x) ∝ b (1 − xi)
Γ(θ) i=1
y es imposible calcular la constante de inte-
gración, la media . . . sin el uso de integración
numérica.
85
Aplicación 1: Inferencia para colas Marko-
vianas

Consideramos inferencia para el sistema de co-

las M/M/1 siguiendo Armero y Bayarri (1994a).

El sistema de colas M (λ)/M (µ)/1 es un sis-

tema con tiempos entre llegadas X ∼ E(λ),
tiempos de servicio Y ∼ E(µ) y un servidor.
Si llega un cliente y el servidor no está vacio,
entonces espera su torno en la cola.

Tı́picamente se interesa por variables como el

tamaño de la cola o el tiempo de espera de un
cliente. Si, en principio, el sistema está vacia
la distribución de N (t), el número de clientes
en el sistema después de un tiempo t tiene una
fórmula complicada (ver Gross y Harris 1985):

86
  
λ
(n−n0 )/2   
−(λ+µ)t
P (N (t) = n) = e In−n0 2 λµt
µ
 (n−n0 −1)/2   
λ
+ In+n0 +1 2 λµt
µ

   n/2 ∞   
λ λ
+ 1− Ij 2 λµt 
µ µ
j=n+n +1
0

donde Ij (c) es una función Bessel modificada,

es decir
∞
(c/2)j+2k
Ij (c) = .
k=0
k!(j + k)!

Es interesante considerar bajo cuáles condi-

ciones converge está distribución cuando t →
∞. Claramente, si la tasa de llegadas es más
alta que la tasa de servicios, se esperarı́a que
el número de personas esperando crecerı́a con
el tiempo.
87
Condición del equilibrio del sistema

El sistema es estable si y sólo si la intensidad

de tráfico ρ = λ/µ < 1.

En este caso la distribución de N (t) se acerca

a un lı́mite cuando t → ∞:

P (N = n|ρ) = lı́m P (N (t) = n|ρ)

t→∞
= (1 − ρ)ρn.
Ver por ejemplo Gross y Harris (1985).

Igualmente, el tiempo, W , pasado en el sis-

tema por un cliente llegando en equilibrio es
exponencial

W |λ, µ ∼ E(µ − λ).

88
Observación 22 También existen expresiones
explı́citas para el tiempo de espera, periodos
de ocupación etc. Ver Gross y Harris (1985).

Para tomar decisiones sobre la introducción de

nuevos servidores etc. se necesita estudiar la
estabilidad del sistema y las distribuciones en
el equilibrio.

89
Experimento y Inferencia

Se supone que λ y µ son desconocidos.

Un experimento sencillo es el de observar nl

tiempos entre llegadas y ns tiempos de sevicio.
En este caso, la verosimilitud es

l(λ, µ|x, y) ∝ λnl e−λtl µns e−µts

donde tl y ts son las sumas de los tiempos de
llegadas y servicios respectivamente.

Suponiendo unas distribuciones a priori (inde-

pendientes) gamma;

λ ∼ G(αl , βl ) µ ∼ G(αs, βs)

entonces las distribuciones a posteriori son tam-
bién (independientes) gamma:

λ|x ∼ G(αl + nl , βl + tl )
µ|y ∼ G(αs + ns, βs + ts)

90
Estimación de la intensidad de tráfico

Es fácil estimar la intensidad de tráfico ρ. Su

esperanza es:

E[ρ|x, y] = E [λ/µ|x, y]
= E[λ|x]E [1/µ|y]
αl + nl βs + ts
=
βl + tl αs + ns − 1

Además es fácil evaluar la distribución de ρ.

A posteriori, λ y µ son independientes y

2(βl +tl )λ|x ∼ χ2

2(αl +nl ) 2(βs+ts)µ|y ∼ χ2
2(αs +ns)
y entonces
(βl + tl )(αs + ns) 2(α +n )
ρ|x, y ∼ F2(α l +nl ) .
(αl + nl )(βs + ts) s s

La probabilidad a posteriori de que el sistema

sea estable es p = P (ρ < 1|x, y). Se la puede
evaluar mediante integración numérica.
91
Estimación del tamaño de la cola en equi-
librio

Si p es grande, es natural suponer que la cola

sea estable. En este caso, se puede hacer infer-
encia sobre el tamaño de la cola en equilibrio.

P (N = n|x, y, equilibrio) = P (N = n|x, y, ρ < 1)

 1
= (1 − ρ)ρnf (ρ|x, y, ρ < 1) dρ
0
 1
1
= (1 − ρ)ρnf (ρ|x, y) dρ
p 0

Se puede expresar este integral en términos de

la función Gauss hipergeométrica:

92
 α∗l Γ(α∗l + n)
P (N = n|·) = ×
Γ(α∗l + n + 2)
 
∗ ∗ ∗ ∗ βl∗
2 F1 αl + αs , αl + n; αl + n + 2; − β ∗
  s

∗ ∗ ∗ ∗ βl∗
2 F1 αl + αs , αl ; αl + 1; − β ∗
s

donde α∗l = αl + nl , βl∗ = βl + tl , α∗s = αs + ns,

βs∗ = βs +ts y la función Gauss hipergeométrica
es
 1
Γ(c)
F
2 1 (a, b; c; d) = z b−1 (1−z)c−b−1 (1−dz)−a dz.
Γ(b)Γ(c − b) 0

No obstante, se necesitan métodos numéricos

para evaluar la integral. Existen varias posibili-
dades (integración numérica, Monte Carlo ...)

Observación 23 Si αs es un número entero,

se puede evaluar el integral explı́citamente.

93
Resultados explı́citos cuando αs es un número
entero

P (N = n|x, y) = E [(1 − ρ)ρn|x, y, ρ < 1]

= E [ρn|x, y, ρ < 1] −
 
n+1
E ρ |x, y, ρ < 1

Hallamos los momentos de ρ condicionado en

el equilibrio.

 r 
λ
E [ρr |x, y, ρ < 1] = E |x, y, λ < µ
µ
 ∞  ∞  r
1 λ
= f (λ|x)f (µ|y) dµ dλ
p 0 µ
 ∞ λ  ∞
1
= λr f (λ|x) µ−r f (µ|y) dµ dλ
p 0
 ∞ λ
1 Γ(α∗s − r) ∗ r
= λ f (λ|x)
r
(βs )
p 0 Γ(α∗s )
 ∞
g(µ|α∗s − r, βs∗) dµ dλ
λ

94
donde g(µ|α∗s − r, βs∗) representa una función de
densidad gamma (Erlang) con parámetros α∗s −
r y βs∗.

Recordamos un resultado para un proceso Pois-

son. Si el número de sucesos en una unidad de
tiempo es Poisson con tasa θ, entonces el tiem-
po hasta el primer suceso es exponencial con
la misma tasa y el tiempo hasta el suceso r es
Erlang o gamma con parámetros r y θ.


Luego λ∞ g(µ|α∗s − r, βs∗) dµ es la probabilidad
de que el suceso número α∗s − r ocurra después
de un tiempo λ. Esto es igual a decir que el
número de sucesos en [0, λ) es menor de α∗s − r
y luego
 ∞ α∗s −r−1
(βs∗λ)j −βs∗λ
g(µ|α∗s − r, βs∗) dµ = e .
λ j=0 j!

95
Entonces,
 ∞
1 Γ(α∗s − r) ∗ r
E [ρ |x, y, ρ < 1] =
r
λ f (λ|x)
r
(βs )
p 0 Γ(α∗s )
α −r−1
∗

(βs∗λ)j −βs∗ λ
s

e dλ
j=0
j!
1 Γ(α∗s − r) ∗ r
= (βs )
p Γ(α∗s )
α∗s −r−1 ∗ j  ∞
βs
λr e−βs λ f (λ|x) dλ
∗

j=0
j! 0
α∗s −r−1
1 Γ(α∗s − r) ∗ r βs∗j
= (βs )
p Γ(α∗s ) j=0
j!
∗
Γ(α∗l + j) βl∗αl
Γ(α∗l ) (βl∗ + βs∗)αl +j
∗

1 Γ(α∗s − r) ∗ r
= ∗
(βs ) P (T ≤ α∗s − r − 1)
p Γ(αs )
 ∗ 
βl
donde T ∼ BN α∗l , ∗
βl +βs∗
.

Luego se tiene una expresión explı́cita para

P (N = n|x, y) y de la misma manera es posible
hallar la distribución predictiva de W .
96
Estimación de la distribución de N (t) por
simulación Monte Carlo

La idea es simular una muestra de las dis-

tribuciones a posteriori de λ y µ suponiendo la
condición del equilibrio y usar los datos mues-
trales para estimar las probabilidades de in-
terés.

Se usa el siguiente algoritmo para estimar las

distribuciones predictivas de N y W a través de
una muestra Monte Carlo de tamaño S. Para
estimar la distribución de N (t), no es necesario
imponer la condición del equilibrio.

97
Algoritmo

1. Para i = 1, . . . , S

a) Generar λi ∼ G(αl + nl , βl + tl )

b) Generar µi ∼ G(αs + ns, βs + ts)

λi
c) Definir ρi = µ .
i

2.
1 S
P (N (t) = n|data) ≈ P (N (t) = n|λi, µi)
S i=1

Se puede usar el mismo tipo de algoritmo para

estimar la densidad de un periodo de ocupación.

98
Estimación de los momentos

A pesar de que los momentos de N (t) existen

para cualquier t, las variables N y W no tienen
media.

E[N |x, y, ρ < 1] = E[E[N |ρ]|x, y, ρ < 1]

 
ρ
= E |x, y, ρ < 1
1−ρ

1 1 ρ
= f (ρ|x, y) dρ
p 0 1−ρ
= ∞

Observación 24 Es una caracterı́stica general

en inferencia para colas. El mismo resultado
ocurre con el tiempo de espera o el periodo
de ocupación. Ver Armero y Bayarri (1994a),
Wiper (1998).

Se resuelve el problema parcialmente supon-

iendo a priori que P (ρ ≥ 1) = 0. Ver Armero y
Bayarri (1994b) o Ruggeri et al (1996).
99
Ejemplo 25 Hall (1991) presenta datos de un
cajero automático en Berkeley. El siguiente his-
tograma muestra los tiempos entre llegadas.
Parece que se ajusten a una distribución expo-
nencial.
0.7

0.6

0.5

0.4
f

0.3

0.2

0.1

0
0 1 2 3 4 5 6
x
a

100
(Los servicios no se ajustan tan bien pero ...)

Los estadı́sticos suficientes son

na = ns = 98, ta = 119,71 ts = 81,35

Suponiendo distribuciones a priori no-informativas

1 1
f (λ) ∝ f (µ) ∝
λ µ
las distribuciones a posteriori son

λ|x ∼ G(98, 119,71) µ|y ∼ G(98, 81,35)

y entonces E[ρ|x, y] ≈ 0,668 y se puede de-
mostrar que P (ρ < 1|x, y) ≈ 0,997.

Parece razonable suponer equilibrio.

El primer dibujo muestra las probabilidades para

diferentes valores del tamaño de la cola en
equilibrio.
101
 0.35

0.3

0.25

0.2
P(N=n|data)

0.15

0.1

0.05

0
0 2 4 6 8 10 12 14 16 18 20
n

El siguiente diagrama muestra la distribución

predictiva de W

102
 1

0.9

0.8

0.7

0.6
FW(w|data)

0.5

0.4

0.3

0.2

0.1

0
0 2 4 6 8 10 12 14 16 18 20
w

El último dibujo muestra la distribución de N (t)

para varios tiempos t sin suponer el equilibrio.

103
 0.7

0.6 t=1
t=10

t=50
0.5

0.4
P(N(t)|data)

0.3

0.2

0.1

0
0 2 4 6 8 10 12 14 16 18 20
N(t)

Se ve la convergencia a la distribución de N .

104
Aplicación 2: Comprobación de software
con múltiples usuarios

Estudiamos el problema de la detección de er-

rores en software siguiendo Wiper y Wilson
(2007).

¿Porqué comprobar el software?

Se quiere detectar y depurar los errores (bugs

o faults) en el programa y mejorar la calidad
del producto final.

El modelo tipo académico es lo siguiente.

Se generan entradas del perfil de un usuario

tı́pico y se graban los tiempos entre fallos (fail-
ures) sucesivos. Después de cada fallo, se depu-
ra el programa.

105
A lo largo del tiempo, la fiabilidad del soft-
ware se mejora y se lanza el producto al mer-
cado cuando sea suficientemente fiable. Ver,
por ejemplo Singpurwalla y Wilson (1999).

En comprobación beta (beta testing), varios

usuarios comprueban el software durante un
periodo del tiempo comentando los varios fal-
los del programa que encuentran.

Cada usuario puede observar varios fallos

debidos al mismo error.

Puede que un usuario olvide a mandar el

informe sobre un error encontrado.

Al final del periodo, se depura el progra-

ma y se decide entre seguir comprobando
o lanzar el producto.

106
Un modelo para comprobación beta

El programa contiene N errores. Cada fallo

es debido a exactamente un error.

Para un usuario dado, el tiempo Sj hasta

observar el primer fallo debido al error j se
distribuye como

Sj |λj ∼ E(λj )
independientemente de los otros errores.

⇒ el tiempo, T para observar el primer fallo

T = mı́n{S1, . . . , SN }.
N
T |N, λ ∼ E(λ0) donde λ0 = λj .
j=1

107
 Si Z es un indicador de la causa de un fallo,
entonces
λj
P (Z = j|N, λ) = para j = 1, . . . , N
λ0

Z y T son independientes.

Los usuarios

M0 usuarios intercambiables comprueban el

software durante un tiempo T0.

A menudo, si se observa un fallo, se puede

identificar la causa, pero no siempre.

Supongamos: p = P (se identifica la causa

del fallo) para todos los fallos y todos los
usuarios.

108
Función de verosimilitud e inferencia

La verosimilitud es:
K

r−r0 r0 r
l(N, λ, p|data) = p (1 − p) λkk
k=0
exp (−M0T0λ0) para N ≥ K,
donde se han identificado K errores y rk fallos
debidos a estos errores, para k = 1, . . . , K.

r0 es el número de fallos no identificados y


r= K k=0 rk .

El EMV para N es N = K, lo que no sirve de

mucho. Mejor emplear inferencia bayesiana.

109
Distribuciones a priori

A menudo, se tiene información a priori sobre

N y λ0, pero no sobre los errores individuales.

Se reparameteriza en términos de λ0 y ρ donde

ρj = λj /λ0. Entonces:
N

r−r0 r0 r
l(p, N, λ0, ρ|datos) = p (1 − p) ρkk
k=1
r
λ0 exp (−M0T0λ0)

Distribuciones a priori p ∼ B(v, w) y λ0 ∼ G(a, b)

son conjugadas.

110
Distribuciones a priori de N y ρ

Supongamos que N ∼ P(θ). Consideramos dos

posibilidades para ρ.

1. Modelo fijo (F )

ρj = 1/N para j = 1, . . . , N .
Se supone que todos los errores tienen la
misma importancia, parecido al modelo de
Jelinski y Moranda (1972).

2. Distribución Dirichlet. (D)

ρ|N ∼ D(φ, . . . , φ) para algún φ > 0

Algunos errores pueden ser más grandes

que otros.

Se generaliza este modelo poniendo una

distribución jerárquica sobre φ.

111
Distribuciones a posteriori

p|datos ∼ B(v + r − r0, w + r0)

λ0 ∼ G(a + r, b + M0T0)

1. Modelo fijo.
θN 1
P (N |F , datos) ∝
N ! N r−r0
para N = K, K + 1, . . . Es fácil estimar la
media de errores en el programa etc. trun-
cando los sumatorios.

2. Modelo Dirichlet.

ρ|D, N, datos ∼ D(φ + r1, . . . , φ + rK , φ, . . . , φ)

θN Γ(N φ)
P (N |D, datos) ∝
N ! Γ(N φ + r − r0)

112
Un segundo periodo de comprobación

Después del primer periodo de ocupación, se

depura el software y se debe decidir seguir com-
probando o no.

Consideramos sólo la posibilidad de comprobar

con M1 usuarios durante un periodo T1. ¿Cómo
elegimos M1 y T1?

Definimos una función de costes y minimizamos

el coste esperado.

Observación 25 Elegimos la decisión Bayes,

es decir la decisión que maximiza la utilidad
esperada. Ver el capı́tulo 5.

113
Una función de costes

c1 = coste por cada usuario por unidad de

tiempo. (Sueldo)

Coste c2 por unidad de tiempo. (Coste de

oportunidad)

Coste c3 por cada error encontrado durante

el periodo. (Coste de depurar)

Coste c4 por cada fallo por unidad de tiem-

po después del lanzamiento del producto.

114
La función de costes es

C(M1, T1) = c1M1T1+ c2T1 + c3B 

N
+c4λ0  ρi(1 − I(i))
i=K+1
donde B es el número de errores encontrados
en el periodo de comprobación y I(i) es un
indicador de si se ha descubierto el error i.

Evaluación del coste esperado

N
Observamos que B = i=K+1 I(i) y entonces
N
C(M1, T1) = c1M1T1 + c2T1 + c4λ0 ρi
i=K+1
N
+ (c3 − c4λ0ρi)I(i)
i=K+1

Se necesita hallar el valor de E[I(i)|N, p, λ0, ρ].

¿Cómo hacerlo?
115
 E[I(i)|N, p, λ0, ρ] = P (se descubre el error i)
y condicionamos sobre el número de fallos F
observados en el periodo de comprobación.

Se tiene
 N f
j=K+1,j = i ρj
P (descubierto|F = f ) = 1−1 − p + p N 
j=K+1 ρj
Además
 
N
F |N, λ0, ρ, datos ∼ P λ0M1T1 ρi 
i=K+1
y, combinando estos dos resultados, se puede
hallar el coste esperado dados los parámetros
del modelo.

Ahora es posible sacar el coste predictivo ba-

jo cada modelo integrando sobre las distribu-
ciones a posteriori de p, λ0, ρ y N .
116
Coste esperado para el modelo fijo

= c1 M1T1 + c2 T1 +
∞
a+r N −K
c4 P (N |F, datos)
b + M0T0 N
N =K+1
∞ ∞
 f  
f
+ (N − K)P (N |F, datos) 1−
s
N =K+1 f =1 s=0
 s
B(v + r − r0 + s, w + r0 + f − s) N − K − 1
B(v + r − r0, w + r0) N −K
Γ(a + r + f ) f
pN (1 − pN )a+r
f !Γ(a + r)


a+r+f
c3 − c4 −K
N (b + M0T0 + M1T1 N N )
−K
M1 T1 N N
donde pN = N −K .
b+M0T0 +M1 T1 N

Existe una expresión parecida para el mode-

lo Dirichlet. (Se necesita integración numérica
unidimensional).

117
Ejemplo

Un programa contiene N = 20 errores con

tasa de fallos total λ0 = 10.

ρ generado de Dirichlet(1, . . . , 1).

Un usuario con periodo de comprobación

de 20 dı́as.

p = 0,9.

Fallos debidos a 14 errores distintos.

r0 = 15 fallos no identificados en r = 199.

Media a priori de N = 20. A prioris no in-

formativos para los demás parámetros.

118
Distribución a posteriori de N

Modelo
N F D H
14 1,0000 ,2836 ,8794
15 0 ,2241 ,1096
16 0 ,1676 ,0101
17 0 ,1189 ,0008
18 0 ,0804 ,0001
19 0 ,0518 ,0000
20 0 ,0320 ,0000
> 20 0 ,0416 0
E[N |datos] 14,0000 16,0256 14,1784

H es un modelo jerárquico, más avanzado.

119
Funciones de fiabilidad esperadas

Se supone la depuración anterior del software.

La fiabilidad es P (T > t) donde T = tiempo al
siguiente fallo.

1.00

Dirichlet
Fixed
0.95 Hierarchical

0.90
P(T>t|data)

0.85

0.80

0.75

0 20 40 60 80 100
t

120
Coste esperado de más comprobación

Consideramos hasta 3 usuarios durante un tiem-

po máximo de 20 dı́as. Costes c1 = 1, c2 =
c3 = 0,1, y c4 = 100000

1000
Dirichlet
Hierarchical
1 Fixed

800
Expected Cost

2
600
3

1
400
2

3
200
3
2
1
0

0 10 20 30 40
Test time

121
Extensiones

Aplicaciones a datos del internet. Ver Wiper

y Wislon (2007).

¿Cómo elegimos un modelo? Utilizamos fac-

tores Bayes. Ver el capı́tulo 5. Otra pos-
bilidad es hacer un promedio de las predic-
ciones mediante los distintos modelos (mod-
el averaging).

Usuarios no homogeneos.

¿Que hacer si no podemos fijar el valor de

M1 ?

122