Stata The Stata listserver
[Date Prev][Date Next][Thread Prev][Thread Next][Date index][Thread index]

Re: st: calculation of bivariate normal probabilities


From   Arnold Kester <arnold.kester@stat.unimaas.nl>
To   statalist@hsphsun2.harvard.edu
Subject   Re: st: calculation of bivariate normal probabilities
Date   Wed, 25 May 2005 14:40:07 +0200

Op 05/25/2005 08:18 AM schreef Jenkins S P:
Stata 9 provides the function -binormal(h,k,r)-, which returns the joint cumulative distribution of the bivariate normal with correlation r; cumulative over (-inf,h]x(-inf,k]. [The equivalent function in version 8.2 was -binorm(h,k,r)-.)

I am not 100% confident about calculation of bivariate normal probabilities once we move outside the standard bivariate normal distribution (or indeed standard univariate normal distribution), and seek advice on this.

My understanding is that -binormal()- evaluates the double integral

INT(-inf,h] INT(-inf,k) { f(x1, x2)dx1dx2 }

where f(x1,x2) refers to the standard bivariate normal density with variates X1, X2 each having mean zero, unit variance, and correlation of r. (Put another way, the correlation matrix R has ones on the leading diagonal, and r in the off-diagonal cells)

By contrast, I want to evaluate

INT(-inf,p] INT(-inf,q) { g(x1,x2)dx1 dx2 }

where g(x1,x2) refers to a general bivariate normal density with variates X1, X2 having means m1 and m2, variances v1, v2, and covariance c (the covariance matrix is C).

To use -binormal()- for the calculation, my inclination was to derive the correlation matrix R from covariance matrix C [-matrix R = corr(C)-] and to extract r from this, and also to set the upper integration points as follows:

h = (p - m1)/sqrt(v1)
q = (q - m2)/sqrt(v2)

Is this correct?
Yes, your transformations standardise the distribution. The covariance matrix is standardised by dividing the first row and first column by sd1 and the second row and column by sd2 (sd=sqrt(v)).


In fact, I am also seeking an answer for the general multivariate normal case. (If one had a function available to evaluate standard multivariate normal distribution probabilities analogous to -binormal()-, would evaluation for the general case proceed in a similar way? I.e. derive R from C, 'standardize' the upper integration points, and then substitute these into the function.)
The statistical computing program R (cran.r-project.org) has the package mvnorm which computes multivariate normal (and T) distribution functions, also for higher than two dimensions.

--
Met vriendelijke groet,
Arnold Kester

*
* For searches and help try:
* http://www.stata.com/support/faqs/res/findit.html
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/




© Copyright 1996–2014 StataCorp LP   |   Terms of use   |   Privacy   |   Contact us   |   What's new   |   Site index