9  The Black-Scholes Formula

In the previous chapter, we developed a general framework for derivative pricing based on the absence of arbitrage. We saw how linear pricing leads to state prices, risk-neutral probabilities, and the fundamental pricing formula. We also introduced change of numeraire methods, which allow us to compute option values as expectations under different probability measures. In this chapter, we apply these principles to derive and analyze the celebrated Black-Scholes formula for European option pricing.

A fundamental insight is that standard European options can be expressed as combinations of digital options. Understanding this connection provides both theoretical insight and practical computational advantages.

The Black-Scholes model assumes that the underlying asset pays a constant dividend yield \(q \,\mathrm{d} t\) (in other words the stock pays the dividend \(q S_t \,\mathrm{d}t\)) and has an ex-dividend price \(S\) satisfying

\[\frac{\mathrm{d}S}{S} = (\mu - q) \,\mathrm{d}t + \sigma \,\mathrm{d}B\,.\]

for a Brownian motion \(B\), where \(\sigma\) and \(\mu\) are assumed to be constant. We also assume a constant continuously-compounded risk-free rate \(r\). In this model, \(\mathrm{d}S\) represents the capital gain and \(q S_t \,\mathrm{d}t\) is the dividend portion of the return.

9.1 Digital Options as Building Blocks

Before deriving the Black-Scholes formula, we need to understand how standard options relate to digital options. A digital option (also called binary option) has a discontinuous payoff — it pays a fixed amount if a condition is met and zero otherwise.

Two Types of Digitals

Consider a call option with strike \(K\) and maturity \(T\). Its payoff \((S_T - K)^+\) can be decomposed as

\[(S_T - K)^+ = \mathbf{1}_{\{S_T > K\}} S_T - K \mathbf{1}_{\{S_T > K\}}\,.\]

This reveals that a call option is the difference between two digital options:

1. Share digital: Pays \(S_T\) when \(S_T > K\), zero otherwise
2. Cash digital: Pays \(K\) when \(S_T > K\), zero otherwise

Similarly, a put option with payoff \((K - S_T)^+\) can be written as

\[(K - S_T)^+ = K \mathbf{1}_{\{S_T < K\}} - \mathbf{1}_{\{S_T < K\}} S_T\,.\]

This is \(K\) cash digitals minus share digitals, both paying when \(S_T < K\).

9.2 Girsanov’s Theorem in the Black-Scholes Model

To value digital options, we need to understand changes of probability measures in continuous time. When we change probability measures, a process \(B\) that was a Brownian motion cannot be expected to remain a Brownian motion under the new measure. However, Girsanov’s Theorem tells us exactly how the drift changes.

In the Black-Scholes model, the dividend-reinvested stock price, denoted \(Y_t\), evolves as

\[\mathrm{d}Y_t = \mu Y_t \,\mathrm{d}t + \sigma Y_t \,\mathrm{d}B_t\,.\]

We can rewrite this as

\[\mathrm{d}Y_t = r Y_t \,\mathrm{d}t + \sigma Y_t \,\mathrm{d}\left(B_t + \frac{\mu - r}{\sigma}t\right)\,.\]

Define \(\kappa = \frac{\mu - r}{\sigma}\). This is called the market price of risk or Sharpe ratio.

Important Principle

Girsanov’s Theorem: Define the exponential martingale:

\[Z_t = \exp\left(-\frac{1}{2}\kappa^2 t - \kappa B_t\right)\,.\]

Using the rule for expectations of lognormal variables, we can see that \(\mathbb{E}[Z_T] = 1\). Under the probability measure defined by:

\[\mathbb{E}^\kappa[Y] = \mathbb{E}[Z_T X]\]

for any function \(X\) of \(B_t\) (\(0 \leq t \leq T\)), the process \(B_t^\kappa = B_t + \kappa t\) is a Brownian motion. Moreover, if \(M_t\) is a martingale under the new measure, then \(Z_t M_t\) is a martingale under the original measure.

This theorem tells us how to construct the risk-neutral measure from the original measure. Using Ito’s lemma, we can verify that \(Z_t \frac{Y_t}{R_t}\) is indeed a martingale:

\[Z_t \frac{Y_t}{R_t} = S_0 \exp\left(-\frac{1}{2}(\sigma - \kappa)^2 t + (\sigma - \kappa)B_t\right)\,.\]

Therefore this is a lognormal random variable and taking expectations gives:

\[\mathbb{E}\left[Z_t \frac{Y_t}{R_t}\right] = S_0\,.\]

Dividend Reinvested Stock as Numeraire

When using the stock as numeraire, we need the process \(\frac{R_t}{Y_t}\) to be a martingale. Following similar analysis with Girsanov’s Theorem, we find that under the stock numeraire measure:

\[\mathrm{d}Y_t = (r + \sigma^2) Y_t \,\mathrm{d}t + \sigma V_t \,\mathrm{d}B_t^Y\,.\]

where \(B_t^Y = B_t^R - \sigma t\) is a Brownian motion under the dividend reinvested stock numeraire measure. It then follows \[\frac{R_t}{Y_t} = \frac{R_0}{S_0} \exp\left(-\frac{1}{2} \sigma^2 t - \sigma B_t^Y \right)\,, \]

which is an exponential martigale.

9.3 Derivation of the Black-Scholes Formula

Now we apply the valuation theory from Chapter 3 (Arbitrage Pricing) and the tail probability calculations from Chapter 2 (Geometric Brownian Motion) to derive the Black-Scholes formula.

Cash Digital Options

A cash digital option paying $1 when \(S_T > K\) has value under the risk-neutral measure:

\[e^{-r(T-t)} \mathbb{E}_t^R[\mathbf{1}_{\{S_T > K\}}] = e^{-r(T-t)} \text{prob}^R(S_T > K)\,.\]

