statistics.java
001 /*
002  * Java Genetic Algorithm Library (jenetics-1.6.0).
003  * Copyright (c) 2007-2014 Franz Wilhelmstötter
004  *
005  * Licensed under the Apache License, Version 2.0 (the "License");
006  * you may not use this file except in compliance with the License.
007  * You may obtain a copy of the License at
008  *
009  *      http://www.apache.org/licenses/LICENSE-2.0
010  *
011  * Unless required by applicable law or agreed to in writing, software
012  * distributed under the License is distributed on an "AS IS" BASIS,
013  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014  * See the License for the specific language governing permissions and
015  * limitations under the License.
016  *
017  * Author:
018  *    Franz Wilhelmstötter (franz.wilhelmstoetter@gmx.at)
019  */
020 package org.jenetics.internal.math;
021 
022 import static java.lang.Math.PI;
023 import static java.lang.Math.abs;
024 import static java.lang.Math.exp;
025 import static java.lang.Math.log;
026 import static java.lang.Math.sqrt;
027 
028 import org.jenetics.util.StaticObject;
029 
030 /**
031  * Some statistical special functions.
032  *
033  @author <a href="mailto:franz.wilhelmstoetter@gmx.at">Franz Wilhelmstötter</a>
034  @since 1.0
035  @version 1.4 &mdash; <em>$Date: 2014-02-15 $</em>
036  */
037 public final class statistics extends StaticObject {
038     private statistics() {}
039 
040 
041     /**
042      * Uses Lanczos approximation formula. See Numerical Recipes 6.1.
043      *
044      @param x the value to calculate the log gamma for.
045      @return the log gamma value.
046      */
047     public static double logGamma(final double x) {
048         final double ser = 1.0 76.18009173 /
049                             (x + 086.50532033 /
050                             (x + 124.01409822 /
051                             (x + 21.231739516 /
052                             (x + 30.00120858003 /
053                             (x + 40.00000536382 /
054                             (x + 5);
055 
056         return (x - 0.5)*log(x + 4.5(x + 4.5+ log(ser*sqrt(* PI));
057     }
058 
059     public static double gamma(final double x) {
060         return exp(logGamma(x));
061     }
062 
063     public static double Γ(final double x) {
064         return gamma(x);
065     }
066 
067 
068 
069     /**
070      * Return the <i>error function</i> of {@code z}. The fractional error of
071      * this implementation is less than 1.2E-7.
072      *
073      @param z the value to calculate the error function for.
074      @return the error function for {@code z}.
075      */
076     public static double erf(final double z) {
077         final double t = 1.0/(1.0 0.5*abs(z));
078 
079         // Horner's method
080         final double result = - t*exp(
081                 -z*z - 1.26551223 +
082                 t*1.00002368 +
083                 t*0.37409196 +
084                 t*0.09678418 +
085                 t*(-0.18628806 +
086                 t*0.27886807 +
087                 t*(-1.13520398 +
088                 t*1.48851587 +
089                 t*(-0.82215223 +
090                 t*(0.17087277))))))))));
091 
092         return z >= ? result : -result;
093     }
094 
095     /**
096      * Return φ(x), the standard Gaussian pdf.
097      *
098      @see #φ(double)
099      @param x the value to calculate φ for.
100      @return the φ value for x.
101      */
102     public static double phi(final double x) {
103         return exp(-x*x/2.0/ sqrt(2.0*PI);
104     }
105 
106     /**
107      * Return φ(x), the standard Gaussian pdf.
108      *
109      @see #phi(double)
110      @param x the value to calculate φ for.
111      @return the φ value for x.
112      */
113     public static double φ(final double x) {
114         return phi(x);
115     }
116 
117     /**
118      * Return φ(x, µ, σ), the standard Gaussian pdf with mean µ and stddev σ.
119      *
120      @see #phi(double, double, double)
121      @param x the value to calculate φ for.
122      @param mu the mean value.
123      @param sigma the stddev.
124      @return the φ value for x.
125      */
126     public static double phi(final double x, final double mu, final double sigma) {
127         return phi((x - mu)/sigma)/sigma;
128     }
129 
130     /**
131      * Return φ(x, µ, σ), the standard Gaussian pdf with mean µ and stddev σ.
132      *
133      @see #phi(double, double, double)
134      @param x the value to calculate φ for.
135      @param µ the mean value.
136      @param σ the stddev.
137      @return the φ value for x.
138      */
139     public static double φ(final double x, final double µ, final double σ) {
140         return phi(x, µ, σ);
141     }
142 
143     /**
144      * Return Φ(z), the standard Gaussian cdf using Taylor approximation.
145      *
146      @param z the value to calculate Φ for.
147      @return the Φ for value z.
148      */
149     public static double Phi(final double z) {
150         if (z < -8.0) {
151             return 0.0;
152         }
153         if (z >  8.0) {
154             return 1.0;
155         }
156 
157         double s = 0.0;
158         double t = z;
159         for (int i = 3; s + t != s; i += 2) {
160             s  = s + t;
161             t = t*z*z/i;
162         }
163         return 0.5 + s*phi(z);
164     }
165 
166     /**
167      * Return Φ(z), the standard Gaussian cdf using Taylor approximation.
168      *
169      @param z the value to calculate Φ for.
170      @return the Φ for value z.
171      */
172     public static double Φ(final double z) {
173         return Phi(z);
174     }
175 
176     /**
177      * Return Φ(z, µ, σ), the standard Gaussian cdf with mean µ and stddev σ.
178      *
179      @see #phi(double, double, double)
180      @param z the value to calculate Φ for.
181      @param mu the mean value.
182      @param sigma the stddev.
183      @return the φ value for x.
184      */
185     public static double Phi(final double z, final double mu, final double sigma) {
186         return Phi((z - mu)/sigma);
187     }
188 
189     /**
190      * Return Φ(z, µ, σ), the standard Gaussian cdf with mean µ and stddev σ.
191      *
192      @see #phi(double, double, double)
193      @param z the value to calculate Φ for.
194      @param µ the mean value.
195      @param σ the stddev.
196      @return the φ value for x.
197      */
198     public static double Φ(final double z, final double µ, final double σ) {
199         return Phi(z, µ, σ);
200     }
201 
202 }