Multivariate normal CDF

General
#1

Hi all,

I need to use the 10 dimensional normal CDF inside Stan. Given that Stan only officially supports the one-dimensional normal CDF, I wonder what my options are?

E.g., is numerical integration possible for my purpose? Is it possible to implement some algorithm that underlies the multivariate normal CDF packages in R or other languages (like the mvtnorm package) (which I’m willing to do myself)?

I have used a 2-dimensional normal CDF inside Stan with the help of a code circulated on this forum. But according to some post on this forum dated several years ago, one reason Stan doesn’t have higher-dimensional normal CDF is the lack of a deterministic way to calculate normal CDF. However, my current project does require me to use a 10-dimensional normal CDF, and I don’t think there is a way to reduce the dimension (because all elements of the 10-dimensional covariance matrix are not zero). If I don’t need the CDF to be super precise or deterministic, would that open up some possibilities for implementing the CDF?

Any suggestion is greatly appreciated! Thank you!

1 Like
#2

If I recall correctly, @spinkney Recently worked a copula-based method for the truncated MVN CDF. Sean did you happen to have a non-truncated CDF implemented as well?

1 Like
#3

I haven’t gotten around to testing the implementation in the TruncatedNormal R package referenced in Multivariate normal cdf. The neat thing is that it allows the non truncated and truncated versions.

I think it would be really interesting to see if it works with the algebra solver. You can try coding it in stan and we can help out if you run into any issues. If it works well then it’s a great candidate to be added to stanmath.

3 Likes
#4

Thanks! That sounds exciting. I’ll look into that.

#5

I’ve had a reasonable success in using nuisance parameters to have Stan integrate the CDF for me during sampling. I feared the models would be super slow, but the slowdown is actually quite manageable so far. The best way to do this turned out to be (for my case - multivariate probit) via the “bgoodri parametrization” (example-models/probit-multi-good.stan at master · stan-dev/example-models · GitHub) which uses the GHK algorithm to build the the univariate distributions.

How directly that’s applicable to your case would depend on how does the CDF fit in with the rest of your model…

Based on advice by @CerulloE , it also turned out to be useful to use logit approximation to Phi and inv_Phi

real approx_Phi(real x) {
        return inv_logit(x * 1.702);
      }

      real approx_inv_Phi(real x) {
        return logit(x) / 1.702;
      }
2 Likes