Serving the Quantitative Finance Community

 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

The "best" method to compute bivariate cumulative normal distribution

January 14th, 2016, 10:04 am

Is there an 'owner' of the GW code?
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

The "best" method to compute bivariate cumulative normal distribution

January 14th, 2016, 1:15 pm

Very good.Some updates1. Even with NX = NY ~ 200/300 I get 6 figures accuracy which is nice.2. My schemes is O(hx^2 + hy^2) and I used Richardson extrapolation to get O(hx^4 + hy^4).3. I can reproduce GenzWest performance; when abs(rho) < 0.925 and rho < 0 the code calls Asin(), Sin() and Exp(). This slows things down big time. In that case, the code is 20 times slower than mine; if rho >= 0.925 it is 7 times slower. That's the bottleneck for sure. Graeme' VBA code is in Espen Haug's book (page 479). See BVND The C and VBA were translated from the Fortran.. Fortran != C.
Last edited by Cuchulainn on January 13th, 2016, 11:00 pm, edited 1 time in total.
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

The "best" method to compute bivariate cumulative normal distribution

January 15th, 2016, 8:30 pm

QuoteOriginally posted by: EBalQuoteOriginally posted by: CuchulainnWhat is the best method forTaking what list1 suggested one step further, you can transform this double integral to an ordinary integral, and use your favorite method to compute it[$]\frac{1}{\sqrt{2\pi}} \int_{-\infty}^b dy\, e^{-y^2/2} N\left(\frac{a-\rho y}{\sqrt{1-\rho^2}}\right)[$]where [$]N[$] is cumulative normal function. Most likely, there is no closed form solution to the ordinary integral, except for special values of [$]a[$], [$]\rho[$].Is it possible to differentiate this double integral twice with respect to a and b to get a PDE u_ab = Known right hand side function The tricky part is the lower limit of inner integral is not constant.?
Last edited by Cuchulainn on January 14th, 2016, 11:00 pm, edited 1 time in total.
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

The "best" method to compute bivariate cumulative normal distribution

January 17th, 2016, 9:54 am

[$]\frac{1}{\sqrt{2\pi}} \int_{-\infty}^b dy\, e^{-y^2/2} N\left(\frac{a-\rho y}{\sqrt{1-\rho^2}}\right)[$]This is essentially Abramowitz/Stegun 26.3.3 (2nd line) in separated form as the product of two integrals. I used the erf(x) to approximate the inner integral and Tanh Rule for the outer integral (the usual integrators take forever). Very good approximation and desired error can be tuned.Also: the formula works for rho = {-1, 1}, which seems to be an issue for other methods.So, 3 methods1. Goursat2. Genz3. Erf + Tanh
Last edited by Cuchulainn on January 16th, 2016, 11:00 pm, edited 1 time in total.
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

The "best" method to compute bivariate cumulative normal distribution

January 18th, 2016, 3:19 pm

QuoteSo, 3 methods1. Goursat2. Genz3. Erf + TanhMoving up a gear, I tested each of these methods on a number of 2-asset problems in Collector's book and the results are consistent across the board. Goursat is still the fastest.It is a good test of C++11 erf as pointed out by outrun.So, for me the 2d case is finished.
Last edited by Cuchulainn on January 17th, 2016, 11:00 pm, edited 1 time in total.
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 10:54 am

Taking what list1 suggested one step further, you can transform this double integral to an ordinary integral, and use your favorite method to compute it[$]\frac{1}{\sqrt{2\pi}} \int_{-\infty}^b dy\, e^{-y^2/2} N\left(\frac{a-\rho y}{\sqrt{1-\rho^2}}\right)[$]where [$]N[$] is cumulative normal function. Most likely, there is no closed form solution to the ordinary integral, except for special values of [$]a[$], [$]\rho[$].
This is certainly the most elegant solution, mathematically IMO and takes the least number lines of code, and the most accurate if  you use the Tanh Integration rule. The latter is a bit slower than Genz and PDE approach.
Is this formula applicable to trivariate?
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 12:49 pm

Isn't the fact that the integration limits are no longer rectangular a problem for most methods?
That makes A&S 26,3,3 tricky, in principle. Lambda functions are useful

auto f = [&](double d) { return Z(d)*CdfNormal((y-rho*d)*rho2); };



// Tanh Rule
result = DEIntegrator::Integrate(f, AL, x, 1.0e-16, evals, error);
 
User avatar
outrun
Posts: 4573
Joined: January 1st, 1970, 12:00 am

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 12:57 pm

The latter is a bit slower than Genz and PDE approach.
If we require to compute something like M(x<=0.4, y<=0.12, rho=0.9) with 7 digits precision how long wil Tanh, Genz and PDE take? The PDE is *a lot* slower than Tanh or Genz in this case. You need to explain how you measure PDE speed and what it's precision is and why my simple example conflicts with that.
I think Genz pretty much does what's stated, translate it to a 1d integral.
 
User avatar
Cuchulainn
Topic Author
Posts: 20250
Joined: July 16th, 2004, 7:38 am
Location: 20, 000

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 2:33 pm

