SSJ
V. 2.6.

umontreal.iro.lecuyer.gof
Class GofStat

java.lang.Object
  extended by umontreal.iro.lecuyer.gof.GofStat

public class GofStat
extends Object

This class provides methods to compute several types of EDF goodness-of-fit test statistics and to apply certain transformations to a set of observations. This includes the probability integral transformation Ui = F(Xi), as well as the power ratio and iterated spacings transformations. Here, U(0),..., U(n-1) stand for n observations U0,..., Un-1 sorted by increasing order, where 0 <= Ui <= 1.

Note: This class uses the Colt library.


Nested Class Summary
static class GofStat.OutcomeCategoriesChi2
          This class helps managing the partitions of possible outcomes into categories for applying chi-square tests.
 
Field Summary
static double EPSILONAD
           
 
Method Summary
static double andersonDarling(double[] sortedData)
          Computes and returns the Anderson-Darling statistic An2.
static double[] andersonDarling(double[] data, ContinuousDistribution dist)
          Computes the Anderson-Darling statistic An2 and the corresponding p-value p.
static double andersonDarling(cern.colt.list.DoubleArrayList sortedData)
          Computes and returns the Anderson-Darling statistic An2 (see method andersonDarling).
static double chi2(double[] nbExp, int[] count, int smin, int smax)
          Computes and returns the chi-square statistic for the observations oi in count[smin...smax], for which the corresponding expected values ei are in nbExp[smin...smax].
static double chi2(GofStat.OutcomeCategoriesChi2 cat, int[] count)
          Computes and returns the chi-square statistic for the observations oi in count, for which the corresponding expected values ei are in cat.
static double chi2(cern.colt.list.IntArrayList data, DiscreteDistributionInt dist, int smin, int smax, double minExp, int[] numCat)
          Computes and returns the chi-square statistic for the observations stored in data, assuming that these observations follow the discrete distribution dist.
static double chi2Equal(cern.colt.list.DoubleArrayList data)
          Equivalent to chi2Equal (data, 10).
static double chi2Equal(cern.colt.list.DoubleArrayList data, double minExp)
          Computes the chi-square statistic for a continuous distribution.
static double chi2Equal(double nbExp, int[] count, int smin, int smax)
          Similar to chi2, except that the expected number of observations per category is assumed to be the same for all categories, and equal to nbExp.
static double cramerVonMises(cern.colt.list.DoubleArrayList sortedData)
          Computes and returns the Cramér-von Mises statistic Wn2.
static void diff(cern.colt.list.DoubleArrayList sortedData, cern.colt.list.DoubleArrayList spacings, int n1, int n2, double a, double b)
          Same as method diff(IntArrayList,IntArrayList,int,int,int,int), but for the continuous case.
static void diff(cern.colt.list.IntArrayList sortedData, cern.colt.list.IntArrayList spacings, int n1, int n2, int a, int b)
          Assumes that the real-valued observations U0,..., Un-1 contained in sortedData are already sorted in increasing order and computes the differences between the successive observations.
static void iterateSpacings(cern.colt.list.DoubleArrayList data, cern.colt.list.DoubleArrayList spacings)
          Applies one iteration of the iterated spacings transformation.
static double[] kolmogorovSmirnov(double[] sortedData)
          Computes the Kolmogorov-Smirnov (KS) test statistics Dn+, Dn-, and Dn (see method kolmogorovSmirnov).
static void kolmogorovSmirnov(double[] data, ContinuousDistribution dist, double[] sval, double[] pval)
          Computes the KolmogorovSmirnov (KS) test statistics and their p-values.
static double[] kolmogorovSmirnov(cern.colt.list.DoubleArrayList sortedData)
          Computes the Kolmogorov-Smirnov (KS) test statistics Dn+, Dn-, and Dn.
static double[] kolmogorovSmirnovJumpOne(cern.colt.list.DoubleArrayList sortedData, double a)
          Compute the KS statistics Dn+(a) and Dn-(a) defined in the description of the method FDist.kolmogorovSmirnovPlusJumpOne, assuming that F is the uniform distribution over [0, 1] and that U(1),..., U(n) are in sortedData.
static double pDisc(double pL, double pR)
          Computes a variant of the p-value p whenever a test statistic has a discrete probability distribution.
