001: package JSci.maths.statistics;
002:
003: import JSci.maths.*;
004:
005: /**
006: * The NormalDistribution class provides an object for encapsulating normal distributions.
007: * @version 1.1
008: * @author Jaco van Kooten
009: */
010: public final class NormalDistribution extends ProbabilityDistribution
011: implements NumericalConstants {
012: private double mean, variance;
013: private double pdfDenominator, cdfDenominator;
014:
015: /**
016: * Constructs the standard normal distribution (zero mean and unity variance).
017: */
018: public NormalDistribution() {
019: this (0.0, 1.0);
020: }
021:
022: /**
023: * Constructs a normal distribution.
024: * @param mu the mean.
025: * @param var the variance.
026: */
027: public NormalDistribution(double mu, double var) {
028: mean = mu;
029: if (var <= 0.0)
030: throw new OutOfRangeException(
031: "The variance should be (strictly) positive.");
032: variance = var;
033: pdfDenominator = SQRT2PI * Math.sqrt(variance);
034: cdfDenominator = SQRT2 * Math.sqrt(variance);
035: }
036:
037: /**
038: * Constructs a normal distribution from a data set.
039: * @param array a sample.
040: * @author Mark Hale
041: */
042: public NormalDistribution(double array[]) {
043: double sumX = array[0];
044: double sumX2 = array[0] * array[0];
045: for (int i = 1; i < array.length; i++) {
046: sumX += array[i];
047: sumX2 += array[i] * array[i];
048: }
049: mean = sumX / array.length;
050: variance = (sumX2 - array.length * mean * mean)
051: / (array.length - 1);
052: pdfDenominator = SQRT2PI * Math.sqrt(variance);
053: cdfDenominator = SQRT2 * Math.sqrt(variance);
054: }
055:
056: /**
057: * Returns the mean.
058: */
059: public double getMean() {
060: return mean;
061: }
062:
063: /**
064: * Returns the variance.
065: */
066: public double getVariance() {
067: return variance;
068: }
069:
070: /**
071: * Probability density function of a normal (Gaussian) distribution.
072: * @return the probability that a stochastic variable x has the value X, i.e. P(x=X).
073: */
074: public double probability(double X) {
075: return Math.exp(-(X - mean) * (X - mean) / (2 * variance))
076: / pdfDenominator;
077: }
078:
079: /**
080: * Cumulative normal distribution function.
081: * @return the probability that a stochastic variable x is less then X, i.e. P(x<X).
082: */
083: public double cumulative(double X) {
084: return SpecialMath.complementaryError(-(X - mean)
085: / cdfDenominator) / 2;
086: }
087:
088: /**
089: * Inverse of the cumulative normal distribution function.
090: * @return the value X for which P(x<X).
091: */
092: public double inverse(double probability) {
093: checkRange(probability);
094: if (probability == 0.0)
095: return -Double.MAX_VALUE;
096: if (probability == 1.0)
097: return Double.MAX_VALUE;
098: if (probability == 0.5)
099: return mean;
100: // To ensure numerical stability we need to rescale the distribution
101: double meanSave = mean, varSave = variance;
102: double pdfDSave = pdfDenominator, cdfDSave = cdfDenominator;
103: mean = 0.0;
104: variance = 1.0;
105: pdfDenominator = Math.sqrt(TWO_PI);
106: cdfDenominator = SQRT2;
107: double X = findRoot(probability, 0.0, -100.0, 100.0);
108: // Scale back
109: mean = meanSave;
110: variance = varSave;
111: pdfDenominator = pdfDSave;
112: cdfDenominator = cdfDSave;
113: return X * Math.sqrt(variance) + mean;
114: }
115: }
|