Trinomial tree option pricing python
|
I've just started picking up Python and have built binomial and trinomial models just to test my understanding, especially about arrays. Show
For terminal stock price array, I have stepped down from InitialStock * u**(iSteps - i) for i going from zero to 2*iSteps+ 1. ie. I start at top-most terminal stock price InitialStock * u**(iSteps) and travel down until I get to InitalStock * u**(-iSteps) and populate the TerminalStock array from 0 to 2*iSteps+1 (ie. last array element is InitialStock * u**(-iSteps). Because I come from the land of VBA, I haven't yet picked up the eloquent (and hence, brief and sometimes difficult-to-read Python style), my loops look like this: I have initialised dblTrinOpt array as dblTrinOpt = np.ndarray( ( 2*iSteps+1, iSteps+1), float) so that it has 2*iSteps price elements at terminal stock price and iStep time steps, where iSteps = OptionTerm / NumberOfSteps. Yes, apology for my clunky python code! (but I have to read it). after that, my option calc is as follows (yes, again the code is inelegant): after that, my Trinomial function passes dblTrinOpt[0][0] as final discounted option price. I wanted to write this in a more readable form so that it became easier to add two components: (a) an American option facility, hence tracking stock prices at each tree node (and so then things like barriers and Bermudans would be an easy addition, and (b) add two additional terminal prices one at each end and so that at time zero, I'd have three stock and option values so that I could compute delta and gamma without needing to recompute the option tree (i.e.. avoid the bump-and-recompute approach). I hope that helps, even though I didn't fix your specific code, but explained the logic. Trinomial trees in options pricingIn the binomial tree, each node leads to two other nodes in the next time step. Similarly in a trinomial tree, each node leads to three other nodes in the next time step. Besides having up and down states, the middle node of the trinomial tree indicates no change in state. When extended over more than two time steps, the trinomial tree can be thought of as a recombining tree, where the middle nodes always retain the same values as the previous time step. Let's consider the Boyle trinomial tree, where the tree is calibrated such that the probability of up, down, and flat movements, Let's create a TrinomialTreeOption class, inheriting from the BinomialTreeOption class. The methods for the TrinomialTreeOption are provided as follows:
def setup_parameters(self): """ Required calculations for the model """ self.u = math.exp(self.sigma*math.sqrt(2.*self.dt)) self.d = 1/self.u self.m = 1 self.qu = ((math.exp((self.r-self.div) * self.dt/2.) - math.exp(-self.sigma * math.sqrt(self.dt/2.))) / (math.exp(self.sigma * math.sqrt(self.dt/2.)) - math.exp(-self.sigma * math.sqrt(self.dt/2.))))**2 self.qd = ((math.exp(self.sigma * math.sqrt(self.dt/2.)) - math.exp((self.r-self.div) ... Trinomial Model: the Auld TriangleA two-jump process for the asset price over each discrete time step was developed in the binomial lattice. Boyle expanded this frame of reference and explored the feasibility of option valuation by allowing for an extra jump in the stochastic process. In keeping with Black Scholes, Boyle examined an asset (S) with a lognormal distribution of returns. Over a small time interval, this distribution can be approximated by a three-point jump process in such a way that the expected return on the asset is the riskless rate, and the variance of the discrete distribution is equal to the variance of the corresponding lognormal distribution. The three point jump process was introduced by Phelim Boyle(1986) as a trinomial tree to price options and the effect has been momentous in the finance literature. Perhaps shamrock mythology or the well-known ballad associated with Brendan Behan inspired the Boyle insight to include a third jump in lattice valuation. His trinomial paper has spawned a huge amount of ground breaking research. In the trinomial model, the asset price S is assumed to jump uS or mS or dS after one time period (dt = T/n), where u > m > d. Joshi (2008) point out that the trinomial model is characterized by the following five parameters: (1) the probability of an up move pu, (2) the probability of an down move pd, (3) the multiplier on the stock price for an up move u, (4) the multiplier on the stock price for a middle move m, (5) the multiplier on the stock price for a down move d. A recombining tree is computationally more efficient so we require ud = m*m M = exp (r∆t), V = exp (σ 2∆t), dt or ∆t = T/N where where N is the total number of steps of a trinomial tree. For a tree to be risk-neutral, the mean and variance across each time steps must be asymptotically correct. Boyle (1986) chose the parameters to be: m = 1, u = exp(λσ√ ∆t), d = 1/u pu =( md − M(m + d) + (M^2)*V )/ (u − d)(u − m) , pd =( um − M(u + m) + (M^2)*V )/ (u − d)(m − d) Boyle suggested that the choice of value for λ should exceed 1 and the best results were obtained when λ is approximately 1.20. One approach to constructing trinomial trees is to develop two steps of a binomial in combination as a single step of a trinomial tree. This can be engineered with many binomials CRR(1979), JR(1979) and Tian (1993) where the volatility is constant. For example, a two-step presentation of the CRR(1979) is presented below in C++: b = r - q; dt = T/n; u = exp(v*sqrt(2.0*dt)); d = 1.0/u; pu = pow((exp(0.5*b*dt) - exp(-v*sqrt(0.5*dt))) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); pd = pow((exp(v*sqrt(0.5*dt)) - exp(0.5*b*dt)) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); pm = 1.0 - pu - pd; or in R using Google Colab # Time Interval: dt = Time/n # Up-and-down jump sizes: u = exp(+sigma * sqrt(2*dt)) d = exp(-sigma * sqrt(2*dt)) # Probabilities of going up and down: pu = ((exp(b * dt/2) - exp( -sigma * sqrt(dt/2))) / (exp(sigma * sqrt(dt/2)) - exp(-sigma * sqrt(dt/2)))) ^ 2 pd = (( exp(sigma * sqrt(dt/2)) - exp( b * dt/2)) / (exp(sigma * sqrt(dt/2)) - exp(-sigma * sqrt(dt/2)))) ^ 2 # Probability of staying at the same asset price level: pm = 1 - pu - pd We initially present below two snippets of C++ code which are investigated for efficiency. Given the numerically intensive nature of lattices we propose making use of dynamic memory techniques in the lattice estimation. We compare the performance of a static vs dynamic memory set up in pricing models not unlike the approach adapted for ESOs in the previous page and consistent with techniques previously outlined for the binomial estimation where memory usage was optimized. The static code was extracted from Volopta.com developed by Fabrice Rouah. The Dynamic representation was transposed from from Espen Haug's Excel VBA trinomial tree. Espen's dynamic tree design appears to be more efficient and can reach higher step size on the online C++ compilers. We also present how to manually set up the trinomial tree in an excel workbook so readers/viewers can compare the basic visual representation against the actual C++ code presented below. Static Trinomial Tree in C++// Static Trinomial Tree C++ // by Fabrice Douglas Rouah, FRouah.com and Volopta.com //#include "StdAfx.h" #include #include #include #include #include #include #include #include using namespace std; double Trinomial(double Spot,double K,double r,double q,double v,double T,char PutCall,char EuroAmer,int n) { int i,j; double b = r - q; double dt = T/n; double u = exp(v*sqrt(2.0*dt)); double d = 1.0/u; double pu = pow((exp(0.5*b*dt) - exp(-v*sqrt(0.5*dt))) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); double pd = pow((exp(v*sqrt(0.5*dt)) - exp(0.5*b*dt)) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); double pm = 1.0 - pu - pd; // Build the trinomial tree for the stock price evolution vector for (j=0; j<=n; j++) for (i=0; i<=2*j; i++) S[i][j] = Spot*pow(u,double(j-i)); // Not Boost Matrices for the Euro and American calls and puts vector // Compute terminal payoffs for (i=0; i<=2*n; i++) { if (PutCall == 'C') V[i][n] = max(S[i][n] - K, 0.0); else if (PutCall == 'P') V[i][n] = max(K - S[i][n], 0.0); } // Backward recursion through the tree for (j=n-1; j>=0; j--) { for (i=0; i<=2*j; i++) { if (EuroAmer == 'E') V[i][j] = exp(-r*dt)*(pu*V[i][j+1] + pm*V[i+1][j+1] + pd*V[i+2][j+1]); else if (EuroAmer == 'A') { if (PutCall == 'C') V[i][j] = max(S[i][j] - K, exp(-r*dt)*(pu*V[i][j+1] + pm*V[i+1][j+1] + pd*V[i+2][j+1])); else if (PutCall == 'P') V[i][j] = max(K - S[i][j], exp(-r*dt)*(pu*V[i][j+1] + pm*V[i+1][j+1] + pd*V[i+2][j+1])); } } } // Output the results return V[0][0]; } int main() { clock_t start_time, end_time; start_time = clock(); // Stock price, strike, maturity, volatility, risk free rate, dividend yield //double S = 40.0; //double K = 40.0; //double T = 0.25; //double v = 0.20; //double r = 0.05; //double q = 0.0; //char PutCall = 'P'; //int n = 1000; double Spot = 100.0; // Spot Price double K = 100.0; // Strike Price double r = 0.03; // Risk Free Rate double q = 0.07; // Dividend double v = 0.20; // Volatility double T = 3; // Years to maturity char PutCall = 'C'; char EuroAmer = 'E'; int n = 3000; cout << setprecision(10); cout << "The Trinomial Output " << Trinomial(Spot,K,r,q,v,T,PutCall,EuroAmer,n) << endl; end_time = clock(); cout << " " << start_time << " " << end_time << " " << (end_time - start_time) << endl; cout << " " << start_time << " " << end_time << " " << (end_time - start_time) / (double)CLOCKS_PER_SEC << " seconds" << endl; } In this video clip, we further outline a manual Trinomial tree and investigate the robustness of both the static and dynamic C++ code. Using a 9-step manual tree in Excel, we verify that results obtained in C++. Dynamic Trinomial Tree in C++// Dynamic Trinomial Tree // by Fabrice Douglas Rouah, Haug, Shang, Byrne //#include "StdAfx.h" #include #include #include #include #include #include #include #include using namespace std; double Trinomial(double S, double K, double r, double q, double v, double T, char PutCall, char OpStyle, int n) { int i, j, z; double b = r - q; double dt = T / n; double u = exp(v*sqrt(2.0*dt)); double d = 1.0 / u; double pu = pow((exp(0.5*b*dt) - exp(-v*sqrt(0.5*dt))) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); double pd = pow((exp(v*sqrt(0.5*dt)) - exp(0.5*b*dt)) / (exp(v*sqrt(0.5*dt)) - exp(-v*sqrt(0.5*dt))), 2); double pm = 1.0 - pu - pd; if (PutCall == 'C') { z = 1; } else if (PutCall == 'P') { z = -1; } vector //resize the column OptionValue.resize(2*(n+1)); for (i = 0; i <= (2 * n); i++) { OptionValue[i] = max(z*(S*pow(u, (i - n)) - K), 0.0); //This code is like Haug but different because I don't use double max; } // Backward recursion through the tree for (j = n - 1; j >= 0; j--) for (i = 0; i <= (2 * j); i++) { if (OpStyle == 'E') OptionValue[i] = exp(-r*dt)*((pu*OptionValue[i + 2]) + (pm*OptionValue[i + 1]) + (pd*OptionValue[i])); else { OptionValue[i] = max(z*(S*pow(u, (i - j)) - K), exp(-r*dt)*((pu*OptionValue[i + 2]) + (pm*OptionValue[i + 1]) + (pd*OptionValue[i]))); //OptionValue[i] = max(z*(S*pow(u, (2 * i - j)) - K), exp(-r*dt)*(p*(OptionValue[i + 1]) + (1.0 - p)*(OptionValue[i]))); } } // Return the option price return OptionValue[0]; } int main() { clock_t start_time, end_time; start_time = clock(); double S = 100.0; // Spot Price double K = 100.0; // Strike Price double r = 0.03; // Risk Free Rate double q = 0.07; // Dividend double v = 0.20; // Volatility double T = 3; // Years to maturity char PutCall = 'C'; char OpStyle = 'E'; int n = 3000; cout << setprecision(10); cout << "The Trinomial Output " << Trinomial(S, K, r, q, v, T, PutCall, OpStyle, n) << endl; end_time = clock(); cout << " " << start_time << " " << end_time << " " << (end_time - start_time) << endl; cout << " " << start_time << " " << end_time << " " << (end_time - start_time) / (double)CLOCKS_PER_SEC << " seconds" << endl; //system("PAUSE"); } In the video below, we work out the value of an American Put option using a 9 step manual tree. Again the results are consistent with Espen Haug's vba code and the C++ code obtained from Fabrice Rouah. Espen's Excel VBA code was translated to C++ with some modifications. The Put Option is developed by using the previous worksheet on the American Call option and modifying simply the intrinsic estimation at each node. Estimation speed using Colab for Static and Dynamic C++ CodeBelow we implement the C++ code in Google Colab. A key advantage of using Colab relates to the ability to share work and collaborate. The processing power provided by google is also quite generous and access to estimation is through your browser. Access to same can be achieved via tablet, smart phone and computer - expanding opportunities to run algorithms. Estimation of the Trinomial Option Pricing model in Google Colab using Static and Dynamic C++ design and R CodeBelow we explore more fully Google Colab, On the relationship between Binomial and Trinomial Option Pricing ModelsRubinstein (2000) p. 3 noted that "Several years ago, when I studied the Brennan/Schwartz paper carefully, it was again natural to ask the question of how binomial trees were related to the trinomial tree produced by explicit finite differences. I had a sense that if the world worked elegantly, a binomial tree ought to be a special case of a trinomial tree. Trinomial trees ought to be composed of binomial trees, somewhat like molecules are composed of atoms. That is, one could hope that the trinomial tree could be parameterized so that it could be interpreted more fundamentally as a binomial tree in which every other period were simply skipped." For an interesting overview and comprehensive discussion of links between trinomial and binomial models please see Haahtela (2010). In the video and colab below we set out estimations that show how to set up the trinomial as a special case of the binomial. Execute Trinomial Model using Matlab code in Google ColabInstall Octave and Execute Matlab code in Google Colab. In the video below, we provide a simple example of how to add two numbers using matlab code in Google Colab. Then we take a Black Scholes example that makes use of matlab code and estimate the value of a European Call Option. Then we make use of the Trinomial Code made available by Mathworks. The major benefits of matlab/mathworks product relate to the MATLAB mathematical function library. This is an expansive collection of computational algorithms ranging from elementary functions like sum, sine, cosine, and complex arithmetic, to more sophisticated functions like matrix inverse, matrix eigenvalues, Bessel functions, and fast Fourier transforms. BERT Excel add in for implementing a R user defined Function in an Excel SpreadsheetThe Trinomial model is implemented here in R with Google Colab and then re-purposed as a user-defined function in excel. We implement R code obtained from R Forge RMetrics and execute in RStudio. We then save R Script file to the BERT2 folder which excel can read. The R function then can be called into the spreadsheet we verify estimation using a manual trinomial tree described above. |