Under the risk-neutral measure, the stock price follows:

\[\frac{\mathrm{d}S}{S} = (r - q) \,\mathrm{d}t + \sigma \,\mathrm{d}B^R\,.\]

This means:

\[\mathrm{d}\log S = \left(r - q - \frac{\sigma^2}{2}\right) \,\mathrm{d}t + \sigma \,\mathrm{d}B^R\,.\]

From Chapter 2, we know this has the solution:

\[S_T = S_t \exp\left(\left(r - q - \frac{\sigma^2}{2}\right)(T-t) + \sigma B_{T-t}^R\right) \,.\]

Applying the tail probability formula from Section 7.2 with \(\alpha = r - q - \sigma^2/2\) gives us

\[\text{prob}^R(S_T > K) = N(d_2)\,,\]

where

\[d_2 = \frac{\log(S_t/K) + (r - q - \frac{\sigma^2}{2})(T-t)}{\sigma\sqrt{T-t}}\,.\]

Share Digital Options

For the share digital paying \(S_T \mathbf{1}_{\{S_T > K\}}\), we use the dividend reinvested stock as numeraire. We have

\[S_T = e^{-q T} Y_T\,; \] therefore, \[\frac{S_T}{Y_T} \mathbf{1}_{\{S_T > K\}} = e^{-q T} \mathbf{1}_{\{S_T > K\}}\,.\]

From the change of numeraire theory in Chapter 3:

\[\text{Value} = Y_t \mathbb{E}_t^V[e^{-q T} \mathbf{1}_{\{S_T > K\}}] = e^{-q (T-t) } S_t \cdot \text{prob}^Y(S_T > K)\,.\]

Under the stock numeraire measure,

\[\mathrm{d}\log S = \left(r - q + \frac{\sigma^2}{2}\right) \,\mathrm{d}t + \sigma \,\mathrm{d}B^Y\,.\]

Applying the tail probability formula from Section 7.2 with \(\alpha = r - q + \sigma^2/2\):

\[\text{prob}^Y(S_T > K) = N(d_1)\,.\]

where:

\[d_1 = \frac{\log(S_t/K) + (r - q + \frac{\sigma^2}{2})(T-t)}{\sigma\sqrt{T-t}} = d_2 + \sigma\sqrt{T-t}\,.\]

Combining the Results

Now we can derive the Black-Scholes formula by combining our digital option results. A European call option has payoff:

\[(S_T - K)^+ = \mathbf{1}_{\{S_T > K\}} S_T - K \mathbf{1}_{\{S_T > K\}}\,.\]

The first term is a share digital worth \(e^{-q (T-t)}S_t \cdot \text{prob}^S(S_T > K) = e^{- q (T-t)} S_t N(d_1)\).

The second term is \(K\) cash digitals worth \(K \cdot e^{-r(T-t)} \cdot \text{prob}^R(S_T > K) = K e^{-r(T-t)} N(d_2)\).

We need \(e^{-q(T-t)}\) in front of the share digital to account for dividends not paid to the call holder upon exercise, the call value is:

The Black-Scholes Formula

This derivation shows how the fundamental arbitrage pricing theory leads directly to the Black-Scholes formula. The key steps were:

1. Linear pricing: From Chapter 3, any security’s value is a linear combination of Arrow securities (digitals)
2. Change of numeraire: Different measures give different probabilities for the same event
   
3. Tail probabilities: From Chapter 2, we can compute \(\text{prob}(S_T > K)\) under any measure
4. Girsanov’s Theorem: Provides the mathematical foundation for changing measures

Key Result

The value of a European call option is:

\[C(t, S_t) = e^{-q(T-t)} S_t N(d_1) - e^{-r(T-t)} K N(d_2)\,.\]

where \(d_1\) and \(d_2\) are defined above and \(N(\cdot)\) is the standard normal cumulative distribution function.

The value of a European put option is similarily derived as an exercise :

Key Result

\[P(t, S_t) = e^{-r(T-t)} K N(-d_2) - e^{-q(T-t)} S_t N(-d_1)\,.\]

The call and put values satisfy a put-call parity relationship:

\[e^{-r(T-t)} K + C(t, S_t) = e^{-q(T-t)} S_t + P(t, S_t)\,.\]

Interactive Black-Scholes Explorer

The following interactive tool allows you to explore how the Black-Scholes formula behaves as you vary the input parameters. You can see how option prices change with stock price, volatility, time to expiration, interest rates, and dividend yields.

Figure 9.1: Black-Scholes Explorer. This interactive calculator shows how call and put option prices vary with the underlying parameters. Notice how the option prices change smoothly with the stock price, and observe the effects of time-to-maturity and volatility on option values.

9.4 Replication and Delta Hedging

The Black Scholes formula can also be derived from a no-arbitrage argument since we can replicate option payoffs through dynamic trading. Following Merton’s argument, we construct a portfolio holding \(\delta_t\) shares of stock and \(\alpha_t\) units of the risk-free asset, with value:

\[W_t = \alpha_t e^{rt} + \delta_t S_t\,.\]

For continuous trading with no cash inflows or outflows:

\[\mathrm{d}W_t = \delta_t \,\mathrm{d}S_t + \delta_t q S_t \,\mathrm{d}t + \alpha_t r e^{rt} \,\mathrm{d}t\,.\]

If the option value is \(C(t, S_t)\), then by Ito’s lemma,

\[\mathrm{d}C = \frac{\partial C}{\partial t} \,\mathrm{d}t + \frac{\partial C}{\partial S} \,\mathrm{d}S + \frac{1}{2} \frac{\partial^2 C}{\partial S^2} \sigma^2 S^2 \,\mathrm{d}t\,.\]

Matching the stochastic terms requires

\[\delta_t = \frac{\partial C}{\partial S}\,.\]

