com.hfg.math.StandardNormalDistribution Maven / Gradle / Ivy
Go to download
Show more of this group Show more artifacts with this name
Show all versions of com_hfg Show documentation
Show all versions of com_hfg Show documentation
com.hfg xml, html, svg, and bioinformatics utility library
package com.hfg.math;
//------------------------------------------------------------------------------
/**
For calculations related to standard normal distributions.
@author J. Alex Taylor, hairyfatguy.com
*/
//------------------------------------------------------------------------------
// com.hfg XML/HTML Coding Library
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
//
// This library 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
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
//
// J. Alex Taylor, President, Founder, CEO, COO, CFO, OOPS hairyfatguy.com
// [email protected]
//------------------------------------------------------------------------------
public class StandardNormalDistribution
{
private Double mMean;
private Double mSampleStandardDeviation;
// Coefficients in rational approximations for the inverse normal cumulative distribution function
private static final double[] sA_coeff = new double[] {
-3.969683028665376e+01, 2.209460984245205e+02,
-2.759285104469687e+02, 1.383577518672690e+02,
-3.066479806614716e+01, 2.506628277459239e+00 };
private static final double[] sB_coeff = new double[] {
-5.447609879822406e+01, 1.615858368580409e+02,
-1.556989798598866e+02, 6.680131188771972e+01,
-1.328068155288572e+01 };
private static final double[] sC_coeff = new double[] {
-7.784894002430293e-03, -3.223964580411365e-01,
-2.400758277161838e+00, -2.549732539343734e+00,
4.374664141464968e+00, 2.938163982698783e+00 };
private static final double[] sD_coeff = new double[] {
7.784695709041462e-03, 3.224671290700398e-01,
2.445134137142996e+00, 3.754408661907416e+00 };
//--------------------------------------------------------------------------
public StandardNormalDistribution()
{
}
//--------------------------------------------------------------------------
public StandardNormalDistribution setMean(Double inValue)
{
mMean = inValue;
return this;
}
//--------------------------------------------------------------------------
public Double getMean()
{
return mMean;
}
//--------------------------------------------------------------------------
public StandardNormalDistribution setSampleStandardDeviation(Double inValue)
{
mSampleStandardDeviation = inValue;
return this;
}
//--------------------------------------------------------------------------
public Double getSampleStandardDeviation()
{
return mSampleStandardDeviation;
}
//--------------------------------------------------------------------------
public double getZScore(double inValue)
{
// Z = (X - μ) / σ
return (inValue - mMean) / mSampleStandardDeviation;
}
//--------------------------------------------------------------------------
public double getProbability(double inValue)
{
double zScore = getZScore(inValue);
double probability = 0;
if (zScore < -7)
{
probability = 0.0;
}
else if (zScore > 7)
{
probability = 1.0;
}
else
{
boolean invert = (zScore > 0.0);
zScore = Math.abs(zScore);
double b = 0.0;
double s = Math.sqrt(2) / 3 * zScore;
double HH = .5;
for (int i = 0; i < 12; i++)
{
double a = Math.exp(-HH * HH / 9) * Math.sin(HH * s) / HH;
b = b + a;
HH = HH + 1.0;
}
probability = .5 - b / Math.PI;
if (invert)
{
probability = 1.0 - probability;
}
}
return probability;
}
//--------------------------------------------------------------------------
/**
Given a probability P, it returns an approximation to the X satisfying P = Pr{Z <= X} where Z is a
random variable from the standard normal distribution. Uses the lower tail quantile for standard normal distribution function from
http://home.online.no/~pjacklam/notes/invnorm/.
*/
public double getValueThresholdForProbability(float inProbability)
{
double z = 0;
// Define break-points.
double pLow = 0.02425;
double pHigh = 1 - pLow;
if ( inProbability < pLow )
{
// Rational approximation for lower region:
double q = Math.sqrt(-2*Math.log(inProbability));
z = (((((sC_coeff[0]*q+sC_coeff[1])*q+sC_coeff[2])*q+sC_coeff[3])*q+sC_coeff[4])*q+sC_coeff[5]) /
((((sD_coeff[0]*q+sD_coeff[1])*q+sD_coeff[2])*q+sD_coeff[3])*q+1);
}
else if ( pHigh < inProbability )
{
// Rational approximation for upper region:
double q = Math.sqrt(-2*Math.log(1 - inProbability));
z = -(((((sC_coeff[0]*q+sC_coeff[1])*q+sC_coeff[2])*q+sC_coeff[3])*q+sC_coeff[4])*q+sC_coeff[5]) /
((((sD_coeff[0]*q+sD_coeff[1])*q+sD_coeff[2])*q+sD_coeff[3])*q+1);
}
else
{
// Rational approximation for central region:
double q = inProbability - 0.5;
double r = q * q;
z = (((((sA_coeff[0] * r + sA_coeff[1]) * r + sA_coeff[2]) * r + sA_coeff[3]) * r + sA_coeff[4]) * r + sA_coeff[5]) * q /
(((((sB_coeff[0] * r + sB_coeff[1]) * r + sB_coeff[2]) * r + sB_coeff[3]) * r + sB_coeff[4]) * r + 1);
}
// X = μ + Zσ
return getMean() + (z * getSampleStandardDeviation());
}
}