15
$\begingroup$

I've been using low discrepancy sequences for a while for Uniform Distributions, as I've found their properties useful (mainly in computer graphics for their random appearance and their ability to densely cover [0,1] in an incremental fashion).

For example, random values above, Halton sequence values below:

https://i.stack.imgur.com/gCXrh.png

I was considering using them for some financial analysis planning, but I need different distributions than just uniform. I started off trying to generate a normal distribution from my uniform distributions via the Marsaglia polar algorithm, but the results don't seem as good as with the uniform distribution.

Another example, again random above, Halton below:

https://i.stack.imgur.com/NNcdH.png

My question is: What is the best method for getting a normal distribution with the properties I get from a uniform low discrepancy sequence - coverage, incremental fill-in, non-correlation across multiple dimensions? Am I on the right track, or should I be taking a completely different approach?

(Python code for uniform and normal distributions I use above: Gist 2566569)

Improve this question
$\endgroup$
7
  • 3
    $\begingroup$ Surely you have looked at the natural and obvious answer, which is to apply the probability transform to the uniform values. In what respect(s) does this approach fail to meet your needs? $\endgroup$
    – whuber ♦
    May 1, 2012 at 15:55
  • $\begingroup$ That was what I had assumed was the normal way of creating such distributions from a uniform input, but I had trouble finding any such algorithms (I'm a programmer, not a statistician, so I was starting from scratch in my search). The methods I'd found were Box Muller, Marsaglia Polar, Ziggarat and Ratio of Uniforms, all of which require mixing more than one input sample which appears to destroy the characteristics I was trying to preserve in the Low Discrepancy sequences. I've just found the following page: home.online.no/~pjacklam/notes/invnorm which has given me the info I want. $\endgroup$
    – Edouard Poor
    May 2, 2012 at 1:37
  • 3
    $\begingroup$ Yes, that site implements the probability transform for the Normal distribution. You may have had trouble researching this because it is so extremely well known that virtually any alternative method of generating random values implicitly assumes you could use this method and then seeks more efficient ways in special cases (including the four methods you list). $\endgroup$
    – whuber ♦
    May 2, 2012 at 13:47
  • 3
    $\begingroup$ Inversion will be the only valid method for this case! $\endgroup$
    – kjetil b halvorsen ♦
    Oct 15, 2012 at 14:48
  • $\begingroup$ I wonder if we can also get some other distribution's quasi-random values, e.g. beta distribution's quasi-random value? $\endgroup$
    – user315648
    Apr 14, 2014 at 3:15

4 Answers 4

Reset to default
7
$\begingroup$

You can transform from $\mathcal U(0, 1)$ random variables to any other distribution using the inverse of the CDF, also called the percent point function or quantile function. It's implemented in scipy as scipy.stats.norm.ppf.

Improve this answer
$\endgroup$
6
$\begingroup$

I recently stumbled upon this problem. Naively I thought that any transformation from uniform would work, so I plugged in a 1D Sobol (and Halton) sequence as if the sequence where a random number generator into an std::normal_distribution<> variate. To my surprise it didn't work, it obviously generated a non normal distribution.

Ok, then I took the Numerical Recipes Third Edition Chapter 7.3.9 Normal_dev function to generate normal numbers from the Sobol or Halton sequences by the method of "Ratio-of-Uniforms" and it failed in the same way. Then I though, ok, if you look at the code, it takes two uniform random numbers to generate two normally distributed random numbers. Perhaps if I used a Sobol (or Halton) 2D sequence it will work. Well, it failed again.

The I remembered about the "Box-Muller method" (mentioned in the comments) and since it has a more geometric interpretation then I though it could work. Well, it did work! I was very excited an starting doing other test, the distribution looks normal.

The problem I saw was that the distribution was no better than random, it terms of filling, so I was a bit disappointed, but ready to publish the result.

Then I did a deeper search (now that I knew what to look for), and it turn out that there is already a paper on this subject: http://www.sciencedirect.com/science/article/pii/S0895717710005935

In this paper it is actually claimed

Two well known methods used with pseudorandom numbers are the Box–Muller and the inverse transformation methods. Some researchers and financial engineers have claimed that it is incorrect to use the Box–Muller method with low-discrepancy sequences, and instead, the inverse transformation method should be used. In this paper we prove that the Box–Muller method can be used with low-discrepancy sequences, and discuss when its use could actually be advantageous.

So the overall conclusion is this:

1) You can use the Box-Muller on 2D low discrepancy sequences to obtain normally distributed sequences. But my few experiments seem to show that the low discrepancy/space, e.g. filling properties are lost in the normal-transformed sequence.

