Apache
/
commons-math
mirror of https://github.com/apache/commons-math.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

PR: http://nagoya.apache.org/bugzilla/show_bug.cgi?id=20570 

Submitted by:	brent@worden.org


git-svn-id: https://svn.apache.org/repos/asf/jakarta/commons/proper/math/trunk@140898 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Mark R. Diggory 2003-06-07 13:57:54 +00:00
parent a6c324e669
commit 8d1482c461
7 changed files with 611 additions and 26 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

7
src/java/org/apache/commons/math/distribution/DistributionFactory.java
View File

@ -110,4 +110,11 @@ public abstract class DistributionFactory {
*/
public abstract GammaDistribution createGammaDistribution(
double alpha, double beta);
/**
* Create a new t distribution with the given degrees of freedom.
* @param degreesOfFreedom degrees of freedom.
* @return a new t distribution.
*/
public abstract TDistribution createTDistribution(double degreesOfFreedom);
}

9
src/java/org/apache/commons/math/distribution/DistributionFactoryImpl.java
View File

@ -90,4 +90,13 @@ public class DistributionFactoryImpl extends DistributionFactory {
return new GammaDistributionImpl(alpha, beta);
}
/**
* Create a new t distribution with the given degrees of freedom.
* @param degreesOfFreedom degrees of freedom.
* @return a new t distribution.
*/
public TDistribution createTDistribution(double degreesOfFreedom) {
return new TDistributionImpl(degreesOfFreedom);
}
}

86
src/java/org/apache/commons/math/distribution/TDistribution.java Normal file
View File

@ -0,0 +1,86 @@
/* ====================================================================
* The Apache Software License, Version 1.1
*
* Copyright (c) 2003 The Apache Software Foundation. All rights
* reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
* are met:
*
* 1. Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in
* the documentation and/or other materials provided with the
* distribution.
*
* 3. The end-user documentation included with the redistribution, if
* any, must include the following acknowlegement:
* "This product includes software developed by the
* Apache Software Foundation (http://www.apache.org/)."
* Alternately, this acknowlegement may appear in the software itself,
* if and wherever such third-party acknowlegements normally appear.
*
* 4. The names "The Jakarta Project", "Commons", and "Apache Software
* Foundation" must not be used to endorse or promote products derived
* from this software without prior written permission. For written
* permission, please contact apache@apache.org.
*
* 5. Products derived from this software may not be called "Apache"
* nor may "Apache" appear in their names without prior written
* permission of the Apache Software Foundation.
*
* THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
* WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
* DISCLAIMED. IN NO EVENT SHALL THE APACHE SOFTWARE FOUNDATION OR
* ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
* USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
* OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
* ====================================================================
*
* This software consists of voluntary contributions made by many
* individuals on behalf of the Apache Software Foundation. For more
* information on the Apache Software Foundation, please see
* <http://www.apache.org/>.
*/
package org.apache.commons.math.stat.distribution;
/**
* <p>
* Student's t-Distribution.
* </p>
*
* <p>
* Instances of TDistribution objects should be created using
* {@link DistributionFactory#createTDistribution(double)}
* </p>
*
* <p>
* Reference:<br/>
* <a href="http://mathworld.wolfram.com/Studentst-Distribution.html">
* Student's t-Distribution</a>
* </p>
*
* @author Brent Worden
*/
public interface TDistribution extends ContinuousDistribution {
/**
* Modify the degrees of freedom.
* @param degreesOfFreedom the new degrees of freedom.
*/
void setDegreesOfFreedom(double degreesOfFreedom);
/**
* Access the degrees of freedom.
* @return the degrees of freedom.
*/
double getDegreesOfFreedom();
}

158
src/java/org/apache/commons/math/distribution/TDistributionImpl.java Normal file
View File

