Let $\mu$ and $\sigma > 0$ and $\beta_1 \ge 0 $ and $\beta_2 \ge 0$ be real numbers. Consider a stochastic process $X_t$ that satisfies the following stochastic differential equation:
\begin{equation} d X_t = \mu X_t^{\beta_1} dt + \sigma X_t^{\beta_2} dB_t \tag{1} \end{equation} where $B_t$ is the Brownian motion.
My objective is to derive the transition probability density function $P_x\left( X_t \in dz\right) := P\left( X_t \in (z,z+dz) | X_0 = x \right)$ for the process in question. Here we will be using a result known from literature (see section 4.11 pages 149-161 in K Ito& H P McKean Jr, "Diffusion Processes and their Sample Paths") as below.
Theorem:
The Laplace transform (with respect to time) of the transition probability density in question reads:
\begin{equation} \int\limits_0^\infty e^{-\lambda t} P_x\left( X_t \in dz\right) dt = G_\lambda(z,x) \cdot m(dz) \tag{2} \end{equation}
Here $m(z):= \left( 1/2 \sigma^2(z) S(z) \right)^{-1}$ is the speed measure with $S(y) := \exp\left(- \int\limits_1^y 2 \mu(\xi)/\sigma^2(\xi) d\xi \right)$ being the scale measure.
Then $G_\lambda(z,x)$ is the Green's function of the operator $\lambda - {\mathfrak G}^{*}$ with ${\mathfrak G}^{*}$ being the infinitesimal generator, i.e. such an operator that its action on a trial function $f$ reads $\left({\mathfrak G}^{*} f\right)(x) = \lim\limits_{\epsilon \rightarrow 0} E \left[ f(X_\epsilon) - f(x) \right]/\epsilon$ subject to $X_0=x$.
Note that the Green's function is easily determined by the eigen-functions and it reads $G_\lambda(z,x) = S(x)/{\mathfrak W}_\lambda(x) \left( G_\lambda(x) F_\lambda(z) 1_{z < x} + F_\lambda(x) G_\lambda(z) 1_{z > x} \right)$ where $F_\lambda(z)$ and $G_\lambda(z)$ is a strictly increasing and a strictly decreasing eigen-function of the infinitesimal generator to an eigenvalue $\lambda$ and ${\mathfrak W}_\lambda(x) := F_\lambda^{'}(x) G_\lambda(x) - G_\lambda^{'}(x) F_\lambda(x)$ is the Wronskian of the generator in question.
What have I managed to achieve?:
I have found the result for $(\beta_1,\beta_2) = (1, \beta)$ with $\beta >1$ and I present the derivation below.
Here the infinitesimal generator reads ${\mathfrak G}_z := \mu z \cdot d/dz + 1/2 \sigma^2 z^{2 \beta} d^2/ d z^2$ and its eigenfunctions (compare my previous question on a similar topic) read:
\begin{eqnarray} F(z) &=& U\left(\frac{\lambda}{2(-1+\beta) \mu}, 1+ \frac{1}{2(-1+\beta)}, \frac{\mu z^{2-2\beta} }{(-1+\beta) \sigma^2} \right) \tag{3a} \\ G(z) &=& F_{1,1}\left(\frac{\lambda}{2(-1+\beta) \mu}, 1+ \frac{1}{2(-1+\beta)}, \frac{\mu z^{2-2\beta} }{(-1+\beta) \sigma^2} \right) \tag{3b} \end{eqnarray} where $U$ is the confluent hypergeometric function.
Note that $F^{'}(z) = \frac{\lambda z^{1-2 \beta } U\left(\frac{\lambda }{2 (\beta -1) \mu }+1,2+\frac{1}{2 (\beta -1)},\frac{z^{2-2 \beta } \mu }{s^2 (\beta -1)}\right)}{(\beta -1) s^2} >0$ and $G^{'}(z) = -\frac{2 \lambda z^{1-2 \beta } \, _1F_1\left(\frac{\lambda }{2 (\beta -1) \mu }+1;2+\frac{1}{2 (\beta -1)};\frac{z^{2-2 \beta } \mu }{s^2 (\beta -1)}\right)}{(2 \beta -1) s^2} < 0$ (note that from the integral representations we see that both $U$ and $F_{1,1}$ are positive --well up to a multiplicative constant in the second case) and as such $F(z)$ and $G(z)$ are strictly increasing and strictly decreasing as it should be.
Now the Wronskian is computed as follows:
\begin{eqnarray} {\mathfrak W}_\lambda(x) &=& {\mathfrak W}_\lambda(1) \cdot \exp\left( -\int\limits_1^x \frac{2\mu(\xi)}{\sigma^2(\xi)} d\xi \right) \tag{4a} \\ &=& {\mathfrak W}_\lambda(1) \cdot \exp\left( - \frac{2 \mu}{\sigma^2} \frac{x^{2-2\beta} -1}{2-2\beta}\right) \tag{4b} \\ &=& -\frac{2 (\beta -1) \Gamma \left(1+\frac{1}{2 (\beta -1)}\right) e^{\frac{\mu }{(\beta -1) \sigma ^2}} \left(\frac{\mu }{(\beta -1) \sigma ^2}\right)^{\frac{1}{2-2 \beta }}}{\Gamma \left(\frac{\lambda }{(2 \beta -2) \mu }\right)} \cdot \exp\left( - \frac{2 \mu}{\sigma^2} \frac{x^{2-2\beta} -1}{2-2\beta}\right) \tag{4c} \end{eqnarray} The first line above simply follows form the first order differential equation being satisfied by the Wronskian. The second line follows form the fact that $\mu(\xi), \sigma(\xi) = \mu \cdot \xi, \sigma \cdot \xi^\beta$. The step to the final line is the most laborious, but not really that mind boggling. Here one employs the integral representations in order to express the Wronskian at unity as a double integral in which one changes variables to obtain the neat final result.
Now we compute the inverse Laplace transform as the Bromwich integral and we use Cauchy theorem (being applied to the usual contour in the negative complex half-plane-- the second & the third quadrants) to evaluate that integral. As seen from $(4c)$ the Laplace transform has infinitely many simple poles at $\lambda_n = (2\beta-2) \mu (-n)$ for $n=0,1,2,\cdots$. Then the final result reads:
\begin{eqnarray} &&\left. P_x \left( X_t \in dz \right)/dz = \frac{2 \mu \left(\frac{\mu }{(\beta -1) \sigma ^2}\right)^{\frac{1}{2 (\beta -1)}}}{\sigma ^2 \Gamma \left(1+\frac{1}{2 (\beta -1)}\right)} \cdot z^{-2 \beta } e^{-\frac{\mu z^{2-2 \beta }}{(\beta -1) \sigma ^2}} \cdot \right.\\ && \left. \sum\limits_{n=0}^\infty \frac{(-1)^n}{n!} % \, _1F_1\left(-n;1+\frac{1}{2 (\beta -1)};\frac{\mu (x \vee z)^{2-2 \beta }}{(\beta -1) \sigma ^2}\right) U\left(-n,1+\frac{1}{2 (\beta -1)},\frac{\mu (x \wedge z)^{2-2 \beta }}{(\beta -1) \sigma ^2}\right) % \cdot e^{-2(-1+\beta) \cdot \mu \cdot n \cdot t} \right. \tag{5} \end{eqnarray}
Now we have checked numerically that the result above matches with the solution known from literature (see equations (4) and (5) in V Linetsky, R Mendoza, The Constant Elasticity of Variance Model). Here we go:
In[1991]:= (*Verify whether our solution matches that known from \
literature?*)
Clear[mu, s, b, z, lmb, phi, v, W, W1, W0, x];
{mu, s, b} = {1/2, 1/3 + 3/10, 5/2};
{x, z} = RandomReal[{1, 3}, 2, WorkingPrecision -> 50];
t = 3/2;
(*Our solution:*)
res1 = (2 mu ((mu/((-1 + b) s^2))^((1/(2 (-1 + b))))) )/(
s^2 Gamma[1 + 1/(2 (-1 + b))])
E^(-((mu z^(2 - 2 b))/((-1 + b) s^2))) z^(-2 b)
Table[(-1)^n/
n! HypergeometricU[-n, 1 + 1/(2 (-1 + b)),
mu/((-1 + b) s^2) Min[x, z]^(2 - 2 b)] Hypergeometric1F1[-n,
1 + 1/(2 (-1 + b)),
mu/((-1 + b) s^2) Max[x, z]^(2 - 2 b)] E^(-
2 (-1 + b) mu n t), {n, 0, M}];
Total[res1] // N
(Solution known from literature: Equations (4) & (5) page 2 in V
Linetsky & R Mendoza The Constant Elasticity of Variance Model.Here
the solution was obtained,using different techniques,meaning by
re-scaling the current price and by changing time in the solution of
a simpler problem,i.e.a problem with the drift parameter being set to
zero.)
Clear[tau]; tau[t_] := (Exp[2 mu (b - 1) t] - 1)/(2 mu (b - 1));
res2 = Exp[-mu t] ((Exp[-mu t] z)^(-2 b + 1/2) x^(1/2))/(
s^2 Abs[b - 1] tau[t])
BesselI[1/(2 Abs[ b - 1]), (x^(-b + 1) (Exp[-mu t] z)^(-b + 1))/(
s^2 (b - 1)^2 tau[
t])] Exp[-((x^(-2 b + 2) + (Exp[-mu t] z)^(-2 b + 2))/(
2 s^2 (b - 1)^2 tau[t]))];
res2 // N
Out[1996]= 0.0481819
Out[1999]= 0.0481819
Having said all this my question is how do we find the solution in the generic case, i.e. when $\beta_1 \neq 1$? Another obvious question would be to provide a literature reference, if it exists, on that problem solution.
