Question
I am looking to find a simplification of the expression below. There seems to be a mistake in my attempt (steps shown below) using the Taylor (Maclaurin) series.
Could someone please verify if this approach is correct / provide additional pointers / alternative solutions?
Expression for Approximation
$$ \beta\left(x\right) = \left[x-\int_{0}^{x}\left[\frac{\Phi\left(\frac{lny-\mu}{\sigma}\right)}{\Phi\left(\frac{lnx-\mu}{\sigma}\right)}\right]^{M-1}dy\right] $$
Here, $\Phi(u)=\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{u}e^{-t^{2}/2}dt$ , the standard normal cumulative distribution and $X=e^{W}$where, $W\sim N\left(\mu,\sigma\right)$. $M$ is a non-zero positive integer. $\beta(x)=0$, when $x=0$
Steps Tried
\begin{eqnarray*} \beta\left(x\right) & = & \left[x-\frac{\int_{0}^{x}\left[\Phi\left(\frac{lny-\mu}{\sigma}\right)\right]^{M-1}dy}{\left[\Phi\left(\frac{lnx-\mu}{\sigma}\right)\right]^{M-1}}\right]\\ & = & \left[x-\frac{\left\{ \int_{0}^{x}\left[\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\left(\frac{lny-\mu}{\sigma}\right)}e^{-t^{2}/2}dt\right]^{M-1}dy\right\} }{\left[\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\left(\frac{lnx-\mu}{\sigma}\right)}e^{-t^{2}/2}dt\right]^{M-1}}\right] \end{eqnarray*} Let, \begin{eqnarray*} h\left(y\right) & = & \left[\int_{-\infty}^{\left(\frac{lny-\mu}{\sigma}\right)}e^{-t^{2}/2}dt\right]^{M-1}\\ j\left(y\right) & = & \int h\left(y\right)dy \end{eqnarray*} We then have, \begin{eqnarray*} \beta\left(x\right) & = & \left[x-\frac{\left\{ \int_{0}^{x}h\left(y\right)dy\right\} }{h\left(x\right)}\right]\\ & = & \left[x-\frac{\left|j\left(y\right)\right|_{0}^{x}}{h\left(x\right)}\right]=\left[x-\left\{ \frac{j\left(x\right)-j\left(0\right)}{h\left(x\right)}\right\} \right]\\ & \approx & \left[x-\frac{j'\left(0\right)x}{h\left(x\right)}\right]\;\left\{ \because j\left(x\right)-j\left(0\right)\simeq j'\left(0\right)x\;,\quad Maclaurin\; Series\right\} \\ & = & \left[x-\frac{h\left(0\right)x}{h\left(x\right)}\right]=x\left[1-\frac{h\left(0\right)}{h\left(x\right)}\right] \end{eqnarray*} \begin{eqnarray*} \Rightarrow\beta\left(x\right) & = & x\;\left[\;\because h\left(0\right)=\left[\int_{-\infty}^{\left(\frac{ln0-\mu}{\sigma}\right)}e^{-t^{2}/2}dt\right]^{M-1}=\left[\int_{-\infty}^{-\infty}e^{-t^{2}/2}dt\right]^{M-1}=0\right] \end{eqnarray*} We could include additional terms, for greater precision, using the subsequent terms of the Maclaurin series, as follows, \begin{eqnarray*} \beta\left(x\right) & \approx & x\left[1-\frac{h\left(0\right)}{h\left(x\right)}-\frac{x}{2}\frac{h'\left(0\right)}{h\left(x\right)}\right]\;\left\{ \because j\left(x\right)-j\left(0\right)\simeq j'\left(0\right)x+\frac{j''\left(0\right)x^{2}}{2!}\right\} \\ \beta\left(x\right) & = & x\left[1-\frac{h\left(0\right)}{h\left(x\right)}-\frac{1}{2}\left\{ 1-\frac{h\left(0\right)}{h\left(x\right)}\right\} \right]\;\left\{ \because\frac{h\left(x\right)-h\left(0\right)}{h\left(x\right)}\simeq\frac{h'\left(0\right)x}{h\left(x\right)}\right\} \\ \Rightarrow\beta\left(x\right) & = & \frac{x}{2}\;\left[\;\because h\left(0\right)=0\right] \end{eqnarray*} \begin{eqnarray*} \beta\left(x\right) & \approx & x\left[1-\frac{h\left(0\right)}{h\left(x\right)}-\frac{x}{2}\frac{h'\left(0\right)}{h\left(x\right)}-\frac{x^{2}}{6}\frac{h''\left(0\right)}{h\left(x\right)}\right]\\ & & \;\left\{ \because j\left(x\right)-j\left(0\right)\simeq j'\left(0\right)x+\frac{j''\left(0\right)x^{2}}{2!}+\frac{j'''\left(0\right)x^{3}}{3!}\right\} \end{eqnarray*} \begin{eqnarray*} \beta\left(x\right) & \approx & x\left[1-\frac{h\left(0\right)}{h\left(x\right)}-\frac{x}{2}\frac{h'\left(0\right)}{h\left(x\right)}-\frac{1}{3}\left\{ 1-\frac{h\left(0\right)}{h\left(x\right)}-x\frac{h'\left(0\right)}{h\left(x\right)}\right\} \right]\\ & & \;\left\{ \because\frac{h\left(x\right)-h\left(0\right)-h'\left(0\right)x}{h\left(x\right)}\simeq\frac{1}{2}\frac{h''\left(0\right)x^{2}}{h\left(x\right)}\right\} \end{eqnarray*} \begin{eqnarray*} \beta\left(x\right) & \approx & x\left[1-\frac{h\left(0\right)}{h\left(x\right)}-\frac{x}{2}\frac{h'\left(0\right)}{h\left(x\right)}-\frac{1}{3}+\frac{1}{3}\frac{h\left(0\right)}{h\left(x\right)}+\frac{x}{3}\frac{h'\left(0\right)}{h\left(x\right)}\right]\\ & = & x\left[\frac{2}{3}-\frac{2}{3}\frac{h\left(0\right)}{h\left(x\right)}-\frac{1}{6}\left\{ 1-\frac{h\left(0\right)}{h\left(x\right)}\right\} \right]\\ \Rightarrow\beta\left(x\right) & = & \frac{x}{2}\;\left[\;\because h\left(0\right)=0\right] \end{eqnarray*}
Reason for Concern
I checked this for correctness using this website: http://www.integral-calculator.com and Wolfram Alpha Definite Integral Calculator
The input in the text box under "Calculate the Integral of …" ((1+erf((ln(y)-5)/(sqrt(2)*10)))/(1+erf((ln(4)-5)/(sqrt(2)*10))))^24 with lower limit of integration "0" and upper limit "4" under options on the right with "Integrate numerically only?" checked.
I am getting different values and hence concerned that this theoretical approximation is incorrect. Also I tried to write R code, that gives a different value from the above theoretical approximation but matches the website.
So unsure what is the real problem with the above theoretical approximation.