@ -0,0 +1,158 @@
/* ====================================================================
* The Apache Software License, Version 1.1
*
* Copyright (c) 2003 The Apache Software Foundation. All rights
* reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
* are met:
*
* 1. Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in
* the documentation and/or other materials provided with the
* distribution.
*
* 3. The end-user documentation included with the redistribution, if
* any, must include the following acknowlegement:
* "This product includes software developed by the
* Apache Software Foundation (http://www.apache.org/)."
* Alternately, this acknowlegement may appear in the software itself,
* if and wherever such third-party acknowlegements normally appear.
*
* 4. The names "The Jakarta Project", "Commons", and "Apache Software
* Foundation" must not be used to endorse or promote products derived
* from this software without prior written permission. For written
* permission, please contact apache@apache.org.
*
* 5. Products derived from this software may not be called "Apache"
* nor may "Apache" appear in their names without prior written
* permission of the Apache Software Foundation.
*
* THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
* WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
* DISCLAIMED. IN NO EVENT SHALL THE APACHE SOFTWARE FOUNDATION OR
* ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
* USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
* OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
* ====================================================================
*
* This software consists of voluntary contributions made by many
* individuals on behalf of the Apache Software Foundation. For more
* information on the Apache Software Foundation, please see
* <http://www.apache.org/>.
*/
package org.apache.commons.math.stat.distribution;
import org.apache.commons.math.special.Beta;
/**
* Default implementation of
* {@link org.apache.commons.math.stat.distribution.TDistribution}.
*
* @author Brent Worden
*/
public class TDistributionImpl
extends AbstractContinuousDistribution
implements TDistribution {
/** The degrees of freedom*/
private double degreesOfFreedom;
/**
* Create a t distribution using the given degrees of freedom.
* @param degreesOfFreedom the degrees of freedom.
*/
public TDistributionImpl(double degreesOfFreedom){
super();
setDegreesOfFreedom(degreesOfFreedom);
}
/**
* Modify the degrees of freedom.
* @param degreesOfFreedom the new degrees of freedom.
*/
public void setDegreesOfFreedom(double degreesOfFreedom) {
if(degreesOfFreedom <= 0.0){
throw new IllegalArgumentException(
"degrees of freedom must be positive.");
}
this.degreesOfFreedom = degreesOfFreedom;
}
/**
* Access the degrees of freedom.
* @return the degrees of freedom.
*/
public double getDegreesOfFreedom() {
return degreesOfFreedom;
}
/**
* For this disbution, X, this method returns P(X &lt; <code>x</code>).
* @param x the value at which the CDF is evaluated.
* @return CDF evaluted at <code>x</code>.
*/
public double cummulativeProbability(double x) {
double t = Beta.regularizedBeta(
getDegreesOfFreedom() / (getDegreesOfFreedom() + (x * x)),
0.5 * getDegreesOfFreedom(), 0.5);
double ret;
if(x < 0.0){
ret = 0.5 * t;
} else if(x > 0.0){
ret = 1.0 - 0.5 * t;
} else {
ret = 0.5;
}
return ret;
}
/**
* Access the domain value lower bound, based on <code>p</code>, used to
* bracket a CDF root. This method is used by
* {@link #inverseCummulativeProbability(double)} to find critical values.
*
* @param p the desired probability for the critical value
* @return domain value lower bound, i.e.
* P(X &lt; <i>lower bound</i>) &lt; <code>p</code>
*/
protected double getDomainLowerBound(double p){
return -Double.MAX_VALUE;
}
/**
* Access the domain value upper bound, based on <code>p</code>, used to
* bracket a CDF root. This method is used by
* {@link #inverseCummulativeProbability(double)} to find critical values.
*
* @param p the desired probability for the critical value
* @return domain value upper bound, i.e.
* P(X &lt; <i>upper bound</i>) &gt; <code>p</code>
*/
protected double getDomainUpperBound(double p){
return Double.MAX_VALUE;
}
/**
* Access the initial domain value, based on <code>p</code>, used to
* bracket a CDF root. This method is used by
* {@link #inverseCummulativeProbability(double)} to find critical values.
*
* @param p the desired probability for the critical value
* @return initial domain value
*/
protected double getInitialDomain(double p){
return 0.0;
}
}

196
src/java/org/apache/commons/math/special/Beta.java Normal file
View File