2) You can use the inverse method, presumably the low discrepancy/space filling properties will be preserved.

3) Ratio-of-Uniforms cannot be used.

EDIT: This https://mathoverflow.net/a/144234 points to the same conclusions.

I made an illustration (the first figure (Ratio-of-uniforms on Sobol) shows that the distribution obtained is not normal but the ohters (Box-Muller and random for comparison) are ):

figures

EDIT2:

The main point is that, even if you find a method that can transform the "distribution" of a low discrepancy sequence, it is not obvious that you will preserve the good filling properties. So you are not better than with a truly random (standard) normal distribution. I have yet to find a method that is low discrepancy and yet it fills nicely with a non uniform distribution. I bet such method is very non-obvious and perhaps an open problem.

Improve this answer
$\endgroup$
2
  • $\begingroup$ In the above mentioned paper (sciencedirect.com/science/article/pii/S0895717710005935) the last sentence reads: "Based on these results, we conclude that Box–Muller is a good alternative to the inverse transformation method to generate low-discrepancy sequences from the normal distribution." Doesn't this contradict point 1 of the stated conclusions? $\endgroup$
    – ngiann
    Sep 13, 2019 at 20:58
  • $\begingroup$ @ngiann, Not contradictory. The paper reaches my same conclusion that the distribution is correctly obtained by the BM transformation but not by the other methods. The second part of the sentence ("But the low discrepancy/space filling properties seem to be lost.") is related to a secondary property one looks for in low-discrepancy sequences, such as efficient filling of the space. My experiment didn't show that property translated to the normal-transformed sequence. If I remember correctly the paper doesn't get into that specific problem. I think this needs more study to be conclusive. $\endgroup$
    – alfC
    Sep 13, 2019 at 21:24
1
$\begingroup$

There are two good methods. First, as noted above, an accurate approximation to the inverse of the Gaussian distribution can be used. Then one can transform any low-discrepancy sequence into a Gaussian.

The second method is the Box-Muller. This method requires two in put numbers (R and A) and generates two outputs. A two-dimensional low-discrepancy sequence is needed. One takes (for example in the Halton Sequence), pairs of primes are used, one for the radial component (R) and one for the angular component (A). One gets Sqrt(-2*Log(R)) for the radial component and Sin(2*Pi*A) and Cos(2*Pi*A) for the angular components. Multiplying the radial by the two angular components (separately) gives two Gaussians. The efficiency is the same as above; two quasi-random inputs and two Gaussian outputs.

Any multi-dimensional low-discrepancy sequence can be used, depending on the dimensionality of the problem.

Improve this answer
$\endgroup$
0
$\begingroup$

The most native method would be indeed using the inverse CDF to transform into normal Gaussian, but there are also a problems with this. If you have e.g. a LDS point set create by rank-1 lattices, then it would be that the starting point is always (0,0), so to transform it you need a little shift, best to have the same gap as for the corner (1,1).

So far no problem still, but for an ideal Gaussian distribution N(0,1)+N(0,1) should give the same distribution as the difference. However, this would be not the case using rank-1 lattice LDS and iCDF on each variable, because the starting point in each variable would give a certain iCDF, like $-3\sigma$ (depending on N), so the difference would be $-6\sigma$.

And that is a too extreme value, leading really to a systematic error (e.g. you will not get $+6\sigma$ at the other side). Best inspect your transformed LDS also for sum and differences, check for such extreme points and also for skew and kurtosis.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.