1. Introduction

Pricing of Electricity commodity is very well known to exhibit mean reversion characteristics. This is because energy prices are often driven by supply and demand. Prices of electricity fluctuated about the equilibrium, therefore, electricity pricing models must have some mean reverting property to capture the mean reverting behaviour of electricity prices (^{[3, 4, 16, 17]}. The demand, supply and equilibrium price of electricity are also associated with spikes (^{[5, 6, 18]}).

There are many models in the literature for pricing electricity from market data, example are the Ornstein – Uhlenbeck, jump diffusion, Box-Cox models and other time series models like Garch and Arima for calibrating models (^{[4, 17, 28]}). Traditional financial models start with the Black-Scholes model with the assumption that prices are log-normal or obeys Geometric Brownian Motion (^{[24, 30]}, and ^[31]]). This assumption does not make sense in the context of electricity pricing for many reasons including the non-predictability of electricity prices (^[16]). A model which has been used in practice is known as the Ornstein-Uhlenbeck process (OUP). This is a continuous time model which permits autocorrelation. It is necessary to incorporate mean reversion when modelling electricity prices because of the electricity price jumps due to unexpected events. Hence OUP gives information on modelling of electricity prices.

When electric power is produced, there is the need for transportation from the power stations to where it is consumed. Basically, this is done through power lines and these power lines constitute the power grid. The main element of the electrical transmission and distribution grid is the power line which enables power flow from one point to the other in the grid ^[26]. When current is flowing through the metallic conductors, electric and magnetic fields are created and the electric power is carried by these fields. In power systems, the primary task of a line is to secure the transportation of power in the grid. Depending on the length of a power line and the required power capacity, an economic optimization results in a certain voltage level. Higher voltages lower the transmission losses but leads to increased costs of the transmission. Transmission of electricity over long distances results in power loss in the form of heat (^[19]). This can be explained by the telegraph equation and the Cobb-Douglas equation.

Electrical power losses are economically important because they constitute the electrical power that needs to be produced by generators along with electric power delivered to customers not minding if they are dissipated as heat in transmission and distribution systems. Since they grow with the square of electric current transmission, voltage is usually increased so as to reduce the current and losses to economic levels when transmitting large quantities of electric power (^[26]).

There are so many literatures from researches all over the world on the efficient transmission and distribution of electricity, electricity pricing models and demand and supply of electricity. Julio and Eduardo ^[27] examined the importance of the regular patterns in the behaviour of electricity prices and its implications for the purpose of derivative pricing. They analyzed the Nordic power exchange’s spot, future and forward prices.

Montero et. al. ^[27] looked at the Spanish electricity market which has suddenly gained an increase in its complexity behaviour caused by deregulation. They looked at asymmetrical patterns of the volatility of the Spanish electricity spot prices with much attention towards the direct or inverse leverage effect. Avishka ^[7] looked at the increased uncertainty in electricity prices (spikes) as a result of the deregulation of electricity markets. He used the Barlow ^[8] and Aid et. al. ^[2] models to price electricity and discovered that the Barlow model can be viewed as a continuous version of the Aid et al. ^[2].

Hogan and Brendan ^[21] examined alternative pricing frameworks for dealing with non-convexities and other complications in models used in day-ahead electricity markets (see ^[20]). The work by O’Neil et. al. ^[29] provided the building block into emphasizing an application for pricing and settlement systems needed for electricity markets. Buzoianu et. al. ^[12] introduced a new model for electricity prices based on the principle of supply and demand equilibrium. The model was applied to study the Californian wholesale electricity prices over a three year period. Oke and Bamgbola ^[26] developed a mathematical model of losses along electric power transmission lines using a combination of Ohmic and Corona losses. Lavaei and Low ^[24] studied how minimum power loss in a power system is related to its network topology. Existing algorithms in the literature all exploit non-linear, heuristic or local search algorithms to find the minimum power loss which makes them blind to the network topology.

Oyelami and Adedoyin ^[17] used The Harvey logistic model to predict the demand and supply of electricity in the Nigeria from 2005 to 2026. Obtained results on the demand and supply of electricity showed that electricity has mixture of spikes because of mean reverse behaviour about the mean values and the demand outweighed the supply. Oyelami ^[16] carried out research on models for pricing electricity from market data using the Ornstein – Uhlenbeck process and other time series models like Garch and Arima to calibrate model ^[23].

