package umontreal.iro.lecuyer.probdist;
import java.util.Formatter;
import java.util.Locale;
import java.util.Arrays;
import java.io.IOException;
import java.io.Reader;
import java.io.BufferedReader;
import umontreal.iro.lecuyer.util.Misc;
import umontreal.iro.lecuyer.util.PrintfFormat;
/**
* Extends {@link DiscreteDistribution} to an empirical
* distribution function,
* based on the observations
* X(1),..., X(n) (sorted by increasing order).
* The distribution is uniform over the n observations, so the
* distribution function has a jump of 1/n at each of the n observations.
*
*/
public class EmpiricalDist extends DiscreteDistribution {
private double[] sortedObs;
private int n = 0;
private double sampleMean;
private double sampleVariance;
private double sampleStandardDeviation;
/**
* Constructs a new empirical distribution using
* all the observations stored in obs, and which are assumed
* to have been sorted in increasing numerical order. 1 These observations are copied into an internal array.
*
*/
public EmpiricalDist (double[] obs) {
if (obs.length <= 1)
throw new IllegalArgumentException
("Two or more observations are needed");
n = obs.length;
sortedObs = new double[n];
System.arraycopy (obs, 0, sortedObs, 0, n);
init();
}
/**
* Constructs a new empirical distribution using
* the observations read from the reader in. This constructor
* will read the first double of each line in the stream.
* Any line that does not start with a +, -, or a decimal digit,
* is ignored. One must be careful about lines starting with a blank.
* This format is the same as in UNURAN. The observations read are assumed
* to have been sorted in increasing numerical order.
*
*/
public EmpiricalDist (Reader in) throws IOException {
BufferedReader inb = new BufferedReader (in);
double[] data = new double[5];
String li;
while ((li = inb.readLine()) != null) {
li = li.trim();
// look for the first non-digit character on the read line
int index = 0;
while (index < li.length() &&
(li.charAt (index) == '+' || li.charAt (index) == '-' ||
li.charAt (index) == 'e' || li.charAt (index) == 'E' ||
li.charAt (index) == '.' || Character.isDigit (li.charAt (index))))
++index;
// truncate the line
li = li.substring (0, index);
if (!li.equals ("")) {
try {
data[n++] = Double.parseDouble (li);
if (n >= data.length) {
double[] newData = new double[2*n];
System.arraycopy (data, 0, newData, 0, data.length);
data = newData;
}
}
catch (NumberFormatException nfe) {}
}
}
sortedObs = new double[n];
System.arraycopy (data, 0, sortedObs, 0, n);
init();
}
public double prob (int i) {
if (i >= 0 && i < n)
return 1.0 / n;
throw new IllegalStateException();
}
public double cdf (double x) {
if (x < sortedObs[0])
return 0;
if (x >= sortedObs[n-1])
return 1;
for (int i = 0; i < (n-1); i++) {
if (x >= sortedObs[i] && x < sortedObs[i+1])
return (double)(i + 1)/n;
}
throw new IllegalStateException();
}
public double barF (double x) {
if (x <= sortedObs[0])
return 1;
if (x > sortedObs[n-1])
return 0;
for (int i = 0; i < (n-1); i++) {
if (x > sortedObs[i] && x <= sortedObs[i+1])
return ((double)n-1-i)/n;
}
throw new IllegalStateException();
}
public double inverseF (double u) {
if (u < 0 || u > 1)
throw new IllegalArgumentException ("u is not in [0,1]");
if (u == 1.0)
return sortedObs[n-1];
int i = (int)Math.floor ((double)n * u);
return sortedObs[i];
}
private void init() {
// Arrays.sort (sortedObs);
// n = sortedObs.length;
double sum = 0.0;
for (int i = 0; i < sortedObs.length; i++) {
sum += sortedObs[i];
}
sampleMean = sum / n;
sum = 0.0;
for (int i = 0; i < n; i++) {
double coeff = (sortedObs[i] - sampleMean);
sum += coeff * coeff;
}
sampleVariance = sum / (n-1);
sampleStandardDeviation = Math.sqrt (sampleVariance);
supportA = sortedObs[0];
supportB = sortedObs[n-1];
}
public double getMean() {
return sampleMean;
}
public double getVariance() {
return sampleVariance;
}
public double getStandardDeviation() {
return sampleStandardDeviation;
}
/**
* Returns the
* n/2th item of the sorted observations when the number
* of items is odd, and the mean of the
* n/2th and the
*
* (n/2 + 1)th items when the number of items is even.
*
*/
public double getMedian () {
if ((n % 2) == 0)
return ((obs[n / 2 - 1] + obs[n / 2]) / 2.0);
else
return obs[(n - 1) / 2];
}
/**
* Returns the
* n/2th item of the array obs when the number
* of items is odd, and the mean of the
* n/2th and the
*
* (n/2 + 1)th items when the number of items is even.
* The array does not have to be sorted.
*
* @param obs the array of observations
*
* @param n the number of observations
*
* @return return the median of the observations
*
*/
public static double getMedian (double obs[], int n) {
if ((n % 2) == 0)
return ((Misc.quickSelect (obs, n, (n / 2 - 1)) +
Misc.quickSelect (obs, n, (n / 2))) / 2.0);
else
return Misc.quickSelect (obs, n, ((n - 1) / 2));
}
/**
* Returns n, the number of observations.
*
*/
public int getN() {
return n;
}
/**
* Returns the value of X(i).
*
*/
public double getObs (int i) {
return sortedObs[i];
}
/**
* Returns the sample mean of the observations.
*
*/
public double getSampleMean() {
return sampleMean;
}
/**
* Returns the sample variance of the observations.
*
*/
public double getSampleVariance() {
return sampleVariance;
}
/**
* Returns the sample standard deviation of the observations.
*
*/
public double getSampleStandardDeviation() {
return sampleStandardDeviation;
}
/**
* Returns the interquartile range of the observations,
* defined as the difference between the third and first quartiles.
*
*/
public double getInterQuartileRange() {
int j = n/2;
double lowerqrt=0, upperqrt=0;
if (j % 2 == 1) {
lowerqrt = sortedObs[(j+1)/2-1];
upperqrt = sortedObs[n-(j+1)/2];
}
else {
lowerqrt = 0.5 * (sortedObs[j/2-1] + sortedObs[j/2+1-1]);
upperqrt = 0.5 * (sortedObs[n-j/2] + sortedObs[n-j/2-1]);
}
double h =upperqrt - lowerqrt;
if (h < 0)
throw new IllegalStateException("Observations MUST be sorted");
return h;
}
/**
* Return a table containing parameters of the current distribution.
*
*/
public double[] getParams () {
double[] retour = new double[n];
System.arraycopy (sortedObs, 0, retour, 0, n);
return retour;
}
/**
* Returns a String containing information about the current distribution.
*
*/
public String toString () {
StringBuilder sb = new StringBuilder();
Formatter formatter = new Formatter(sb, Locale.US);
formatter.format(getClass().getSimpleName() + PrintfFormat.NEWLINE);
for(int i = 0; i
Ads help maintain this website.