static void powerRatios(cern.colt.list.DoubleArrayList sortedData)
          Applies the power ratios transformation W.
static int scan(cern.colt.list.DoubleArrayList sortedData, double d)
          Computes and returns the scan statistic Sn(d ), defined in FBar.scan.
static cern.colt.list.DoubleArrayList unifTransform(cern.colt.list.DoubleArrayList data, ContinuousDistribution dist)
          Applies the probability integral transformation Ui = F(Vi) for i = 0, 1,…, n - 1, where F is a continuous distribution function, and returns the result as an array of length n.
static cern.colt.list.DoubleArrayList unifTransform(cern.colt.list.DoubleArrayList data, DiscreteDistribution dist)
          Applies the transformation Ui = F(Vi) for i = 0, 1,…, n - 1, where F is a discrete distribution function, and returns the result as an array of length n.
static double watsonG(cern.colt.list.DoubleArrayList sortedData)
          Computes and returns the Watson statistic Gn.
static double watsonU(cern.colt.list.DoubleArrayList sortedData)
          Computes and returns the Watson statistic Un2.
 
Methods inherited from class java.lang.Object
equals, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait
 

Field Detail

EPSILONAD

public static double EPSILONAD
Method Detail

unifTransform

public static cern.colt.list.DoubleArrayList unifTransform(cern.colt.list.DoubleArrayList data,
                                                           ContinuousDistribution dist)
Applies the probability integral transformation Ui = F(Vi) for i = 0, 1,…, n - 1, where F is a continuous distribution function, and returns the result as an array of length n. V represents the n observations contained in data, and U, the returned transformed observations. If data contains random variables from the distribution function dist, then the result will contain uniform random variables over [0, 1].

Parameters:
data - array of observations to be transformed
dist - assumed distribution of the observations
Returns:
the array of transformed observations

unifTransform

public static cern.colt.list.DoubleArrayList unifTransform(cern.colt.list.DoubleArrayList data,
                                                           DiscreteDistribution dist)
Applies the transformation Ui = F(Vi) for i = 0, 1,…, n - 1, where F is a discrete distribution function, and returns the result as an array of length n. V represents the n observations contained in data, and U, the returned transformed observations.

Note: If V are the values of random variables with distribution function dist, then the result will contain the values of discrete random variables distributed over the set of values taken by dist, not uniform random variables over [0, 1].

Parameters:
data - array of observations to be transformed
dist - assumed distribution of the observations
Returns:
the array of transformed observations

diff

public static void diff(cern.colt.list.IntArrayList sortedData,
                        cern.colt.list.IntArrayList spacings,
                        int n1,
                        int n2,
                        int a,
                        int b)
Assumes that the real-valued observations U0,..., Un-1 contained in sortedData are already sorted in increasing order and computes the differences between the successive observations. Let D be the differences returned in spacings. The difference Ui - Ui-1 is put in Di for n1 < i <= n2, whereas Un1 - a is put into Dn1 and b - Un2 is put into Dn2+1. The number of observations must be greater or equal than n2, we must have n1 < n2, and n1 and n2 are greater than 0. The size of spacings will be at least n + 1 after the call returns.

Parameters:
sortedData - array of sorted observations
spacings - pointer to an array object that will be filled with spacings
n1 - starting index, in sortedData, of the processed observations
n2 - ending index, in sortedData of the processed observations
a - minimum value of the observations
b - maximum value of the observations

diff

public static void diff(cern.colt.list.DoubleArrayList sortedData,
                        cern.colt.list.DoubleArrayList spacings,
                        int n1,
                        int n2,
                        double a,
                        double b)
Same as method diff(IntArrayList,IntArrayList,int,int,int,int), but for the continuous case.

Parameters:
sortedData - array of sorted observations
spacings - pointer to an array object that will be filled with spacings
n1 - starting index of the processed observations in sortedData
n2 - ending index, in sortedData of the processed observations
a - minimum value of the observations
b - maximum value of the observations

iterateSpacings

public static void iterateSpacings(cern.colt.list.DoubleArrayList data,
                                   cern.colt.list.DoubleArrayList spacings)
