BinomialTest.java

  1. /*
  2.  * Licensed to the Apache Software Foundation (ASF) under one or more
  3.  * contributor license agreements.  See the NOTICE file distributed with
  4.  * this work for additional information regarding copyright ownership.
  5.  * The ASF licenses this file to You under the Apache License, Version 2.0
  6.  * (the "License"); you may not use this file except in compliance with
  7.  * the License.  You may obtain a copy of the License at
  8.  *
  9.  *      http://www.apache.org/licenses/LICENSE-2.0
  10.  *
  11.  * Unless required by applicable law or agreed to in writing, software
  12.  * distributed under the License is distributed on an "AS IS" BASIS,
  13.  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  14.  * See the License for the specific language governing permissions and
  15.  * limitations under the License.
  16.  */
  17. package org.apache.commons.math4.legacy.stat.inference;

  18. import org.apache.commons.statistics.distribution.BinomialDistribution;
  19. import org.apache.commons.math4.legacy.exception.MathIllegalArgumentException;
  20. import org.apache.commons.math4.legacy.exception.MathInternalError;
  21. import org.apache.commons.math4.legacy.exception.NotPositiveException;
  22. import org.apache.commons.math4.legacy.exception.NullArgumentException;
  23. import org.apache.commons.math4.legacy.exception.OutOfRangeException;
  24. import org.apache.commons.math4.legacy.exception.util.LocalizedFormats;

  25. /**
  26.  * Implements binomial test statistics.
  27.  * <p>
  28.  * Exact test for the statistical significance of deviations from a
  29.  * theoretically expected distribution of observations into two categories.
  30.  *
  31.  * @see <a href="http://en.wikipedia.org/wiki/Binomial_test">Binomial test (Wikipedia)</a>
  32.  * @since 3.3
  33.  */
  34. public class BinomialTest {

  35.     /**
  36.      * Returns whether the null hypothesis can be rejected with the given confidence level.
  37.      * <p>
  38.      * <strong>Preconditions</strong>:
  39.      * <ul>
  40.      * <li>Number of trials must be &ge; 0.</li>
  41.      * <li>Number of successes must be &ge; 0.</li>
  42.      * <li>Number of successes must be &le; number of trials.</li>
  43.      * <li>Probability must be &ge; 0 and &le; 1.</li>
  44.      * </ul>
  45.      *
  46.      * @param numberOfTrials number of trials performed
  47.      * @param numberOfSuccesses number of successes observed
  48.      * @param probability assumed probability of a single trial under the null hypothesis
  49.      * @param alternativeHypothesis type of hypothesis being evaluated (one- or two-sided)
  50.      * @param alpha significance level of the test
  51.      * @return true if the null hypothesis can be rejected with confidence {@code 1 - alpha}
  52.      * @throws NotPositiveException if {@code numberOfTrials} or {@code numberOfSuccesses} is negative
  53.      * @throws OutOfRangeException if {@code probability} is not between 0 and 1
  54.      * @throws MathIllegalArgumentException if {@code numberOfTrials} &lt; {@code numberOfSuccesses} or
  55.      * if {@code alternateHypothesis} is null.
  56.      * @see AlternativeHypothesis
  57.      */
  58.     public boolean binomialTest(int numberOfTrials, int numberOfSuccesses, double probability,
  59.                                 AlternativeHypothesis alternativeHypothesis, double alpha) {
  60.         double pValue = binomialTest(numberOfTrials, numberOfSuccesses, probability, alternativeHypothesis);
  61.         return pValue < alpha;
  62.     }

  63.     /**
  64.      * Returns the <i>observed significance level</i>, or
  65.      * <a href="http://www.cas.lancs.ac.uk/glossary_v1.1/hyptest.html#pvalue">p-value</a>,
  66.      * associated with a <a href="http://en.wikipedia.org/wiki/Binomial_test"> Binomial test</a>.
  67.      * <p>
  68.      * The number returned is the smallest significance level at which one can reject the null hypothesis.
  69.      * The form of the hypothesis depends on {@code alternativeHypothesis}.</p>
  70.      * <p>
  71.      * The p-Value represents the likelihood of getting a result at least as extreme as the sample,
  72.      * given the provided {@code probability} of success on a single trial. For single-sided tests,
  73.      * this value can be directly derived from the Binomial distribution. For the two-sided test,
  74.      * the implementation works as follows: we start by looking at the most extreme cases
  75.      * (0 success and n success where n is the number of trials from the sample) and determine their likelihood.
  76.      * The lower value is added to the p-Value (if both values are equal, both are added). Then we continue with
  77.      * the next extreme value, until we added the value for the actual observed sample.</p>
  78.      * <p>
  79.      * <strong>Preconditions</strong>:
  80.      * <ul>
  81.      * <li>Number of trials must be &ge; 0.</li>
  82.      * <li>Number of successes must be &ge; 0.</li>
  83.      * <li>Number of successes must be &le; number of trials.</li>
  84.      * <li>Probability must be &ge; 0 and &le; 1.</li>
  85.      * </ul>
  86.      *
  87.      * @param numberOfTrials number of trials performed
  88.      * @param numberOfSuccesses number of successes observed
  89.      * @param probability assumed probability of a single trial under the null hypothesis
  90.      * @param alternativeHypothesis type of hypothesis being evaluated (one- or two-sided)
  91.      * @return p-value
  92.      * @throws NotPositiveException if {@code numberOfTrials} or {@code numberOfSuccesses} is negative
  93.      * @throws OutOfRangeException if {@code probability} is not between 0 and 1
  94.      * @throws MathIllegalArgumentException if {@code numberOfTrials} &lt; {@code numberOfSuccesses} or
  95.      * if {@code alternateHypothesis} is null.
  96.      * @see AlternativeHypothesis
  97.      */
  98.     public double binomialTest(int numberOfTrials, int numberOfSuccesses, double probability,
  99.                                AlternativeHypothesis alternativeHypothesis) {
  100.         if (numberOfTrials < 0) {
  101.             throw new NotPositiveException(numberOfTrials);
  102.         }
  103.         if (numberOfSuccesses < 0) {
  104.             throw new NotPositiveException(numberOfSuccesses);
  105.         }
  106.         if (probability < 0 || probability > 1) {
  107.             throw new OutOfRangeException(probability, 0, 1);
  108.         }
  109.         if (numberOfTrials < numberOfSuccesses) {
  110.             throw new MathIllegalArgumentException(
  111.                 LocalizedFormats.BINOMIAL_INVALID_PARAMETERS_ORDER,
  112.                 numberOfTrials, numberOfSuccesses);
  113.         }
  114.         if (alternativeHypothesis == null) {
  115.             throw new NullArgumentException();
  116.         }

  117.         final BinomialDistribution distribution = BinomialDistribution.of(numberOfTrials, probability);
  118.         switch (alternativeHypothesis) {
  119.         case GREATER_THAN:
  120.             return distribution.survivalProbability(numberOfSuccesses - 1);
  121.         case LESS_THAN:
  122.             return distribution.cumulativeProbability(numberOfSuccesses);
  123.         case TWO_SIDED:
  124.             int criticalValueLow = 0;
  125.             int criticalValueHigh = numberOfTrials;
  126.             double pTotal = 0;

  127.             while (true) {
  128.                 double pLow = distribution.probability(criticalValueLow);
  129.                 double pHigh = distribution.probability(criticalValueHigh);

  130.                 if (pLow == pHigh) {
  131.                     if (criticalValueLow == criticalValueHigh) {
  132.                         pTotal += pLow;
  133.                     } else {
  134.                         pTotal += 2 * Math.nextDown(pLow);
  135.                     }
  136.                     criticalValueLow++;
  137.                     criticalValueHigh--;
  138.                 } else if (pLow < pHigh) {
  139.                     pTotal += pLow;
  140.                     criticalValueLow++;
  141.                 } else {
  142.                     pTotal += pHigh;
  143.                     criticalValueHigh--;
  144.                 }

  145.                 if (criticalValueLow > numberOfSuccesses || criticalValueHigh < numberOfSuccesses) {
  146.                     break;
  147.                 }
  148.             }
  149.             return pTotal;
  150.         default:
  151.             throw new MathInternalError(LocalizedFormats. OUT_OF_RANGE_SIMPLE, alternativeHypothesis,
  152.                       AlternativeHypothesis.TWO_SIDED, AlternativeHypothesis.LESS_THAN);
  153.         }
  154.     }
  155. }