# Evaluating an integral of geometric Brownian motion

I am a Mathematica novice. For an operations research application, I am trying to work with the following wealth process in Mathematica, $$d W_t = (\rho W_t + P_t) dt\ ,$$ where $$P_t$$ is a stream of payoffs that itself follows a geometric Brownian motion: $$d P_t = \mu P_t dt + \sigma P_t dZ_t\ .$$ Unfortunately, I cannot represent this problem as a standard Ito process and thus take advantage of Mathematica's rich in-built features. Just as in that case, I would like to be able to perform simulations, compute expected values of functions of terminal wealth, etc.

Here is what I've tried so far to compute, for example, the expected value.

WealthEvolve :=
Simplify[DSolve[{W'[t] == (ρ W[t] + P[t]), W[0] == w - k}, W[t],
t]]

FinalWealth[t_] := W[t] /. WealthEvolve[[1]]


This gives the expected output

E^(t ρ) (-k + w + Inactive[Integrate][E^(-ρ K[1]) P[K[1]], {K[1], 0, t}])


If I now go ahead and evaluate

Simplify[Activate[FinalWealth[t] /.
P -> Mean[GeometricBrownianMotionProcess[μ, σ, p]]]]


It does not evaluate, even though the expression looks correct.

E^(t ρ) (-k + w + \!$$\*SubsuperscriptBox[\(∫$$, $$0$$, $$t$$]$$\( \*SuperscriptBox[\(E$$, $$\(-ρ$$\ K[1]\)]\ $$Mean[ GeometricBrownianMotionProcess[μ, σ, p]]$$[
K[1]]\) \[DifferentialD]K[1]\)\))


Does anyone have any insight as to why?

I am also stuck when it comes to simulations. My intuition tells me that for any $$t$$, I need to need to sample from LogNormalDistribution and then evaluate the following integral for that $$t$$, but I am not sure how to do this. I would be grateful for any insight!

## 1 Answer

You are not replacing P correctly because you are forgeting the "time" argument. Namely, Mean[GeometricBrownianMotionProcess[μ, σ, p]] does not evaluate because it should be Mean[GeometricBrownianMotionProcess[μ, σ, p][t]].

FinalWealth[t] /.
P[t_] :> Mean[GeometricBrownianMotionProcess[μ, σ, p][t]] // Activate

(* E^(t ρ) (-k + w + ((-1 + E^(t (μ - ρ))) p)/(μ - ρ)) *)

• Thanks a lot for this! Would you have any idea how to simulate with this process? I am still unable to add code to the question, but I would like to wrap the code below in a "Simulate" function to generate paths over a specified horizon. Commented Apr 25, 2023 at 10:18
• FinalWealth[t] /. P -> GeometricBrownianMotionProcess[[Mu], [Sigma], p] /. {[Rho] -> 0.02, [Mu] -> 0.05, [Sigma] -> 0.2, w -> 100, k -> 10, p -> 0.5} Commented Apr 25, 2023 at 10:18