Applies one iteration of the iterated spacings transformation. Let U be the n observations contained into data, and let S be the spacings contained into spacings, Assumes that S[0..n] contains the spacings between n real numbers U0,..., Un-1 in the interval [0, 1]. These spacings are defined by

Si = U(i) - U(i-1),        1 <= i < n,

where U(0) = 0, U(n-1) = 1, and U(0),..., U(n-1), are the Ui sorted in increasing order. These spacings may have been obtained by calling diff. This method transforms the spacings into new spacings: it sorts S0,..., Sn to obtain S(0) <= S(1) <= S(2) <=  ...  <= S(n), computes the weighted differences
S0 = (n + 1)S(0),  
S1 = n(S(1) - S(0)),  
S2 = (n - 1)(S(2) - S(1)),  
  ...    
Sn = S(n) - S(n-1),  

and computes Vi = S0 + S1 + ... + Si for 0 <= i < n. It then returns S0,..., Sn in S[0..n] and V1,..., Vn in V[1..n].

Under the assumption that the Ui are i.i.d. U(0, 1), the new Si can be considered as a new set of spacings having the same distribution as the original spacings, and the Vi are a new sample of i.i.d. U(0, 1) random variables, sorted by increasing order.

This transformation is useful to detect clustering in a data set: A pair of observations that are close to each other is transformed into an observation close to zero. A data set with unusually clustered observations is thus transformed to a data set with an accumulation of observations near zero, which is easily detected by the Anderson-Darling GOF test.

Parameters:
data - array of observations
spacings - spacings between the observations, will be filled with the new spacings

powerRatios

public static void powerRatios(cern.colt.list.DoubleArrayList sortedData)
Applies the power ratios transformation W. Let U be the n observations contained into sortedData. Assumes that U contains n real numbers U(0),..., U(n-1) from the interval [0, 1], already sorted in increasing order, and computes the transformations:

U'i = (U(i)/U(i+1))i+1,        i = 0,..., n - 1,

with U(n) = 1. These U'i are sorted in increasing order and put back in U[1...n]. If the U(i) are i.i.d. U(0, 1) sorted by increasing order, then the U'i are also i.i.d. U(0, 1).

This transformation is useful to detect clustering, as explained in iterateSpacings, except that here a pair of observations close to each other is transformed into an observation close to 1. An accumulation of observations near 1 is also easily detected by the Anderson-Darling GOF test.

Parameters:
sortedData - sorted array of real-valued observations in the interval [0, 1] that will be overwritten with the transformed observations

chi2

public static double chi2(double[] nbExp,
                          int[] count,
                          int smin,
                          int smax)
Computes and returns the chi-square statistic for the observations oi in count[smin...smax], for which the corresponding expected values ei are in nbExp[smin...smax]. Assuming that i goes from 1 to k, where k = smax-smin+1 is the number of categories, the chi-square statistic is defined as

X2 = ∑i=1k(oi - ei)2/ei.

Under the hypothesis that the ei are the correct expectations and if these ei are large enough, X2 follows approximately the chi-square distribution with k - 1 degrees of freedom. If some of the ei are too small, one can use OutcomeCategoriesChi2 to regroup categories.

Parameters:
nbExp - numbers expected in each category
count - numbers observed in each category
smin - index of the first valid data in count and nbExp
smax - index of the last valid data in count and nbExp
Returns:
the X2 statistic

chi2

public static double chi2(GofStat.OutcomeCategoriesChi2 cat,
                          int[] count)
Computes and returns the chi-square statistic for the observations oi in count, for which the corresponding expected values ei are in cat. This assumes that cat.regroupCategories has been called before to regroup categories in order to make sure that the expected numbers in each category are large enough for the chi-square test.

Parameters:
cat - numbers expected in each category
count - numbers observed in each category
Returns:
the X2 statistic

chi2

public static double chi2(cern.colt.list.IntArrayList data,
                          DiscreteDistributionInt dist,
                          int smin,
                          int smax,
                          double minExp,
                          int[] numCat)
Computes and returns the chi-square statistic for the observations stored in data, assuming that these observations follow the discrete distribution dist. For dist, we assume that there is one set S = {a, a + 1,..., b - 1, b}, where a < b and a >=  0, for which p(s) > 0 if sS and p(s) = 0 otherwise.