In year 2013 there was total deregulation in the power sector in Nigeria, therefore, there are now three stages involved in the business of electricity in Nigeria; these stages are the generation, transmission and distribution of electricity (see ^[17]). The basic problems facing the power sector are how to use scientific method to compute fair price for selling of electricity to consumers and how the electricity transmission companies and distribution companies will hedge out the cost of electricity from the electricity generating companies considering the market volatility of demand and supply of electricity in the country as reported in (^[17]). There may also be some other problems that have to do with technology associated with the generation and distribution of electricity (^{[1, 4, 6, 12, 15]}).

In this paper, we intend to develop simulation for computing the fair price for electric tariff for electric consumption for domestic users in Nigeria. We will make of plain vanilla European option pricing mechanism to hedge out the costs of generation and distribution electricity from power stations and use Ornstein-Uhlenbeck model with the aid of Monte Carlo simulation to obtain the fair electricity price for domestic consumers. Price sensitivity analyses are made by the use of Greeks in the Matlab environment.

2. Preliminaries and Notation

2.1. Brownian Motion ([30] & [31])

A scalar standard Brownian motion or standard Wiener process over [0, T] is a random variable) that depends continuously on and satisfies the following conditions:

1. w(0)=0 (with probability 1)

2. For the random variable given by the increments is normally distributed with mean zero and variance t-s, Equivalently, is the normal distribution with zero mean and unit variance.

3. For the increments w (t)-w(s) and W (u) – w (v) is independent.

3. Statement of Problem

The prime motivation for the paper is to develop simulation for pricing electricity in Nigeria. Moreover, formulate the strategy on how to hedge out the transmission cost and distribution cost of electricity to consumers. We will use models with mean reverse property devoid of complexities to achieve our goal.

4. Methods

4.1 Ornstein-Uhlenbeck Model

Consider stochastic differential equation given as the Ornstein-Uhlenbeck model as follows:

(1)

Where

is the long run equilibrium price or mean reversion level.

is the volatility

is the mean reversion rate

is the standard Brownian motion

is the spot price

The above model is known for its mean-reversion property in the sense that the process is always pulled back towards a long term level of which stands for the mean reversion rate or the speed of adjustment (^{[13, 31]}). The mean reversion parameter controls the size of the expected adjustment toward the long term level. The drift component is governed by the distance between the spot price and the mean reversion level as well as by the mean reversion rate. As long as the spot price stays above the mean reversion level, the mean reversion component will be negative resulting in a decrease of the spot price. If we have the spot price below the mean reversion level, the mean reversion component will be positive thus resulting into an increment for the spot price. Therefore, the price will always be pulled towards the mean reversion level at a speed determined by the mean reversion rate. We will find analytic and numerical solution for the above model.

4.2. Discretization of the Ornstein-Uhlenbeck Model Using Monte Carlo method

Let us consider the interval and the partition such that . If the intervals are evenly spaced such that we consider a Brownian process which is a real valued stochastic process such that are independent increments. We will to simulate, let be a constant time increment such that . Hence can be written as the sum of increments is follows:

Let Z be a standard normally distributed random number, that is, for each k. we can therefore compute Brownian process as . A forward simulation of over T is given by . Where are independent Gaussian standard variables. Since we have evenly spaced intervals then .

From the properties of stochastic integrals then given any integral of the form

Then the above integral becomes a normal random variable with zero mean and variance .

Calibration can be achieved using the relationship between consecutive observations. The relationship between the linear fit and the model parameter is given as where . We can rewrite as

4.3. Option Pricing of Electricity Commodity

For European Call and Put options driven by the Ornstein-Uhlenbeck process, we express the price of the option as the expected value of the discounted payoff at maturity. Therefore, three types of pricing we will be made as follows:

(1) The fair price for domestic consumption of electricity.

(2) Fair price for buying Call options for electricity commodity by transmission companies from power generating companies.

(3) Fair price for selling Put options on electricity commodity by transmission companies to distribution companies.

Therefore, for European call option

(2)

For the Put option

(3)

Where is the solution to the Ornstein-Uhlenbeck model at maturity period .

Hence for every obtained from the discretized price process, we can then obtain a respective Call Option

(4)

And for Put option

(5)

4.4. Calibration of the Ornstein-Uhlenbeck Model

The price process is a function of normalized demand which follows an Ornstein-Uhlenbeck process (OUP). Let the exact solution to OUP for any .

