Suma de variables aleatorias truncadas normales

Supongamos que tengo variables aleatorias normales independientes $n$

X_{1} \sim norte (μ_{1}, σ_{1}^{2}) X_{2} \sim norte (μ_{2}, σ_{2}^{2}) ⋮ X_{norte} \sim norte (μ_{norte}, σ_{norte}^{2})

$X_1 \sim \mathrm{N}(\mu_1, \sigma_1^2)\\X_2 \sim \mathrm{N}(\mu_2, \sigma_2^2)\\\vdots\\X_n \sim \mathrm{N}(\mu_n, \sigma_n^2)$

e . ¿Cómo caracterizaría la densidad de si la distribución de cada se trunca dentro de ? En otras palabras, estoy muestras de distribuciones normales independientes, descartando muestras que no dentro de de cada media y sumando. $Y=X_1+X_2+\dotsm+X_n$ $Y$ $X_i$ $(\mu_i - 2\sigma_i, \mu_i + 2\sigma_i)$ $n$ $2\sigma_i$

En este momento, estoy haciendo esto con el código R a continuación:

x_mu <- c(12, 18, 7)
x_sd <- c(1.5, 2, 0.8)
a <- x_mu - 2 * x_sd
b <- x_mu + 2 * x_sd

samples <- sapply(1:3, function(i) {
  return(rtruncnorm(100000, a[i], b[i], x_mu[i], x_sd[i]))
})

y <- rowSums(samples)

¿Hay algún método para generar la densidad de directamente? $Y$

r normal-distribution pdf truncation saddlepoint-approximation Devin
fuente

Su pregunta implica que conoce todos los . ¿Es ese realmente el caso o los estás estimando ? ¡Hay una gran diferencia! Por curiosidad, ¿por qué estás tirando esos datos? Dependiendo de sus objetivos, sospecho que existen (mucho) mejores procedimientos.

σ_{i}

$\sigma_i$

whuber

Sí conozco todos los medios y SD para mis datos, sí.

Devin

Creo que podrías caracterizarlo como "un desastre". Este artículo, jstor.org/stable/2236545 , examina el asunto con más rigor científico.

Alecos Papadopoulos

Aproximación externa a través de CLT, esto es relativamente complicado. Supongo que si es lo suficientemente pequeño, podría intentar la convolución numérica.

n

$n$

Glen_b -Reinstate Monica

@Silverfish Dependiendo de la implementación, plataforma y cuán fina es tolerable una grilla, cientos deberían estar bien (quizás más); Además de la velocidad, sin embargo, con suficientes términos, debe ser mucho más cuidadoso con los detalles de implementación o pueden surgir una serie de problemas numéricos.

Glen_b -Reinstate Monica

Respuestas:

Puede usar la aproximación por métodos de punto de silla, para la suma de las normales truncadas. No daré los detalles ahora, puede ver mi respuesta a la suma general de distribuciones gamma para obtener pistas. Lo que necesitamos es encontrar la función generadora de momento para una normal truncada, lo cual es fácil. Lo haré aquí para un estándar normal truncado en , que tiene densidad donde aquí son densidad y cdf para una normal estándar, respectivamente. $\pm 2$

F (X) = {\begin{cases} \frac{1}{C} ϕ (X), & El | X El | \leq 2 \\ 0 0, & El | X El | > 2 \end{cases}

$f(x) =\begin{cases} \frac1{C} \phi(x), & |x| \le 2 \\ 0, & |x| > 2 \end{cases}$

C = Φ (2) - Φ (- 2)

$C=\Phi(2) - \Phi(-2)$

ϕ (x), Φ (x)

$\phi(x), \Phi(x)$

La función de generación de momentos se puede calcular como y luego podemos usar aproximaciones de punto de silla de montar.

METRO (t) = mi {mi}^{t X} = \frac{1}{C} \int_{- 2}^{2} {mi}^{t X} ϕ (X) re X = \frac{1}{C} {mi}^{\frac{1}{2} t^{2}} [Φ (2 - t) - Φ (- 2 - t)]

$\DeclareMathOperator{\E}{\mathbb{E}} M(t) = \E e^{tX}=\frac1{C}\int_{-2}^2 e^{tx} \phi(x)\; dx=\frac1{C}e^{\frac12 t^2} [\Phi(2-t)-\Phi(-2-t) ]$

kjetil b halvorsen
fuente

-3

Tengo curiosidad por qué, pero sí, hay una manera simple de generar el pdf de esta suma de distribuciones:

## install.packages("truncnorm")
## install.packages("caTools")
library(truncnorm)

x.mu <- c(12, 18, 7)
x.sd <- c(1.5, 2, 0.8)
x.a <- x.mu - 2*x.sd
x.b <- x.mu + 2*x.sd

dmulti <- function(x, a, b, mu, sd)
  rowSums(
    sapply(1:length(mu),
           function(idx)
             dtruncnorm(x, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)
pmulti <- function(q, a, b, mu, sd)
  rowSums(
    sapply(1:length(mu),
           function(idx)
             ptruncnorm(q, a=a[idx], b=b[idx], mean=mu[idx], sd=sd[idx])))/length(mu)

pointrange <- range(c(x.a, x.b))
pointseq <- seq(pointrange[1], pointrange[2], length.out=100)
## Plot the probability density function
plot(pointseq, dmulti(pointseq, x.a, x.b, x.mu, x.sd),
     type="l")

## Plot the cumulative distribution function
plot(pointseq, pmulti(pointseq, x.a, x.b, x.mu, x.sd),
     type="l")

Bill Denney
fuente

Si leo este código correctamente, parece que estás implementando algo como una mezcla en lugar de una suma. La trama que produce este código es lamentablemente incorrecta. ¡Ni siquiera es una función de densidad de probabilidad válida!

whuber

@whuber, gracias por la captura. Normalicé el pdf y agregué el cdf.

Bill Denney

Gracias. Sin embargo, el error básico persiste: está calculando una distribución de mezcla en lugar de la suma.

whuber