Numerical methods for American Option Pricing

Rancy Chepchirchir
8 min readJan 9, 2023

Definition:

An American option is a type of financial derivative that gives the holder the right, but not the obligation, to buy or sell an underlying asset at a specified price on or before a certain date. The key difference between an American option and a European option is that an American option can be exercised at any time before the expiration date, while a European option can only be exercised at the expiration date.

Approaches:

There are several numerical methods that can be used to price American options. One popular method is the finite difference method, which involves discretizing the underlying partial differential equation that describes the option price and solving the resulting system of equations using numerical techniques. Other methods include the binomial options model and the Monte Carlo simulation method.

  1. ) Binomial option model

The Binomial options pricing model is a widely used method for valuing American options. It is based on the idea of breaking the time to expiration of the option into a series of discrete periods and modeling the price of the underlying asset at the end of each period. The model assumes that the underlying asset can either increase or decrease in price over each time period by a fixed percentage. By working backward from the expiration date, it is possible to calculate the value of the option at each time period, taking into account the probability of the option being exercised at each point in time. The final value of the option is then obtained by taking the present value of all of these future option values.

To solve for the value of an American option using the binomial model in Python, you can use the following steps:

  1. Define the parameters of the option, including the current price of the underlying asset, the strike price, the risk-free interest rate, the time to expiration, the volatility, and the number of steps in the binomial tree.
  2. Create a binomial tree using these parameters. The binomial tree is a discrete model that represents the possible prices of the underlying asset at different points in time. It consists of a series of nodes, with each node representing a possible price of the underlying asset at a given point in time.
  3. Starting at the final node of the tree, work backward to calculate the value of the option at each node. To do this, you will need to use the Black-Scholes formula to calculate the value of the option at each node.
  4. Once the value of the option has been calculated at each node of the tree, you can use these values to determine the value of the option at the current time.

Here is some sample Python code that demonstrates how to use the binomial model to value an American option:

import numpy as np
import matplotlib.pyplot as plt

def binomial_option_values(S, X, r, T, N, sigma, option_type):
"""
S: current price of the underlying asset
X: exercise price of the option
r: risk-free interest rate
T: time to expiration of the option (in years)
N: number of time periods in the binomial model
option_type: "call" or "put"
"""
dt = T/N
u = np.exp(sigma * np.sqrt(dt))
d = 1/u
p = (np.exp(r*dt) - d)/(u - d)
values = np.zeros((N+1, N+1))

for j in range(N+1):
values[N][j] = max(0, (S*u**j*d**(N-j) - X))

for i in range(N-1, -1, -1):
for j in range(i+1):
values[i][j] = np.exp(-r*dt)*(p*values[i+1][j+1] + (1-p)*values[i+1][j])

return values

# Example usage:
S = 50
X = 45
r = 0.05
T = 0.5
N = 20
sigma = 0.7
option_type = "call"

option_values = binomial_option_values(S, X, r, T, N, sigma, option_type)

# Plot the option values over time
plt.plot(option_values)
plt.xlabel("Time period")
plt.ylabel("Option value")
plt.show()

The plot shows the option values at each time period, with the time periods on the x-axis and the option values on the y-axis.

2.) Finite difference method

To use the finite difference method to price an American option, you first need to define the grid for the option’s underlying variables, such as the stock price and time. You can then use the partial differential equation for the option price, known as the Black-Scholes equation, to derive a set of finite difference equations for the option price at each point on the grid.

Once you have the finite difference equations, you can use numerical techniques, such as the explicit or implicit method, to solve for the option price at each point on the grid. The option price at the final time step will give you the price of the American option.

There are several benefits to using the finite difference method to price American options. It is relatively simple to implement, and it can handle a wide range of option types and underlying assets. However, it can be computationally intensive, and it may not always be the most accurate method for pricing options.

Here is a simple example of how the finite difference method can be implemented in Python to price an American call option on a non-dividend paying stock:

import numpy as np
import matplotlib.pyplot as plt