This is the option’s delta - the number of shares needed to hedge the option. Matching the drift terms in both expressions \[ \frac{\partial C}{\partial S} \mu S_t + \alpha_t r R_t = \frac{\partial C}{\partial t} + \frac{\partial C}{\partial S} (\mu-q) S_t + \frac{1}{2} \frac{\partial^2 C}{\partial S^2} \sigma^2 S_t^2\,.\] which can be solved to give \[ \alpha_t r R_t=r \left(W_t -\frac{\partial C}{\partial S} S_t\right) = \frac{\partial C}{\partial t}-\frac{\partial C}{\partial S} q S_t + \frac{1}{2} \frac{\partial^2 C}{\partial S^2} \sigma^2 S_t^2\,,\] which gives the equation \[r W_t = \frac{\partial C}{\partial t} + \frac{\partial C}{\partial S}(r-q) S_t+\frac{1}{2} \frac{\partial^2 C}{\partial S^2} \sigma^2 S_t^2\] with a boundary condition \(W_T = (S_T-K)^{+}\). However, no-arbitrage suggests \(W_t = C(t,S_t)\), which gives us the partial differential equation \[r C = \frac{\partial C}{\partial t} + \frac{\partial C}{\partial S} (r-q) S+\frac{1}{2} \frac{\partial^2 C}{\partial S^2} \sigma^2 S^2\] with a boundary condition \(C(T,S_T)= (S_T - K)^{+}\). This is a partial differential equation and a fairly tedious set of calcuations show the Black Scholes formula is a solution (in fact it is the only positive solution). Close observation of the right hand side we see this is the drift term of Ito expansion for \(C\) if we work in the risk-neutral probability. The right hand side then says in the risk-neutral probability, the call option earns the risk free return.

However, there is nothing special about a call option. The same argument will apply for any European style option. The only difference is the boundary condition.

Important Principle

The price \(V(S_t,t)\) of any European derivative with payoff \(V(S_T)\) at maturity \(T\) must satisfy:

\[r V = \frac{\partial V}{\partial t} + \frac{\partial V }{\partial S}(r - q)S + \frac{1}{2} \frac{\partial^2 V}{\partial S^2} \sigma^2 S^2\]

with boundary condition \(V(S_T, T) = V(S_T)\).

This PDE has several important interpretations:

1. Risk-neutral expectation: The solution is \(V(t, S_t) = \mathbb{E}_t^R[e^{-r(T-t)} V(S_T)]\). This sometimes called the Feynman Kac Theorem.¹
2. Replication: The portfolio holds \(\frac{\partial V}{\partial S}\) shares and \(V - \frac{\partial V}{\partial S}S\) in cash
3. Hedging condition: The drift term equals the risk-free return when perfectly hedged

The Black-Scholes formula satisfies this PDE with the appropriate boundary conditions for calls and puts.

Simulating Portfolio Insurance

Recall, that a protective put position buys a put and buys a share and the payoff at the expiration of the put is given by \(\max(K,S_T)\). The reason for the name protective put is apparent since the position can pay off no less than \(K\). The cost of this insurance is the price of the put. However, if the put is not traded, we can synthetically replicate this payoff using the prodedure above assuming we can trade continuously. The basic recipe is to start with inital wealth equal to that for a protective put position: \(W_0 = P(0,S_0)+S_0\). The delta of the protective put position can be calcuated to be the delta of the put plus 1 which is \(N(d_1)\), where \(d_1\) is calculated at each point in time. However, in practice we cannot trade continuously. A simple discrete strategy would rebalance at intervals \(\Delta t\). The strategy calculates \(N(d_1)\) at time 0 and holds \(P(0,S_0)+S_0 - N(d_1)S_0\) dollars in the risk free asset and \(N(d_1)\) shares of the asset. Thereafter these holdings are adjusted. The change in portfolio value over the interval \(\Delta t\) is \[\Delta W= W_{i \Delta t}- W_{(i-1)\Delta t} \] \[ = (P((i-1)\Delta t,S_{(i-1)\Delta t})+S_{(i-1)\Delta t} - N(d_1-)S_{(i-1)\Delta t})(R_{i\Delta t}-R_{(i-1)\Delta t}) + N(d1-)(S_{i\Delta t}-S_{(i-1)\Delta t}) \,, \] where \(N(d1-)\) is the delta chosen at time \((i-1)\Delta t\). The question is if the Black Scoles model is correct, how accurate can a discrete rebalancing scheme be? This is simulated in the following code:

import numpy as np
# from bsfunctions import *
import matplotlib.pyplot as plt
import time
from math import pow, exp, sqrt
from scipy import stats
# incs = np.genfromtxt('incs.csv',delimiter=",",skip_header=1)
def blackscholes(S0, K, r, q, sig, T, call = True):
    '''Calculate option price using B-S formula.

    Args:
        S0 (num): initial price of underlying asset.
        K (num): strick price.
        r (num): risk free rate.
        q (num): dividend yield
        sig (num): Black-Scholes volatility.
        T (num): maturity.
        call (bool): True returns call price, False returns put price.

    Returns:
        num
    '''
    d1 = (np.log(S0/K) + (r -q + sig**2/2) * T)/(sig*np.sqrt(T))
    d2 = d1 - sig*np.sqrt(T)
#     norm = sp.stats.norm
    norm = stats.norm
    if call:
        return np.exp(-q*T)*S0 * norm.cdf(d1,0,1) - K * np.exp(-r * T) * norm.cdf(d2,0, 1)
    else:
        return np.exp(-q*T)*S0 * -norm.cdf(-d1,0,1) + K * np.exp(-r * T) * norm.cdf(-d2,0, 1)

