Many resources are needed to download a project. Please understand that we have to compensate our server costs. Thank you in advance.
Project price only 1 $
You can buy this project and download/modify it how often you want.
net.librec.math.algorithm.Randoms Maven / Gradle / Ivy
// Copyright (C) 2014-2015 Guibing Guo
//
// This file is part of LibRec.
//
// LibRec is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// LibRec is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with LibRec. If not, see _tempList = new ArrayList<>();
/**
* Random generate an integer in [0, range)
*
* @param range range of the interval
* @return an integer random generated in [0, range)
*/
public static int uniform(int range) {
return uniform(0, range);
}
public static void seed(long seed) {
r = new Random(seed);
}
/**
* Random generate an integer in [min, max)
*
* @param min minimum of the range
* @param max maximum of the range
* @return an integer random generated in [min, max)
*/
public static int uniform(int min, int max) {
return min + r.nextInt(max - min);
}
/**
* Return real number uniformly in [0, 1).
*
* @return real number uniformly in [0, 1).
*/
public static double random() {
return uniform();
}
/**
* Return a random number from a given list of numbers.
*
* @param type parameter
* @param data a given list of numbers
* @return a random number from a given list of numbers
*/
public static T random(List data) {
int idx = uniform(data.size());
return data.get(idx);
}
/**
* A random double array with values in [0, 1)
*
* @param size the size of random array
* @return a random double array with values in [0, 1)
*/
public static double[] doubles(int size) {
double[] array = new double[size];
for (int i = 0; i < size; i++)
array[i] = random();
return array;
}
/**
* A random double array with values in [min, max).
*
* @param min minimum
* @param max maximum
* @param size the size of random array
* @return a random double array with values in [min, max)
*/
public static double[] doubles(double min, double max, int size) {
double[] array = new double[size];
for (int i = 0; i < size; i++)
array[i] = uniform(min, max);
return array;
}
/**
* Random (uniformly distributed) double in [0, 1)
*
* @return Random (uniformly distributed) double in [0, 1)
*/
public static double uniform() {
return uniform(0.0, 1.0);
}
/**
* Random (uniformly distributed) double in [min, max)
*
* @param max max of the range
* @param min min of the range
* @return Random (uniformly distributed) double in [min, max)
*/
public static double uniform(double min, double max) {
return min + (max - min) * r.nextDouble();
}
/**
* Return a boolean, which is true with probability p, and false otherwise.
*
* @param p probability p
* @return a boolean, which is true with probability p, and false otherwise.
*/
public static boolean bernoulli(double p) {
return uniform() < p;
}
/**
* Return a boolean, which is true with probability .5, and false otherwise.
*
* @return a boolean, which is true with probability .5, and false otherwise.
*/
public static boolean bernoulli() {
return bernoulli(0.5);
}
/**
* Return a real number from a Gaussian distribution with given mean and stddev.
*
* @param mu mean
* @param sigma stddev
* @return a real number from a Gaussian distribution with given mean and stddev
*/
public static double gaussian(double mu, double sigma) {
return mu + sigma * r.nextGaussian();
}
/**
* Randomly sample 1 point from Gamma Distribution with the given parameters. The code is from Mahout
* (http://mahout.apache.org/), available under Apache 2 license.
*
* @param alpha alpha parameter for Gamma Distribution.
* @param scale scale parameter for Gamma Distribution.
* @return a sample point randomly drawn from the given distribution.
*/
public static double gamma(double alpha, double scale) {
double rate = 1 / scale;
if (alpha <= 0.0 || rate <= 0.0) {
throw new IllegalArgumentException();
}
double gds;
double b = 0.0;
// CASE A: Acceptance rejection algorithm gs
if (alpha < 1.0) {
b = 1.0 + 0.36788794412 * alpha; // Step 1
while (true) {
double p = b * r.nextDouble();
// Step 2. Case gds <= 1
if (p <= 1.0) {
gds = Math.exp(Math.log(p) / alpha);
if (Math.log(r.nextDouble()) <= -gds) {
return gds / rate;
}
}
// Step 3. Case gds > 1
else {
gds = -Math.log((b - p) / alpha);
if (Math.log(r.nextDouble()) <= ((alpha - 1.0) * Math.log(gds))) {
return gds / rate;
}
}
}
}
// CASE B: Acceptance complement algorithm gd (gaussian distribution,
// box muller transformation)
else {
double ss = 0.0;
double s = 0.0;
double d = 0.0;
// Step 1. Preparations
if (alpha != -1.0) {
ss = alpha - 0.5;
s = Math.sqrt(ss);
d = 5.656854249 - 12.0 * s;
}
// Step 2. Normal deviate
double v12;
double v1;
do {
v1 = 2.0 * r.nextDouble() - 1.0;
double v2 = 2.0 * r.nextDouble() - 1.0;
v12 = v1 * v1 + v2 * v2;
} while (v12 > 1.0);
double t = v1 * Math.sqrt(-2.0 * Math.log(v12) / v12);
double x = s + 0.5 * t;
gds = x * x;
if (t >= 0.0) { // Immediate acceptance
return gds / rate;
}
double u = r.nextDouble();
if (d * u <= t * t * t) { // Squeeze acceptance
return gds / rate;
}
double q0 = 0.0;
double si = 0.0;
double c = 0.0;
// Step 4. Set-up for hat case
if (alpha != -1.0) {
double rr = 1.0 / alpha;
double q9 = 0.0001710320;
double q8 = -0.0004701849;
double q7 = 0.0006053049;
double q6 = 0.0003340332;
double q5 = -0.0003349403;
double q4 = 0.0015746717;
double q3 = 0.0079849875;
double q2 = 0.0208333723;
double q1 = 0.0416666664;
q0 = ((((((((q9 * rr + q8) * rr + q7) * rr + q6) * rr + q5) * rr + q4) * rr + q3) * rr + q2) * rr + q1)
* rr;
if (alpha > 3.686) {
if (alpha > 13.022) {
b = 1.77;
si = 0.75;
c = 0.1515 / s;
} else {
b = 1.654 + 0.0076 * ss;
si = 1.68 / s + 0.275;
c = 0.062 / s + 0.024;
}
} else {
b = 0.463 + s - 0.178 * ss;
si = 1.235;
c = 0.195 / s - 0.079 + 0.016 * s;
}
}
double v, q;
double a9 = 0.104089866;
double a8 = -0.112750886;
double a7 = 0.110368310;
double a6 = -0.124385581;
double a5 = 0.142873973;
double a4 = -0.166677482;
double a3 = 0.199999867;
double a2 = -0.249999949;
double a1 = 0.333333333;
// Step 5. Calculation of q
if (x > 0.0) {
// Step 6.
v = t / (s + s);
if (Math.abs(v) > 0.25) {
q = q0 - s * t + 0.25 * t * t + (ss + ss) * Math.log(1.0 + v);
}
// Step 7. Quotient acceptance
else {
q = q0
+ 0.5
* t
* t
* ((((((((a9 * v + a8) * v + a7) * v + a6) * v + a5) * v + a4) * v + a3) * v + a2) * v + a1)
* v;
}
if (Math.log(1.0 - u) <= q) {
return gds / rate;
}
}
double e7 = 0.000247453;
double e6 = 0.001353826;
double e5 = 0.008345522;
double e4 = 0.041664508;
double e3 = 0.166666848;
double e2 = 0.499999994;
double e1 = 1.000000000;
// Step 8. Double exponential deviate t
while (true) {
double sign_u;
double e;
do { // Step 9. Rejection of t
e = -Math.log(r.nextDouble());
u = r.nextDouble();
u = u + u - 1.0;
sign_u = (u > 0) ? 1.0 : -1.0;
t = b + (e * si) * sign_u;
} while (t <= -0.71874483771719);
// Step 10. New q(t)
v = t / (s + s);
if (Math.abs(v) > 0.25) {
q = q0 - s * t + 0.25 * t * t + (ss + ss) * Math.log(1.0 + v);
} else {
q = q0
+ 0.5
* t
* t
* ((((((((a9 * v + a8) * v + a7) * v + a6) * v + a5) * v + a4) * v + a3) * v + a2) * v + a1)
* v;
}
// Step 11.
if (q <= 0.0) {
continue;
}
// Step 12. Hat acceptance
double w;
if (q > 0.5) {
w = Math.exp(q) - 1.0;
} else {
w = ((((((e7 * q + e6) * q + e5) * q + e4) * q + e3) * q + e2) * q + e1) * q;
}
if (c * u * sign_u <= w * Math.exp(e - 0.5 * t * t)) {
x = s + 0.5 * t;
return x * x / rate;
}
}
}
}
/**
* Randomly sample a matrix from Wishart Distribution with the given parameters.
*
* @param scale scale parameter for Wishart Distribution.
* @param df degree of freedom for Wishart Distribution.
* @return the sample randomly drawn from the given distribution.
* @throws LibrecException if error occurs
*/
public static DenseMatrix wishart(DenseMatrix scale, double df) throws LibrecException {
DenseMatrix A = scale.cholesky();
if (A == null)
return null;
int p = scale.numRows();
DenseMatrix z = new DenseMatrix(p, p);
for (int i = 0; i < p; i++) {
for (int j = 0; j < p; j++) {
z.set(i, j, Randoms.gaussian(0, 1));
}
}
SparseVector y = new SparseVector(p);
for (int i = 0; i < p; i++)
y.set(i, Randoms.gamma((df - (i + 1)) / 2, 2));
DenseMatrix B = new DenseMatrix(p, p);
B.set(0, 0, y.get(0));
if (p > 1) {
// rest of diagonal:
for (int j = 1; j < p; j++) {
SparseVector zz = new SparseVector(j);
for (int k = 0; k < j; k++)
zz.set(k, z.get(k, j));
B.set(j, j, y.get(j) + zz.inner(zz));
}
// first row and column:
for (int j = 1; j < p; j++) {
B.set(0, j, z.get(0, j) * Math.sqrt(y.get(0)));
B.set(j, 0, B.get(0, j)); // mirror
}
}
if (p > 2) {
for (int j = 2; j < p; j++) {
for (int i = 1; i <= j - 1; i++) {
SparseVector zki = new SparseVector(i);
SparseVector zkj = new SparseVector(i);
for (int k = 0; k <= i - 1; k++) {
zki.set(k, z.get(k, i));
zkj.set(k, z.get(k, j));
}
B.set(i, j, z.get(i, j) * Math.sqrt(y.get(i)) + zki.inner(zkj));
B.set(j, i, B.get(i, j)); // mirror
}
}
}
return A.transpose().mult(B).mult(A);
}
/**
* Return an integer with a Poisson distribution with mean lambda.
*
* @param lambda mean lambda
* @return an integer with a Poisson distribution with mean lambda
*/
public static int poisson(double lambda) {
// using algorithm given by Knuth
// see http://en.wikipedia.org/wiki/Poisson_distribution
int k = 0;
double p = 1.0;
double L = Math.exp(-lambda);
do {
k++;
p *= uniform();
} while (p >= L);
return k - 1;
}
/**
* Return a real number with a Pareto distribution with parameter alpha.
*
* @param alpha parameter alpha
* @return a real number with a Pareto distribution with parameter alpha.
*/
public static double pareto(double alpha) {
return Math.pow(1 - uniform(), -1.0 / alpha) - 1.0;
}
/**
* Return a real number with a Cauchy distribution.
*
* @return a real number with a Cauchy distribution
*/
public static double cauchy() {
return Math.tan(Math.PI * (uniform() - 0.5));
}
/**
* Return a number from a discrete distribution: i with probability a[i]. Precondition: array entries are
* nonnegative and their sum (very nearly) equals 1.0.
*
* @param a probability a[i]
* @return a number from the discrete distribution
*/
public static int discrete(double[] a) {
double EPSILON = 1E-6;
double sum = 0.0;
for (int i = 0; i < a.length; i++) {
if (a[i] < 0.0)
throw new IllegalArgumentException("array entry " + i + " is negative: " + a[i]);
sum = sum + a[i];
}
if (sum > 1.0 + EPSILON || sum < 1.0 - EPSILON)
throw new IllegalArgumentException("sum of array entries not equal to one: " + sum);
// the for loop may not return a value when both r is (nearly) 1.0 and when the cumulative sum is less than 1.0 (as a result of floating-point roundoff error)
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;
}
}
}
/**
* Return a real number from an exponential distribution with rate lambda.
*
* @param lambda rate lambda
* @return a real number from an exponential distribution with rate lambda.
*/
public static double exp(double lambda) {
return -Math.log(1 - uniform()) / lambda;
}
/**
* generate next random integer in a range besides exceptions
*
* @param range the range located of next integer
* @param exceptions the exception values when generating integers, sorted first
* @return next no-repeated random integer
*/
public static int nextInt(int range, int... exceptions) {
return nextInt(0, range, exceptions);
}
/**
* generate next random integer in a range [min, max) besides exceptions
*
* @param min the minimum of range
* @param max the maximum of range
* @param exceptions the exception values when generating integers, sorted first
* @return next no-repeated random integer
*/
public static int nextInt(int min, int max, int... exceptions) {
int next;
while (true) {
next = min + r.nextInt(max - min);
if (exceptions != null && exceptions.length > 0 && Arrays.binarySearch(exceptions, next) >= 0) {
continue;
}
if (_tempList.contains(next))
continue;
else {
_tempList.add(next);
break;
}
}
return next;
}
/**
* Generate no repeat {@code size} indexes from {@code min} to {@code max}
*
* @param min min of the range
* @param max max of the range
* @param size size of the index array
* @return no repeat {@code size} indexes from {@code min} to {@code max}
*/
public static int[] indexs(int size, int min, int max) {
Set used = new HashSet();
int[] index = new int[size];
for (int i = 0; i < index.length; i++) {
while (true) {
int ind = uniform(min, max);
if (!used.contains(ind)) {
index[i] = ind;
used.add(ind);
break;
}
}
}
return index;
}
public static void clearCache() {
_tempList.clear();
}
/**
* generate next integers array with no repeated elements
*
* @param length the length of the array
* @param range the index range of the array, default [0, range)
* @return ascending sorted integer array
* @throws Exception if the range is less than length, an exception will be thrown
*/
public static int[] nextIntArray(int length, int range) throws Exception {
return nextIntArray(length, 0, range, null);
}
public static int[] nextIntArray(int length, int range, int... exceptions) throws Exception {
return nextIntArray(length, 0, range, exceptions);
}
public static int[] nextIntArray(int length, int min, int max) throws Exception {
return nextIntArray(length, min, max, null);
}
public static int[] nextIntArray(int length, int min, int max, int... exceptions) throws Exception {
int maxLen = max - min; // because max itself is not counted
if (maxLen < length)
throw new Exception("The range is less than legth");
int[] index = new int[length];
if (maxLen == length) {
for (int i = 0; i < length; i++)
index[i] = min + i;
} else {
Randoms.clearCache();
for (int i = 0; i < index.length; i++)
index[i] = Randoms.nextInt(min, max, exceptions);
Arrays.sort(index);
}
return index;
}
/**
* Generate a set of random (unique) integers in the range [min, max) with length {@code length}
*
* @param max max of the range
* @param min min of the range
* @param length length of the List
* @return a set of unique integers
* @throws Exception if error occurs
*/
public static List randInts(int length, int min, int max) throws Exception {
int len = max - min;
if (len < length)
throw new Exception("The range is less than legth");
Set ints = new HashSet();
while (true) {
int rand = min + r.nextInt(max - min);
ints.add(rand);
if (ints.size() >= length)
break;
}
List res = new ArrayList(ints);
Collections.sort(res);
return res;
}
/**
* Get a normalize array of probabilities
*
* @param size array size
* @return a normalize array of probabilities
*/
public static double[] randProbs(int size) {
if (size < 1)
throw new IllegalArgumentException("The size param must be greate than zero");
double[] pros = new double[size];
int sum = 0;
for (int i = 0; i < pros.length; i++) {
//avoid zero
pros[i] = r.nextInt(size) + 1;
sum += pros[i];
}
//normalize
for (int i = 0; i < pros.length; i++) {
pros[i] = pros[i] / sum;
}
return pros;
}
public static int[] ints(int range, int size) {
int[] data = new int[size];
for (int i = 0; i < size; i++)
data[i] = uniform(range);
return data;
}
public static int[] ints(int min, int max, int size) {
int[] data = new int[size];
for (int i = 0; i < size; i++)
data[i] = uniform(min, max);
return data;
}
public static List list(int size) {
return list(size, 0, 1, false);
}
public static List list(int size, int min, int max) {
return list(size, min, max, false);
}
public static List list(int size, int min, int max, boolean isInteger) {
List list = new ArrayList(size);
for (int i = 0; i < size; i++) {
if (isInteger)
list.add(uniform(min, max) + 0.0);
else
list.add(uniform(min + 0.0, max + 0.0));
}
return list;
}
/**
* Generate a permutation from min to max
*
* @param min the minimum value
* @param max the maximum value
* @return a permutation
*/
public static List permute(int min, int max) {
List list = new ArrayList();
int len = max - min + 1;
for (int i = 0; i < len; i++) {
while (true) {
int index = uniform(min, max + 1);
if (!list.contains(index)) {
list.add(index);
break;
}
}
}
return list;
}
}