The mean .

The variance=.

The equation for Simulation of at time is given as

(6)

Where the standard normal distribution with zero mean and unit standard deviation. Hence, for Call we use Matlab to simulate the expression

(7)

For Put option we use

(8)

if is a sequence of random variables each of which has same probability information as , and the series is convergent, then by Kolmogorov criterion a sufficient condition for the strong law of large numbers to be applied to is guaranteed.

5. Results and Discussions

The computation of market volatility for the demand and supply is very useful to our study. We make use of the result we obtained in ^[16] and ^[17]. The Figure 1 and Figure 2 are graphs of electricity demand and supplied in Nigeria in the given period as obtained from the National Bureau of Statistics.

The solution to the equation (1) can be found by employing the Ito’s formula.

Let then the corresponding Ito’s expression for the Ornstein-Uhlenbeck model is given as

(9)

Hence by partial differentiations of give .

Hence,

(10)

Integrating through the equation (10), we get

(11)

Therefore, from our solution of the Ornstein-Uhlenbeck equation, we have

The stochastic integral has zero mean and the variance is

(12)

Hence the volatility is .

Discretizing the Ornstein-Uhlenbeck model,

(13)

Where are independent and identically distributed and.

Consider the case where the OUP model, is purely deterministic, then we have

(14)

Therefore, by the use of the financial tools from Maple 17 software using

(15)

And for for with 100 replications, we have the solution to the deterministic OUP model as

(16)

The graph of the solutions to OUB is shown in the figures below:

We will therefore calibrate the model using the in-built maximum likelihood tool in the Financial toolbox of Matlab 7 for price sensitivity analysis we also employed the demand and supply of electricity data obtained from the Nigerian National Bureau of Statistics (see Oyelami and Adedoyin (2014)).We therefore obtained the following result using the Matlab (see ^[17]):

We choose , because 13.96 naira (about US$0.087) is the current flat rate being charged in Nigeria for domestic consumption of electricity per unit per kilowatt of electricity used. Figure 3 shows the graph of OUM for the deterministic case (i.e. ). If we make use of the equation (6) we can generate the Monte Carlo simulation for the pricing of the electricity at the volatility of by the Maple 17 code in the Appendix II and graph is the Figure 4. In Figure 5 prices of electricity increases with time.

Figure 6 is the graph which shows the price path for the Ornstein –Uhlenbeck model and there appears to be no long run average level about which the series evolves which is an example of a non-stationary time series.

Figure 6 is a graph showing us the price path of the Ornstein –Uhlenbeck process and there appears to be no long run average level about which the series evolves which is an example of a non-stationary time series. We can see the different paths for the different volatilities which is a good explanation for the reason why higher volatilities always come with a higher price since the risk involved is also high. The continuously compounded return series associated with the same price series is shown in Figure 7. The return series appear to be quite stable over time and produce a stationary time series.

Figure 7 also describes the price processes of the Ornstein-Uhlenbeck model at different volatilities ranging from the market volatility obtained from the market at 0.0812 up to 0.7812. Basically, we looked at .

5.1. Option Pricing of the Electricity Commodity

The call and put prices are generated from the equation (7) and equation (8) respectively using the Monte Carlo simulation in the Matlab environment see the Appendix I for the code used and the results are displayed in the Table 1 and Table 2 respectively. From the Figure 8 the call price increases with the time to maturity time while the put price decreases the time to maturity time.

Suppose that the transmission companies are to buy electricity from generating companies at a price fixed ahead of time, there would be the need to set up derivative contracts between the two parties. As shown in Figure 8, we simulated a fair price for a call option using the Monte Carlo simulation and obtained 0.3100 as the simulated fair price. In the same vein, we also simulated a fair price for a put option contract written by the transmission company to the distribution company, we obtained 0.3071 as the fair price.

5.2. Sensitivity Analyses for the Call and Put Options

Figure 9 is the simulation of put price as time approaches maturity for selling of electricity by transmitting companies to distribution companies. The portfolio created is accompanied by very high risk. This risk is as a result price uncertainty. We therefore the risks must be managed by providing effective sensitivity analysis through a proper description of how the Greeks affect the value of the portfolio described above.

From Figure 10, the value of gamma rose around the low electricity price region which shows that around that price region, the delta is highly sensitive to any electricity price change. It would be risky to leave a delta-neutral portfolio unchanged without frequent adjustment. (Adjustment could be in terms of rebalancing the hedge ratio).