@ -0,0 +1,196 @@
/* ====================================================================
* The Apache Software License, Version 1.1
*
* Copyright (c) 2003 The Apache Software Foundation. All rights
* reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
* are met:
*
* 1. Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in
* the documentation and/or other materials provided with the
* distribution.
*
* 3. The end-user documentation included with the redistribution, if
* any, must include the following acknowlegement:
* "This product includes software developed by the
* Apache Software Foundation (http://www.apache.org/)."
* Alternately, this acknowlegement may appear in the software itself,
* if and wherever such third-party acknowlegements normally appear.
*
* 4. The names "The Jakarta Project", "Commons", and "Apache Software
* Foundation" must not be used to endorse or promote products derived
* from this software without prior written permission. For written
* permission, please contact apache@apache.org.
*
* 5. Products derived from this software may not be called "Apache"
* nor may "Apache" appear in their names without prior written
* permission of the Apache Software Foundation.
*
* THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
* WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
* DISCLAIMED. IN NO EVENT SHALL THE APACHE SOFTWARE FOUNDATION OR
* ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
* USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
* OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
* ====================================================================
*
* This software consists of voluntary contributions made by many
* individuals on behalf of the Apache Software Foundation. For more
* information on the Apache Software Foundation, please see
* <http://www.apache.org/>.
*/
package org.apache.commons.math.special;
/**
* This is a utility class that provides computation methods related to the
* Gamma family of functions.
*
* @author Brent Worden
*/
public class Beta {
/** Maximum allowed numerical error. */
private static final double DEFAULT_EPSILON = 10e-9;
/**
* Default constructor. Prohibit instantiation.
*/
private Beta() {
super();
}
/**
* <p>
* Returns the regularized beta function I(x, a, b).
* </p>
*
* @param x ???
* @param a ???
* @param b ???
* @return the regularized beta function I(x, a, b)
*/
public static double regularizedBeta(double x, double a, double b) {
return regularizedBeta(x, a, b, DEFAULT_EPSILON, Integer.MAX_VALUE);
}
/**
* <p>
* Returns the regularized beta function I(x, a, b).
* </p>
*
* @param x ???
* @param a ???
* @param b ???
* @return the regularized beta function I(x, a, b)
*/
public static double regularizedBeta(double x, double a, double b, double epsilon) {
return regularizedBeta(x, a, b, epsilon, Integer.MAX_VALUE);
}
/**
* <p>
* Returns the regularized beta function I(x, a, b).
* </p>
*
* @param x ???
* @param a ???
* @param b ???
* @return the regularized beta function I(x, a, b)
*/
public static double regularizedBeta(double x, double a, double b, int maxIterations) {
return regularizedBeta(x, a, b, DEFAULT_EPSILON, maxIterations);
}
/**
* <p>
* Returns the regularized beta function I(x, a, b).
* </p>
*
* <p>
* The implementation of this method is based on:
* <ul>
* <li>
* <a href="http://mathworld.wolfram.com/RegularizedBetaFunction.html">
* Regularized Beta Function</a>.</li>
* <li>
* <a href="http://mathworld.wolfram.com/IncompleteBetaFunction.html">
* Incomplete Beta Function</a>.</li>
* </ul>
* </p>
*
* @param x ???
* @param a ???
* @param b ???
* @return the regularized beta function I(x, a, b)
*/
public static double regularizedBeta(double x, final double a, final double b, double epsilon, int maxIterations) {
double ret;
if (a <= 0.0) {
throw new IllegalArgumentException("a must be positive");
} else if (b <= 0.0) {
throw new IllegalArgumentException("b must be positive");
} else if (x < 0.0 || x > 1.0) {
throw new IllegalArgumentException(
"x must be between 0.0 and 1.0, inclusive");
} else if(x == 0.0){
ret = 0.0;
} else if(x == 1.0){
ret = 1.0;
} else {
double n = 0.0;
double an = 1.0 / a;
double s = an;
while(Math.abs(an) > epsilon && n < maxIterations){
n = n + 1.0;
an = an * (n - b) / n * x / (a + n) * (a + n - 1);
s = s + an;
}
ret = Math.exp(a * Math.log(x) - logBeta(a, b)) * s;
}
return ret;
}
/**
* <p>
* Returns the natural logarithm of the beta function B(a, b).
* </p>
*
* <p>
* The implementation of this method is based on:
* <ul>
* <li><a href="http://mathworld.wolfram.com/BetaFunction.html">
* Beta Function</a>, equation (1).</li>
* </ul>
* </p>
*
* @param x ???
* @return log(B(a, b))
*/
public static double logBeta(double a, double b) {
double ret;
if (a <= 0.0) {
throw new IllegalArgumentException("a must be positive");
} else if (b <= 0.0) {
throw new IllegalArgumentException("b must be positive");
} else {
ret = Gamma.logGamma(a) + Gamma.logGamma(b)
- Gamma.logGamma(a + b);
}
return ret;
}
}

52
src/test/org/apache/commons/math/distribution/DistributionFactoryImplTest.java
View File

