Supongamos que tenemos un fijo, conocido, $n$ y cada $x_1 \ldots x_n$ es un número aleatorio generado uniformemente sobre $[0,1]$ . ¿Cuál es la probabilidad de que $x_n$ es el mayor valor de la secuencia?
Respuestas
¿Demasiados anuncios?Mayúsculas $X$ denota variables aleatorias, en minúsculas $x$ denota realizaciones de variables aleatorias.
Si el realización de la secuencia se conoce hasta $n-1$ entonces sabemos cuál es el valor del máximo hasta entonces, digamos $x_{(n-1)}$ . Entonces queremos la probabilidad
$$\Pr\big(X_n > x_{(n-1)}) = 1-x_{(n-1)}$$
Si no conocemos las realizaciones reales, entonces la probabilidad que buscamos es la probabilidad de la diferencia de dos variables aleatorias independientes:
$$\Pr\big(X_n > X_{(n-1)} \big) = \Pr\big(X_{(n-1)} -X_n \leq 0 \big)$$
(como tenemos RVs continuas la aparición de la desigualdad débil para ajustarse a las expresiones estándar es neutral para los resultados).
La distribución de $X_n$ es $U(0,1)$ . $X_{(n-1)}$ es el estadístico de orden máximo de una muestra de $n-1$ independiente $U(0,1)$ variables aleatorias por lo que su distribución y funciones de densidad de probabilidad son
$$F_{X_{(n-1)}} (x_{(n-1)}) = x_{(n-1)}^{n-1},\;\;\; f_{X_{(n-1)}} (x_{(n-1)}) = (n-1)x_{(n-1)}^{n-2}$$
Derivando la distribución de su diferencia $Z = X_{(n-1)} -X_n$ requiere cierto cuidado. $Z$ rangos en $(-1,1)$ mientras que sus componentes definitorios son no negativos. En cambio, podemos optar por la probabilidad específica que buscamos: para que $Z\leq 0$ debe ser el caso que $X_{(n-1)} \leq X_n$ . La fórmula general para calcular dicha probabilidad es
$$P(Y < W) = \int_{S_w}\int_{\{y<w\}}f_{YW}(y,w) {\rm d}y{\rm d}w$$ que en nuestro caso se convierte (descomponiendo la densidad conjunta en el producto de los marginales debido a la independencia) en
$$\Pr(Z\leq 0) = \int_0^1\int_0^{x_n}(n-1)x_{(n-1)}^{n-2}\cdot 1 {\rm d}x_{(n-1)}{\rm d}x_n$$
$$= \int_0^1x_n^{n-1} {\rm d}x_n = \frac 1n$$
Que esta probabilidad disminuye con $n$ (y así con $n-1$ ) es intuitivo: cuantas más posibilidades tenga la secuencia de obtener el máximo teórico, más probable es que se acerque cada vez más a él, y por tanto menos probable es que el siguiente valor sea aún mayor.
Las simulaciones verifican el resultado anterior.
No estoy seguro de la configuración. Aquí hay algunas respuestas.
Si observa $x_1, x_2, \ldots, x_{n-1}$ y se pregunta si la próxima observación superará la serie, entonces la respuesta es $1-F(M)$ , donde $F$ es la función de distribución acumulativa y $M$ es el máximo de la primera $n-1$ observaciones.
Por otro lado, si sólo observa $x_n$ y se preguntan por la probabilidad de que este fuera el máximo de la serie, entonces la respuesta es $$F(x_n)^{n-1}$$ .
Y por otro lado, si no se observa nada en absoluto, pero se pregunta por la probabilidad de que el máximo se produzca en la última observación, entonces la probabilidad es $1/n$ .
descargo de responsabilidad: no soy estadístico ni algebrista.
Así que descubrí que hay que trabajar mucho para conseguir un gran número de muestras que sean secuencialmente mayores a partir de una distribución uniforme. Quizás Chebychev hable de esto.
Mi resultado fue esta trama:
Este es el resumen para un ajuste decente (la línea roja):
Call:
lm(formula = I(log10(y3)) ~ x)
Residuals:
Min 1Q Median 3Q Max
-0.13084 -0.06066 0.01368 0.03426 0.11696
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) -0.199568 0.056497 -3.532 0.00771 **
x 0.383633 0.009105 42.133 1.11e-10 ***
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.0827 on 8 degrees of freedom
Multiple R-squared: 0.9955, Adjusted R-squared: 0.995
F-statistic: 1775 on 1 and 8 DF, p-value: 1.11e-10
Eso significa que esta es la relación:
$ y = 10^{0.383\cdot x-0.2}$
Y aquí está mi código:
#my CV 132949
#housekeeping
rm(list = ls())
#I am going to draw a random number from a uniform distribution.
# I am then going to count how many other numbers I have to draw until one is bigger than it.
# I am going to do this 10 times in a row.
# I am going to repeat this 300 times
nbig <- 300
nsmall <- 10
nmax <- 300000
#pre-declare storage
data <- matrix(0,nrow=nsmall,ncol=nbig)
#big for
for (i in 1:nbig){
print(i)
#little for
for (j in 1:nsmall){
#if this is the first time
#draw a random number
if (j==1){
y <- 0.5
}
#if this is not the first time
#update the random number
if (j>1){
y <- z
}
#initialize z and n
z <- y
n <- 0
#draw random numbers until the next one is bigger
# recording how long it takes
while ((z <= y)&(n < nmax)) {
z <- runif(1);
n <- n+1;
}
#store in a variable
if (n < nmax){
data[j,i] <- n
} else {
#this is to retain the median value below
data[j,i] <- 1e12
}
}
}
y3 <- (rowMeans(data))
for (i in 1:nsmall){
y3[i] <- median(data[i,])
}
x <- 1:nsmall
plot(x,log10(y3),
xlab=c("number of terms in sequence"),
ylab=c("log_{10} median number new samples"))
est <- lm( I(log10(y3))~x )
summary(est)
lines(x,predict(lm( I(log10(y3))~x ) ),col="Red")
Observaciones:
Este era el valor medio, no un percentil superior. Si estuviera apostando, me inclinaría por el resultado más probable. Si estuviera construyendo un enfoque de seguridad empresarial, entonces me gustaría envolver todas las cosas malas y miraría el percentil 99 (o más) y no la mediana.