Figure 13 show that theta is negative which is in line with the theoretical property of the European call option Greeks. Since the theta is negative, the gamma for this call option must be positive. This theoretical fact supports why the gamma discussed above is positive. The implication of a negative theta is that as the contract expires, option becomes less valuable.

Figure 14 describes how the change in interest rate affects the option value. From the Figure, standing at a point where rho is 20 means that a 1% change in the risk-free interest rate says from 5% to 6% will lead to an increase in the value of the contract by 0.01 x 20 = 0.2.

For the put option, Theta is observed to be positive from Figure 15 which is responsible for the reason why gamma is negative in Figure 12. This is in line with the theoretical background of the Greeks. The implication of this is that as gamma appears to reduce in value (negative), the changes in delta would be very slow. Hence transmission companies would have to worry less about adjusting their portfolio in order to keep it delta-neutral. If they have to adjust, it should be made but not frequently.

From Figure 16, it is noted that Vega stayed strongly positive all through the region. This implies that the portfolio is going to be very sensitive to small changes in volatility. Rho is observed negative in

Figure 17. If for example, we take the point where rho is -10, this means that a 1% change in the risk free interest rate say from 4% to 5% would reduce the value of the option contract by 0.01 x (-10) = -0.2.

6. Conclusion

Determination of fair price for electricity for industry and pricing according to occupational sector, busy hour periods, holidays etc are not considered in this paper. We did not consider the use of exotic options pricing too. However, this does underscore this study in this paper. The Nigerian power sector has just been fully deregulated. This is expected to open up investments in the power sector in Nigeria and followed trading of financial instruments based on the current market situation. This work is more or less introducing use of mean reversing models and opting pricing for trading in the power sector in Nigeria. Our future research would be on complex system modelling of electricity commodity in Nigeria and simulations would be carrying out subject to available data.

7. Acknowledgements

The authors are grateful to the National Mathematical Centre, Abuja and The University of Abuja, Nigeria.

References

[1]	Ackerman (1985) Short-Term Load Prediction for Electric-Utility Control of Generating Units. Comparative models for electrical Load forecasting (pp. 32-42) New York.
	In article
[2]	Aid, R., Campi, L., Nguyenhuu, A. & Touzi, N. (2009) ‘A Structural Risk Neutral Model of Electricity Prices’ International Journal of Theoretical and Applied Finance 12 (7), 925-947
	In article		CrossRef
[3]	Al-Alawi, S.M and Syed, M.Islam (1996) Principles of Electricity Demand Forecasting: Part 1 Methodologies. Power Engineering Journal 10 (3), 139-143.
	In article		CrossRef
[4]	Al-Garni, A.Z., Syed, M.Z. and Javeed, S.N. (1994) A Regression Model for Electric-Energy-consumption forecasting in Eastern Saudi Arabia. Energy, 19 (10), 1043-1049.
	In article		CrossRef
[5]	Amarawickrama, H.A. and Hunt, L.C (2008) Electricity Demand for Sri Lanka: A Time series Analysis.
	In article
[6]	Anton, E., Michael, L., Lefevre, T., O’Leary, D., Peters, R., Svensson, B., and Wilkinson, R. (2000) “Electricity Supply and Demand Side managent Options” Prepared for the World commission on Dams.
	In article
[7]	Avishka Raghoonundun (2010) ‘Models for Electricity Prices’ Msc dissertation submitted to the University of British Colombia (Vancouver).
	In article
[8]	Barlow, M.T. (2002), ‘A diffusion model for electricity prices’ Mathematical Finance 12 (4), 287-298.
	In article		CrossRef
[9]	Barton, B (2003) “Does electricity Market Liberalization Contribute to energy Sustainability?” in Energy Law and Sustainable Development, Ed. By Bradbrook, A.J and Ottinger, R.I., pp. 2117-231 Cambridge: IUCN Publications.
	In article
[10]	Batstone, S. (2000) “Risk Management for Deregulated electricity Market: Simulation Results from a Hydro Management model” in Proceedings of the 33rd Annual Conference, Ed. By Henderson, A. and Orsnz, P.
	In article
[11]	Bodger P.S and Tay, H.S (1987) ‘Logistic and energy substitution Models for Electricity Forecasting: A Comparison Using New Zealand Consumption Data’ Technological Forecasting and Social Change, Vol. 31, pp. 27-48.
	In article		CrossRef