def blackscholes_delta(S0, K, r, q, sig, T, call = True):
    '''Calculate option price using B-S formula.

    Args:
        S0 (num): initial price of underlying asset.
        K (num): strick price.
        r (num): risk free rate.
        q (num): dividend yield
        sig (num): Black-Scholes volatility.
        T (num): maturity.
        call (bool): True returns call price, False returns put price.

    Returns:
        num
    '''
    d1 = (np.log(S0/K) + (r -q + sig**2/2) * T)/(sig*np.sqrt(T))
    d2 = d1 - sig*np.sqrt(T)
#     norm = sp.stats.norm
    norm = stats.norm
    if type(call) == bool:
        if call:
            return np.exp(-q*T)*norm.cdf(d1,0,1)
        else:
            return np.exp(-q*T)*norm.cdf(-d1,0,1)
    else:
        print("Not a valid value for call")

# parameters
# number of paths
# n = incs.shape[1]
n = 100000
# number of divisions
# m = incs.shape[0]
m = 100
# interest rate
r = .1
# dividend yield
q=0.0
# true drift
mu = .15
# volatility
sig = .2
# Initial Stock Price
S0 = 42
# Strike Price
K = 42
# Maturity
T = 0.5


# seed for random generator
seed= 1234
# define a random generator
rg = np.random.RandomState(seed)
# initialize


# generate normal random vairables
dt= T/m
vol=sig*np.sqrt(dt)
incs = rg.normal(0,vol,[m,n])


tline = np.linspace(0,T,m+1)


St = np.zeros((m+1,n))
#St1 = np.zeros((m+1,n))

V_vec = np.zeros((m+1,n))

delta = np.zeros((m,n))

put= blackscholes(S0,K,r, q, sig,T,call=False)

incs_cumsum =  np.concatenate((np.zeros((1,n)),incs),axis=0).cumsum(axis=0)
V_vec = np.zeros((m+1,n))
t_mat =  np.repeat(tline.reshape((m+1,1)), n, axis=1)
drift_cumsum = (mu -q -0.5*sig**2) * t_mat

St = S0 * np.exp(incs_cumsum + drift_cumsum)

delta = blackscholes_delta(St[:-1,:],K,r, q, sig,T-t_mat[:-1,:])

V_vec[0,:] = S0 + put

for i in range(1,m+1):
    V_vec[i,:] = V_vec[i-1,:] + (np.exp(r*dt)-1) * (V_vec[i-1,:] - delta[i-1,:] * St[i-1,:])+ delta[i-1,:] * (St[i,:]-St[i-1,:])

# Uses actual simulated changes in riskfree and stock price not the dt and dB approximations 
# plot ST versus VT
plt.scatter(St[m,:],V_vec[m,:])
plt.xlabel('Stock Price at Maturity')
plt.ylabel('Value of Portfolio Insurance')
plt.show()

With \(m=100\) rebalancing dates over \(T=0.5\) for the parameters chosen the repliacting strategy does a pretty good job. The hedging errors occur when the stock price is close to the strike price. This is not surprising since the delta changes (measured by the gamma) fastest around this point. A gamma hedge would potentially improve the performance.

The portfolio insurance rebalancing scheme involves sell stock and buying bonds when the stock price goes down and buying stocks and selling bonds when the stock price goes up. This can be destabilizing and was identified as a contributor to the 1987 stock market crash.

9.5 Greeks

The sensitivities of option values to various inputs are called Greeks. From the Black-Scholes formula:

Delta (\(\delta\))

The sensitivity to stock price changes:

\[\delta_{call} = e^{-q(T-t)} N(d_1)\,.\] \[\delta_{put} = -e^{-q(T-t)} N(-d_1)\,.\]

Gamma (\(\Gamma\))

The rate of change of delta:

\[\Gamma = \frac{e^{-q(T-t)} n(d_1)}{S \sigma \sqrt{T-t}}\,,\]

where \(n(\cdot)\) is the standard normal density. Gamma is the same for calls and puts.

Theta (\(\Theta\))

The time decay (using \(\partial/\partial(-T)\) so positive theta means value increases as time passes):

\[\Theta_{call} = -\frac{e^{-q(T-t)} S n(d_1) \sigma}{2\sqrt{T-t}} + q e^{-q(T-t)} S N(d_1) - r e^{-r(T-t)} K N(d_2)\,.\]

Vega (\(\mathcal{V}\))

The sensitivity to volatility:

\[\mathcal{V} = e^{-q(T-t)} S n(d_1) \sqrt{T-t}\,.\]

Rho (\(\rho\))

The sensitivity to interest rates:

\[\rho_{call} = (T-t) e^{-r(T-t)} K N(d_2)\,.\]

Interactive Greeks Explorer

The Greeks measure how option prices change with respect to various parameters. The following interactive tool lets you visualize how the Greeks behave across different stock prices and parameter values.

Figure 9.2: Black-Scholes Greeks Explorer. This tool shows how delta, gamma, theta, vega, and rho vary with the stock price and other parameters. Pay particular attention to how gamma peaks near the strike price and how theta becomes more negative as expiration approaches.

9.6 Theta and Gamma in Delta Hedges

The relationship between theta and gamma provides deep insight into option hedging. Consider a delta-hedged portfolio (short call, long \(\delta\) shares). From the fundamental PDE and our Greeks:

\[\Theta + \frac{1}{2} \Gamma \sigma^2 S^2 = r C - (r - q) S \delta\,.\]

In continuous time, the portfolio change is

\[-\Theta \,\mathrm{d}t - \frac{1}{2} \Gamma \sigma^2 S^2 \,\mathrm{d}t + q \delta S \,\mathrm{d}t + (C - \delta S) r \,\mathrm{d}t\,.\]

These terms exactly cancel—the time decay and dividends received offset losses from being short gamma and interest payments. This perfect cancellation only holds with continuous rebalancing.

