# Mathematica 8 – FinancialDerivative American options issue

Again while testing my routines (Java code to implement the binomial tree to price American options) I came across an issue with the FinancialDerivative command of Mathematica 8 (I am working with the linux version). The default algorithm chosen by Mathematica for American style options seem to converge to a wrong result. To get the correct result  one must manually set the option  “Method” -> “Binomial”. As a simple example of the issue I post here the code to compute the price of a European and American call option with zero dividends. The two prices should be the same in this case. Here what I get instead:

In:= \$Version
Out= "8.0 for Linux x86 (32-bit) (October 10, 2011)"

In:= FinancialDerivative[{"American", "Call"}, {"StrikePrice" -> 110.00, "Expiration" -> 1},  {"InterestRate" -> 0.08, "Volatility" -> 0.,"CurrentPrice" -> 100, "Dividend" -> 0.0}]

Out= 7.10517

In:= FinancialDerivative[{"European","Call"}, {"StrikePrice" -> 110.00,"Expiration" -> 1},  {"InterestRate" -> 0.08, "Volatility" -> 0.2,"CurrentPrice" -> 100, "Dividend" -> 0.0}]

Out= 7.27904



the two results are different. To get the correct price also for the American option try this

In:= FinancialDerivative[{"American","Call"}, {"StrikePrice" -> 110.00,Expiration" -> 1},  {"InterestRate" -> 0.08, "Volatility" -> 0.2,"CurrentPrice" -> 100, "Dividend" -> 0.0}, "Method" -> "Binomial"]

Out= 7.27904

I reported the issue to Wolfram support also asking what is the default method for American options and how to get the full list of available methods. By now I just received a short answer saying

[…] I have
forwarded your example to our developers so that they can take a look into
this and resolve the issue for a future version of Mathematica.

so let’s hope that this will be fixed in Mathematica v9. If I will get more information I’ll modify the post, by now just use the Binomial method to get the correct price for American options.

# Understanding the Higgs Understanding the Higgs –  28’th November 2012. A workshop focused on the discovery of the new boson at LHC and on the current theoretical implications. A great opportunity to understand more on Higgs discovery and future perspectives from some of the theorists and experimentalists which are playing a key role in Higgs physics research. Also an opportunity to visit Turin, a beautiful city in the northern Italy.

To have an overview on Higgs discovery look at this nice video. For more technical information try this.

# Wanna be a quant: QuantNet 2012-2013 International Guide to Programs in Financial Engineering The QuantNet 2012-2013 International Guide to Programs in Financial Engineering have been released! Awesome guide, download it. Check out the table of contents: What do Financial Engineers Do? – Efficient Ways to Set Up a Successful Career – Which Master’s Degree is Right for You? – Do’s and Don’ts of MFE Applications – C++ Programming for Financial Engineering – Hunting the Quant Job: Questions and Answers – Questions Asked at Quant Interviews – Reading List: Books to Help Prepare for Quant Interviews, and more.

Don’t miss the chance to have a look at all the information you can find on QuantNet. For example,  if you deem yourself strong at C++ try their good (and free) C++ on-line tests.

# Hats and capricious secretary: a MonteCarlo trick to compute “e”

A popular question on basic probability  appearing in quant job interviews is the classical matching problem described in amusing forms: hats in a check-room are mixed and randomly distributed to the guests; a secretary has mixed a set of letters and their envelopes. What is the probability that no hat will match the right guest? What is the probability that exactly 3 letters will reach the right destination?

The more classical formulation of the same problem is about seeking for a match among two identical decks of N cards randomly mixed. Being $A_i$ the probability that the $i^{th}$ card will match (regardless to other possible matchings!), we have $P(A_i) = \frac{1}{N}\,\,, P(A_iA_j) = \frac{1}{N(N-1)}\,\,, P(A_{i_1}...A_{i_{i_r}}) = \frac{(N-r)!}{N!}$
The probability of at least one match is given by $P(A_1\cup A_2 \cup ... \cup A_N)$ $= \sum_{i=1}^{N} P(A_i)- \sum_{i $= 1- \frac{1}{2!}+\frac{1}{3!}-...+(-1)^{N-1}\frac{1}{N!}$
Therefore, in the limit of large N the probability approaches $1-1/e$, which is a quite surprisingly result. This result can be used to compute the probability of zero matches which is simply 1 minus the probability of at least one match.

From this simple result the amusing idea to try to compute $e$ with a MonteCarlo, just as it is usually done for $\pi$. Here my Mathematica code which does the job:

MonteCarloNepero[NN_, L_] := Module[{mean, check, indicatorCount},
indicatorCount = 0;
Do[
check = (# - RandomSample[#]) &@Table[l, {l, 1, L}];
If[FreeQ[check, 0],
indicatorCount = indicatorCount + 1;
];
, {NN}];
mean = indicatorCount/NN;
Return[N[1/mean]];
];

Where NN is the number of MonteCarlo runs and L is the length of the card deck. The convergence in the number of cards is very fast so with L=15 the result is already fair. Far more slow is the convergence in NN… For example we may get

In:= MonteCarloNepero[1000, 15]
Out:= 2.75482

In:= E // N
Out= 2.71828

Not that accurate, but fun.

# Wanna be a quant: job interview books reviews

Just updated my Books page with a few lines of review about some popular “Quant job interviews” books.

Well written and useful, as the other Joshi’s books. Many non-trivial questions and original solutions. Interesting and unusual the section about C++ programming.
The classic reference for quant job interviews. On top of the probability-statistics and math finance problems, there is big section on “brainteasers” (I don’t like them…) and a chapter on non-technical questions.
A more mathematical tone with respect the other interview books. Even the brainteasers are non-irritating. Each section of questions is introduced by a quick and clever review of the math that should be used to solve the problems.
Different from all the other interviews books: the FAQ are a good opportunity to jump through different interesting topics. The answers and the references given in the book are the starting point for a further study.

# Mathematica v8 – Linux 32-bit – FinancialDerivative issues

– fixed in Mathematica 9 –

Some problems with FinancialDerivative in Mathematica 8? Wait for v9 (should be released in a couple of months) or try the the Finance Platform.

While testing my app European Options, I used as a cross-check the excellent FinancialDerivative function of Mathematica 8. FinancialDerivative gives price and greeks for a wide choice of financial products, European and American style. Nevertheless at some point I came across some problems:

• For the “Out” BarrierOptions the “Rebate” term is not implemented. So you get the following wrong result by typing
In:= FinancialDerivative[{"BarrierUpOut", "European",
"Call"}, {"StrikePrice" -> 110.00, "Expiration" -> 0.5,
"Barrier" -> 105, "Rebate" -> 3}, {"InterestRate" -> 0.08,
"Volatility" -> 0.25, "CurrentPrice" -> 100,
"Dividend" -> 0.04}, {"Value", "Greeks"}]

Out= {0., {"Delta" -> 0., "Gamma" -> 0., "Rho" -> 0.,
"Theta" -> 0., "Vega" -> 0.}}

Which is correct with a zero rebate, while in this case the correct result should be

 Out=
{2.34535,{Delta->0.130242,Gamma->0.000861361,Rho->1.57311,Theta->-0.61502,Vega->2.02904}}
• The “Gamma” (the sensitivity of Delta with respect the spot price) for all the financial products I tested so far is surprisingly inaccurate. As an example I show here what we get for simple vanilla call and put for which it easily derived by call-put paritythat the Gamma should be the same.
In:= FinancialDerivative[{"Vanilla", "European",
"Put"}, {"StrikePrice" -> 110.00,
"Expiration" -> 0.5}, {"InterestRate" -> 0.08, "Volatility" -> 0.25,
"CurrentPrice" -> 100, "Dividend" -> 0.04}, {"Value", "Greeks"}]

Out= {11.6465, {"Delta" -> -0.619661, "Gamma" -> 0.020876,
"Rho" -> -36.8063, "Theta" -> -3.11937, "Vega" -> 26.1189}}

In:= FinancialDerivative[{"Vanilla", "European",
"Call"}, {"StrikePrice" -> 110.00,
"Expiration" -> 0.5}, {"InterestRate" -> 0.08, "Volatility" -> 0.25,
"CurrentPrice" -> 100, "Dividend" -> 0.04}, {"Value", "Greeks"}]

Out= {3.97952, {"Delta" -> 0.360538, "Gamma" -> 0.0208805,
"Rho" -> 16.0371, "Theta" -> -7.65352, "Vega" -> 26.1189}}

The two Gamma are different from each other and different from the correct result which is in this case Gamma = 0.0208952.

I reported the bug to the Wolphram Support and after few days I received a first answer. The guy at support was very kind, but told me that the problem was not showing up on his version (Windows and Finance Platform). After a few try the final answer has been the following:
Hi Stefano

I have carried out the same evaluations as you mentioned for FinancialDerivative for BarrierDownOut and BarrierUpOut options under Linux and have been able to replicate your problems. Thank you for bringing this to our attention. It has now been fixed for version 9 of Mathematica which will be released in a couple of months. Also the problems with the inaccuracy of Gamma were not showing up in my version because I was running it on the Finance Platform, which had fixed that bug. All these errors have been resolved in version 9 of Mathematica. Thanks again for sending in your emails about this as they help to improve our product.

So… keep waiting for v9.

– fixed in Mathematica 9 –

# Free Open-Source Statistics Cookbook Ryan Swanstrom

Matthias Vallentin, a computer science PhD student at UC Berkeley, has published a Probability and Statistics Cookbook. The book can be freely downloaded in PDF format via the website. Also, the latex source is available on Github. Matthias states that others are free to fork the source and make changes.

The book is not a textbook. It is more of a cheatsheet. It contains many of the common probability and statistics techniques and the associated formula. I would consider this book to be an excellent resource to have around.

View original post

# Sampling from Cauchy (and Gaussian) distribution

How we construct  a random draw  from a Cauchy distribution, if we start with a draw from uniform distribution? That’s what asks the Exercise 7.3 of the book The Concepts and Practice of Mathematical Finance by Mark Joshi. We post here some C++ code as an answer. This is easily done by explicitly computing the inverse cumulative distribution function $F^{-1}(x)$ of the Cauchy distribution and plugging then random numbers generated with a uniform distribution (this is done with the function GetOneRandomByICD). $f(x) = \frac{1}{\pi}\frac{1}{1+x^2},\; F(x) = \frac{1}{\pi}\arctan (x) + \frac{1}{2},$ $F^{-1}(x) = \tan [\pi (x-1/2)]$

More tricky is the case in which we don’t have the $F^{-1}(x)$. In that case we can obtain our sample by using two random uniforms – GetOneCauchyByUniforms –  or using two random numbers with gaussian distribution – GetOneCauchyByGaussians. In the latter case of course we need also a routine  to compute gaussian numbers. We present here the routine GetOneGaussianByRatio, which is similar to GetOneCauchyByUniforms (both are based on ratio of uniforms – see Probability and Random Processes by Grimmet and Stirzaker). The mathematical proofs are in the probability section of my notes.

//Random.h
#ifndef RANDOM_H
#define RANDOM_H

double CauchyICD(double uniform);
double GetOneRandomByICD(double (*func_ptr)(double));

double GetOneCauchyByUniforms();
double GetOneCauchyByGaussian();

double GetOneGaussianByRatio();

#endif
//Random.cpp
#include "Random.h"
#include <cmath>
#include <cstdlib>

double CauchyICD(double x){
return std::tan(M_PI*(x-0.5));
}
double GetOneRandomByICD(double (*func_ptr)(double)){
return (*func_ptr)(rand()/static_cast<double>(RAND_MAX));
}

double GetOneGaussianByRatio(){

double x;
double y;

double ratio;

do {
x = 2.*rand()/static_cast(RAND_MAX)-1.;
y = 2.*rand()/static_cast(RAND_MAX)-1.;
ratio = x/y;
}
while((y > std::exp(-ratio*ratio*0.25))
|| (y < -std::exp(-ratio*ratio*0.25)) );

return ratio;
}

double GetOneCauchyByUniforms(){

double x;
double y;

double sizeSquared;

do {
x = 2.*rand()/static_cast<double>(RAND_MAX)-1.;
y = 2.*rand()/static_cast<double>(RAND_MAX)-1.;
sizeSquared = x*x + y*y;
}
while(sizeSquared >= 4.0);

return y/x;
}

double GetOneCauchyByGaussian(){
double x = GetOneGaussianByRatio();
double y = GetOneGaussianByRatio();

return x/y;
}

A simple test program which write to files samples of Cauchy distributed numbers:

//testRandom.cpp
#include "Random.h"
#include <iostream>
#include <fstream>
#include <cmath>
#include <cstdlib>

int main(){

long sampleLength = 10000;

std::ofstream myfile1;
myfile1.open ("inverse_cauchy.txt");
std::ofstream myfile2;
myfile2.open ("ratio_cauchy.txt");
std::ofstream myfile3;
myfile3.open ("gauss_cauchy.txt");

for(long n=0; n< sampleLength; n++){
myfile1 << GetOneRandomByICD(&CauchyICD) << "\n";
myfile2 << GetOneCauchyByRatio() << "\n";
myfile3 << GetOneCauchyByGaussian() << "\n";
}

myfile1.close();
myfile2.close();
myfile3.close();

return 0;
}

Now just compile and try it.

# Pricing on Trees notes

Just added to the Notes page a few pages about option pricing on trees. We explicitly show the equivalence of replicating, hedging or using a risk-neutral approach to price options on a binomial tree. A Mathematica implementation of binomial and trinomial tree will be added soon, with discussion on convergence issues.

# European Options on Google Play

Easily compute price and greeks of european-style option derivatives. European Options is an handy app to compute price and greeks of Vanilla, Digital and Barrier european-style options (for definitions we refer to the standard book “Hull – Options, futures and other derivatives”). It is planned to add more exotic options. Price and greeks are computed via analytical formulas.

The source code is available at my GitHub repo. 