[12]	Buzoianu, M., Seppi, D.J. & Brockwell, E. Anthony (2011) “A Dynamic Supply-Demand Model for Electricity Prices” Technical Report submitted to the department of statistics, Carnegie Mellon University, Pittsburgh, PA.
	In article
[13]	Conelo, A.J., Contreras, R.E. and Plazas, M. (2005) “Forecasting Electricity Prices for a day-ahead pool-based electric energy market” international Journal of Forecasting, 21, pp. 435-462.
	In article		CrossRef
[14]	Dapice, D. “Electricity demand and supply in Myanmar”. Electric Power Research Institute, (1977) “Forecasting and Modelling Time-of-Day and Seasonal Electricity Demands” EPRI EA-578-SR, Palo Alto, CA. 2012.
	In article
[15]	Eydeland, A & Wolyniec, K (2002): Energy and Power Risk Management. New Development in Modelling, Pricing and hedging. John Wiley and Sons, New York, NY.
	In article
[16]	Oyelami Benjamin Oyediran (2013) “Models for pricing of electricity commodity”. Mathematical theory and modelling volume 3 No. 1.
	In article
[17]	Oyelami Benjamin Oyediran, and Adedoyin Adewumi Adedamola (2004) “Models for Forecasting the Demand and Supply of Electricity in Nigeria.” American Journal of Modelling and Optimization, vol. 2, no. 1: 25-33.
	In article		CrossRef
