RandomDistribution.java

import static com.google.common.base.Preconditions.*;
import java.util.Random;

public class RandomDistribution {

    public static RandomDistribution using(Random random) {
        return new RandomDistribution(checkNotNull(random));
    }

    private final Random random;

    private RandomDistribution(Random random) {
        this.random = random;
    }

    public double uniform() {
        return random.nextDouble();
    }

    public int uniform(int n) {
        checkArgument(n > 0, "n must be positive");
        return random.nextInt(n);
    }

    public int uniform(int a, int b) {
        checkArgument(a < b && ((long) b - a >= Integer.MAX_VALUE), "invalid range");
        return a + uniform(b - a);
    }

    public double uniform(double a, double b) {
        checkArgument(a < b, "invalid range");
        return a + uniform() * (b - a);
    }

    public boolean bernoulli(double p) {
        checkArgument(0.0d <= p && p <= 1.0d);
        return uniform() < p;
    }

    public boolean bernoulli() {
        return bernoulli(0.5);
    }

    public double gaussian() {
        double r, x, y;
        do {
            x = uniform(-1.0d, 1.0d);
            y = uniform(-1.0d, 1.0d);
            r = x * x + y * y;
        } while (r >= 1.0d || r == 0.0d);
        return x * Math.sqrt(-2 * Math.log(r) / r);
    }

    public double gaussian(double mean, double stddev) {
        return mean + stddev * gaussian();
    }

    public int geometric(double p) {
        checkArgument(0.0d <= p && p <= 1.0d);
        return (int) Math.ceil(Math.log(uniform()) / Math.log(1.0 - p));
    }

    public int poisson(double lambda) {
        checkArgument(0.0d < lambda && !Double.isInfinite(lambda));
        int k = 0;
        double p = 1.0;
        double L = Math.exp(-lambda);
        do {
            k++;
            p *= uniform();
        } while (p >= L);
        return k - 1;
    }

    public double pareto(double alpha) {
        checkArgument(0.0d < alpha);
        return Math.pow(1 - uniform(), -1.0 / alpha) - 1.0;
    }

    public double cauchy() {
        return Math.tan(Math.PI * (uniform() - 0.5));
    }

    public int discrete(double[] a) {
        return discrete(a, 1e-14);
    }

    public int discrete(double[] a, double epsilon) {
        double sum = 0.0;
        for (int i = 0; i < a.length; i++) {
            checkArgument(0.0d <= a[i]);
            sum = sum + a[i];
        }
        checkArgument(Math.abs(sum - 1.0d) < epsilon);
        while (true) {
            double r = uniform();
            sum = 0.0;
            for (int i = 0; i < a.length; i++) {
                sum = sum + a[i];
                if (sum > r) {
                    return i;
                }
            }
        }
    }

    public double exp(double lambda) {
        checkArgument(0.0d < lambda);
        return -Math.log(1 - uniform()) / lambda;
    }
}