Exercise: Realised vs implied — hedging error simulation

The continuous-time hedging-error formula states

P\&L = \tfrac{1}{2}\int_0^T \Gamma_t S_t^2 (\sigma_{r,t}^2 - \sigma_{h}^2)\,dt.

Implement a discrete-time simulation and verify this relation numerically.

Parameters: $S_0 = 100$ , $K = 100$ , $r = 0$ , $T = 0.25$ (3 months), hedge vol $\sigma_h = 0.20$ , $n = 63$ daily rebalances. Consider realised vols $\sigma_r \in \{0.15, 0.20, 0.25\}$ , constant across the horizon.

Tasks

Implement a simulator: generate underlying paths under GBM with realised vol $\sigma_r$ ; at each rebalance time, set the hedge portfolio to $\Delta^{BS}(S_t, t; \sigma_h)$ shares; aggregate P&L over the path and across $M = 5000$ paths.
For each $\sigma_r$ , compute (a) mean P&L per path, and (b) standard deviation of P&L per path.
Compute the theoretical expected P&L using the formula. For $\sigma_r = 0.20$ (same as hedge vol) the expected P&L is $0$ . For $\sigma_r = 0.25$ , estimate the integral $\tfrac12 \int \mathbb{E}[\Gamma_t S_t^2] (0.25^2 - 0.20^2)\, dt$ (use the fact that under risk-neutral GBM with $r = 0$ , $\mathbb{E}[\Gamma_t S_t^2]$ has a closed-form or compute it numerically).
Compare simulated mean P&L with theoretical. They should agree within Monte Carlo error.
Why is the standard deviation of P&L much larger than the absolute mean when $|\sigma_r - \sigma_h| = 0.05$ ? Interpret the result as "options are a noisy vol play."

Hint

$\mathbb{E}[\Gamma_t S_t^2]$ can be computed from $\Gamma_t = \varphi(d_1(t))/(S_t\sigma_h\sqrt{T-t})$ and Monte Carlo averaging.