Generally, it is not possible to divide the integers in intervals satisfying nP(a0 <= s < a1) = nP(a1 <= s < a2) = ... = nP(aj-1 <= s < aj) for a discrete distribution, where n is the sample size, i.e., the number of observations stored into data. To perform a general chi-square test, the method starts from smin and finds the first non-negligible probability p(s) >= ε, where ε = DiscreteDistributionInt.EPSILON. It uses smax to allocate an array storing the number of expected observations (np(s)) for each s >=  smin. Starting from s = smin, the np(s) terms are computed and the allocated array grows if required until a negligible probability term is found. This gives the number of expected elements for each category, where an outcome category corresponds here to an interval in which sample observations could lie. The categories are regrouped to have at least minExp observations per category. The method then counts the number of samples in each categories and calls chi2 to get the chi-square test statistic. If numCat is not null, the number of categories after regrouping is returned in numCat[0]. The number of degrees of freedom is equal to numCat[0]-1. We usually choose minExp = 10.

Parameters:
data - observations, not necessarily sorted
dist - assumed probability distribution
smin - estimated minimum value of s for which p(s) > 0
smax - estimated maximum value of s for which p(s) > 0
minExp - minimum number of expected observations in each interval
numCat - one-element array that will be filled with the number of categories after regrouping
Returns:
the chi-square statistic for a discrete distribution

chi2Equal

public static double chi2Equal(double nbExp,
                               int[] count,
                               int smin,
                               int smax)
Similar to chi2, except that the expected number of observations per category is assumed to be the same for all categories, and equal to nbExp.

Parameters:
nbExp - number of expected observations in each category (or interval)
count - number of counted observations in each category
smin - index of the first valid data in count and nbExp
smax - index of the last valid data in count and nbExp
Returns:
the X2 statistic

chi2Equal

public static double chi2Equal(cern.colt.list.DoubleArrayList data,
                               double minExp)
Computes the chi-square statistic for a continuous distribution. Here, the equiprobable case can be used. Assuming that data contains observations coming from the uniform distribution, the [0, 1] interval is divided into 1/p subintervals, where p = minExp/n, n being the sample size, i.e., the number of observations stored in data. For each subinterval, the method counts the number of contained observations and the chi-square statistic is computed using chi2Equal. We usually choose minExp = 10.

Parameters:
data - array of observations in [0, 1)
minExp - minimum number of expected observations in each subintervals
Returns:
the chi-square statistic for a continuous distribution

chi2Equal

public static double chi2Equal(cern.colt.list.DoubleArrayList data)
Equivalent to chi2Equal (data, 10).

Parameters:
data - array of observations in [0, 1)
Returns:
the chi-square statistic for a continuous distribution

scan

public static int scan(cern.colt.list.DoubleArrayList sortedData,
                       double d)
Computes and returns the scan statistic Sn(d ), defined in FBar.scan. Let U be the n observations contained into sortedData. The n observations in U[0..n - 1] must be real numbers in the interval [0, 1], sorted in increasing order. (See FBar.scan for the distribution function of Sn(d )).

Parameters:
sortedData - sorted array of real-valued observations in the interval [0, 1]
d - length of the test interval (∈(0, 1))
Returns:
the scan statistic

cramerVonMises

public static double cramerVonMises(cern.colt.list.DoubleArrayList sortedData)
Computes and returns the Cramér-von Mises statistic Wn2. It is defined by

Wn2 = 1/(12n) + ∑j=0n-1(U(j) - (j + 0.5)/n)2,

assuming that sortedData contains U(0),..., U(n-1) sorted in increasing order.

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the Cramér-von Mises statistic

watsonG

public static double watsonG(cern.colt.list.DoubleArrayList sortedData)
Computes and returns the Watson statistic Gn. It is defined by
Gn = (n)1/2max0 <= j <= n-1{(j + 1)/n - U(j) + bar(U)n -1/2}  
  = (n)1/2(Dn+ + bar(U)n - 1/2),  

where bar(U)n is the average of the observations U(j), assuming that sortedData contains the sorted U(0),..., U(n-1).

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the Watson statistic Gn

watsonU

public static double watsonU(cern.colt.list.DoubleArrayList sortedData)
Computes and returns the Watson statistic Un2. It is defined by
Wn2 = 1/(12n) + ∑j=0n-1{U(j) - (j + 0.5)/n}2,  
Un2 = Wn2 - n(bar(U)n -1/2)2.

 

