umontreal.iro.lecuyer.probdist.ChiDist Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of ssj Show documentation
Show all versions of ssj Show documentation
SSJ is a Java library for stochastic simulation, developed under the direction of Pierre L'Ecuyer,
in the Département d'Informatique et de Recherche Opérationnelle (DIRO), at the Université de Montréal.
It provides facilities for generating uniform and nonuniform random variates, computing different
measures related to probability distributions, performing goodness-of-fit tests, applying quasi-Monte
Carlo methods, collecting (elementary) statistics, and programming discrete-event simulations with both
events and processes.
The newest version!
/*
* Class: ChiDist
* Description: chi distribution
* Environment: Java
* Software: SSJ
* Copyright (C) 2001 Pierre L'Ecuyer and Université de Montréal
* Organization: DIRO, Université de Montréal
* @author
* @since
* SSJ is free software: you can redistribute it and/or modify it under
* the terms of the GNU General Public License (GPL) as published by the
* Free Software Foundation, either version 3 of the License, or
* any later version.
* SSJ 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.
* A copy of the GNU General Public License is available at
GPL licence site.
*/
package umontreal.iro.lecuyer.probdist;
import umontreal.iro.lecuyer.util.Num;
import umontreal.iro.lecuyer.functions.MathFunction;
/**
* Extends the class {@link ContinuousDistribution} for the chi
* distribution with shape parameter
* v > 0, where the number of degrees of freedom
* v is a positive integer.
* The density function is given by
*
*
*
* f (x) = e-x2/2xv-1/(2(v/2)-1Γ(v/2)) for x > 0,
*
* where
* Γ(x) is the gamma function defined in
* {@link GammaDist}.
* The distribution function is
*
*
*
* F(x) = 1/Γ(v/2)∫0x2/2tv/2-1e-t dt.
*
* It is equivalent to the gamma distribution function with parameters
*
* α = v/2 and λ = 1, evaluated at x2/2.
*
*/
public class ChiDist extends ContinuousDistribution {
private int nu;
private double C1;
private static class Function implements MathFunction {
protected int n;
protected double sum;
public Function (double s, int n)
{
this.n = n;
this.sum = s;
}
public double evaluate (double k)
{
if (k < 1.0) return 1.0e200;
return (sum + n * (Num.lnGamma (k / 2.0) - 0.5*(Num.LN2) - Num.lnGamma ((k + 1.0) / 2.0)));
}
}
/**
* Constructs a ChiDist object.
*
*/
public ChiDist (int nu) {
setNu (nu);
}
public double density (double x) {
if (x <= 0.0)
return 0.0;
return Math.exp ((nu - 1)*Math.log (x) - x*x/2.0 - C1);
}
public double cdf (double x) {
return cdf (nu, x);
}
public double barF (double x) {
return barF (nu, x);
}
public double inverseF (double u) {
return inverseF (nu, u);
}
public double getMean() {
return ChiDist.getMean (nu);
}
public double getVariance() {
return ChiDist.getVariance (nu);
}
public double getStandardDeviation() {
return ChiDist.getStandardDeviation (nu);
}
/**
* Computes the density function.
*
*/
public static double density (int nu, double x) {
if (nu <= 0)
throw new IllegalArgumentException ("nu <= 0");
if (x <= 0.0)
return 0.0;
return Math.exp ((nu - 1)*Math.log (x) - x*x/2.0
- (nu/2.0 - 1.0)*Num.LN2 - Num.lnGamma (nu/2.0));
}
/**
* Computes the distribution function by using the
* gamma distribution function.
*
*/
public static double cdf (int nu, double x) {
if (x <= 0.0)
return 0.0;
return GammaDist.cdf (nu/2.0, 15, x*x/2.0);
}
/**
* Computes the complementary distribution.
*
*/
public static double barF (int nu, double x) {
if (x <= 0.0)
return 1.0;
return GammaDist.barF (nu/2.0, 15, x*x/2.0);
}
/**
* Returns the inverse distribution function computed
* using the gamma inversion.
*
*/
public static double inverseF (int nu, double u) {
double res = GammaDist.inverseF (nu/2.0, 15, u);
return Math.sqrt (2*res);
}
/**
* Estimates the parameter ν of the chi distribution
* using the maximum likelihood method, from the n observations
* x[i],
* i = 0, 1,…, n - 1. The estimate is returned in element 0
* of the returned array.
*
* @param x the list of observations to use to evaluate parameters
*
* @param n the number of observations to use to evaluate parameters
*
* @return returns the parameter [hat(ν)]
*
*/
public static double[] getMLE (double[] x, int n) {
double[] parameters = new double[1];
double mean = 0.0;
for (int i = 0; i < n; i++)
mean += x[i];
mean /= (double) n;
double var = 0.0;
for (int i = 0; i < n; i++)
var += ((x[i] - mean) * (x[i] - mean));
var /= (double) n;
double k = Math.round (var + mean * mean) - 5.0;
if (k < 1.0)
k = 1.0;
double sum = 0.0;
for (int i = 0; i < n; i++) {
if (x[i] > 0.0)
sum += Math.log (x[i]);
else
sum -= 709.0;
}
Function f = new Function (sum, n);
while (f.evaluate(k) > 0.0)
k++;
parameters[0] = k;
return parameters;
}
/**
* Creates a new instance of a chi distribution with parameter ν estimated using
* the maximum likelihood method based on the n observations x[i],
*
* i = 0, 1,…, n - 1.
*
* @param x the list of observations to use to evaluate parameters
*
* @param n the number of observations to use to evaluate parameters
*
*
*/
public static ChiDist getInstanceFromMLE (double[] x, int n) {
double parameters[] = getMLE (x, n);
return new ChiDist ((int) parameters[0]);
}
/**
* Computes and returns the mean
* of the chi distribution with parameter ν.
*
* @return the mean of the chi distribution
*
* E[X] = (2)1/2Γ((ν +1)/2)/Γ(ν/2)
*
*/
public static double getMean (int nu) {
if (nu <= 0)
throw new IllegalArgumentException ("nu <= 0");
return Num.RAC2 * Num.gammaRatioHalf(nu / 2.0);
}
/**
* Computes and returns the variance
* of the chi distribution with parameter ν.
*
* @return the variance of the chi distribution
*
* Var[X] = 2[Γ(ν/2)Γ(1 + ν/2) - Γ2(1/2(ν +1))]/Γ(ν/2)
*
*/
public static double getVariance (int nu) {
if (nu <= 0)
throw new IllegalArgumentException ("nu <= 0");
double mean = ChiDist.getMean(nu);
return (nu - (mean * mean));
}
/**
* Computes and returns the standard deviation of the chi distribution
* with parameter ν.
*
* @return the standard deviation of the chi distribution
*
*/
public static double getStandardDeviation (int nu) {
return Math.sqrt (ChiDist.getVariance (nu));
}
/**
* Returns the value of ν for this object.
*
*/
public int getNu() {
return nu;
}
/**
* Sets the value of ν for this object.
*
*/
public void setNu (int nu) {
if (nu <= 0)
throw new IllegalArgumentException ("nu <= 0");
this.nu = nu;
supportA = 0.0;
C1 = (nu/2.0 - 1.0)*Num.LN2 + Num.lnGamma (nu/2.0);
}
/**
* Return a table containing parameters of the current distribution.
*
*
*/
public double[] getParams () {
double[] retour = {nu};
return retour;
}
public String toString () {
return getClass().getSimpleName() + " : nu = " + nu;
}
}