With discrete rebalancing, the hedge is imperfect. The portfolio is short gamma, meaning it loses money when the stock moves significantly between rebalances. This loss is approximately:

\[\frac{1}{2} \Gamma (\Delta S)^2\,,\]

where \(\Delta S\) is the stock price change. The hedging error increases with:

• Larger gamma (near the strike at expiration)
• Higher volatility
• Longer time between rebalances

Discretely Rebalanced Delta Hedges

The theoretical analysis above assumes continuous rebalancing, but in practice, hedging must be done at discrete intervals. This introduces hedging error that we can analyze through simulation.

The Discrete Hedging Problem

Consider a delta hedge that is rebalanced at discrete times \(0 = t_0 < t_1 < \cdots < t_N = T\) with \(\Delta t = T/N\). Between rebalancing dates, the hedge portfolio consists of:

• Short position in the option
• Long \(\delta_{t_i}\) shares of stock (where \(\delta_{t_i}\) is computed at the last rebalancing date)
• Cash position to finance the hedge

The key insight from the analysis in ?sec-s_deltahedging is that perfect continuous hedging relies on the relationship:

\[-\Theta \,\mathrm{d} t - \frac{1}{2}\Gamma \sigma^2S^2\,\mathrm{d} t+ q \delta S\,\mathrm{d} t+(C-\delta S)r\,\mathrm{d} t = 0\,.\]

With discrete rebalancing, this balance is disrupted. The portfolio gains and losses come from:

1. Time decay: \(-\Theta \Delta t\) (predictable, typically positive for short option)
2. Gamma exposure: \(-\frac{1}{2}\Gamma (\Delta S)^2\) (stochastic losses when short gamma)
3. Interest and dividends: Financing costs and dividend income
4. Delta drift: Changes in delta between rebalancing dates

Simulation of Discrete Hedging

The following code simulates the performance of discretely rebalanced delta hedges. It runs 10,000 simulations of a stock price path, with trading to an updated delta hedge each period.²

import numpy as np
from scipy.stats import norm
from numba import jit

def blackscholes(S0, K, r, q, sig, T, call = True):
    '''Calculate option price using Black-Scholes formula.'''
    d1 = (np.log(S0/K) + (r -q + sig**2/2) * T)/(sig*np.sqrt(T))
    d2 = d1 - sig*np.sqrt(T)
    if call:
        return np.exp(-q*T)*S0 * norm.cdf(d1,0,1) - K * np.exp(-r * T) * norm.cdf(d2,0, 1)
    else:
        return np.exp(-q*T)*S0 * (-norm.cdf(-d1,0,1)) + K * np.exp(-r * T) * norm.cdf(-d2,0, 1)

def blackscholes_delta(S0, K, r, q, sig, T, call = True):
    '''Calculate option delta using Black-Scholes formula.'''
    d1 = (np.log(S0/K) + (r -q + sig**2/2) * T)/(sig*np.sqrt(T))
    if call:
        return np.exp(-q*T)*norm.cdf(d1,0,1)
    else:
        return np.exp(-q*T)*norm.cdf(-d1,0,1)

@jit(nopython=True)
def norm_cdf_fast(x):
    """
    Fast approximation of cumulative normal distribution.
    Accurate to about 7 decimal places using Abramowitz & Stegun approximation.
    """
    # Constants for the approximation
    a1 =  0.254829592
    a2 = -0.284496736
    a3 =  1.421413741
    a4 = -1.453152027
    a5 =  1.061405429
    p  =  0.3275911

    # Save the sign of x
    sign = 1.0 if x >= 0.0 else -1.0
    x = abs(x) / np.sqrt(2.0)

    # Abramowitz & Stegun formula
    t = 1.0 / (1.0 + p * x)
    y = 1.0 - (((((a5 * t + a4) * t) + a3) * t + a2) * t + a1) * t * np.exp(-x * x)

    return 0.5 * (1.0 + sign * y)

@jit(nopython=True)
def blackscholes_delta_fast(S, K, r, q, sigma, T):
    """Fast Black-Scholes delta for call option using Numba."""
    if T <= 0:
        return 1.0 if S > K else 0.0
    d1 = (np.log(S / K) + (r - q + sigma**2 / 2) * T) / (sigma * np.sqrt(T))
    return np.exp(-q * T) * norm_cdf_fast(d1)

@jit(nopython=True)
def simulate_single_path(S0, K, r, sigma, q, T, mu, N, drift, vol_dt,
                        discount, div_factor, delta_0, cash_0, random_normals):
    """
    Simulate a single delta hedge path.
    Optimized with Numba JIT compilation.
    """
    dt = T / N
    log_s = np.log(S0)
    cash = cash_0
    s = S0
    delta = delta_0

    # Rebalancing loop
    for j in range(1, N):
        # Stock price evolution
        log_s += drift + vol_dt * random_normals[j]
        new_s = np.exp(log_s)

        # New delta for rebalancing
        time_remaining = T - j * dt
        new_delta = blackscholes_delta_fast(new_s, K, r, q, sigma, time_remaining)

        # Update cash position: interest + dividends - rebalancing cost
        cash = (discount * cash +
               delta * s * div_factor -
               (new_delta - delta) * new_s)

        s = new_s
        delta = new_delta

    # Final period
    log_s += drift + vol_dt * random_normals[N]
    final_s = np.exp(log_s)

    # Final hedge value
    final_hedge_value = (discount * cash +
                       delta * s * div_factor +
                       delta * final_s)

    # Profit = hedge value - option payoff
    option_payoff = max(final_s - K, 0.0)
    return final_hedge_value - option_payoff