[18]	O’Neill Richard P., Paul M. Sotkiewicz, Benjamin F. Hobbs, Michael H. Rothkopf, William R. Stewart, Jr., “Efficient Market-Clearing Prices in Markets with Nonconvexities,” December 9, 2002, available on the web page of the Harvard Electricity Policy Group. (https://ksgwww.harvard.edu/hepg/).
	In article
[19]	Haque, R. (2006) Transmission loss allocation using artificial neural networks. Master of Science thesis submitted to the Department of Electrical Engineering, University of Saskatchewan, Saskatoon, Saskatchewan.
	In article
[20]	Hogan W W (2003) Transmission Market Design. John F. Kennedy School of Government or Harvard University Working paper series, RWP03-040.
	In article
[21]	Hogan, William & Brendan J. Ring (2003) “On Minimum Uplift Pricing for Electricity Markets” available on the web page of Harvard Kennedy group. (https://www.hks.havard.edu)
	In article
[22]	Hull J. (2002): Options, Futures and Other derivatives 6th ed. Pearson, Prentice Hall, NJ.
	In article
[23]	Kwok Y.K Mathematical Model of Financial Derivatives (1998) Springer Verlag. Berlin.
	In article
[24]	Lavaei, J & Low, S.H. (2010) ‘Relationship Between Power Loss and Network Topology in Power Systems’ Technical report, California Institute of Technology. https://caltechcdstr.library.caltech.edu
	In article
[25]	North America Electric Reliability Council (1990) “Electricity Supply and Demand for 1990-1999”. Princeton New jersey
	In article
[26]	Oke, M. & Bambgola, M. (2013). ‘Minimization of Losses on Electric Power Transmission Line’ Mathematical Theory and Modelling. Vol. 3 (7), 28-31.
	In article
[27]	Julio, J.L &Eduardo, S.S (2002) “Electricity Prices and Power Derivatives: Evidence from the Nordic Power Exchange” Review of Derivative Research Vol. 5, issue1, pp 5-50, Springer.
	In article
[28]	Montero, Jose M. Maria C. Garcia –Centeno, Gema Fernandez-Aviles (2011) Modelling the Volatility of Spanish Wholesale Electricity Spot Market: Asymmetric Garch models Vs Threshold ARSV model, Estudios de Economia Aplicada, vol. 29, no 2, pp 597-616 Association International de Economia Aplicada Espana.
	In article
[29]	Richard P. O’Neill, Paul M. Sotkiewicz, Benjamin F. Hobbs, Michael H. Rothkopf, William R. Stewart, Jr., “Efficient Market-Clearing Prices in Markets with Nonconvexities,” December 9, 2002, available on the web page of the Harvard Electricity Policy Group. (https://ksgwww.harvard.edu/hepg/).
	In article
[30]	Paulo Glasserman Monte Carlo Methods in Financial engineering. Application of Mathematics Stochastic modelling and applied probability. Springer Science, USA 2004.
	In article
[31]	Steven E. Shreve Stochastic Calculus for finance II. Continuous Time Model USA, 2004 Springer Finance
	In article

Appendix I

>>% matlab code for price, call and put option simulation begins here

>> a=0.0948;

>>mu=46.8371;

>>sigma=0.0812;

>>n=7000;

>>alpha=1.100;

>>t=1;

>>epsilon=0.008;

>>h=t/n;

>>tvec= [1: n];

>>x (1) =13.96;

>> for t=2: n

x(t)=miu+((x(t-1)-miu)*exp(-a*h))+sigma*sqrt((1-exp(-2*a*h)))/(2*a)*randn(1);

end

>> for t=1: n

if 1+alpha*x(t)>epsilon

S (t) = (1+alpha*x (t)) ^ (1/alpha);

else

S (t) =epsilon^(1/alpha);

end

end

>> s=mean(S)

s = 13.9629

>> K=13.96;

>> r=0.05;

>> Ct=max (0, S (1: n)-K)*exp (-r*h);

>> Pt=max (0, K-S (1: n))*exp (-r*h);

>> CALL=mean (Ct);

>> PUT=mean (Pt);

>> CALL

CALL= 0.3100

>> PUT

PUT = 0.3071

>> % Matlab code for parameter estimation using maximum likelihood begins here

>>n = length(S)-1;

>>Sx = sum(S (1: end-1));

>>Sy = sum(S (2: end));

>>Sxx = sum(S (1: end-1).^2);

>>Sxy = sum(S (1: end-1).*S (2: end));

>>Syy = sum(S (2: end).^2);

>>mu = (Sy*Sxx - Sx*Sxy) / (n*(Sxx - Sxy) - (Sx^2 - Sx*Sy));

>>lambda = - (1/deltat)*log ((Sxy - mu*Sx - mu*Sy + n*mu^2) / (Sxx -2*mu*Sx + n*mu^2));

>>alpha = exp (- lambda*deltat);

>>alpha2 = exp (-2*lambda*deltat);

>>sigmahat2 = (1/n)*(Syy - 2*alpha*Sxy + alpha2*Sxx -...

2*mu*(1-alpha)*(Sy - alpha*Sx) + n*mu^2*(1-alpha)^2);

>>sigma = sqrt (sigmahat2*2*lambda/ (1-alpha2));

End

>>% matlab code for call option portfolio sensitivity begins here

if ~exist ('normpdf')

Msgbox ('the statistics toolbox is required to run this’, ‘product dependency')

return

end

>> range=10:20;

>> plen=length(range);

>> K= [10.80 11 11.60 11.60 12 12.60 12.80 12.80 13.96 14.70 13.00 16.60];

>> r=0.05*ones (12, 1);

>> time= [36 36 36 27 18 18 18 9 9 9 6 6];

>> sigma=0.0812*ones(12, 1);

>> numopt= [1.1704 1.5605 1.3054 1.9506 1.9144 1.3292 1.2679 1.9506 1.5020 1.2679 3.9012 1.3654];

>> zval=zeros (36, plen);

>> for i=1:12

Pad=ones (time (i), plen);

newr=range (ones (time (i), 1), :);

t= (1: time (i))';

newt=t (:, ones (plen, 1));

zval(36-time(i)+1:36,:)=zval(36-time(i)+1:36,:)+numopt(i)*blsgamma(newr,K(i)*pad,r(i)*pad,newt/36,sigma(i)*pad);

color(36-time(i)+1:36,:)=zval(36-time(i)+1:36,:)+numopt(i)*blsdelta(newr,K(i)*pad,r(i)*pad,newt/36,sigma(i)*pad);

end

>> mesh (range, 1:36, zval, color);

>> View (50, 50)

>> set (gca, 'xdir', 'reverse', 'tag', 'mesh_axes_3');

>> set (gca, 'box', 'on');

>> Colorbar ('horiz');

MATLAB CODE FOR PARAMETER ESTIMATION USING MAXIMUM LIKELIHOOD

>>n = length(S)-1;

>>Sx = sum(S (1: end-1));

>>Sy = sum(S (2: end));

>>Sxx = sum(S (1: end-1).^2 );

>>Sxy = sum(S (1: end-1).*S (2: end));

>>Syy = sum(S (2: end).^2);

>>mu = (Sy*Sxx - Sx*Sxy) / ( n*(Sxx - Sxy) - (Sx^2 - Sx*Sy) );

>>lambda = -(1/deltat)*log((Sxy - mu*Sx - mu*Sy + n*mu^2) / (Sxx -2*mu*Sx + n*mu^2));

>>alpha = exp (- lambda*deltat);

>>alpha2 = exp (-2*lambda*deltat);

>>sigmahat2 = (1/n)*(Syy - 2*alpha*Sxy + alpha2*Sxx -...

2*mu*(1-alpha)*(Sy - alpha*Sx) + n*mu^2*(1-alpha)^2);

>>sigma = sqrt (sigmahat2*2*lambda/ (1-alpha2));

End

MATLAB CODE FOR CALL OPTION RHO SENSITIVITY

if ~exist ('normpdf')

msgbox ('the statistics toolbox is required to run this’, ‘product dependency')

return

end

>> range=10:20;

>> plen=length (range);

>> K= [10.80 11 11.60 11.60 12 12.60 12.80 12.80 13.96 14.70 13.00 16.60];

>> r=0.05*ones (12, 1);

>> time= [36 36 36 27 18 18 18 9 9 9 6 6];

>> Sigma=0.0812*ones (12, 1);

>> numopt= [1.1704 1.5605 1.3054 1.9506 1.9144 1.3292 1.2679 1.9506 1.5020 1.2679 3.9012 1.3654];

>> zval=zeros (36, plen);

>> for i=1:12

pad=ones (time (i), plen);

newr=range (ones (time (i), 1), :);

t= (1: time (i))';

newt=t (: ones (plen, 1));

color(36-time(i)+1:36,:)=zval(36-time(i)+1:36,:)+numopt(i)*blsrho(newr,K(i)*pad,r(i)*pad,newt/36,sigma(i)*pad);

end

>> mesh (range, 1:36, zval, color);

>> view (50, 50)

>> set (gca, 'xdir', 'reverse', 'tag', 'mesh_axes_3');

>> set (gca, 'box', 'on');

>> colorbar ('horiz');

MATLAB CODE FOR CALL OPTION THETA SENSITIVITY

if ~exist ('normpdf')

msgbox ('the statistics toolbox is required to run this’, ‘product dependency')

return

end

>> range=10:20;

>> plen=length (range);

>> K= [10.80 11 11.60 11.60 12 12.60 12.80 12.80 13.96 14.70 13.00 16.60];

>> r=0.05*ones (12, 1);

>> time= [36 36 36 27 18 18 18 9 9 9 6 6];

>> sigma=0.0812*ones (12, 1);

>> numopt= [1.1704 1.5605 1.3054 1.9506 1.9144 1.3292 1.2679 1.9506 1.5020 1.2679 3.9012 1.3654];

>> zval=zeros (36, plen);

>> for i=1:12

pad=ones (time (i), plen);

newr=range (ones (time (i), 1), :);

t= (1: time (i))';

newt=t (: ones (plen, 1));

color(36-time(i)+1:36,:)=zval(36-time(i)+1:36,:)+numopt(i)*blstheta(newr,K(i)*pad,r(i)*pad,newt/36,sigma(i)*pad);

end

>> mesh (range, 1:36, zval, color);

>> view (50, 50)

>> set (gca, 'xdir', 'reverse', 'tag', 'mesh_axes_3');

>> set (gca, 'box', 'on');

>> colorbar ('horiz');

Appendix II

with (Finance):

> r:= OrnsteinUhlenbeckProcess (r₀, μ, θ, σ)

r:= _X

> Drift (r(t))

θ(μ-_X(t))

> Diffusion (r(t))

σ

r₀:= 13.96;

> sigma:= 0.0812;

σ:= 0.0812

> theta:= 0.0948

θ:= 0.0948

> mu:= 46.8371

μ:= 46.8371

> PathPlot (r(t), t=0..5, timesteps = 100, replications = 2, color =red .. blue, axes = BOXED)

q: OrnsteinUhlenbeckProcess (r₀, μ, θ, σ, scheme = unbiased)

q: _X1

>p: dsolve

p:=

>P1:= PathPlot (r(t), t=0..5, timesteps=5, replication=1, color= red)

>P2:= PathPlot (q(t), t=0..5, timesteps=5, replication=1, color= blue)

>P3:= plot (qrhs(t), t=0..5, color= green)

>plots [display] (P1, P2, P3, thickness = 3, axes =BOXELD)