where bar(U)n is the average of the observations U(j), assuming that sortedData contains the sorted U(0),..., U(n-1).

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the Watson statistic Un2

andersonDarling

public static double andersonDarling(cern.colt.list.DoubleArrayList sortedData)
Computes and returns the Anderson-Darling statistic An2 (see method andersonDarling).

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the Anderson-Darling statistic

andersonDarling

public static double andersonDarling(double[] sortedData)
Computes and returns the Anderson-Darling statistic An2. It is defined by
An2 = - n - 1/n    ∑j=0n-1{(2j + 1)ln(U(j)) + (2n - 1 - 2j)ln(1 - U(j))},  

assuming that sortedData contains U(0),..., U(n-1) sorted in increasing order.

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the Anderson-Darling statistic

andersonDarling

public static double[] andersonDarling(double[] data,
                                       ContinuousDistribution dist)
Computes the Anderson-Darling statistic An2 and the corresponding p-value p. The n (unsorted) observations in data are assumed to be independent and to come from the continuous distribution dist. Returns the 2-elements array [An2, p].

Parameters:
data - array of observations
dist - assumed distribution of the observations
Returns:
the array [An2, p].

kolmogorovSmirnov

public static double[] kolmogorovSmirnov(double[] sortedData)
Computes the Kolmogorov-Smirnov (KS) test statistics Dn+, Dn-, and Dn (see method kolmogorovSmirnov). Returns the array [Dn+, Dn-, Dn].

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the array [Dn+, Dn-, Dn]

kolmogorovSmirnov

public static double[] kolmogorovSmirnov(cern.colt.list.DoubleArrayList sortedData)
Computes the Kolmogorov-Smirnov (KS) test statistics Dn+, Dn-, and Dn. It is defined by
Dn+ = max0 <= j <= n-1((j + 1)/n - U(j)),  
Dn- = max0 <= j <= n-1(U(j) - j/n),  
Dn = max (Dn+, Dn-).  

and returns an array of length 3 that contains [Dn+, Dn-, Dn]. These statistics compare the empirical distribution of U(1),..., U(n), which are assumed to be in sortedData, with the uniform distribution over [0, 1].

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
Returns:
the array [Dn+, Dn-, Dn]

kolmogorovSmirnov

public static void kolmogorovSmirnov(double[] data,
                                     ContinuousDistribution dist,
                                     double[] sval,
                                     double[] pval)
Computes the KolmogorovSmirnov (KS) test statistics and their p-values. This is to compare the empirical distribution of the (unsorted) observations in data with the theoretical distribution dist. The KS statistics Dn+, Dn- and Dn are returned in sval[0], sval[1], and sval[2] respectively, and their corresponding p-values are returned in pval[0], pval[1], and pval[2].

Parameters:
data - array of observations to be tested
dist - assumed distribution of the observations
sval - values of the 3 KS statistics
pval - p-values for the 3 KS statistics

kolmogorovSmirnovJumpOne

public static double[] kolmogorovSmirnovJumpOne(cern.colt.list.DoubleArrayList sortedData,
                                                double a)
Compute the KS statistics Dn+(a) and Dn-(a) defined in the description of the method FDist.kolmogorovSmirnovPlusJumpOne, assuming that F is the uniform distribution over [0, 1] and that U(1),..., U(n) are in sortedData. Returns the array [Dn+, Dn-].

Parameters:
sortedData - array of sorted real-valued observations in the interval [0, 1]
a - size of the jump
Returns:
the array [Dn+, Dn-]

pDisc

public static double pDisc(double pL,
                           double pR)
Computes a variant of the p-value p whenever a test statistic has a discrete probability distribution. This p-value is defined as follows:
pL = P[Y <= y]  
pR = P[Y >= y]  

p = pR,         if pR < pL,
p = 1 - pL,         if pR >= pL and pL < 0.5,
p = 0.5         otherwise.

The function takes pL and pR as input and returns p.

Parameters:
pL - left p-value
pR - right p-value
Returns:
the p-value for a test on a discrete distribution

SSJ
V. 2.6.

To submit a bug or ask questions, send an e-mail to Pierre L'Ecuyer.