¿Es posible integrar analíticamente

En primer lugar, por integración analítica, quiero decir, ¿hay una regla de integración para resolver esto en oposición a los análisis numéricos (como las reglas trapezoidales, Gauss-Legendre o Simpson)?

Tengo una función donde es la función de densidad de probabilidad de una distribución lognormal con parámetros y . A continuación, abreviaré la notación a y usaré para la función de distribución acumulativa.$\newcommand{\rd}{\mathrm{d}}f(x) = x g(x; \mu, \sigma)$

$$g(x; \mu, \sigma) = \frac{1}{\sigma x \sqrt{2\pi}} e^{-\frac{1}{2\sigma^2}(\log(x) - \mu)^2}$$ $\mu$$\sigma$$g(x)$$G(x)$

Necesito calcular la integral

$$\int_{a}^{b} f(x) \,\rd x \>.$$

Actualmente, estoy haciendo esto con integración numérica usando el método Gauss-Legendre. Como necesito ejecutar esto muchas veces, el rendimiento es importante. Antes de analizar la optimización de los análisis numéricos / otras piezas, me gustaría saber si hay alguna regla de integración para resolver esto.

Intenté aplicar la regla de integración por partes, y llegué a esto, donde estoy atascado nuevamente,

1. $\int u \,\mathrm{d}v = u v - \int v \mathrm{d}u$ .

2. $u=x \implies \rd u = \rd x$

3. $\rd v = g(x) \rd x \implies v = G(x)$

4. $u v - \int v \rd x = x G(x) - \int G(x) \rd x$

Estoy atascado, ya que no puedo evaluar el .$\int G(x) \rd x$

Esto es para un paquete de software que estoy construyendo.

Respuesta corta : No, no es posible, al menos en términos de funciones elementales. Sin embargo, existen algoritmos numéricos muy buenos (¡y razonablemente rápidos!) Para calcular dicha cantidad y, en este caso, deberían preferirse a cualquier técnica de integración numérica.

Cantidad de interés en términos de cdf normal

La cantidad que le interesa en realidad está estrechamente relacionada con la media condicional de una variable aleatoria lognormal. Es decir, si se distribuye como lognormal con los parámetros y , entonces, usando su notación, $X$$\mu$$\sigma$

$$\newcommand{\e}{\mathbb{E}}\renewcommand{\Pr}{\mathbb{P}}\newcommand{\rd}{\mathrm{d}} \int_a^b f(x) \rd x = \int_a^b \frac{1}{\sigma \sqrt{2\pi}} e^{-\frac{1}{2\sigma^2}(\log(x) - \mu)^2} \rd x = \Pr(a \leq X \leq b) \e(X \mid a \leq X \leq b) \>.$$

Para obtener una expresión para esta integral, realice la sustitución . Al principio, esto puede parecer un poco desmotivado. Pero, tenga en cuenta que usando esta sustitución, y simplemente cambiando las variables, obtenemos donde y .$z = (\log(x) - (\mu + \sigma^2))/\sigma$$x = e^{\mu + \sigma^2} e^{\sigma z}$

$$\int_a^b f(x) \rd x = e^{\mu + \frac{1}{2}\sigma^2} \int_{\alpha}^{\beta} \frac{1}{\sqrt{2\pi}} e^{-\frac{1}{2} z^2} \rd z \> ,$$ $\alpha = (\log(a) - (\mu + \sigma^2))/\sigma$$\beta = (\log(b) - (\mu + \sigma^2))/\sigma$

Por lo tanto, donde es el estándar Función de distribución acumulativa normal.

$$\int_a^b f(x) \rd x = e^{\mu + \frac{1}{2}\sigma^2} \big( \Phi(\beta) - \Phi(\alpha) \big) \>,$$ $\Phi(x) = \int_{-\infty}^x \frac{1}{\sqrt{2\pi}} e^{-z^2/2} \rd z$

Aproximación numérica

A menudo se afirma que no existe una expresión de forma cerrada conocida para . Sin embargo, un teorema de Liouville de principios de 1800 afirma algo más fuerte: no existe una expresión de forma cerrada para esta función . (Para la prueba en este caso particular, vea el artículo de Brian Conrad ).$\Phi(x)$

Por lo tanto, nos queda usar un algoritmo numérico para aproximar la cantidad deseada. Esto se puede hacer dentro del punto flotante de doble precisión IEEE a través de un algoritmo de WJ Cody's. Es el algoritmo estándar para este problema y, al utilizar expresiones racionales de un orden bastante bajo, también es bastante eficiente.

Aquí hay una referencia que analiza la aproximación:

WJ Cody, aproximaciones racionales de Chebyshev para la función de error , matemática. Comp. , 1969, pp.631-637.

También es la implementación utilizada tanto en MATLAB como en , entre otros, en caso de que sea más fácil obtener un código de ejemplo.$R$

Aquí hay una pregunta relacionada, en caso de que esté interesado.