You are viewing a plain text version of this content. The canonical link for it is here.
Posted to dev@commons.apache.org by md...@apache.org on 2003/06/14 06:17:49 UTC
cvs commit: jakarta-commons-sandbox/math/src/test/org/apache/commons/math/special BetaTest.java
mdiggory 2003/06/13 21:17:49
Modified: math/src/java/org/apache/commons/math/special Beta.java
math/src/java/org/apache/commons/math/stat
UnivariateImpl.java
math/src/java/org/apache/commons/math/stat/distribution
TDistributionImpl.java
math/src/test/org/apache/commons/math/stat/distribution
TDistributionTest.java
Added: math/src/test/org/apache/commons/math/special BetaTest.java
Log:
PR: http://nagoya.apache.org/bugzilla/show_bug.cgi?id=20766
Submitted by: brent@worden.org
Revision Changes Path
1.3 +5 -15 jakarta-commons-sandbox/math/src/java/org/apache/commons/math/special/Beta.java
Index: Beta.java
===================================================================
RCS file: /home/cvs/jakarta-commons-sandbox/math/src/java/org/apache/commons/math/special/Beta.java,v
retrieving revision 1.2
retrieving revision 1.3
diff -u -r1.2 -r1.3
--- Beta.java 11 Jun 2003 01:19:18 -0000 1.2
+++ Beta.java 14 Jun 2003 04:17:48 -0000 1.3
@@ -139,17 +139,9 @@
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;
+ if (Double.isNaN(x) || Double.isNaN(a) || Double.isNaN(b) || (x < 0)
+ || (x > 1) || (a <= 0.0) || (b <= 0.0)) {
+ ret = Double.NaN;
} else {
ContinuedFraction fraction = new ContinuedFraction() {
protected double getB(int n, double x) {
@@ -228,10 +220,8 @@
public static double logBeta(double a, 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");
+ if (Double.isNaN(a) || Double.isNaN(b) || (a <= 0.0) || (b <= 0.0)) {
+ ret = Double.NaN;
} else {
ret = Gamma.logGamma(a, epsilon, maxIterations) + Gamma.logGamma(b, epsilon, maxIterations)
- Gamma.logGamma(a + b, epsilon, maxIterations);
1.4 +125 -119 jakarta-commons-sandbox/math/src/java/org/apache/commons/math/stat/UnivariateImpl.java
Index: UnivariateImpl.java
===================================================================
RCS file: /home/cvs/jakarta-commons-sandbox/math/src/java/org/apache/commons/math/stat/UnivariateImpl.java,v
retrieving revision 1.3
retrieving revision 1.4
diff -u -r1.3 -r1.4
--- UnivariateImpl.java 4 Jun 2003 12:23:44 -0000 1.3
+++ UnivariateImpl.java 14 Jun 2003 04:17:49 -0000 1.4
@@ -89,12 +89,12 @@
/** running sum of squares that have been added */
private double sumsq = 0.0;
- /** running sum of 3rd powers that have been added */
- private double sumCube = 0.0;
-
- /** running sum of 4th powers that have been added */
- private double sumQuad = 0.0;
-
+ /** running sum of 3rd powers that have been added */
+ private double sumCube = 0.0;
+
+ /** running sum of 4th powers that have been added */
+ private double sumQuad = 0.0;
+
/** count of values that have been added */
private int n = 0;
@@ -120,17 +120,17 @@
/**
- * @see org.apache.commons.math.stat.Univariate#addValue(double)
- */
- public void addValue(double v) {
+ * @see org.apache.commons.math.stat.Univariate#addValue(double)
+ */
+ public void addValue(double v) {
insertValue(v);
}
/**
- * @see org.apache.commons.math.stat.Univariate#getMean()
- */
- public double getMean() {
+ * @see org.apache.commons.math.stat.Univariate#getMean()
+ */
+ public double getMean() {
if (n == 0) {
return Double.NaN;
} else {
@@ -140,97 +140,97 @@
/**
- * @see org.apache.commons.math.stat.Univariate#getGeometricMean()
- */
- public double getGeometricMean() {
+ * @see org.apache.commons.math.stat.Univariate#getGeometricMean()
+ */
+ public double getGeometricMean() {
if ((product <= 0.0) || (n == 0)) {
return Double.NaN;
} else {
- return Math.pow(product,( 1.0/(double)n ) );
+ return Math.pow(product,( 1.0 / (double) n ) );
}
}
/**
- * @see org.apache.commons.math.stat.Univariate#getProduct()
- */
- public double getProduct() {
+ * @see org.apache.commons.math.stat.Univariate#getProduct()
+ */
+ public double getProduct() {
return product;
}
- /**
- * @see org.apache.commons.math.stat.Univariate#getStandardDeviation()
- */
- public double getStandardDeviation() {
- double variance = getVariance();
- if ((variance == 0.0) || (variance == Double.NaN)) {
- return variance;
- } else {
- return Math.sqrt(variance);
- }
- }
-
- /**
- * Returns the variance of the values that have been added as described by
- * <a href="http://mathworld.wolfram.com/k-Statistic.html">Equation (5) for k-Statistics</a>.
- *
- * @return The variance of a set of values. Double.NaN is returned for
- * an empty set of values and 0.0 is returned for a <= 1 value set.
- */
- public double getVariance() {
- double variance = Double.NaN;
-
- if( n == 1 ) {
- variance = 0.0;
- } else if( n > 1 ) {
- variance = (((double)n)*sumsq - (sum * sum)) / (double) (n * (n - 1));
- }
+ /**
+ * @see org.apache.commons.math.stat.Univariate#getStandardDeviation()
+ */
+ public double getStandardDeviation() {
+ double variance = getVariance();
+ if ((variance == 0.0) || (variance == Double.NaN)) {
+ return variance;
+ } else {
+ return Math.sqrt(variance);
+ }
+ }
+
+ /**
+ * Returns the variance of the values that have been added as described by
+ * <a href="http://mathworld.wolfram.com/k-Statistic.html">Equation (5) for k-Statistics</a>.
+ *
+ * @return The variance of a set of values. Double.NaN is returned for
+ * an empty set of values and 0.0 is returned for a <= 1 value set.
+ */
+ public double getVariance() {
+ double variance = Double.NaN;
- return variance < 0 ? 0.0 : variance;
- }
+ if( n == 1 ) {
+ variance = 0.0;
+ } else if( n > 1 ) {
+ variance = (((double) n) * sumsq - (sum * sum)) / (double) (n * (n - 1));
+ }
+
+ return variance < 0 ? 0.0 : variance;
+ }
- /**
- * Returns the skewness of the values that have been added as described by
+ /**
+ * Returns the skewness of the values that have been added as described by
* <a href="http://mathworld.wolfram.com/k-Statistic.html">Equation (6) for k-Statistics</a>.
*
- * @return The skew of a set of values. Double.NaN is returned for
- * an empty set of values and 0.0 is returned for a <= 2 value set.
- */
- public double getSkewness() {
-
- if( n < 1) return Double.NaN;
- if( n <= 2 ) return 0.0;
-
- return ( 2*Math.pow(sum,3) - 3*sum*sumsq + ((double)(n*n))*sumCube ) /
- ( (double)(n*(n-1)*(n-2)) ) ;
- }
-
- /**
- * Returns the kurtosis of the values that have been added as described by
+ * @return The skew of a set of values. Double.NaN is returned for
+ * an empty set of values and 0.0 is returned for a <= 2 value set.
+ */
+ public double getSkewness() {
+
+ if( n < 1) return Double.NaN;
+ if( n <= 2 ) return 0.0;
+
+ return ( 2 * Math.pow(sum, 3) - 3 * sum * sumsq + ((double) (n * n)) * sumCube ) /
+ ( (double) (n * (n - 1) * (n - 2)) ) ;
+ }
+
+ /**
+ * Returns the kurtosis of the values that have been added as described by
* <a href="http://mathworld.wolfram.com/k-Statistic.html">Equation (7) for k-Statistics</a>.
*
- * @return The kurtosis of a set of values. Double.NaN is returned for
- * an empty set of values and 0.0 is returned for a <= 3 value set.
- */
- public double getKurtosis() {
-
- if( n < 1) return Double.NaN;
- if( n <= 3 ) return 0.0;
-
- double x1 = -6*Math.pow(sum,4);
- double x2 = 12*((double)n)*Math.pow(sum,2)*sumsq;
- double x3 = -3*((double)(n*(n-1)))*Math.pow(sumsq,2);
- double x4 = -4*((double)(n*(n+1)))*sum*sumCube;
- double x5 = Math.pow(((double)n),2)*((double)(n+1))*sumQuad;
-
- return (x1 + x2 + x3 + x4 + x5) /
- ( (double)(n*(n-1)*(n-2)*(n-3)) );
- }
-
+ * @return The kurtosis of a set of values. Double.NaN is returned for
+ * an empty set of values and 0.0 is returned for a <= 3 value set.
+ */
+ public double getKurtosis() {
+
+ if( n < 1) return Double.NaN;
+ if( n <= 3 ) return 0.0;
+
+ double x1 = -6 * Math.pow(sum, 4);
+ double x2 = 12 * ((double) n) * Math.pow(sum, 2) * sumsq;
+ double x3 = -3 * ((double) (n * (n - 1))) * Math.pow(sumsq,2);
+ double x4 = -4 * ((double) (n * (n + 1))) * sum * sumCube;
+ double x5 = Math.pow(((double) n),2) * ((double) (n+1)) * sumQuad;
+
+ return (x1 + x2 + x3 + x4 + x5) /
+ ( (double) (n * (n - 1) * (n - 2) * (n - 3)) );
+ }
+
/**
* Called in "addValue" to insert a new value into the statistic.
- * @param v The value to be added.
- */
- private void insertValue(double v) {
+ * @param v The value to be added.
+ */
+ private void insertValue(double v) {
// The default value of product is NaN, if you
// try to retrieve the product for a univariate with
@@ -250,34 +250,36 @@
// Remove the influence of the discarded
sum -= discarded;
sumsq -= discarded * discarded;
- sumCube -= Math.pow(discarded,3);
- sumQuad -= Math.pow(discarded,4);
-
+ sumCube -= Math.pow(discarded, 3);
+ sumQuad -= Math.pow(discarded, 4);
+
if(discarded == min) {
min = doubleArray.getMin();
- } else {
- if(discarded == max){
+ } else if(discarded == max){
max = doubleArray.getMax();
- }
}
if(product != 0.0){
// can safely remove discarded value
- product *= v/discarded;
+ product *= v / discarded;
} else if(discarded == 0.0){
// need to recompute product
product = 1.0;
double[] elements = doubleArray.getElements();
for( int i = 0; i < elements.length; i++ ) {
- product *= elements[i];
+ product *= elements[i];
}
} // else product = 0 and will still be 0 after discard
} else {
- doubleArray.addElement( v );
+ doubleArray.addElement( v );
n += 1.0;
- if (v < min) min = v;
- if (v > max) max = v;
+ if (v < min) {
+ min = v;
+ }
+ if (v > max) {
+ max = v;
+ }
product *= v;
}
} else {
@@ -285,15 +287,19 @@
// worry about storing any values. We don't need to discard the
// influence of any single item.
n += 1.0;
- if (v < min) min = v;
- if (v > max) max = v;
+ if (v < min) {
+ min = v;
+ }
+ if (v > max) {
+ max = v;
+ }
product *= v;
}
- sum += v;
- sumsq += v*v;
- sumCube += Math.pow(v,3);
- sumQuad += Math.pow(v,4);
+ sum += v;
+ sumsq += v * v;
+ sumCube += Math.pow(v,3);
+ sumQuad += Math.pow(v,4);
}
/** Getter for property max.
@@ -339,20 +345,20 @@
return sumsq;
}
- /** Getter for property sumCube.
- * @return Value of property sumCube.
- */
- public double getSumCube() {
- return sumCube;
- }
-
- /** Getter for property sumQuad.
- * @return Value of property sumQuad.
- */
- public double getSumQuad() {
- return sumQuad;
- }
-
+ /** Getter for property sumCube.
+ * @return Value of property sumCube.
+ */
+ public double getSumCube() {
+ return sumCube;
+ }
+
+ /** Getter for property sumQuad.
+ * @return Value of property sumQuad.
+ */
+ public double getSumQuad() {
+ return sumQuad;
+ }
+
/**
* Generates a text report displaying
* univariate statistics from values that
@@ -367,8 +373,8 @@
outBuffer.append("max: " + max + "\n");
outBuffer.append("mean: " + getMean() + "\n");
outBuffer.append("std dev: " + getStandardDeviation() + "\n");
- outBuffer.append("skewness: " + getSkewness() + "\n");
- outBuffer.append("kurtosis: " + getKurtosis() + "\n");
+ outBuffer.append("skewness: " + getSkewness() + "\n");
+ outBuffer.append("kurtosis: " + getKurtosis() + "\n");
return outBuffer.toString();
}
1.2 +12 -9 jakarta-commons-sandbox/math/src/java/org/apache/commons/math/stat/distribution/TDistributionImpl.java
Index: TDistributionImpl.java
===================================================================
RCS file: /home/cvs/jakarta-commons-sandbox/math/src/java/org/apache/commons/math/stat/distribution/TDistributionImpl.java,v
retrieving revision 1.1
retrieving revision 1.2
diff -u -r1.1 -r1.2
--- TDistributionImpl.java 7 Jun 2003 13:57:54 -0000 1.1
+++ TDistributionImpl.java 14 Jun 2003 04:17:49 -0000 1.2
@@ -103,18 +103,21 @@
* @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 {
+ if(x == 0.0){
ret = 0.5;
+ } else {
+ double t = Beta.regularizedBeta(
+ getDegreesOfFreedom() / (getDegreesOfFreedom() + (x * x)),
+ 0.5 * getDegreesOfFreedom(), 0.5);
+
+ if(x < 0.0){
+ ret = 0.5 * t;
+ } else {
+ ret = 1.0 - 0.5 * t;
+ }
}
+
return ret;
}
1.2 +67 -19 jakarta-commons-sandbox/math/src/test/org/apache/commons/math/stat/distribution/TDistributionTest.java
Index: TDistributionTest.java
===================================================================
RCS file: /home/cvs/jakarta-commons-sandbox/math/src/test/org/apache/commons/math/stat/distribution/TDistributionTest.java,v
retrieving revision 1.1
retrieving revision 1.2
diff -u -r1.1 -r1.2
--- TDistributionTest.java 7 Jun 2003 13:57:54 -0000 1.1
+++ TDistributionTest.java 14 Jun 2003 04:17:49 -0000 1.2
@@ -85,45 +85,93 @@
super.tearDown();
}
- public void testLowerTailProbability(){
+ public void testInverseCummulativeProbability001() {
+ testValue(-5.893, .001);
+ }
+
+ public void testInverseCumulativeProbability010() {
+ testValue(-3.365, .010);
+ }
+
+ public void testInverseCumulativeProbability025() {
+ testValue(-2.571, .025);
+ }
+
+ public void testInverseCumulativeProbability050() {
+ testValue(-2.015, .050);
+ }
+
+ public void testInverseCumulativeProbability100() {
+ testValue(-1.476, .100);
+ }
+
+ public void testInverseCummulativeProbability999() {
+ testValue(5.893, .999);
+ }
+
+ public void testInverseCumulativeProbability990() {
+ testValue(3.365, .990);
+ }
+
+ public void testInverseCumulativeProbability975() {
+ testValue(2.571, .975);
+ }
+
+ public void testInverseCumulativeProbability950() {
+ testValue(2.015, .950);
+ }
+
+ public void testInverseCumulativeProbability900() {
+ testValue(1.476, .900);
+ }
+
+ public void testCummulativeProbability001() {
testProbability(-5.893, .001);
+ }
+
+ public void testCumulativeProbability010() {
testProbability(-3.365, .010);
+ }
+
+ public void testCumulativeProbability025() {
testProbability(-2.571, .025);
+ }
+
+ public void testCumulativeProbability050() {
testProbability(-2.015, .050);
+ }
+
+ public void testCumulativeProbability100() {
testProbability(-1.476, .100);
}
- public void testUpperTailProbability(){
+ public void testCummulativeProbability999() {
testProbability(5.893, .999);
+ }
+
+ public void testCumulativeProbability990() {
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 testCumulativeProbability975() {
+ testProbability(2.571, .975);
+ }
+
+ public void testCumulativeProbability950() {
+ testProbability(2.015, .950);
}
- public void testUpperTailValues(){
- testValue(5.893, .999);
- testValue(3.365, .990);
- testValue(2.571, .975);
- testValue(2.015, .950);
- testValue(1.476, .900);
+ public void testCumulativeProbability900() {
+ testProbability(1.476, .900);
}
private void testProbability(double x, double expected){
double actual = t.cummulativeProbability(x);
- assertEquals("probability for " + x, expected, actual, 10e-4);
+ assertEquals(expected, actual, 10e-4);
}
private void testValue(double expected, double p){
double actual = t.inverseCummulativeProbability(p);
- assertEquals("value for " + p, expected, actual, 10e-4);
+ assertEquals(expected, actual, 10e-4);
}
}
1.1 jakarta-commons-sandbox/math/src/test/org/apache/commons/math/special/BetaTest.java
Index: BetaTest.java
===================================================================
/* ====================================================================
* 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;
import junit.framework.TestCase;
/**
* @author Brent Worden
*/
public class BetaTest extends TestCase {
/**
* Constructor for BetaTest.
* @param name
*/
public BetaTest(String name) {
super(name);
}
private void testRegularizedBeta(double expected, double x, double a, double b) {
double actual = Beta.regularizedBeta(x, a, b);
if (Double.isNaN(expected)) {
assertTrue(Double.isNaN(actual));
} else {
assertEquals(expected, actual, 10e-5);
}
}
private void testLogBeta(double expected, double a, double b) {
double actual = Beta.logBeta(a, b);
if (Double.isNaN(expected)) {
assertTrue(Double.isNaN(actual));
} else {
assertEquals(expected, actual, 10e-5);
}
}
public void testRegularizedBetaNanPositivePositive() {
testRegularizedBeta(Double.NaN, Double.NaN, 1.0, 1.0);
}
public void testRegularizedBetaPositiveNanPositive() {
testRegularizedBeta(Double.NaN, 0.5, Double.NaN, 1.0);
}
public void testRegularizedBetaPositivePositiveNan() {
testRegularizedBeta(Double.NaN, 0.5, 1.0, Double.NaN);
}
public void testRegularizedBetaNegativePositivePositive() {
testRegularizedBeta(Double.NaN, -0.5, 1.0, 2.0);
}
public void testRegularizedBetaPositiveNegativePositive() {
testRegularizedBeta(Double.NaN, 0.5, -1.0, 2.0);
}
public void testRegularizedBetaPositivePositiveNegative() {
testRegularizedBeta(Double.NaN, 0.5, 1.0, -2.0);
}
public void testRegularizedBetaZeroPositivePositive() {
testRegularizedBeta(0.0, 0.0, 1.0, 2.0);
}
public void testRegularizedBetaPositiveZeroPositive() {
testRegularizedBeta(Double.NaN, 0.5, 0.0, 2.0);
}
public void testRegularizedBetaPositivePositiveZero() {
testRegularizedBeta(Double.NaN, 0.5, 1.0, 0.0);
}
public void testRegularizedBetaPositivePositivePositive() {
testRegularizedBeta(0.75, 0.5, 1.0, 2.0);
}
public void testLogBetaNanPositive() {
testLogBeta(Double.NaN, Double.NaN, 2.0);
}
public void testLogBetaPositiveNan() {
testLogBeta(Double.NaN, 1.0, Double.NaN);
}
public void testLogBetaNegativePositive() {
testLogBeta(Double.NaN, -1.0, 2.0);
}
public void testLogBetaPositiveNegative() {
testLogBeta(Double.NaN, 1.0, -2.0);
}
public void testLogBetaZeroPositive() {
testLogBeta(Double.NaN, 0.0, 2.0);
}
public void testLogBetaPositiveZero() {
testLogBeta(Double.NaN, 1.0, 0.0);
}
public void testLogBetaPositivePositive() {
testLogBeta(-0.693147, 1.0, 2.0);
}
}
---------------------------------------------------------------------
To unsubscribe, e-mail: commons-dev-unsubscribe@jakarta.apache.org
For additional commands, e-mail: commons-dev-help@jakarta.apache.org