Monte Carlo Simulation of Merton’s Portfolio Problem
Summary
The document describes an attempt to check the optimal consumption and investment policy in Merton’s portfolio problem by simulating wealth paths and comparing average terminal utility with an analytical value function. It gives a discretized wealth update, recalculates consumption from current wealth at each step, and lists parameters for a sample run.
The author reports a large difference between the simulation and analytical result and asks whether the SDE implementation, consumption updates, and terminal condition are correct. The excerpt does not include the cited slides’ formulas or a resolution, so it does not establish which implementation detail causes the discrepancy. The code also uses both W0 and X0 for initial wealth, making it important to check that the simulated and analytical calculations start from the same value. Results from this Euler-style simulation may also depend on discretization and sampling error; the document does not quantify either.
Key ideas
- The simulation evolves wealth using investment returns and consumption over discrete time steps.
- Consumption is recalculated from current wealth and the time-dependent value function factor.
- The author compares mean terminal utility across simulated paths with an analytical value at time zero.
- The excerpt reports a mismatch but does not provide an answer or diagnose its source.
- The code uses different symbols for simulated and analytical initial wealth, which merits checking.
Tags
Full text
# Merton's Portfolio Problem Monte Carlo Simulation
# Merton's Portfolio Problem Monte Carlo Simulation
I am currently studying Merton's Problem Portfolio from the slides here: https://stanford.edu/~ashlearn/RLForFinanceBook/MertonPortfolio.pdf
Where they stated the optimal solution are:
I tried to do a simple Monte Carlo simulation using the optimal solutions and the SDE:
and compare it with the analytical formula for V* at time 0 as below
```
def Sim(mu, r, sigma, gamma, rho, k, X0, T, N, M): # N = number of timestep, M = sample size
dt = T / N
sqrt_dt = sqrt(dt)
normals = sqrt_dt * np.random.standard_normal(size=(N, M)) # generate a NxM matrix of (0,dt) normals
epsilon = k ** (1 / gamma)
pi_star = (mu - r) / (gamma * sigma * sigma)
alpha_factor = (((mu - r) ** 2) / (2 * sigma * sigma * gamma)) + r
alpha = (rho - (1 - gamma) * alpha_factor) / gamma
def f(t):
if alpha != 0.:
f_num = 1 + (alpha * epsilon - 1) * exp(-(alpha * (T - t)))
f = f_num / alpha
else:
f = T - t + epsilon
return f
W = np.ones(M) * W0
for n in range(N):
c_star = W / f(n * dt)
W += ((r + (mu - r) * pi_star) * W - c_star) * dt + sigma * pi_star * W * normals[n]
u_star = k * np.sign(W) * (np.abs(W) ** (1 - gamma)) / (1 - gamma) # terminal wealth
V0_star_MC = np.mean(u_star)
V0_star = (f(0) ** gamma) * (X0 ** (1 - gamma)) / (1 - gamma)
print(V0_star_MC)
print(V0_star)
```
However, I got a big difference between the Monte Carlo result and the analytical formula. I checked all the numbers and constants by hand and they seem to work fine. Is there anything wrong with the way I simulate the SDE and recompute C* in each loop? I also assume the terminal wealth to be calculated at each path is the boundary condition:
Is this correct? Thanks a lot!
I used the following parameters for my simulation:
```
mu = 0.15
r = 0.05
sigma = 0.2
gamma = 0.8
k = 0.8
rho = 0.1
W0 = 100
T = 1
N = 1000
M = 10000
```Shown in full with attribution under the source's licence. Licence: CC BY-SA 4.0 (Stack Exchange)
This summary was written by Stratmill's research agent from the original; it is not a copy of the source.