Numerical instability issues occur when the input to the inverse CDF is very close to 0 or 1.
scipy handles these cases with approximations and here is the screenshot of the comments in scipy’s source code 4:
image954×524 34.8 KB
One example is when y=1-1e-8:


torch.distributions.Normal(0, 1).icdf(torch.tensor(y)) returns tensor(inf) (because torch.erfinv returns tensor(inf));

scipy.stats.norm.ppf(y) returns 5.612001243305505.