@jit(nopython=True)
def run_all_simulations(S0, K, r, sigma, q, T, mu, M, N, drift, vol_dt,
                       discount, div_factor, delta_0, cash_0, all_random_normals):
    """Run all M simulations using Numba JIT."""
    profits = np.zeros(M)

    for i in range(M):
        profits[i] = simulate_single_path(
            S0, K, r, sigma, q, T, mu, N, drift, vol_dt,
            discount, div_factor, delta_0, cash_0, all_random_normals[i, :]
        )

    return profits

def simulated_delta_hedge_profit(S0, K, r, sigma, q, T, mu, M, N, pct):
    """
    Simulate profits/losses from discretely rebalanced delta hedge.
    Optimized version using Numba JIT compilation.

    Parameters:
    S0: initial stock price
    K: strike price
    r: risk-free rate
    sigma: volatility
    q: dividend yield
    T: time to maturity
    mu: expected stock return (actual, not risk-neutral)
    M: number of simulations
    N: number of rebalancing periods
    pct: percentile to return
    """
    dt = T / N
    vol_dt = sigma * np.sqrt(dt)
    drift = (mu - q - 0.5 * sigma**2) * dt
    discount = np.exp(r * dt)
    div_factor = np.exp(q * dt) - 1

    # Initial setup (using scipy for initial values)
    call_0 = blackscholes(S0, K, r, q, sigma, T, True)
    delta_0 = blackscholes_delta(S0, K, r, q, sigma, T, True)
    cash_0 = call_0 - delta_0 * S0

    # Generate all random numbers upfront for better performance
    np.random.seed(42)
    all_random_normals = np.random.randn(M, N + 1)

    # Run all simulations using JIT-compiled function
    profits = run_all_simulations(
        S0, K, r, sigma, q, T, mu, M, N, drift, vol_dt,
        discount, div_factor, delta_0, cash_0, all_random_normals
    )

    return np.percentile(profits, pct * 100)

# Example: Hedge a call option with discrete rebalancing
S0 = 100      # Initial stock price
K = 100       # Strike price (at-the-money)
r = 0.05      # Risk-free rate
sigma = 0.2   # Volatility
q = 0.02      # Dividend yield
T = 0.25      # 3 months to maturity
mu = 0.12     # Expected stock return

print("Discrete Delta Hedging Analysis")
print("==============================")
print(f"Parameters: S0={S0}, K={K}, r={r:.2f}, σ={sigma:.2f}, q={q:.2f}, T={T}")
print(f"Expected stock return μ={mu:.2f}\n")

# Test different rebalancing frequencies
rebalancing_frequencies = [3, 13, 63]  
M = 10000  # Number of simulations

print("Hedging Error Analysis (95th percentile of absolute profits):")
print("Rebalancing Frequency | Periods | 95th Percentile Error")
print("-" * 55)

for N in rebalancing_frequencies:
    # Get both tails of the distribution
    p95 = simulated_delta_hedge_profit(S0, K, r, sigma, q, T, mu, M, N, 0.95)
    p5 = simulated_delta_hedge_profit(S0, K, r, sigma, q, T, mu, M, N, 0.05)
    
    # Report the larger absolute error
    max_error = max(abs(p95), abs(p5))
    
    freq_name = {3: "Monthly", 13: "Weekly", 63: "Daily"}[N]
    print(f"{freq_name:<20} | {N:>7} | {max_error:>18.4f}")

print(f"\nNote: Errors should decrease as rebalancing frequency increases.")

Discrete Delta Hedging Analysis
==============================
Parameters: S0=100, K=100, r=0.05, σ=0.20, q=0.02, T=0.25
Expected stock return μ=0.12

Hedging Error Analysis (95th percentile of absolute profits):
Rebalancing Frequency | Periods | 95th Percentile Error
-------------------------------------------------------
Monthly              |       3 |             3.2002
Weekly               |      13 |             1.5593
Daily                |      63 |             0.7157

Note: Errors should decrease as rebalancing frequency increases.

9.7 Delta-Gamma Hedging

To attempt to improve the performance of a discretely rebalanced delta hedge, one can use another option to create a portfolio that is both delta and gamma neutral. Being delta neutral means hedged as in the previous section—the portfolio value has no exposure to changes in the underlying asset price. In other words, it means that the derivative of the portfolio value with respect to the price of the underlying (the portfolio delta) is zero. Being gamma neutral means that the delta of the portfolio has no exposure to changes in the underlying price, which is equivalent to the second derivative of the portfolio value with respect to the price of the underlying (the portfolio gamma) being zero. If the delta truly did not change, then there would be no need to rebalance continuously, and hence no hedging error introduced by only adjusting the portfolio at discrete times rather than continuously. However, there is certainly no guarantee that a discretely rebalanced delta/gamma hedge will perform better than a discretely rebalanced delta hedge.

A delta/gamma hedge can be constructed as follows. Suppose we have written (shorted) a call option and we want to hedge both the delta and gamma using the underlying asset and another option, for example, another call option with a different strike. In practice, one would want to use a liquid option for this purpose, which typically means that the strike of the option will be near the current value of the underlying (i.e., the option used to hedge would be approximately at the money).

Let \(\delta\) and \(\Gamma\) denote the delta and gamma of the written option and let \(\delta^*\) and \(\Gamma^*\) denote the delta and gamma of the option used to hedge. Consider holding \(a\) shares of the stock and \(b\) units of the option used to hedge in conjunction with the short option. The delta of the stock is one (\(\mathrm{d} S/\mathrm{d} S = 1\)), so to obtain a zero portfolio delta we need

\[ 0 = - \delta + a + b\delta^*. \qquad(9.1)\]

The gamma of the stock is zero (\(\mathrm{d}^2 S/\mathrm{d} S^2 = d\,1/\mathrm{d} S = 0\)), so to obtain a zero portfolio gamma we need

\[ 0 = - \Gamma + b\Gamma^*\;. \qquad(9.2)\]