@ -101,30 +101,30 @@ public class DistributionFactoryImplTest extends TestCase {
fail("positive alpah and beta. IllegalArgumentException is not expected");
}
}
//
// public void testCreateTDistributionNegative(){
// try {
// factory.createTDistribution(-1.0);
// fail("negative degrees of freedom. IllegalArgumentException expected");
// } catch (IllegalArgumentException ex) {
// ;
// }
// }
//
// public void testCreateTDistributionZero(){
// try {
// factory.createTDistribution(0.0);
// fail("zero degrees of freedom. IllegalArgumentException expected");
// } catch (IllegalArgumentException ex) {
// ;
// }
// }
//
// public void testCreateTDistributionPositive(){
// try {
// factory.createTDistribution(1.0);
// } catch (IllegalArgumentException ex) {
// fail("positive degrees of freedom. IllegalArgumentException is not expected");
// }
// }
public void testCreateTDistributionNegative(){
try {
factory.createTDistribution(-1.0);
fail("negative degrees of freedom. IllegalArgumentException expected");
} catch (IllegalArgumentException ex) {
;
}
}
public void testCreateTDistributionZero(){
try {
factory.createTDistribution(0.0);
fail("zero degrees of freedom. IllegalArgumentException expected");
} catch (IllegalArgumentException ex) {
;
}
}
public void testCreateTDistributionPositive(){
try {
factory.createTDistribution(1.0);
} catch (IllegalArgumentException ex) {
fail("positive degrees of freedom. IllegalArgumentException is not expected");
}
}
}

129
src/test/org/apache/commons/math/distribution/TDistributionTest.java Normal file
View File

@ -0,0 +1,129 @@
/* ====================================================================
* The Apache Software License, Version 1.1
*
* Copyright (c) 2003 The Apache Software Foundation. All rights
* reserved.
*
* Redistribution and use in source and binary forms, with or without
* modification, are permitted provided that the following conditions
* are met:
*
* 1. Redistributions of source code must retain the above copyright
* notice, this list of conditions and the following disclaimer.
*
* 2. Redistributions in binary form must reproduce the above copyright
* notice, this list of conditions and the following disclaimer in
* the documentation and/or other materials provided with the
* distribution.
*
* 3. The end-user documentation included with the redistribution, if
* any, must include the following acknowlegement:
* "This product includes software developed by the
* Apache Software Foundation (http://www.apache.org/)."
* Alternately, this acknowlegement may appear in the software itself,
* if and wherever such third-party acknowlegements normally appear.
*
* 4. The names "The Jakarta Project", "Commons", and "Apache Software
* Foundation" must not be used to endorse or promote products derived
* from this software without prior written permission. For written
* permission, please contact apache@apache.org.
*
* 5. Products derived from this software may not be called "Apache"
* nor may "Apache" appear in their names without prior written
* permission of the Apache Software Foundation.
*
* THIS SOFTWARE IS PROVIDED ``AS IS'' AND ANY EXPRESSED OR IMPLIED
* WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
* OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
* DISCLAIMED. IN NO EVENT SHALL THE APACHE SOFTWARE FOUNDATION OR
* ITS CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
* SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
* LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF
* USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
* ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
* OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT
* OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
* SUCH DAMAGE.
* ====================================================================
*
* This software consists of voluntary contributions made by many
* individuals on behalf of the Apache Software Foundation. For more
* information on the Apache Software Foundation, please see
* <http://www.apache.org/>.
*/
package org.apache.commons.math.stat.distribution;
import junit.framework.TestCase;
/**
* @author Brent Worden
*/
public class TDistributionTest extends TestCase {
private TDistribution t;
/**
* Constructor for ChiSquareDistributionTest.
* @param name
*/
public TDistributionTest(String name) {
super(name);
}
/*
* @see TestCase#setUp()
*/
protected void setUp() throws Exception {
super.setUp();
t = DistributionFactory.newInstance().createTDistribution(5.0);
}
/*
* @see TestCase#tearDown()
*/
protected void tearDown() throws Exception {
t = null;
super.tearDown();
}
public void testLowerTailProbability(){
testProbability(-5.893, .001);
testProbability(-3.365, .010);
testProbability(-2.571, .025);
testProbability(-2.015, .050);
testProbability(-1.476, .100);
}
public void testUpperTailProbability(){
testProbability(5.893, .999);
testProbability(3.365, .990);
testProbability(2.571, .975);
testProbability(2.015, .950);
testProbability(1.476, .900);
}
public void testLowerTailValues(){
testValue(-5.893, .001);
testValue(-3.365, .010);
testValue(-2.571, .025);
testValue(-2.015, .050);
testValue(-1.476, .100);
}
public void testUpperTailValues(){
testValue(5.893, .999);
testValue(3.365, .990);
testValue(2.571, .975);
testValue(2.015, .950);
testValue(1.476, .900);
}
private void testProbability(double x, double expected){
double actual = t.cummulativeProbability(x);
assertEquals("probability for " + x, expected, actual, 10e-4);
}
private void testValue(double expected, double p){
double actual = t.inverseCummulativeProbability(p);
assertEquals("value for " + p, expected, actual, 10e-4);
}
}