Comparing One-Period and Reoptimized Black–Scholes Strategies

# Comparing optimal strategies in Black-Scholes model with python


# Comparing optimal strategies in Black-Scholes model with python












Consider the 1-dimensional Black-Scholes model $$dS_t = S_t(\mu dt + \sigma dW_t)$$ $$dB_t = rB_t dt$$ Given a maturity time $T$, by martingale method the optimal discounted portfolio at the maturity over the interval $[0,T]$ is given by $$V_T = x + \frac{1}{\gamma} \mathbb{E}[Z_T \log Z_T] - \frac{1}{\gamma}\log Z_T$$ where $x$ is the initial wealth, $$Z_t = e^{\theta W_t - \frac{\theta^2}{2}t}$$ is the density of the unique risk neutral measure and $$\theta = - \frac{\mu-r}{\sigma}.$$ The optimization has taken is account the utility function $$u(x) = -\frac{1}{\gamma} e^{-\gamma x}.$$

Given two time $T_{mid} < T$, I'd like to compare two strategies:

- Simply optimizing on $[0,T]$. Let's call $V_1$ the portfolio at $T$.

- Optimizing between $[0,T_{mid}]$ and then re-optimizing on $[T_{mid}, T]$, with initial wealth the portfolio obtained at $T_{mid}$ with the first optimization. Let's call $V_2$ the final portfolio at $T$.

I expect that in general $$\mathbb{E}[u(V_2)] \leq \mathbb{E}[u(V_1)]$$ but this is not happening.

- Is my expecation wrong?

- Is my code wrong?

I put here my code

> `import numpy as np import matplotlib.pyplot as plt #BS parameters sigma = 0.2 mu = 1.5 r = 1 theta = -(mu-r)/sigma #investor parameters x = 1 gamma_0 = 1 #time T = 1 #final time N = 100 delta = T/N Nmid = 80 Tmid = delta*Nmid #mid time steps = np.arange(0, N+1) t = steps*delta ns = 10000 #brownian motion W = np.zeros((ns, N+1)) W[:,0] = 0 #log of density logZ = np.zeros((ns,N+1)) logZ[:,0] = 0 #compute density for n in steps[:-1]: logZ_n = logZ[:,n] W_n = W[:,n] dW = np.random.normal(loc=0, scale=np.sqrt(delta), size=ns) W[:,n+1] = W_n + dW logZ[:,n+1] = logZ_n - 0.5*theta**2*delta +theta*dW #optimal 2-step portfolio, first step mean_value = np.mean(logZ[:,Nmid]*np.exp(logZ[:,Nmid])) Vmid = x + 1/gamma_0*mean_value - 1/gamma_0*logZ[:,Nmid] #optimal 1-step portfolio mean_value = np.mean(logZ[:,-1]*np.exp(logZ[:,-1])) V1step = x + 1/gamma_0*mean_value - 1/gamma_0*logZ[:,-1] #re-compute logZ for the second step over [Tmid,T] logZ_mid = np.zeros((ns,N-Nmid+1)) logZ_mid[:,0] = 0 for n in range(0, N-Nmid): logZ_mid_n = logZ_mid[:,n] logZ_mid[:,n+1] = logZ_mid_n - 0.5*theta**2*delta +theta*(W[:,Nmid+n+1]-W[:,Nmid+n]) #optimal 2-step portfolio, ending step mean_value = np.mean(logZ_mid[:,-1]*np.exp(logZ_mid[:,-1])) V2step = Vmid + 1/gamma_0*mean_value - 1/gamma_0*logZ_mid[:,-1] #comparing def utility(x): return -1/gamma_0*np.exp(-gamma_0*x) print(f"2-step: {np.mean(utility(V2step))}") print(f"1-step: {np.mean(utility(V1step))}") `

```
import numpy as np
import matplotlib.pyplot as plt

#BS parameters
sigma = 0.2
mu = 1.5
r = 1
theta = -(mu-r)/sigma

#investor parameters
x = 1
gamma_0 = 1

#time
T = 1 #final time
N = 100
delta = T/N
Nmid = 80
Tmid = delta*Nmid #mid time
steps = np.arange(0, N+1)

t = steps*delta
ns = 10000

#brownian motion
W = np.zeros((ns, N+1))
W[:,0] = 0

#log of density
logZ = np.zeros((ns,N+1))
logZ[:,0] = 0

#compute density
for n in steps[:-1]:

   logZ_n = logZ[:,n]
   W_n = W[:,n]

   dW = np.random.normal(loc=0, scale=np.sqrt(delta), size=ns)
   W[:,n+1] = W_n + dW

   logZ[:,n+1] = logZ_n - 0.5*theta**2*delta +theta*dW

#optimal 2-step portfolio, first step
mean_value = np.mean(logZ[:,Nmid]*np.exp(logZ[:,Nmid]))
Vmid = x + 1/gamma_0*mean_value - 1/gamma_0*logZ[:,Nmid]

#optimal 1-step portfolio
mean_value = np.mean(logZ[:,-1]*np.exp(logZ[:,-1]))
V1step = x + 1/gamma_0*mean_value - 1/gamma_0*logZ[:,-1]

#re-compute logZ for the second step over [Tmid,T]
logZ_mid = np.zeros((ns,N-Nmid+1))
logZ_mid[:,0] = 0

for n in range(0, N-Nmid):

   logZ_mid_n = logZ_mid[:,n]
   logZ_mid[:,n+1] = logZ_mid_n - 0.5*theta**2*delta +theta*(W[:,Nmid+n+1]-W[:,Nmid+n])

#optimal 2-step portfolio, ending step
mean_value = np.mean(logZ_mid[:,-1]*np.exp(logZ_mid[:,-1]))
V2step = Vmid + 1/gamma_0*mean_value - 1/gamma_0*logZ_mid[:,-1]

#comparing
def utility(x):
   return -1/gamma_0*np.exp(-gamma_0*x)

print(f"2-step: {np.mean(utility(V2step))}")
print(f"1-step: {np.mean(utility(V1step))}")
```

Thank you for your help!