4

In order to generate normal random numbers, one usually generates random numbers following a uniform distribution $Z \sim \mathcal{U}(0,1)$ and then applies the reverse CDF function on them $X=\Phi^{-1}(Z) \sim \mathcal{N}(0,1)$.

However, I encountered a problematic case when one of the generated $Z$ turns out to give exactly 0. Then, you have $X=\Phi^{-1}(Z)=- \infty$.

This is pretty problematic when you generated random samples because it will usually break all your variance/covariance measure basically returning nan or inf when the samples contain infinite number.

How do you usually handle this? Do you check after each generated random number whether the value is 0 or 1 and shift it slightly (or simply dicard it)?

SRKX
  • 11,126
  • 4
  • 42
  • 83

4 Answers4

2

Indeed for computational purposes, best you can do is use a uniform distribution on another interval $[10^{-10},1-10^{-10}]$, or just discard all occurences of $Z=0,1$.

Discarding $Z=0,1$ is justified, since for continuous distributions $P(Z=0)=P(Z=1)=0$.

emcor
  • 5,795
  • 3
  • 31
  • 57
  • Well for continuous distribution $\mathbb{P}[Z=z]=0 ~ \forall z$... – SRKX Jul 22 '14 at 16:00
  • @SRKX Yes, you can leave out any single values because their probability is $0$. Only you cannot leave out intervals of values. If you leave out $0.5$, the generator will still create infinitely many values around $[0.5+-\epsilon]$. – emcor Jul 22 '14 at 16:06
  • While that is theoretically true, a real RNG in code can only return a discrete set of values, so all values removed represent a range. So you can remove a few problematic values, but it will leave (small) holes in your resulting distribution. – Phil H Jul 28 '14 at 09:29
2

I agree with @encor that it isn't an issue to include some logic to avoid errors. I imagine that most non-uniform random number generators already include that.

I don't think I've used a pseudo-random number generator that's given a 0. I'm not an expert on the topic, but a common implementation (linear congruential generator) relies on modular arithmetic. For instance, the generator $X_{n+1}=f(X_{n})mod\: m$ might produce numbers between 0 and $m-1$ where $m$ is some big number. You could then divide by $m$ and ensure all the numbers are between 0 and 1. The only way that you should get a zero is if $f(X_{n})$ equals $m$. If $m$ is large enough (and most seem to be using $2^{31}$ or $2^{32}$), then that should be very rare.

The only other thing I could think of was an issue with types. For instance, you might need an unsigned int or long to generate the random integers. The division to get between 0 and 1 suggests a float or double type. It might be rounding things down to 0 instead of 0.xxx1.

If you're not happy with the above, there are other techniques for sampling distributions. For instance, you can represent the distribution as a grid as in Meucci's Fully Flexible Extreme Views (though without the Entropy Pooling part of it). You can then use linear interpolation as he describes to simulate from the grid distribution.

John
  • 5,391
  • 22
  • 28
2

This depends on your method to generate the normal random numbers. The problem with normal cdf is that the direct inverse $\phi^{-1}(Z)$ is hard to solve for directly. There are some other methods to generate $N(\mu, \sigma^2)$ from $U(0,1)$. Two notable methods are:

  1. Box-Muller method
  2. Marsaglia polar method

For most purposes you can use the above methods but you need to ensure the uniform random numbers are good. You will not get the problem of $-\inf$ or $\inf$. I prefer Box-Muller since I don't have to discard any generated numbers.

Taran
  • 311
  • 2
  • 6
  • 1
    These are other simulation methods, but they share the same problem still: If $U$ or $s$ respectively become $0$, the value becomes infinite. – emcor Jul 22 '14 at 19:16
  • 1
    So the uniform random generator is generating 0 values? I personally haven't come across exact 0 value from a uniform pseudo-RNG. Can you share which algorithm are you using? It could be very very small where you have to deal with problems of limits on storing them. – Taran Jul 22 '14 at 19:31
0

Most uniform RNGs initially output integers in [0,2^32). When such a value is then being converted to double, which have more than 32 bits of mantissa precision, and then scaled to [0,1), the best route is to just add 2^-33 at the end to center the outputs around 1/2.

This way no values are wasted, no checks are needed and the (0,1) interval is best covered reaching farther away in the tails. Furthermore on modern processors and compilers this can come completely for free.

Discarding/shifting 0 or 1, or scaling to other intervals are both costly and produce a suboptimal output.

Quartz
  • 1,553
  • 10
  • 15