Pricing Asian Options Using Maple
Overview
Solution using PDE
Monte Carlo Simulation
This example demonstrates the use of Maple for computing the price of an Asian option, a derivative security that has gained popularity in financial markets in recent years.
The payoff of an Asian option is based on the difference between an asset's average price over a given time period, and a fixed price called the strike price. Asian options are popular because they tend to have lower volatility than options whose payoffs are based purely on a single price point. It is also harder for big traders to manipulate an average price over an extended period than a single price, so Asian options offer further protection against risk.
One disadvantage of Asian options is that their prices are very hard to compute using standard techniques. Unlike European options, which can be priced using the classic Black-Scholes formula, there is no analytical formula for pricing an Asian option when the underlying asset is assumed to have a lognormal distribution, which is par for the course in financial modeling.
In this application, we use two different approaches for computing these prices numerically: (1) Solving a partial differential equation, and (2) Monte Carlo simulation. We will see that the two numerical solutions that Maple derives are the same, providing strong validation for both the techniques and Maple's numerics.
Arithmetic average Asian options are securities whose payoff depends on the average of the underlying stock price over a certain period of time. To be more precise, the value of the continuous arithmetic Asian call option at time t≤T is given by
Ct,T⁡K=ⅇ−r⁢T−t⁢E⁢1T⋅∫0TSuⅆu−K`+`
where S⁡t is the stock price at time t, T is the expiration date, and K is the strike price.
Since no general closed form solution for the price of the arithmetic average Asian option is known, a variety of numerical methods have been developed. These include formulation as a PDE, Monte Carlo simulation, and numeric inversion of the Laplace transform.
This worksheet demonstrates how these methods can be easily and quickly implemented in Maple.
restart
It was shown by Vecer ("A new PDE approach for pricing arithmetic average Asian options", Journal of Computational Finance, vol. 4, No. 4, 105-113) that the value of the Asian option at time 0 is given by
S0⁢u⁢0,S0−KS0
where the function u⁡t,z satisfies the following PDE
ut+r⁢1−tT−z⁢uz+11−tT−z2⁢σ2⁢uz,z2=0
with the boundary condition
u⁡T,z=z`+`
First set up the equation:
q:=t→1−tT
q≔t↦1−tT
pde:=∂∂t⁢u⁡t,z+r⁢1−tT−z⁢∂∂z⁢u⁡t,z+1⁢1−tT−z2⁢σ2⁢∂2∂z2⁢u⁡t,z2=0
pde≔∂∂tu⁡t,z+r⁢1−tT−z⁢∂∂zu⁡t,z+1−tT−z2⁢σ2⁢∂2∂z2u⁡t,z2=0
State the boundary conditions:
ibc1:=u⁡T,z=max⁡z,0;
ibc1≔u⁡T,z=max⁡0,z
ibc2:=u⁡t,−1=0;
ibc2≔u⁡t,−1=0
ibc3≔D2ut,1=1
ibc3≔D2⁡u⁡t,1=1
Assume a fixed time period of 1 year, a starting price of $100, a strike price of $100, and σ=0.2, r=0.15.
T:=1;
T≔1
S0≔100.;
S0≔100.
K:=100.;
K≔100.
σ:=0.2;
σ≔0.2
r:=0.15;
r≔0.15
Now solve the PDE numerically using Maple's built-in solver.
pds:=pdsolve⁡pde,ibc1,ibc2,ibc3,numeric,spacestep=0.001,timestep=0.001;
pds≔module...end module
val:=pds:-value⁡t=0,output=listprocedure:
v:=eval⁡u⁡t,z,val:
Finally, evaluate S0⁢u⁢0,S0−KS0 at Maple's numerical solution of uz,t.
S0 vS0−KS0;
8.40878796725886
Alternatively, you can use Monte Carlo methods to price the option.
restart:
withFinance:
Use the same parameters as before.
In the classical Black-Scholes setting the underlying process follows a geometric Brownian motion with drift r and volatility σ.
X:=GeometricBrownianMotion⁡S0,r,σ:
Option payoff is given by the following formula:
P:=t→max⁡1⁢∫0tX⁡uⅆut−K,0
P≔t↦max⁡∫0tX⁡uⅆut−K,0
You can compute the risk-neutral price of the option as expected payoff discounted using the risk-free rate.
DiscountFactor⁡r,T⁢ExpectedValue⁡P⁡T,timesteps=100,replications=104,output=value
8.342020802
You can also compute the various statistics related to this simulation.
for t in ExpectedValuePT,'timesteps'=100,'replications'=104,'output'='all' do printf%20s = %a\n, opt end:
value = 9.569798629 standarderror = .1022532908 standarddeviation = 10.22532908 skewness = 1.172246041 kurtosis = 1.209793258 minimum = 0. maximum = 66.07060386 standarderror = .1022532908
Alternatively, you can generate sample data for the option payoff and perform statistical analysis of this data.
A:=SampleValues⁡P⁡T,'timesteps'=100,'replications'=104;
eval(NextAfter);
procoptionbuiltin=NextAfter;end proc
A:=StatisticsRemoveInRange⁡A,−∞ ..NextAfter0, ∞;
StatisticsHistogram⁡A,bincount=20,averageshifted=4, axes = BOXED, gridlines = true;
You can consider a more complicated process for modeling the underlying asset. For example, you can incorporate jumps or regime switching into the model.
X1:=GeometricBrownianMotion⁡S0,r,σ;
X1≔_X4
X2:=GeometricBrownianMotion⁡S0,0.1,0.5;
X2≔_X5
P:=0.5,0.5|0.5,0.5;
P≔0.50.50.50.5
X3:=RegimeSwitchingProcess⁡P,X1,X2,1,12;
X3≔_X6
P:=X,t→max⁡∫0tX⁡uⅆut−K,0;
P≔X,t↦max⁡∫0tX⁡uⅆut−K,0
DiscountFactor⁡r,T⁢ExpectedValue⁡P⁡X1,T,timesteps=120,replications=104,output=value;
8.497423429
DiscountFactor⁡r,T⁢ExpectedValue⁡P⁡X2,T,timesteps=120,replications=104,output=value;
12.45400511
DiscountFactor⁡r,T⁢ExpectedValue⁡P⁡X3,T,timesteps=120,replications=104,output=value;
10.85502698
Here is another example. You can consider another asset, which also follows a geometric Brownian motion but with different parameters. Consider the following contract: if the average price of X is higher than the average price of Y, we get the difference, otherwise we get nothing.
Y:=GeometricBrownianMotion⁡S0,0.1,0.5:
P:=t→max1⁢∫0tX⁡uⅆut−1⁢∫0tY⁡uⅆut,0:
ExpectedValue⁡P⁡T,timesteps=100,replications=104;
value=14.85155399,standarderror=0.1795756877
Again, in order to get the price of the contract, you have to discount this expected value using the risk-free rate.
DiscountFactor⁡r,T⁢ExpectedValue⁡P⁡T,timesteps=100,replications=104,output=value;
12.83559553
Download Help Document