Equation Equation 9.2 shows that we should hold enough of the second option to neutralize the gamma of the option we have shorted; i.e., \[ b= \frac{\Gamma}{\Gamma^*} \] Equation Equation 9.1 shows that we should use the stock to delta hedge the portfolio of options; i.e., \[ a=\delta - \frac{\Gamma}{\Gamma^*}\delta^*\;. \]

9.8 Implied Volatilities

While all other Black-Scholes inputs are observable, volatility must be estimated. Given an option’s market price, we can invert the Black-Scholes formula to find the implied volatility—the \(\sigma\) that produces the observed price.

Implied volatilities serve several purposes:

1. Price quotes: Options are often quoted in terms of implied volatility
2. Relative value: Comparing implied volatilities helps identify expensive or cheap options
3. Market views: Implied volatilities reflect market expectations of future volatility

The Volatility Smile

If the Black-Scholes model were perfect, all options with the same maturity would have the same implied volatility. In practice, plotting implied volatility against strike typically shows:

• Higher implied volatilities for low strikes (out-of-the-money puts)
• Lower implied volatilities near the at-the-money strike
• Slightly increasing implied volatilities for high strikes

This “volatility smile” or “smirk” indicates that market prices reflect: - Fat tails: Higher probability of extreme moves than lognormal - Negative skewness: Larger probability of extreme downward moves

The smile has been particularly pronounced for equity index options since the 1987 crash, suggesting market participants price in crash risk.

9.9 Exercises

Digital Options and Building Blocks

Exercise 9.1 Consider a cash digital option that pays $1 when \(S_T > K\) and a share digital that pays \(S_T\) when \(S_T > K\). Using Black-Scholes parameters \(S_0 = 100\), \(K = 105\), \(r = 5\%\), \(q = 2\%\), \(\sigma = 20\%\), \(T = 0.25\).

1. Calculate the analytical values of both digital options using the formulas from the chapter
2. Verify your results using Monte Carlo simulation with 100,000 paths
3. Show that a call option can be decomposed as: \(C = S_0 N(d_1) e^{-qT} - K N(d_2) e^{-rT}\)
4. Implement this decomposition and verify it matches the standard Black-Scholes formula

Exercise 9.2 The delta of a cash digital option that pays $1 when \(S_T > K\) is \[\Delta_{\text{digital}} = \frac{e^{-rT} n(d_2)}{\sigma S \sqrt{T}}\,,\]

where \(n(\cdot)\) is the standard normal density.

1. Plot this delta against stock price for \(S \in [80, 120]\) with the parameters from the previous exercise
2. Compare the magnitude of digital delta to call option delta near the strike
3. Explain why delta hedging a short digital position near expiration and near the strike is problematic
4. What happens to the digital delta as time to expiration approaches zero?

