Pricing Barrier Options Using Monte Carlo Methods

080000231 080014353 080005791 080014940 040004438

Barrier Options Computational Finance EC5704 Spring 2009

Table of Contents Introduction ........................................................................................................................ 3 Background ......................................................................................................................... 3 In + Out = Vanilla............................................................................................................. 3 Volatility .......................................................................................................................... 5 Analytical Solutions......................................................................................................... 6 Monte Carlo Simulation...................................................................................................... 7 Simulation and Path Generation..................................................................................... 7 Calculating Payoffs.......................................................................................................... 9 Challenges ..................................................................................................................... 11 Variance Reduction Techniques ....................................................................................... 11 Antithetic Variates ........................................................................................................ 12 Control Variates ............................................................................................................ 13 Results............................................................................................................................... 14 Error Analysis ................................................................................................................ 15 Conclusion......................................................................................................................... 17 Appendix A: Results .......................................................................................................... 18 Appendix B: Code.............................................................................................................. 19

Introduction This paper analyses the pricing of barrier options using Monte Carlo methods. Variance reduction techniques are also analysed and implemented in the pricing of various barrier options. The corresponding error is also investigated to reveal any possible relationships between the simulation parameters and the resulting error.

Background A barrier option is a type of path-dependent option where the payoff is determined by whether or not the price of the stock hits a certain level during the life of the option. There are several different types of barrier options. We consider two general types of barrier options, ‘in’ and ‘out’ options. An ‘out’ option only pays off if the stock does not hit the barrier level throughout the life of the option. If the stock hits a specified barrier, then it has knocked out and expires worthless. A knock-in option on the other hand only pays out if the barrier is crossed during the life of the option. For each knockout and knock-in option, we can have a down or an up option.

Call Put

Knock-out Up and Out Down and Out Up and Out Down and Out

Knock-in Up and In Down and In Up and In Down and In

We therefore have four different basic barrier options each having the typical payoffs. It is worth pointing out that the relationship between the barrier and spot price indicates whether the option is an up option or down option. If the barrier price is above the spot price, then the option is an up option; if the barrier is below the spot price, then it is a down option.

In + Out = Vanilla For a given set of parameters, a vanilla option can be replicated by combining an in and an out option of the same type. For example, one could combine a down and out call with a down and in call to replicate a plain vanilla call. These are equivalent due to the fact that if one option gets knocked out, the other is knocked in and vise versa. This idea is illustrated in the figure below by graphing the values for an up and out put and an up and in put using various barrier levels. The vanilla put was valued with the same parameters using MATLAB’s built in Black-Scholes formula for pricing options, blsprice.

Barrier Option 4.5

4 Up and Out Up and In Vanilla Put

3.5

3

Option Price

2.5

2

1.5

1

0.5

0

-0.5 50

55

60

65

70

75

Barrier

Figure 1: Spot: $50, Strike: $50, Interest Rate: 10%, Time to Expiry: 5 Months, Volatility: 40% with barriers ranging from $50-$75 It is quite clear that for any given barrier level, adding the value of an up and out put and an up and in put totals to the value of a vanilla put. It is also apparent that both barrier options are cheaper than their vanilla counterpart. This is because barrier options limit the potential of either getting in the money or staying in the money, thereby reducing the value relative to a vanilla option. The figure also illustrates the relationship between the option price and barrier level. For a knock-out option, increasing the absolute difference between the barrier level and the initial spot price has the effect of increasing the value of the option. As the absolute difference increases, the value of the option converges towards the value of a plain vanilla option. This is simply because the probability of the option knocking out approaches zero as the absolute difference increases. For a knock-in option, increasing the absolute difference has the opposite effect. The probability of the option knocking in approaches zero as the absolute difference increase thus making the option worthless.

Volatility Volatility plays an important role in valuing options. Barrier options are especially sensitive to volatility. For knock-out options, increased volatility has the effect of decreasing the option value as knock-out becomes more probable. Knock-in options however increase in price with increased volatility as knock-in becomes more probable.

4.5 4 3.5

Option Value

3 2.5 2 1.5 1 0.5 0 -0.5 60 0.3

55 0.25 50

0.2 0.15

45 0.1 Spot Price

40

0.05 Volatility

Figure 2: Analytical solutions to an up and out call option. Barrier: $60, Strike: $50, Interest Rate: 10%, Time to Expiry: 1 Year. For the up and out call option, the graph seems quite odd. It seems that there is only a very small range for which the option is valuable. The peak of the curve can be interpreted as the spot price that has the best chance of getting into the money without hitting the barrier. However, as the volatility increases the chance of the option getting knocked out increases and consequently the option price declines.

18 16 14

Option Value

12 10 8 6 4 2 0 60 55

0.3 0.25 50

0.2 0.15

45 0.1 Spot Price

40

0.05 Volatility

Figure 3: Analytical solutions to an up and in call option. Barrier: $60, Strike: $50, Interest Rate: 10%, Time to Expiry: 1 Year It is evident from the figure that increasing the volatility for any given spot price has the effect of increasing the value of the up and in call option. It is also clear that as the spot price approaches the barrier level, the value of the option increases due to the increased probability of the option knocking in.

Analytical Solutions The barrier options priced in this paper all have readily available analytical solutions. The formulae for pricing the barrier options have been converted to MATLAB code to check the Monte Carlo simulation results. See Appendix B for the code.

Monte Carlo Simulation The basic principle behind Monte Carlo simulation is to simulate as many possible scenarios and to average those scenarios to get an expectation. To price options using Monte Carlo methods, many price paths are generated using the evolution of the stock according to geometric Brownian motion:

S (t ) = S (0)e (( r −0.5σ

2

) t +σ t Z

Z ~ N (0,1) The option payoff is then calculated for each price path. An accurate estimation of option price is obtained by discounting the average of all the payoffs. Barrier options are path-dependent options - the payoff is determined by whether or not the price of the asset hits a certain level, the barrier, during the life of the option. Due to this path-dependency, simulation of the entire price evolution is necessary. Therefore, the first step in pricing barrier options using Monte Carlo methods is to simulate the price evolution.

Simulation and Path Generation The easiest way to generate multiple price paths is by using a 2-D matrix. In this matrix, each row represents an asset path. Therefore, by simply increasing the number of rows one can increase the number of simulated asset paths. The columns of this matrix represent the time steps of the evolution. Therefore, each column represents a price evolution. A function, AssetPaths, was created using MATLAB to create a matrix of asset paths. The size of this matrix was determined by the parameters passed, namely the desired number of simulations and the number of steps – price evolutions - required over the specified time interval. This function returns the matrix of asset paths to be used for pricing options. The graph below was created using the function written for generating asset paths to illustrate the concept behind Monte Carlo simulation. Only seven paths were generated in this simulation to maintain clarity. In practice, thousands of simulations would be required to accurately price an option.

Simulated Paths 70

65

Asset Price

60

55

50

45

40

0

0.1

0.2

0.3

0.4

0.5 Time

0.6

0.7

0.8

0.9

1

Figure 4: Simulated price evolutions of a stock with the following characteristics: Spot price: $50, Interest Rate: 10%, Volatility: 15%. Simulation Parameters: Time Interval: 1 Year, Number of Steps: 365 The graph simulating seven asset paths also shows a hypothetical barrier level to help illustrate the concept behind barrier options. In this simulation, two of the asset paths hit the barrier. If the matrix in this simulation were used to price a knock-out option, particularly an up and out option, then the payoffs associated with these two asset paths would be zero. The remaining five payoffs would be evaluated as one would evaluate a vanilla option. The path-dependency of a barrier option also necessitates a sufficient number of steps to accurately model the price evolution.

A graphical interpretation of the number of steps used in the simulation is presented below: Simulated Paths 60 10 Steps 100 Steps 500 Steps

Asset Price

55

50

45

0

0.1

0.2

0.3

0.4

0.5 Time

0.6

0.7

0.8

0.9

1

Figure 5: Number of steps illustration It is clear that if there are an insufficient number of steps used over a given time interval, a barrier might not be hit that might have otherwise been hit had the number of steps been increased. Thus, for a given time interval, increasing the number of steps in the simulation increases the number of price evolutions leading to an increase in the accuracy of the simulation.

Calculating Payoffs Once multiple asset paths have been simulated, the next step is to determine the payoff for each asset path. This is done by evaluating each path to see whether it has hit the barrier. The payoff is then dependent on the type of barrier option and the knowledge of whether or not the barrier level has been hit during the life of the option. Some options will expire worthless if the barrier is reached, others will be worthless unless the barrier is reached. Regardless, coding this for the Monte Carlo simulation is rather trivial once the asset path has been generated. The amount of code required can be significantly reduced by noting that whether the option is an up or down option is determined by the relationship between the barrier level and the initial spot price. A barrier above the initial spot price indicates that the option is an up option, while a barrier below indicates a down option.

This information limits the coding to just four functions to price eight different barrier options: • • • •

knock-out call knock-out put knock-in call knock-in put

Whether the option is an up or down option is then evaluated within each function by checking if the barrier is above or below the initial spot price. Code for pricing a knock-out call is provided to show how the payoffs were calculated. Commented MATLAB code for pricing a knock-out call function [price, std, CI] = ... knockout_call(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and out, down and out call % if the barrier is less than the spot price, then down and out call % if the barrier is above the spot price, then up and out call payoff = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path knocked = barrierCrossing(spot,Sb,path); % determine if up or down % knocked = 1 if barrier has been hit if knocked == 0 % path has always been above/below the barrier payoff(i) = max(0,path(NSteps+1)-K); % use the last asset price end end % return mean and confidence interval of present value payoff vector [price, std, CI] = normfit(exp(-r*T)*payoff);

Challenges Computational Costs Monte Carlo simulation, while a relatively easy way to value options, is computationally slow as accuracy requires an increase in the number of simulations.

0.08

0.07

Time Elapsed

0.06 0.05 0.04 0.03 0.02 0.01 100 80 60

800 40 400

20 Number of Simulations

0

100

200

500

600

900

1000

700

300 Number of Steps

Figure 6: Computational cost of increasing simulations and step size Using the tic, toc function in MATLAB, a graphic depiction of the computational cost of accuracy is made possible. The spikes are likely due to interference of running background processes on the computer used to generate the graph.

Variance Reduction Techniques The concept of variance reduction arises from the following relationship for determining a confidence interval:

Var ( X ) M

Where M is the sample size and Var(X) is the variance. There are only two ways to decrease the size of the confidence interval. Either increase the sample size or decrease the variance. In Monte Carlo simulations, increasing the sample size is computationally costly as was illustrated previously. Variance reduction works by replacing X with a random variable Y that has the same expectation of X, but with a smaller variance. We can compute the expectation of X via Monte Carlo methods using Y. Using Y in this instance allows for the same expectation but with a lower variance.

Antithetic Variates Variance reduction using antithetic variates is based on the following lemma: If X is an arbitrary random variable and f is a monotonically increasing or monotonically decreasing function then,

cov( f ( X ), f (− X )) ≤ 0 Let us now consider the case of a random variable which is of the form f (U), where U is a standard normally distributed random variable, i.e. f(U) ~ N(0,1). The standard normal distribution is symmetric, and hence also – U ~ N(0, 1). It then follows obviously that f (U) ~ f (-U). In particular they have the same expectation. Therefore, in case we want to compute the expectation of X = f (U), we can take Z = f (-U) and

Y=

f (U ) + f (− U ) 2

If we now assume that the map f is monotonically increasing, then we conclude from the above lemma that

cov( f (U ), f (−U )) ≤ 0 and we finally obtain X +Z Ε  = Ε( X )  2   1 X +Z  X + Z X + Z  1  Var  , ( X , Z ) + Var (Z ) ≤ Var ( X )  = cov  = Var ( X ) + 2Cov 14243 123  2 2  4  2   2 ≤0 =Var ( X )   where X = f(U) and Z = f(-U).

Basically, in an option pricing context, the payoff of an option is a function of the asset price at expiration which is generated by a standard normally distributed random variable. Following the method outlined above, we calculate another payoff depending upon a final asset price generated by the opposite of our standard normally distributed random variable. Taking the average of the two discounted payoffs gives us the more accurate option value.

Control Variates Using control variates is another good way to reduce the variance and generate more accurate results. Suppose we again wish to estimate E(X). If we can somehow find another random variable Y, which is close to X in some sense, with known expected value E(Y) then we can define a random variable Z via

Z = X + c(Y − Ε(Y )) where c ∈ ℜ and Y is called the control variate. It is obvious that Ε(Z ) = Ε( X ) , so we can apply Monte Carlo simulation to Z in order to compute E(X). On the other hand, we have

Var (Z ) = Var ( X + cY ) = Var ( X ) + 2cCov( X , Y ) + c ²Var (Y ) As we want to minimize this function we look for the value of c that does exactly that. Straightforward calculation gives c min =

− Cov( X , Y ) Var (Y )

In practice, however, Cov(X,Y) and Var(Y) are generally unknown so we have to estimate them via Monte Carlo simulation. To provide more insight into this method, consider an up and out call option. A suitable control variate has to be chosen depending on the particular problem we are dealing with. As we want to price a barrier option, we can assume a plain vanilla call option is a suitable control variate. Its exact value is known as it can be calculated analytically using the Black-Scholes formula. Furthermore, its expected value and variance can easily be estimated. To estimate the covariance between the up and out call price and the plain vanilla call price, we estimate the corresponding payoffs via Monte Carlo simulations. This may be accomplished by using a set of pilot replications. Next, we are able to calculate the value of c that minimizes the variance of Z. In a final step we run a Monte Carlo simulation once more to estimate E(Z) using the previously obtained value of c.

To avoid inducing bias in the estimate of E(Z), pilot replications are used to select c and to generate an estimate of E(Z). We run another set of replications for the final estimation of E(Z).

Results The results for an up and in call barrier option priced using Monte Carlo simulation methods is presented for discussion. Table 1: Monte Carlo simulation results for an up and in call option. Up and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $8.1541 Estimated $7.3903 $7.9892 $8.1701 $8.1639 σ $10.2415 $11.9190 $12.3546 $12.1117 95% Confidence Interval $5.3582 $9.4224 $7.2496 $8.7288 $7.9280 $8.4123 $8.0888 $8.2390 Error $0.7638 $0.1649 $0.0160 $0.0098 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Table 2: Monte Carlo using Variance Reduction: Antithetic Variates Up and In Call Number of Simulations 250 Steps 100 1000 10000 100000 True $8.1541 Estimated $7.4464 $7.9892 $8.0408 $8.1269 σ $5.5213 $6.2478 $6.3965 $6.3681 95% Confidence Interval $6.3508 $8.5419 $7.6015 $8.3769 $7.9154 $8.1662 $8.0874 $8.1664 Error $0.7077 $0.1649 $0.1133 $0.0272 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

It is clear that in comparing the use of antithetic variates to the standard Monte Carlo simulation, the variance of our estimates have been significantly reduced. We halved the number of steps used to make the comparison fair as the use of antithetic pairs doubles the sample used for the simulation.

Table 3: Monte Carlo using Variance Reduction: Control Variates Up and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True $8.1541 Estimated $7.9724 $8.0699 $8.1147 $8.1209 σ $1.2545 $1.1733 $1.0661 $1.0428 95% Confidence Interval $7.7235 $8.2214 $7.9971 $8.1427 $8.0938 $8.1356 $8.1144 $8.1273 Error $0.1817 $0.0841 $0.0394 $0.0332 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

As expected the use of control variates greatly reduces the size of the confidence interval. Generally, the use of control variates leads to more accurate estimations of the option price. The results here however can be somewhat misleading. The time elapsed in performing the simulations is not included. While variance reduction techniques often lead to increased accuracy for a given number of simulations, the computational cost associated with the increased accuracy is still present. Results for all the barrier options considered are provided in Appendix A.

Error Analysis Increasing the accuracy of Monte Carlo simulations requires a reduction of error. There are likely to be many sources of error in a simulation, such as the estimates of the parameters themselves or the assumptions made within the model itself, however, if these errors can be eliminated, then one of the simplest ways to increase accuracy is to increase the number of simulations. The use of variance reduction techniques helps shortcut the need for an increased number of simulations.

0.5 0.45 0.4

Absolute Error

0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 0

0.05 0.1 1000

2000

0.15 3000

4000

5000

0.2 6000

7000

8000

0.25 9000

10000

0.3 Volatility

Number of Simulations

Figure 7: Error associated with an up and out call: Varying volatility and number of simulations Spot: $50, Barrier: $60, Strike: $50, Time to Expiry: 1 Year, Interest Rate: 10%, Number of Steps: 500. As expected, the error decreases as the number of simulations increases. It is interesting to note however, that the absolute error decreases as the volatility increases.

Error in valuing an up and in call option for various Monte Carlo methods 2 Standard Monte Carlo Monte Carlo: Antithetic Variates Monte Carlo: Control Variates

1.8 1.6 1.4

Error

1.2 1 0.8 0.6 0.4 0.2 0

0

1000

2000

3000

4000 5000 6000 Number of Simulations

7000

8000

9000

10000

Figure 8: Up and in call option with Spot: $50, Barrier: $60, Strike: $50, Interest Rate: 10%, Time to Expiry: 1 Year, Volatility: 30%. Errors plotted against number of simulations used to price the option for various Monte Carlo methods. The standard Monte Carlo method benefits greatly when the number of simulations increase. Using Monte Carlo with antithetic variates, the number of simulations required to obtain the same accuracy as the standard Monte Carlo method is greatly reduced. The use of control variates quickly results in accurate estimates of the option price. It appears as if increasing the number of simulations does virtually nothing to reduce the error. This however may not necessarily be the case as the current figure distorts the interpretation for the control variates result. In any case, it is clear that the Monte Carlo method using control variates is much more effective at estimating the option price for a given number of simulations. However, the time elapsed to value the option is not considered and would surely factor into the evaluation of which method is ultimately more efficient.

Conclusion In this paper an overview of how to estimate barrier option prices via Monte Carlo simulation was provided. We discussed the presence of computational errors and described variance reduction techniques to eliminate these as much as possible. In performing Monte Carlo simulations, many of the advantages and disadvantages of the method have been made apparent. Monte Carlo simulations, while rather easy to perform, come at the cost of computational efficiency. In approaching Monte Carlo simulation, one should be cautious of the assumptions made and be aware of the errors that may be present in the parameters used. The use of variance reduction techniques can significantly increase the accuracy of the estimates leading to a reduction in the number of simulations needed for a given confidence interval; however this increased accuracy also imposes a computational cost which erodes some of the benefits gained.

Appendix A: Results The results were exported from MATLAB to Excel and for this reason they are attached on the subsequent page.

Monte Carlo Down and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.5310 Estimated $0.2368 $0.4275 $0.4609 $0.4774 σ $1.2271 $2.2273 $2.4001 $2.4980 95% Confidence Interval -$0.0066 $0.4803 $0.2893 $0.5657 $0.4138 $0.5079 $0.4619 $0.4929 Error $0.2942 $0.1035 $0.0702 $0.0536 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $7.8360 Estimated $8.2700 $8.0048 $7.7761 $7.8692 σ $11.2728 $11.9747 $12.0294 $12.0074 95% Confidence Interval $6.0333 $10.5068 $7.2617 $8.7479 $7.5403 $8.0119 $7.7948 $7.9436 Error $0.4340 $0.1687 $0.0599 $0.0332 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and In Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $3.2686 Estimated $2.9489 $3.5473 $3.1877 $3.2432 σ $5.5979 $5.9763 $5.6809 $5.7170 95% Confidence Interval $1.8381 $4.0596 $3.1765 $3.9182 $3.0763 $3.2990 $3.2078 $3.2787 Error $0.3197 $0.2787 $0.0809 $0.0254 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.3403 Estimated $0.6581 $0.4271 $0.3933 $0.3774 σ $1.7870 $1.3414 $1.3298 $1.3022 95% Confidence Interval $0.3035 $1.0127 $0.3439 $0.5103 $0.3673 $0.4194 $0.3693 $0.3855 Error $0.3177 $0.0868 $0.0530 $0.0371 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Up and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $8.1541 Estimated $7.3903 $7.9892 $8.1701 $8.1639 σ $10.2415 $11.9190 $12.3546 $12.1117 95% Confidence Interval $5.3582 $9.4224 $7.2496 $8.7288 $7.9280 $8.4123 $8.0888 $8.2390 Error $0.7638 $0.1649 $0.0160 $0.0098 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.2130 Estimated $0.3168 $0.2615 $0.2330 $0.2448 σ $1.1738 $1.0774 $1.0110 $1.0479 95% Confidence Interval $0.0839 $0.5497 $0.1947 $0.3284 $0.2132 $0.2528 $0.2383 $0.2513 Error $0.1038 $0.0486 $0.0200 $0.0318 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and In Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.5559 Estimated $0.3377 $0.4702 $0.4900 $0.4964 σ $1.5891 $1.9919 $1.9785 $2.0367 95% Confidence Interval $0.0224 $0.6530 $0.3466 $0.5938 $0.4513 $0.5288 $0.4837 $0.5090 Error $0.2182 $0.0857 $0.0659 $0.0596 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $3.0530 Estimated $3.4964 $3.1870 $3.1694 $3.0997 σ $5.7948 $5.5438 $5.6284 $5.5334 95% Confidence Interval $2.3466 $4.6462 $2.8430 $3.5311 $3.0590 $3.2797 $3.0654 $3.1340 Error $0.4433 $0.1340 $0.1163 $0.0467 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Monte Carlo using Antithetic Variates Down and In Call Number of Simulations 250 Steps 100 1000 10000 100000 True Value $0.5310 Estimated $0.3338 $0.4818 $0.4270 $0.4511 σ $1.0282 $1.6899 $1.6455 $1.7081 95% Confidence Interval $0.1298 $0.5378 $0.3769 $0.5866 $0.3947 $0.4592 $0.4406 $0.4617 Error $0.1973 $0.0493 $0.1041 $0.0799 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Call Number of Simulations 250 Steps 100 1000 10000 100000 True Value $7.8360 Estimated $8.3832 $8.0819 $8.0490 $7.9325 σ $7.0658 $6.5490 $6.4208 $6.4219 95% Confidence Interval $6.9812 $9.7852 $7.6755 $8.4883 $7.9231 $8.1748 $7.8927 $7.9723 Error $0.5472 $0.2459 $0.2130 $0.0965 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and In Put Number of Simulations 250 Steps 100 1000 10000 100000 True Value $3.2686 Estimated $3.1268 $3.0805 $3.1926 $3.2259 σ $3.1688 $3.3764 $3.3158 $3.3352 95% Confidence Interval $2.4980 $3.7555 $2.8710 $3.2901 $3.1276 $3.2575 $3.2053 $3.2466 Error $0.1418 $0.1881 $0.0761 $0.0427 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Put Number of Simulations 250 Steps 100 1000 10000 100000 True Value $0.3403 Estimated $0.4274 $0.4302 $0.4008 $0.3973 σ $0.8506 $0.9168 $0.9147 $0.9050 95% Confidence Interval $0.2586 $0.5961 $0.3733 $0.4871 $0.3829 $0.4188 $0.3917 $0.4029 Error $0.0870 $0.0898 $0.0605 $0.0569 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Up and In Call Number of Simulations 250 Steps 100 1000 10000 100000 True Value $8.1541 Estimated $7.4464 $7.9892 $8.0408 $8.1269 σ $5.5213 $6.2478 $6.3965 $6.3681 95% Confidence Interval $6.3508 $8.5419 $7.6015 $8.3769 $7.9154 $8.1662 $8.0874 $8.1664 Error $0.7077 $0.1649 $0.1133 $0.0272 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Call Number of Simulations 250 Steps 100 1000 10000 100000 True Value $0.2130 Estimated $0.1392 $0.2320 $0.2587 $0.2625 σ $0.4476 $0.7223 $0.7603 $0.7637 95% Confidence Interval $0.0504 $0.2280 $0.1872 $0.2768 $0.2438 $0.2736 $0.2578 $0.2672 Error $0.0737 $0.0190 $0.0457 $0.0495 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and In Put Number of Simulations 250 Steps 100 1000 10000 100000 True Value $0.5559 Estimated $0.3692 $0.5250 $0.4785 $0.4848 σ $1.2150 $1.4640 $1.3787 $1.3828 95% Confidence Interval $0.1281 $0.6103 $0.4342 $0.6159 $0.4514 $0.5055 $0.4762 $0.4933 Error $0.1867 $0.0309 $0.0774 $0.0712 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Put Number of Simulations 250 Steps 100 1000 10000 100000 True Value $3.0530 Estimated $3.6459 $3.1782 $3.1508 $3.1353 σ $3.4791 $3.2969 $3.2711 $3.2365 95% Confidence Interval $2.9556 $4.3362 $2.9736 $3.3828 $3.0867 $3.2149 $3.1152 $3.1553 Error $0.5929 $0.1252 $0.0978 $0.0823 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Monte Carlo using Control Variates Down and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.5310 $0.7455 $0.3590 $0.4658 $0.4741 Estimated σ $3.1011 $2.0914 $2.3979 $2.5271 95% Confidence Interval $0.1301 $1.3608 $0.2293 $0.4888 $0.4188 $0.5128 $0.4584 $0.4897 Error $0.2144 $0.1720 $0.0653 $0.0570 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $7.8360 $7.5204 $7.9893 $7.9138 $7.9027 Estimated $2.8588 $2.1183 $2.4377 $2.4713 σ 95% Confidence Interval $6.9532 $8.0877 $7.8579 $8.1208 $7.8660 $7.9616 $7.8874 $7.9180 Error $0.3156 $0.1533 $0.0777 $0.0667 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and In Put Number of Simulations 500 Steps 100 1000 10000 100000 $3.2686 True Value $3.1211 $3.2183 $3.2358 $3.2338 Estimated $1.6185 $1.3379 $1.2852 $1.2943 σ 95% Confidence Interval $2.7999 $3.4422 $3.1353 $3.3014 $3.2106 $3.2610 $3.2258 $3.2418 $0.1475 $0.0503 $0.0328 $0.0348 Error S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Down and Out Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.3403 Estimated $0.0578 $0.3952 $0.4001 $0.3774 σ $0.4590 $1.3780 $1.3420 $1.2962 95% Confidence Interval -$0.0333 $0.1489 $0.3097 $0.4808 $0.3738 $0.4264 $0.3694 $0.3855 Error $0.2825 $0.0549 $0.0598 $0.0371 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $40

Up and In Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $8.1541 Estimated $7.9724 $8.0699 $8.1147 $8.1209 σ $1.2545 $1.1733 $1.0661 $1.0428 95% Confidence Interval $7.7235 $8.2214 $7.9971 $8.1427 $8.0938 $8.1356 $8.1144 $8.1273 Error $0.1817 $0.0841 $0.0394 $0.0332 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Call Number of Simulations 500 Steps 100 1000 10000 100000 True Value $0.2130 Estimated $0.2600 $0.2207 $0.2546 $0.2441 σ $1.1069 $0.9568 $1.0687 $1.0362 95% Confidence Interval $0.0403 $0.4796 $0.1613 $0.2801 $0.2337 $0.2756 $0.2377 $0.2505 $0.0470 $0.0077 $0.0416 $0.0311 Error S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and In Put Number of Simulations 100 1000 10000 100000 500 Steps True Value $0.5559 Estimated $0.6797 $0.4610 $0.4958 $0.5000 σ $2.2537 $1.8825 $1.9877 $1.9908 95% Confidence Interval $0.2326 $1.1269 $0.3442 $0.5778 $0.4569 $0.5348 $0.4877 $0.5124 $0.1238 $0.0949 $0.0601 $0.0559 Error S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Up and Out Put Number of Simulations 500 Steps 100 1000 10000 100000 True Value $3.0530 Estimated $3.1737 $3.1822 $3.1104 $3.1027 σ $1.6173 $1.7722 $2.0042 $1.9998 95% Confidence Interval $2.8528 $3.4946 $3.0722 $3.2922 $3.0711 $3.1497 $3.0903 $3.1151 Error $0.1207 $0.1292 $0.0574 $0.0497 S = $50, K = $50, r = 10%, σ = 30%, T = 1 Barrier = $60

Appendix B: Code The following m-files have been appended to the file. They follow the same order as listed and begin on the subsequent page. MATLAB Code for Analytical Solutions: A_DICall.m A_DOCall.m A_DIPut.m A_DOPut.m A_UICall.m A_UOCall.m A_UIPut.m A_UOPut.m MATLAB Code for Monte Carlo: knockin_call knockout_call knockin_put knockout_put MATLAB Code for Monte Carlo using Antithetic Variates: AVknockin_call.m AVknockout_call.m AVknockin_put.m AVknockout_put.m MATLAB Code for Monte Carlo using Control Variates: CVknockin_call.m CVknockout_call.m CVknockin_put.m CVknockout_put.m MATLAB Code used to generate data: AssetPaths.m AssetPathsAV.m AVdataGenerator.m barrierCrossing.m dataGenerator.m CVdataGenerator.m graphingFunction.m plottingFunction.m

function price = A_DICall(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% down and in call if E > Sb price = S*exp(-q*(T))*(b*(1-normcdf(d8)))... -E*exp(-r*(T))*(a*(1-normcdf(d7))); else price = S*exp(-q*(T))*(normcdf(d1)-normcdf(d3)+b*(1-normcdf(d6)))... -E*exp(-r*(T))*(normcdf(d2)-normcdf(d4)+a*(1-normcdf(d5))); end

function price = A_DOCall(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% down and out call if E > Sb price = S*exp(-q*(T))*(normcdf(d1)-b*(1-normcdf(d8)))... -E*exp(-r*(T))*(normcdf(d2)-a*(1-normcdf(d7))); else price = S*exp(-q*(T))*(normcdf(d3)-b*(1-normcdf(d6)))... -E*exp(-r*(T))*(normcdf(d4)-a*(1-normcdf(d5))); end

function price = A_DIPut(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); % d1 % d2 d3 = d4 = d5 = d6 = d7 = d8 =

= (log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); = (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% down and in put price = -S*exp(-q*(T))*(1 - normcdf(d3)+... b*(normcdf(d8)-normcdf(d6)))... +E*exp(-r*(T))*(1 - normcdf(d4)+a*(normcdf(d7)-normcdf(d5)));

function price = A_DOPut(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% down and out put price = -S*exp(-q*(T))*(normcdf(d3)-normcdf(d1)-... b*(normcdf(d8)-normcdf(d6)))... +E*exp(-r*(T))*(normcdf(d4)-normcdf(d2)-a*(normcdf(d7)-normcdf(d5)));

function price = A_UICall(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); % d1 % d2 d3 = d4 = d5 = d6 = d7 = d8 =

= (log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); = (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% up and in call price = S*exp(-q*(T))*(normcdf(d3)+b*(normcdf(d6)-normcdf(d8)))... -E*exp(-r*(T))*(normcdf(d4)+a*(normcdf(d5)-normcdf(d7)));

function price = A_UOCall(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% up and out call price = S*exp(-q*(T))*(normcdf(d1)-normcdf(d3)-... b*(normcdf(d6)-normcdf(d8)))... -E*exp(-r*(T))*(normcdf(d2)-normcdf(d4)-a*(normcdf(d5)-normcdf(d7)));

function price = A_UIPut(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% up and in put if E > Sb price = -S*exp(-q*(T))*(normcdf(d3)-normcdf(d1)+b*(normcdf(d6)))... +E*exp(-r*(T))*(normcdf(d4)-normcdf(d2)+a*(normcdf(d5))); else price = -S*exp(-q*(T))*(b*(normcdf(d8)))... +E*exp(-r*(T))*(a*(normcdf(d7))); end

function price = A_UOPut(S,Sb,E,r,q,T,sigma) a = (Sb/S)^(-1+(2*(r-q)/sigma^2)); b = (Sb/S)^(1+(2*(r-q)/sigma^2)); d1 d2 d3 d4 d5 d6 d7 d8

= = = = = = = =

(log(S/E)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/E)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)+(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S/Sb)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q-0.5*sigma^2)*(T))/(sigma*sqrt(T)); (log(S*E/Sb^2)-(r-q+0.5*sigma^2)*(T))/(sigma*sqrt(T));

% up and out put if E > Sb price = -S*exp(-q*(T))*(1-normcdf(d3)-b*(normcdf(d6)))... +E*exp(-r*(T))*(1-normcdf(d4)-a*(normcdf(d5))); else price = -S*exp(-q*(T))*(1-normcdf(d1)-b*(normcdf(d8)))... +E*exp(-r*(T))*(1-normcdf(d2)-a*(normcdf(d7))); end

function [price, std, CI] = knockin_call(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and in, down and in call % if the barrier is less than the spot price, then down and in call % if the barrier is above the spot price, then up and in call payoff = zeros(NRepl,1); % col vector of payoffs for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier payoff(i) = max(0,path(NSteps+1)-K); end end % return mean and confidence interval of present value payoff vector [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price, std, CI] = knockout_call(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and out, down and out call % if the barrier is less than the spot price, then down and out call % if the barrier is above the spot price, then up and out call payoff = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 0 % path has always been above/below the barrier payoff(i) = max(0,path(NSteps+1)-K); % use the last price end end % return mean and confidence interval of present value payoff vector [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price, std, CI] = knockin_put(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and in, down and in put % if the barrier is less than the spot price, then down and in put % if the barrier is above the spot price, then up and in put payoff = zeros(NRepl,1); % col vector of payoffs for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier payoff(i) = max(0,K-path(NSteps+1)); end end % return mean and confidence interval of present value payoff vector [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price, std, CI] = knockout_put(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and out, down and out put % if the barrier is less than the spot price, then down and out put % if the barrier is above the spot price, then up and out put payoff = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 0 % path has always been above/below the barrier payoff(i) = max(0,K-path(NSteps+1)); % use the last price end end % return mean and confidence interval of present value payoff vector [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price,std,CI] = AVknockin_call(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and in, down and in call % if the barrier is less than the spot price, then down and in call % if the barrier is above the spot price, then up and in call payoffA = zeros(NRepl,1); % col vector of payoffs payoffB = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl [pathA, pathB] = AssetPathsAV(spot,r,sigma,T,NSteps,1); knockedA = barrierCrossing(spot,Sb,pathA); % determine if up or down knockedB = barrierCrossing(spot,Sb,pathB); if knockedA == 1 payoffA(i) = max(0,pathA(NSteps+1)-K); end if knockedB == 1 payoffB(i) = max(0,pathB(NSteps+1)-K); end end payoff = (payoffA + payoffB)./2; [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price,std,CI] = AVknockout_call(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and out, down and out call % if the barrier is less than the spot price, then down and out call % if the barrier is above the spot price, then up and out call payoffA = zeros(NRepl,1); % col vector of payoffs payoffB = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl [pathA, pathB] = AssetPathsAV(spot,r,sigma,T,NSteps,1); knockedA = barrierCrossing(spot,Sb,pathA); % determine if up or down knockedB = barrierCrossing(spot,Sb,pathB); if knockedA == 0 % path has not crossed barrier payoffA(i) = max(0,pathA(NSteps+1)-K); end if knockedB == 0 payoffB(i) = max(0,pathB(NSteps+1)-K); end end payoff = (payoffA + payoffB)./2; [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price,std,CI] = AVknockin_put(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and in, down and in put % if the barrier is less than the spot price, then down and in put % if the barrier is above the spot price, then up and in put payoffA = zeros(NRepl,1); % col vector of payoffs payoffB = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl [pathA, pathB] = AssetPathsAV(spot,r,sigma,T,NSteps,1); knockedA = barrierCrossing(spot,Sb,pathA); % determine if up or down knockedB = barrierCrossing(spot,Sb,pathB); if knockedA == 1 payoffA(i) = max(0,K-pathA(NSteps+1)); end if knockedB == 1 payoffB(i) = max(0,K-pathB(NSteps+1)); end end payoff = (payoffA + payoffB)./2; [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price,std,CI] = AVknockout_put(spot,Sb,K,r,T,sigma,NSteps,NRepl) % up and out, down and out put % if the barrier is less than the spot price, then down and out put % if the barrier is above the spot price, then up and out put payoffA = zeros(NRepl,1); % col vector of payoffs payoffB = zeros(NRepl,1); % col vector of payoffs % loop through the paths to determine payoffs for i=1:NRepl [pathA, pathB] = AssetPathsAV(spot,r,sigma,T,NSteps,1); knockedA = barrierCrossing(spot,Sb,pathA); % determine if up or down knockedB = barrierCrossing(spot,Sb,pathB); if knockedA == 0 % path has not crossed barrier payoffA(i) = max(0,K-pathA(NSteps+1)); end if knockedB == 0 payoffB(i) = max(0,K-pathB(NSteps+1)); end end payoff = (payoffA + payoffB)./2; [price, std, CI] = normfit(exp(-r*T)*payoff);

function [price, std, CI] = CVknockin_call(spot,Sb,K,r,T,sigma,NSteps,NRepl,NPilot) % up and in, down and in call % if the barrier is less than the spot price, then down and in call % if the barrier is above the spot price, then up and in call payoff = zeros(NPilot,1); vanillaPayoff = zeros(NPilot,1); muVanilla = blsprice(spot,K,r,T,sigma); for i=1:NPilot path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path vanillaPayoff(i) = max(0,path(NSteps+1)-K); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier payoff(i) = max(0,path(NSteps+1)-K); end end vanillaPayoff = exp(-r*T)*vanillaPayoff; payoff = exp(-r*T)*payoff; covMat = cov(vanillaPayoff, payoff); varVanilla = var(vanillaPayoff); c = -covMat(1,2)/varVanilla; newPayoff = zeros(NRepl,1); newVanillaPayoff = zeros(NRepl,1); for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path newVanillaPayoff(i) = max(0,path(NSteps+1)-K); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier newPayoff(i) = max(0,path(NSteps+1)-K); end end newVanillaPayoff = exp(-r*T)*newVanillaPayoff; newPayoff = exp(-r*T)*newPayoff; CVpayoff = newPayoff + c*(newVanillaPayoff - muVanilla); [price, std, CI] = normfit(CVpayoff);

function [price, std, CI] = CVknockout_call(spot,Sb,K,r,T,sigma,NSteps,NRepl,NPilot) payoff = zeros(NPilot,1); vanillaPayoff = zeros(NPilot,1); muVanilla = blsprice(spot,K,r,T,sigma); for i=1:NPilot path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path vanillaPayoff(i) = max(0,path(NSteps+1)-K); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 0 payoff(i) = max(0,path(NSteps+1)-K); end end vanillaPayoff = exp(-r*T)*vanillaPayoff; payoff = exp(-r*T)*payoff; covMat = cov(vanillaPayoff, payoff); varVanilla = var(vanillaPayoff); c = -covMat(1,2)/varVanilla; newPayoff = zeros(NRepl,1); newVanillaPayoff = zeros(NRepl,1); for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path newVanillaPayoff(i) = max(0,path(NSteps+1)-K); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 0 newPayoff(i) = max(0,path(NSteps+1)-K); end end newVanillaPayoff = exp(-r*T)*newVanillaPayoff; newPayoff = exp(-r*T)*newPayoff; CVpayoff = newPayoff + c*(newVanillaPayoff - muVanilla); [price, std, CI] = normfit(CVpayoff);

function [price, std, CI] = CVknockin_put(spot,Sb,K,r,T,sigma,NSteps,NRepl,NPilot) payoff = zeros(NPilot,1); vanillaPayoff = zeros(NPilot,1); [aux, muVanilla] = blsprice(spot,K,r,T,sigma); for i=1:NPilot path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path vanillaPayoff(i) = max(0,K-path(NSteps+1)); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier payoff(i) = max(0,K-path(NSteps+1)); end end vanillaPayoff = exp(-r*T)*vanillaPayoff; payoff = exp(-r*T)*payoff; covMat = cov(vanillaPayoff, payoff); varVanilla = var(vanillaPayoff); c = -covMat(1,2)/varVanilla; newPayoff = zeros(NRepl,1); newVanillaPayoff = zeros(NRepl,1); for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); % generate one path newVanillaPayoff(i) = max(0,K-path(NSteps+1)); knocked = barrierCrossing(spot,Sb,path); % determine if up or down if knocked == 1 % path hit the barrier newPayoff(i) = max(0,K-path(NSteps+1)); end end newVanillaPayoff = exp(-r*T)*newVanillaPayoff; newPayoff = exp(-r*T)*newPayoff; CVpayoff = newPayoff + c*(newVanillaPayoff - muVanilla); [price, std, CI] = normfit(CVpayoff);

function [price, std, CI] = CVknockout_put(spot,Sb,K,r,T,sigma,NSteps,NRepl,NPilot) payoff = zeros(NPilot,1); vanillaPayoff = zeros(NPilot,1); [aux, muVanilla] = blsprice(spot,K,r,T,sigma); for i=1:NPilot path = AssetPaths(spot,r,sigma,T,NSteps,1); vanillaPayoff(i) = max(0,K-path(NSteps+1)); knocked = barrierCrossing(spot,Sb,path); if knocked == 0 payoff(i) = max(0,K-path(NSteps+1)); end end vanillaPayoff = exp(-r*T)*vanillaPayoff; payoff = exp(-r*T)*payoff; covMat = cov(vanillaPayoff, payoff); varVanilla = var(vanillaPayoff); c = -covMat(1,2)/varVanilla; newPayoff = zeros(NRepl,1); newVanillaPayoff = zeros(NRepl,1); for i=1:NRepl path = AssetPaths(spot,r,sigma,T,NSteps,1); newVanillaPayoff(i) = max(0,K-path(NSteps+1)); knocked = barrierCrossing(spot,Sb,path); if knocked == 0 newPayoff(i) = max(0,K-path(NSteps+1)); end end newVanillaPayoff = exp(-r*T)*newVanillaPayoff; newPayoff = exp(-r*T)*newPayoff; CVpayoff = newPayoff + c*(newVanillaPayoff - muVanilla); [price, std, CI] = normfit(CVpayoff);

function sim_paths = AssetPaths(spot,r,sigma,T,num_steps,num_sims) sim_paths = zeros(num_sims,num_steps+1); % Each row is a simulated path sim_paths(:,1) = spot; % fill first col w/ spot prices dt = T/num_steps; % time step for i=1:num_sims for j=2:num_steps+1 wt = randn; % generate a random number from a normal distribution sim_paths(i,j) = sim_paths(i,j-1)*exp((r-0.5*sigma^2)... *dt+sigma*sqrt(dt)*wt); end end

function [sim_pathsA, sim_pathsB] = AssetPathsAV(spot,r,sigma,T,num_steps,num_sims) sim_pathsA = zeros(num_sims,num_steps+1); % Each row is a simulated path sim_pathsA(:,1) = spot; % fill first col w/ spot prices sim_pathsB = zeros(num_sims,num_steps+1); % Each row is a simulated path sim_pathsB(:,1) = spot; % fill first col w/ spot prices dt = T/num_steps; % time step for i=1:num_sims for j=2:num_steps+1 wt = randn; % generate a random number from a normal distribution sim_pathsA(i,j) = sim_pathsA(i,j-1)*exp((r-0.5*sigma^2)... *dt+sigma*sqrt(dt)*wt); sim_pathsB(i,j) = sim_pathsB(i,j-1)*exp((r-0.5*sigma^2)... *dt-sigma*sqrt(dt)*wt); end end

barrier_type = 1; filename = 'AVPricing.xls'; spot = 50; sb = 40; strike = 50; r = 0.1; T = 1; sigma = 0.3; num_steps = 250; num_sims = [100 1000 10000 100000]; % % % % % % % %

[c,i,j,k,c,m,l,n,o] i = barriers j = strike k = interest rate (r) l = volatility m = time interval (T) n = num_steps o = num_sims

% [aux sb_length] = size(sb); % get size of row vector % [aux strike_length] = size(strike); % [aux r_length] = size(r); % [aux sigma_length] = size(sigma); % [aux time_length] = size(T); % [aux steps_length] = size(num_steps); [aux sims_length] = size(num_sims);

price_matrix = zeros(1,sims_length); ASol = zeros(1,sims_length); error_matrix = zeros(1,sims_length); sigmahat = zeros(1,sims_length); CI = zeros(2,sims_length); % for i=1:sb_length % for j=1:strike_length % for k=1:r_length % for l=1:sigma_length % for m = 1:time_length switch barrier_type case 1 ASol(:) = A_DICall(spot,sb,strike,r,0,T,sigma); range = 'C13:F19'; case 2 ASol(:) = A_DOCall(spot,sb,strike,r,0,T,sigma); range = 'C23:F29'; case 3 ASol(:) = A_DIPut(spot,sb,strike,r,0,T,sigma); range = 'C33:F39'; case 4 ASol(:) = A_DOPut(spot,sb,strike,r,0,T,sigma); range = 'C43:F49'; case 5 ASol(:) = A_UICall(spot,sb,strike,r,0,T,sigma); range = 'C53:F59'; case 6

ASol(:) range = case 7 ASol(:) range = case 8 ASol(:) range =

= A_UOCall(spot,sb,strike,r,0,T,sigma); 'C63:F69'; = A_UIPut(spot,sb,strike,r,0,T,sigma); 'C73:F79'; = A_UOPut(spot,sb,strike,r,0,T,sigma); 'C83:F89';

end % for n = 1:steps_length % issue here? for i = 1:sims_length switch barrier_type case {1,5} [price_matrix(i), sigmahat(i), CI(:,i)] = ... AVknockin_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {2,6} [price_matrix(i), sigmahat(i), CI(:,i)] = ... AVknockout_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {3,7} [price_matrix(i), sigmahat(i), CI(:,i)] = ... AVknockin_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {4,8} [price_matrix(i), sigmahat(i), CI(:,i)] = ... AVknockout_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); end error_matrix(i) = abs(price_matrix(i)-ASol(i)); end % end % end % end % end % end % end M = [num_sims;ASol;price_matrix;sigmahat;CI;error_matrix]; xlswrite(filename, M, range) disp('Completed successfully!'); % xlswrite(filename, M) writes matrix M to the Excel file filename. % The filename input is a string enclosed in single quotes. % The input matrix M is an m-by-n numeric, character, or cell array, % where m < 65536 and n < 256. The matrix data is written to the first % worksheet in the file, starting at cell A1. % Note: To specify sheet or range, but not both, call xlswrite with just % three inputs. If the third input is a string that includes a colon % character (e.g., 'D2:H4'), it specifies range.

function knocked = barrierCrossing(spot,sb,path) if sb < spot % barrier is less than spot price, down and in/out knocked = any(path <= sb); else % barrier is greater than spot, up and in/out knocked = any(path >= sb); end

barrier_type = 5; filename = 'CVPricing.xls'; spot = 50; sb = 60; strike = 50; r = 0.1; T = 1; sigma = 0.3; num_steps = 500; num_sims = [100 1000 10000 100000]; % % % % % % % %

[c,i,j,k,c,m,l,n,o] i = barriers j = strike k = interest rate (r) l = volatility m = time interval (T) n = num_steps o = num_sims

% [aux sb_length] = size(sb); % get size of row vector % [aux strike_length] = size(strike); % [aux r_length] = size(r); % [aux sigma_length] = size(sigma); % [aux time_length] = size(T); % [aux steps_length] = size(num_steps); [aux sims_length] = size(num_sims);

price_matrix = zeros(1,sims_length); ASol = zeros(1,sims_length); error_matrix = zeros(1,sims_length); sigmahat = zeros(1,sims_length); CI = zeros(2,sims_length); % for i=1:sb_length % for j=1:strike_length % for k=1:r_length % for l=1:sigma_length % for m = 1:time_length switch barrier_type case 1 ASol(:) = A_DICall(spot,sb,strike,r,0,T,sigma); range = 'C13:F19'; case 2 ASol(:) = A_DOCall(spot,sb,strike,r,0,T,sigma); range = 'C23:F29'; case 3 ASol(:) = A_DIPut(spot,sb,strike,r,0,T,sigma); range = 'C33:F39'; case 4 ASol(:) = A_DOPut(spot,sb,strike,r,0,T,sigma); range = 'C43:F49'; case 5 ASol(:) = A_UICall(spot,sb,strike,r,0,T,sigma); range = 'C53:F59'; case 6

ASol(:) range = case 7 ASol(:) range = case 8 ASol(:) range =

= A_UOCall(spot,sb,strike,r,0,T,sigma); 'C63:F69'; = A_UIPut(spot,sb,strike,r,0,T,sigma); 'C73:F79'; = A_UOPut(spot,sb,strike,r,0,T,sigma); 'C83:F89';

end % for n = 1:steps_length % issue here? for i = 1:sims_length switch barrier_type case {1,5} [price_matrix(i), sigmahat(i), CI(:,i)] = ... CVknockin_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i),5000); case {2,6} [price_matrix(i), sigmahat(i), CI(:,i)] = ... CVknockout_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i),5000); case {3,7} [price_matrix(i), sigmahat(i), CI(:,i)] = ... CVknockin_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i),5000); case {4,8} [price_matrix(i), sigmahat(i), CI(:,i)] = ... CVknockout_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i),5000); end error_matrix(i) = abs(price_matrix(i)-ASol(i)); end % end % end % end % end % end % end M = [num_sims;ASol;price_matrix;sigmahat;CI;error_matrix]; xlswrite(filename, M, range) disp('Completed successfully!'); % xlswrite(filename, M) writes matrix M to the Excel file filename. % The filename input is a string enclosed in single quotes. % The input matrix M is an m-by-n numeric, character, or cell array, % where m < 65536 and n < 256. The matrix data is written to the first % worksheet in the file, starting at cell A1. % Note: To specify sheet or range, but not both, call xlswrite with just % three inputs. If the third input is a string that includes a colon % character (e.g., 'D2:H4'), it specifies range.

barrier_type = 4; filename = 'DICall.xls'; spot = 50; sb = 40; strike = 50; r = 0.1; T = 1; sigma = 0.3; num_steps = 500; num_sims = [100 1000 10000 100000]; % % % % % % % %

[c,i,j,k,c,m,l,n,o] i = barriers j = strike k = interest rate (r) l = volatility m = time interval (T) n = num_steps o = num_sims

% [aux sb_length] = size(sb); % get size of row vector % [aux strike_length] = size(strike); % [aux r_length] = size(r); % [aux sigma_length] = size(sigma); % [aux time_length] = size(T); % [aux steps_length] = size(num_steps); [aux sims_length] = size(num_sims);

price_matrix = zeros(1,sims_length); ASol = zeros(1,sims_length); error_matrix = zeros(1,sims_length); sigmahat = zeros(1,sims_length); CI = zeros(2,sims_length); % for i=1:sb_length % for j=1:strike_length % for k=1:r_length % for l=1:sigma_length % for m = 1:time_length switch barrier_type case 1 ASol(:) = A_DICall(spot,sb,strike,r,0,T,sigma); range = 'C13:F19'; case 2 ASol(:) = A_DOCall(spot,sb,strike,r,0,T,sigma); range = 'C23:F29'; case 3 ASol(:) = A_DIPut(spot,sb,strike,r,0,T,sigma); range = 'C33:F39'; case 4 ASol(:) = A_DOPut(spot,sb,strike,r,0,T,sigma); range = 'C43:F49'; case 5 ASol(:) = A_UICall(spot,sb,strike,r,0,T,sigma); range = 'C53:F59'; case 6

ASol(:) range = case 7 ASol(:) range = case 8 ASol(:) range =

= A_UOCall(spot,sb,strike,r,0,T,sigma); 'C63:F69'; = A_UIPut(spot,sb,strike,r,0,T,sigma); 'C73:F79'; = A_UOPut(spot,sb,strike,r,0,T,sigma); 'C83:F89';

end % for n = 1:steps_length % issue here? for i = 1:sims_length switch barrier_type case {1,5} [price_matrix(i), sigmahat(i), CI(:,i)] = ... knockin_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {2,6} [price_matrix(i), sigmahat(i), CI(:,i)] = ... knockout_call(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {3,7} [price_matrix(i), sigmahat(i), CI(:,i)] = ... knockin_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); case {4,8} [price_matrix(i), sigmahat(i), CI(:,i)] = ... knockout_put(spot,sb,strike,r,T,sigma,num_steps,num_sims(i)); end error_matrix(i) = abs(price_matrix(i)-ASol(i)); end % end % end % end % end % end % end M = [num_sims;ASol;price_matrix;sigmahat;CI;error_matrix]; xlswrite(filename, M, range) disp('Completed successfully!'); % xlswrite(filename, M) writes matrix M to the Excel file filename. % The filename input is a string enclosed in single quotes. % The input matrix M is an m-by-n numeric, character, or cell array, % where m < 65536 and n < 256. The matrix data is written to the first % worksheet in the file, starting at cell A1. % Note: To specify sheet or range, but not both, call xlswrite with just % three inputs. If the third input is a string that includes a colon % character (e.g., 'D2:H4'), it specifies range.

function [price_matrix,tElapsed,error_matrix,ASol] = ... graphingFunction(spot,sb,strike,r,T,sigma,... num_steps,num_sims,barrier_type,filename) % passing in vector of possible values for ranges of variables % sb = 0:10:100; % vector of barriers % % % % % % % %

Result: Asset Price || Error changing barriers (sb) changing strike price (K) changing num_steps changing num_sims changing volatility (sigma) changing interest rates (r) changing time interval (T)

% % % % % % % %

[c,i,j,k,c,m,l,n,o] i = barriers j = strike k = interest rate (r) l = volatility m = time interval (T) n = num_steps o = num_sims

[aux [aux [aux [aux [aux [aux [aux

sb_length] = size(sb); % get size of row vector strike_length] = size(strike); r_length] = size(r); sigma_length] = size(sigma); time_length] = size(T); steps_length] = size(num_steps); sims_length] = size(num_sims);

price_matrix = zeros([sb_length strike_length r_length sigma_length ... time_length steps_length sims_length]); tElapsed = zeros([sb_length strike_length r_length sigma_length ... time_length steps_length sims_length]); error_matrix = zeros([sb_length strike_length r_length sigma_length ... time_length steps_length sims_length]); ASol = zeros([sb_length strike_length r_length sigma_length ... time_length]); for i=1:sb_length for j=1:strike_length for k=1:r_length for l=1:sigma_length for m = 1:time_length switch barrier_type case 1 ASol(i,j,k,l,m) = A_DICall(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 2 ASol(i,j,k,l,m) = A_DOCall(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 3 ASol(i,j,k,l,m) = A_DIPut(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 4 ASol(i,j,k,l,m) = A_DOPut(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 5

ASol(i,j,k,l,m) = A_UICall(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 6 ASol(i,j,k,l,m) = A_UOCall(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 7 ASol(i,j,k,l,m) = A_UIPut(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); case 8 ASol(i,j,k,l,m) = A_UOPut(spot,sb(i),strike(j),r(k),0,T(m),sigma(l)); end for n = 1:steps_length % issue here? for o = 1:sims_length switch barrier_type case {1,5} tic; [price_matrix(i,j,k,l,m,n,o), aux, aux] = ... knockin_call(spot,sb(i),strike(j),r(k),T(m),sigma(l),... num_steps(n),num_sims(o)); tElapsed(i,j,k,l,m,n,o) = toc; case {2,6} tic; [price_matrix(i,j,k,l,m,n,o), aux, aux] = ... knockout_call(spot,sb(i),strike(j),r(k),T(m),sigma(l),... num_steps(n),num_sims(o)); tElapsed(i,j,k,l,m,n,o) = toc; case {3,7} tic; [price_matrix(i,j,k,l,m,n,o), aux, aux] = ... knockin_put(spot,sb(i),strike(j),r(k),T(m),sigma(l),... num_steps(n),num_sims(o)); tElapsed(i,j,k,l,m,n,o) = toc; case {4,8} tic; [price_matrix(i,j,k,l,m,n,o), aux, aux] = ... knockout_put(spot,sb(i),strike(j),r(k),T(m),sigma(l),... num_steps(n),num_sims(o)); tElapsed(i,j,k,l,m,n,o) = toc; end error_matrix(i,j,k,l,m,n,o) = abs(price_matrix(i,j,k,l,m,n,o)-... ASol(i,j,k,l,m)); end end end end end end end save(filename, 'price_matrix', 'tElapsed','error_matrix','ASol','spot',... 'sb','strike','r','T','sigma','num_steps','num_sims','barrier_type'); % tic; any_statements; tElapsed = toc; makes the same time measurement, % but assigns the elapsed time output to a variable, tElapsed. MATLAB does % not display the elapsed time unless you omit the terminating semicolon.

% % % % % % %

- changing barriers (sb) - changing strike price (K) changing num_steps changing num_sims - changing volatility (sigma) - changing interest rates (r) - changing time interval (T)

% % % % % % % % % % % % % % % % % % % % % % %

Down and In Call barrier_type = 1; filename = 'DICall.mat'; spot = 50; % sb = 0:1:45; sb = 40; % strike = 0:1:50; strike = 50; % r = 0:.001:.5; r = .1; % T = 0:.1:10; T = 1; % sigma = 0:.01:.5; sigma = .4; % num_steps = [100 1000 10000 100000]; num_steps = 1000; num_sims = 0:5000:100000; % num_sims = 10000;

% % % % % % % % % % %

prices vs. analytical prices as num_sims change as num_steps varies i = sb_length j = strike_length k = r_length l = sigma_length m = time_length n = steps_length o = sims_length [i j k l m n o] = size(price_matrix);

[price_matrix,tElapsed,error_matrix,ASol] = ... graphingFunction(spot,sb,strike,r,T,sigma,... num_steps,num_sims,barrier_type,filename);

% prices = zeros(1,o); % for c=1:o % prices(i) = price_matrix(1,1,1,1,1,1,i); % end %======== Changing Strike =========== %======= MC vs Analytical ========== % prices = price_matrix(1,:,1,1,1,1,1); % plot(strike,prices); % hold on; % aprices = ASol(1,:,1,1,1); % plot(strike,aprices);

% hold off; %======== Changing Barrier =========== %======= MC vs Analytical ========== % prices = price_matrix(:,30,1,1,1,1,1); % plot(sb,prices); % hold on; % aprices = ASol(:,30,1,1,1); % plot(sb,aprices); % hold off; %======== Changing Interest =========== %======= MC vs Analytical ========== % prices = zeros(1,k); % aprices = zeros(1,k); % for c=1:k % prices(c) = price_matrix(1,1,c,1,1,1,1); % end % plot(r,prices); % hold on; % for c=1:k % aprices(c) = ASol(1,1,c,1,1); % end % plot(r,aprices,'-r'); % hold off; %======== Changing Time =========== %======= MC vs Analytical ========== % prices = zeros(1,m); % aprices = zeros(1,m); % error = zeros(1,m); % for c=1:m % prices(c) = price_matrix(1,1,1,1,c,1,1); % end % plot(T,prices); % hold on; % for c=1:m % aprices(c) = ASol(1,1,1,1,c); % end % plot(T,aprices,'-r'); % hold off; % for c=1:m % error(c) = error_matrix(1,1,1,1,c,1,1)/prices(c); % end % plot(T,error,'-r'); % for c=1:l % error(c) = error_matrix(1,1,1,c,1,1,1)/prices(c); % end % plot(T,error,'-r'); %======= Changing Volatitility ======= %======== MC vs Analytical ========= % prices = zeros(1,l); % aprices = zeros(1,l); % error = zeros(1,l); % for c=1:l

% % % % % % % % % % % % % %

prices(c) = price_matrix(1,1,1,c,1,1,1); end prices plot(sigma,prices); hold on; for c=1:l aprices(c) = ASol(1,1,1,c,1); end plot(sigma,aprices,'-r'); hold off; for c=1:l error(c) = error_matrix(1,1,1,c,1,1,1)/prices(c); end plot(sigma,error,'-r');

%======== Changing Sims =========== %======= MC vs Analytical ========== % prices = zeros(1,o); % aprices = zeros(1,o); % error = zeros(1,o); % time = zeros(1,o); % for c=1:o % prices(c) = price_matrix(1,1,1,1,1,1,c); % end % plot(num_sims,prices); % hold on; % for c=1:o % aprices(c) = ASol(1,1,1,1,1); % end % plot(num_sims,aprices,'-r'); % hold off; % for c=1:l % time(c) = tElapsed(1,1,1,c,1,1,1); % end % plot(num_sims,tElapsed(c),'-r');

% Up and Out Call barrier_type = 5; filename = 'UOCall.mat'; spot = 50; % sb = 50:.5:75; sb = 60; % strike = 0:1:50; strike = 50; % r = 0:.001:.5; r = .1; T = 0:.1:10; % T = 1; sigma = .1; % sigma = .05:.01:.3; % num_steps = [100 1000 10000 100000]; num_steps = 1:5:500; num_sims = 1000; % num_sims = 10000;

[price_matrix,tElapsed,error_matrix,ASol] = ... graphingFunction(spot,sb,strike,r,T,sigma,... num_steps,num_sims,barrier_type,filename); [i j k l m n o] = size(price_matrix); %======== Changing Barrier and Volatility =========== %======= MC vs Analytical ========== % % % % % % % % % % %

[X,Y]=meshgrid(sigma,sb); [m n] = size(X); prices = zeros(m,n); aprices = zeros(m,n); error = zeros(m,n); for c=1:m for d=1:n prices(c,d) = price_matrix(c,1,1,d,1,1,1); end end surf(X,Y,prices);

% % % % % % % % % % %

[X,Y]=meshgrid(sigma,sb); [m n] = size(X); prices = zeros(m,n); aprices = zeros(m,n); error = zeros(m,n); for c=1:m for d=1:n error(c,d) = error_matrix(c,1,1,d,1,1,1); end end surf(X,Y,error);

[X,Y]=meshgrid(T,num_steps); [m n] = size(X); prices = zeros(m,n); aprices = zeros(m,n); error = zeros(m,n); for c=1:m for d=1:n error(c,d) = error_matrix(1,1,1,1,d,c,1); end end surf(X,Y,error); xlabel('Time Interval'); ylabel('Number of Steps'); zlabel('Absolute Error');

% % % % % %

prices = zeros(1,l); aprices = zeros(1,l); error = zeros(1,l); for c=1:l prices(c) = price_matrix(1,1,1,c,1,1,1); end

% % % % % % % %

prices; plot(sigma,prices); hold on; for c=1:l aprices(c) = ASol(1,1,1,c,1); end plot(sigma,aprices,'-r'); hold off;

% % % % % % % % % % % % % %

prices = zeros(1,i); aprices = zeros(1,i); error = zeros(1,i); for c=1:i prices(c) = price_matrix(c,1,1,1,1,1,1); end prices; plot(sb,prices); hold on; for c=1:i aprices(c) = ASol(c,1,1,1,1); end plot(sb,aprices,'-r'); hold off;

% for c=1:l % error(c) = error_matrix(1,1,1,c,1,1,1)/prices(c); % end % plot(sigma,error,'-r'); % % % % % % %

plot(sb,prices); hold on; for c=1:i aprices(c) = ASol(c,1,1,1,1); end plot(sb,aprices,'-r'); hold off;

% % % % % % %

Down and Out Call Down and Out Put Down and In Put Up and In Call Up and Out Call Up and Out Put Up and In Put

% [price_matrix,tElapsed,error_matrix] = ... % graphingFunction(spot,sb,strike,r,delta,T,sigma,... % num_steps,num_sims,barrier_type,filename)