The latter is a bit slower than Genz and PDE approach.
If we require to compute something like M(x<=0.4, y<=0.12, rho=0.9) with 7 digits precision how long wil Tanh, Genz and PDE take? The PDE is *a lot* slower than Tanh or Genz in this case. You need to explain how you measure PDE speed and what it's precision is and why my simple example conflicts with that.
I think Genz pretty much does what's stated, translate it to a 1d integral.
Your analysis is incorrect. Your intuition tells you it is slower. I have done all the tests and will post on the other thread. 
And this thread on "best method", Efficiency is only one of the metrics.

BTW have you ever used these methods? What's your run time for Tanh for example? Is it 10, 20, 30 times slower than Genz?

// 
Examples on rough meshes

** Bivariate a,b, rho (one-off scenario): 0.4,0.12,0.9
** Bivariate NX, NY:        15,15
Goursat Classico NX:        0.5275646914439414,
Goursat Extrap 2*NX:        0.519400274728015
Goursat Extrap 4*NX:        0.5195176107403655
Genz:                       0.5195173418566086
Tanh +(A&S 26.3.3)       :  0.5195173411668271
Tanh fun evals and error :  385,2.842170943040401e-15
 
ex2

** Bivariate a,b, rho (one-off scenario): 8,1.1,-0.5

** Bivariate NX, NY:        15,15
Goursat Classico NX:        0.8680812732240816,
Goursat Extrap 2*NX:        0.864315708265091
Goursat Extrap 4*NX:        0.8643339688264521
Genz:                       0.8643339390536167
Tanh +(A&S 26.3.3)       :  0.8643339395362047
Tanh fun evals and error :  385,9.71445146547012e-17
 
User avatar
outrun
Posts: 4573
Joined: January 1st, 1970, 12:00 am

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 3:06 pm

The latter is a bit slower than Genz and PDE approach.
If we require to compute something like M(x<=0.4, y<=0.12, rho=0.9) with 7 digits precision how long wil Tanh, Genz and PDE take? The PDE is *a lot* slower than Tanh or Genz in this case. You need to explain how you measure PDE speed and what it's precision is and why my simple example conflicts with that.
I think Genz pretty much does what's stated, translate it to a 1d integral.
Your analysis is incorrect. Your intuition tells you it is slower. I have done all the tests and will post on the other thread. 
And this thread on "best method", Efficiency is only one of the metrics.

BTW have you ever used these methods? What's your run time for Tanh for example? Is it 10, 20, 30 times slower than Genz?
Your results show me it's slower, there is no intuition.You need very large grid to get a couple of digits precision. 
Shall we do a defined bet then? I bet that:
1) to compute  M(x<=0.4, y<=0.12, rho=0.9) with 7 digits precision we have Tahn and Genz being *much* faster than PDE. Probably at least a factor 100. You claim as quoted that PDE is faster and so I bet it's it slower by at leat a factor 10. You are allowed to use Richardson extrapolation but need to of course include that in your timing. The Genz answer to this bet is 5.195173365E-1

To answer your question: depending on the required precision you need to pick a number of quadrature points for the 2d numerical tanh integral . In my experience both Gauss Lengedre and Tanh quatrature converge very fast (as a function of the number of nodes) fro smooth function like the bivariate Guassian. A 2d numerical Tanh integral is more generic than the specialized Genz method and I expect Genz to be at least 5 times fasters for similar precision.  

With the low precision 7 digits  requirements in the bet the PDE method will also be at least factor 1000 less accurate. So its' both slower *and* less accurate by orders of magnitudes!
 
User avatar
outrun
Posts: 4573
Joined: January 1st, 1970, 12:00 am

Re: The "best" method to compute bivariate cumulative normal distribution

March 8th, 2017, 3:16 pm

** Bivariate NX, NY:        15,15
Goursat Classico NX:        0.5275646914439414,
Goursat Extrap 2*NX:        0.519400274728015
Goursat Extrap 4*NX:        0.5195176107403655
Genz:                       0.5195173418566086
Tanh +(A&S 26.3.3)       :  0.5195173411668271
Let's discuss this example since you've posted it.
the 4x extrapolation seems to have 6 significant digits,  right (519 517)? This means you've computed at some point a 60x60 grid -since you started with a 15x15 grid (Nx = 15.. later on it say's 4*Nx)? So how long did that computation take compared to the Genz calculation?
 
User avatar
katastrofa
Posts: 7440
Joined: August 16th, 2007, 5:36 am
Location: Alpha Centauri

Re: The "best" method to compute bivariate cumulative normal distribution

March 9th, 2017, 2:22 pm

Isn't the fact that the integration limits are no longer rectangular a problem for most methods?
That makes A&S 26,3,3 tricky, in principle. Lambda functions are useful

auto f = [&](double d) { return Z(d)*CdfNormal((y-rho*d)*rho2); };



// Tanh Rule
result = DEIntegrator::Integrate(f, AL, x, 1.0e-16, evals, error);
Image
 
User avatar
outrun
Posts: 4573
Joined: January 1st, 1970, 12:00 am

Re: The "best" method to compute bivariate cumulative normal distribution

March 9th, 2017, 2:46 pm

..about the non-rectangular region: indeed, I've run into that. E.g. this paper computes the value of spread options on two GBMs using quadrature. The non-linear transform they do makes the region on integration non-rectangular.

About that 1e-16: you can change the weights of a quadrature and pick different orthogonal polynomials depending on the type of interval: see wikipedia