::: {#digital puts} The following two questions can be answered in two ways.

What is the vaule of a derivative security which pays \[\mathbf{1}_{\{S_T <K\}}\,.\]

What is the value of a derivative security which pays \[ S_T \mathbf{1}_{\{S_T <K\}}\,.\]

:::

Change of Numeraire and Girsanov’s Theorem

Exercise 9.3 Using simulation, verify the change of numeraire formulas for digital options:

1. Under the risk-neutral measure, simulate \(S_T\) following \(dS/S = (r-q)dt + \sigma dB^R\)
2. Under the stock numeraire measure, simulate \(S_T\) following \(dS/S = (r-q+\sigma^2)dt + \sigma dB^S\)
   
3. Calculate \(\mathbb{E}^R[e^{-rT} \mathbf{1}_{\{S_T > K\}}]\) and \(S_0 \mathbb{E}^S[\mathbf{1}_{\{S_T > K\}}/S_T]\)
4. Verify both give the same digital option value up to simulation error
5. Use 50,000 simulation paths and report the standard errors

Exercise 9.4 Verify the exponential martingale property in Girsanov’s Theorem. With \(\kappa = (\mu - r)/\sigma\):

1. Simulate paths of \(Z_t = \exp(-\frac{1}{2}\kappa^2 t - \kappa B_t)\) where \(B_t\) is standard Brownian motion
2. Verify that \(\mathbb{E}[Z_T] = 1\) for various values of \(T\)
3. Show that \(Z_t S_t/R_t\) is a martingale under the original measure
4. Compare the distribution of \(B_t + \kappa t\) to standard Brownian motion under the new measure

Black-Scholes Formula and Properties

Exercise 9.5 Implement the complete Black-Scholes formula with all Greeks:

1. Create a Python class BlackScholesOption that calculates price, delta, gamma, theta, vega, and rho
2. Include both calls and puts with proper dividend yield handling
3. Verify put-call parity: \(C - P = S_0 e^{-qT} - K e^{-rT}\)
4. Test your implementation against the interactive Black-Scholes explorer values

Exercise 9.6 Analyze the limiting behavior of the Black-Scholes formula:

1. Show that as \(T \to 0\), the call value approaches \(\max(S_0 - K, 0)\)
2. Show that as \(\sigma \to 0\), the formula approaches the deterministic payoff present value
3. What happens to call and put values as \(\sigma \to \infty\)? Explain intuitively.
4. Analyze the behavior as \(r \to 0\) and as \(r \to \infty\)

Exercise 9.7 Using the interactive Black-Scholes explorer or your own implementation:

1. For an at-the-money option, plot call value against volatility for \(\sigma \in [0.1, 1.0]\)
2. Repeat for options that are 10% in-the-money and 10% out-of-the-money
3. For which moneyness is the option price most sensitive to volatility changes?
4. Explain this pattern in terms of the option’s probability of finishing in-the-money

Greeks and Risk Management

Exercise 9.8 Create comprehensive plots of the Greeks using the interactive Greeks explorer as reference:

1. Plot delta, gamma, theta, vega, and rho against stock price for an at-the-money option
2. Repeat for different times to expiration: 1 month, 3 months, 6 months, 1 year
3. Identify where gamma is maximized and explain why this matters for hedging
4. Show how theta changes as expiration approaches and explain the time decay acceleration

Exercise 9.9 Consider delta and gamma hedging a short call option using the underlying and a put with the same strike and maturity:

1. Derive the positions in stock and put needed for a delta-gamma neutral portfolio
2. Show that this hedge never needs adjustment (relate to put-call parity)
3. Implement this strategy and compare hedging errors to delta-only hedging
4. What are the practical limitations of this “perfect” hedge?

The Fundamental PDE and Replication

Show that the stock price and the boond price satisfy the fundamental PDE with the appropriate boundary condition.

Show that the share digital and cash digital prices satisfy the fundamenatl PDE with the appropriate boundary condition.

Use risk neutral valuation to find the value of a derivative security which pays \(\ln(S_T)\). Verify your value satisfies the fundamental PDE with the appropriate boundary condition. What is the replicating portfolio? Is this option harder to hedge than a standard option? Why or why not?

Exercise 9.10 Calculate the Greeks for exotic payoffs:

1. For a derivative paying \(S_T^2\), find the value, delta, and gamma using risk-neutral valuation
2. For a derivative paying \(\log(S_T)\), repeat the calculation
3. Verify both satisfy the fundamental PDE: \(rV = \frac{\partial V}{\partial t} + (r-q)S\frac{\partial V}{\partial S} + \frac{1}{2}\sigma^2 S^2 \frac{\partial^2 V}{\partial S^2}\)
4. Compare the hedging difficulty (gamma exposure) of these payoffs to standard options.

:::{#digital option replication} Modify the code which simulates portfolio insurance to simulate replication of a cash digital option. Where do the largest errors tend to occur? :::

:::{#option replication and transaction costs} Modify the code which simulates portfolio insurance to incur a transaction cost. One model has ask price \(1.01 S_t\) and and bid price share \(0.99 S_t\). Another model has ask price \(S_t +.005\) and bid price \(S_t - .005\). Try different rebalancing frequencies. What rebalancing frequency gives the best hedging? Where do the largest errors tend to occur? :::

Discrete Hedging and Practical Implementation

Exercise 9.11 Extend the discrete hedging simulation from the chapter:

1. Compare hedging errors for different rebalancing frequencies: daily, weekly, monthly
2. Analyze how hedging errors scale with volatility and time to expiration
   
3. Study the impact of transaction costs: assume 0.1% bid-ask spread on each rebalance
4. Find the optimal rebalancing frequency that minimizes total cost (hedging error + transaction costs)

Exercise 9.12 Analyze the gamma scalping strategy:

1. Show that a delta-hedged short option position has P&L approximately equal to \(-\frac{1}{2}\Gamma(\Delta S)^2\)
2. Simulate this P&L for different realized volatilities vs. implied volatility
3. Demonstrate that selling options when implied volatility > realized volatility is profitable
4. Account for the theta decay and show the complete P&L attribution

Implied Volatility and Market Practice

Exercise 9.13 Implement an implied volatility calculator:

1. Use numerical methods (bisection or Newton-Raphson) to invert the Black-Scholes formula
2. Test with market-like option prices and verify convergence
3. Handle edge cases: very deep ITM/OTM options, very short/long expirations
4. Compare computational efficiency of different numerical methods

Delta-Gamma Hedging

Exercise 9.14 Consider delta and gamma hedging a short call option, using the underlying and a put with the same strike and maturity as the call. Calculate the position in the underlying and the put that you should take, using the analysis in the Delta-Gamma Hedging section. Will you ever need to adjust this hedge? Relate your result to put-call parity.

Exercise 9.15 Modify the discrete delta hedging simulation to compute percentiles of gains and losses for an investor who writes a call option and constructs a delta and gamma hedge using the underlying asset and another call option. Include the exercise price of the call option used to hedge as an input, and assume it has the same time to maturity as the option that is written.

Hint: In each period, the updated cash position can be calculated as:

Cash = exp(r*dt)*Cash + a*S*(exp(q*dt)-1) - (Newa-a)*NewS - (Newb-b)*PriceHedge

where a denotes the number of shares of the stock held, b denotes the number of units held of the option that is used for hedging, and PriceHedge denotes the price of the option used for hedging (computed from the Black-Scholes formula each period). This expression embodies the interest earned (paid) on the cash position, the dividends received on the shares of stock and the cash inflows (outflows) from adjusting the hedge. At the final date, the value of the hedge is:

exp(r*dt)*Cash + a*S*(exp(q*dt)-1) + a*NewS + b*max(NewS-KHedge,0)

and the value of the overall portfolio is the value of the hedge less max(NewS-KWritten,0), where KHedge denotes the strike price of the option used to hedge and KWritten denotes the strike of the option that was written.

9.10 Summary

The Black-Scholes model combines the theoretical foundations of arbitrage pricing with practical implementation:

1. Change of numeraire: Digital options under different measures give the formula
2. Replication: Delta hedging shows how to replicate option payoffs
3. Fundamental PDE: All Europena style derivatives satisfy the same partial differential equation with diffrent boundary conditions.
4. Greeks: Sensitivities guide risk management and hedging
5. Market practice: Implied volatilities reveal market expectations and model limitations

The model’s elegance lies in showing that, under its assumptions, options can be perfectly hedged through dynamic trading. While real markets violate these assumptions, the Black-Scholes framework remains the foundation for understanding option pricing and hedging.

---

1. See Karatzas and Shreve (1991).

2. The implementation uses Numba JIT compilation for significant performance improvements (typically 20-100x speedup compared to pure Python loops).