def finite_difference_option_price(S, X, r, T, sigma, N, M, option_type):
"""
S: current price of the underlying asset
X: exercise price of the option
r: risk-free interest rate
T: time to expiration of the option (in years)
sigma: volatility of the underlying asset's price
N: number of time steps in the finite difference grid
M: number of asset price steps in the finite difference grid
option_type: "call" or "put"
"""
# Set up the finite difference grid
dt = T/N
dx = sigma * np.sqrt(dt)
edx = np.exp(dx)
edx2 = edx**2
pu = (r - 0.5*sigma**2) * dt / (2*dx)
pm = 1 - dt*sigma**2/(dx**2)
pd = -(r - 0.5*sigma**2) * dt / (2*dx)
u = np.arange(M+1)
u = X * edx**u
u = u[::-1]
vm = np.zeros((N+1, M+1))

# Initialize the boundary conditions
if option_type == "call":
vm[:, 0] = 0.0
vm[:, M] = u[M] - u
else:
vm[:, 0] = u - u[0]
vm[:, M] = 0.0

# Iterate over the time steps
for i in range(N):
for j in range(1, M):
vm[i+1, j] = pu * vm[i, j-1] + pm * vm[i, j] + pd * vm[i, j+1]
if option_type == "call":
vm[i+1, j] = max(vm[i+1, j], u[j] - u[M-i])
else:
vm[i+1, j] = max(vm[i+1, j], u[0] - u[M-i-j])

# Return the option value at the current time
return vm[0, M/2]

# Example usage:
S = 50
X = 45
r = 0.05
T = 0.5
sigma = 0.2
N = 20
M = 50
option_type = "call"

option_price = finite_difference_option_price(S, X, r, T, sigma, N, M, option_type)

Here, we employ the finite difference method to solve for the option price at each point on the grid. It uses the Black-Scholes partial differential equation to derive the finite difference equations, and it uses the boundary conditions for an American call option to set the option price at the edges of the grid.

The plot below is a simple representation of the errors (with respect to the closed-form Black Scholes model) across the three finite difference methods namely; the explicit, implicit, and Crank-Nicolson methods.

3.) Monte Carlo simulation method

The Monte Carlo simulation method is another numerical method that can be used to price American options. It involves using random sampling to generate a large number of possible scenarios for the underlying asset’s price over time and then using these scenarios to estimate the option’s value.

To use the Monte Carlo simulation method to price an American option, you first need to define the option’s underlying variables, such as the stock price and time. You can then use a random number generator to simulate the asset’s price over time, using a model such as the geometric Brownian motion model.

Once you have simulated the asset’s price over time, you can use the option’s payoffs at each time step to determine the option’s value. You can repeat this process a large number of times to generate a distribution of possible option values, and then use statistical techniques to estimate the option’s expected value.

Here is a simple example of how the Monte Carlo simulation method can be implemented in Python to price an American call option on a non-dividend paying stock:

import numpy as np

# Monte Carlo simulation
def monte_carlo(S0, K, r, sigma, T, N):
dt = T / N # Time step
S = np.empty(N+1) # Stock price
S[0] = S0
for i in range(N):
S[i+1] = S[i] * np.exp((r - sigma**2 / 2) * dt + sigma * np.sqrt(dt) * np.random.normal())
return np.maximum(S[-1] - K, 0) # Payoff of call option

# Parameters
S0 = 100 # Initial stock price
K = 90 # Strike price
r = 0.05 # Risk-free interest rate
sigma = 0.2 # Volatility
T = 1 # Time to expiration (in years)

# Run Monte Carlo simulation
M = 10000 # Number of simulations
payoffs = np.empty(M)
for i in range(M):
payoffs[i] = monte_carlo(S0, K, r, sigma, T, N)

# Estimate option price
C = np.exp(-r * T) * payoffs.mean()

print(C)

In the code above, we define a function monte_carlo that uses the Monte Carlo simulation method to simulate the price of the underlying asset over time. It uses the geometric Brownian motion model to generate a random path for the asset's price, and it returns the payoff of the American call option at the end of the simulation.

The code then runs the monte_carlo function a large number of times and stores the payoffs in an array. It uses the mean of the payoffs to estimate the option price, and it adjusts the estimate for the time value of money using the risk-free interest rate. Finally, it prints the estimated option price.

The plot below is a representation of the Monte Carlo option values in comparison with the closed-form Black Scholes model:

Conclusion:

For the finite difference methods, both the implicit and Crank-Nicolson models achieve higher accuracy and stability as opposed to the explicit method. Nonetheless, the Monte Carlo simulation method can handle a wide range of option types and underlying assets, and it is relatively simple to implement. However, it can be computationally intensive, and it may not always be the most accurate method for pricing options.

Thanks again! Reach out for more discussion on the same!

--

--