001/*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements.  See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License.  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 */
017package org.apache.commons.math3.stat.inference;
018
019import org.apache.commons.math3.distribution.BinomialDistribution;
020import org.apache.commons.math3.exception.MathIllegalArgumentException;
021import org.apache.commons.math3.exception.MathInternalError;
022import org.apache.commons.math3.exception.NotPositiveException;
023import org.apache.commons.math3.exception.NullArgumentException;
024import org.apache.commons.math3.exception.OutOfRangeException;
025import org.apache.commons.math3.exception.util.LocalizedFormats;
026
027/**
028 * Implements binomial test statistics.
029 * <p>
030 * Exact test for the statistical significance of deviations from a
031 * theoretically expected distribution of observations into two categories.
032 *
033 * @see <a href="http://en.wikipedia.org/wiki/Binomial_test">Binomial test (Wikipedia)</a>
034 * @version $Id: BinomialTest.html 885258 2013-11-03 02:46:49Z tn $
035 * @since 3.3
036 */
037public class BinomialTest {
038
039    /**
040     * Returns whether the null hypothesis can be rejected with the given confidence level.
041     * <p>
042     * <strong>Preconditions</strong>:
043     * <ul>
044     * <li>Number of trials must be &ge; 0.</li>
045     * <li>Number of successes must be &ge; 0.</li>
046     * <li>Number of successes must be &le; number of trials.</li>
047     * <li>Probability must be &ge; 0 and &le; 1.</li>
048     * </ul>
049     *
050     * @param numberOfTrials number of trials performed
051     * @param numberOfSuccesses number of successes observed
052     * @param probability assumed probability of a single trial under the null hypothesis
053     * @param alternativeHypothesis type of hypothesis being evaluated (one- or two-sided)
054     * @param alpha significance level of the test
055     * @return true if the null hypothesis can be rejected with confidence {@code 1 - alpha}
056     * @throws NotPositiveException if {@code numberOfTrials} or {@code numberOfSuccesses} is negative
057     * @throws OutOfRangeException if {@code probability} is not between 0 and 1
058     * @throws MathIllegalArgumentException if {@code numberOfTrials} &lt; {@code numberOfSuccesses} or
059     * if {@code alternateHypothesis} is null.
060     * @see AlternativeHypothesis
061     */
062    public boolean binomialTest(int numberOfTrials, int numberOfSuccesses, double probability,
063                                AlternativeHypothesis alternativeHypothesis, double alpha) {
064        double pValue = binomialTest(numberOfTrials, numberOfSuccesses, probability, alternativeHypothesis);
065        return pValue < alpha;
066    }
067
068    /**
069     * Returns the <i>observed significance level</i>, or
070     * <a href="http://www.cas.lancs.ac.uk/glossary_v1.1/hyptest.html#pvalue">p-value</a>,
071     * associated with a <a href="http://en.wikipedia.org/wiki/Binomial_test"> Binomial test</a>.
072     * <p>
073     * The number returned is the smallest significance level at which one can reject the null hypothesis.
074     * The form of the hypothesis depends on {@code alternativeHypothesis}.</p>
075     * <p>
076     * The p-Value represents the likelihood of getting a result at least as extreme as the sample,
077     * given the provided {@code probability} of success on a single trial. For single-sided tests,
078     * this value can be directly derived from the Binomial distribution. For the two-sided test,
079     * the implementation works as follows: we start by looking at the most extreme cases
080     * (0 success and n success where n is the number of trials from the sample) and determine their likelihood.
081     * The lower value is added to the p-Value (if both values are equal, both are added). Then we continue with
082     * the next extreme value, until we added the value for the actual observed sample.</p>
083     * <p>
084     * <strong>Preconditions</strong>:
085     * <ul>
086     * <li>Number of trials must be &ge; 0.</li>
087     * <li>Number of successes must be &ge; 0.</li>
088     * <li>Number of successes must be &le; number of trials.</li>
089     * <li>Probability must be &ge; 0 and &le; 1.</li>
090     * </ul></p>
091     *
092     * @param numberOfTrials number of trials performed
093     * @param numberOfSuccesses number of successes observed
094     * @param probability assumed probability of a single trial under the null hypothesis
095     * @param alternativeHypothesis type of hypothesis being evaluated (one- or two-sided)
096     * @return p-value
097     * @throws NotPositiveException if {@code numberOfTrials} or {@code numberOfSuccesses} is negative
098     * @throws OutOfRangeException if {@code probability} is not between 0 and 1
099     * @throws MathIllegalArgumentException if {@code numberOfTrials} &lt; {@code numberOfSuccesses} or
100     * if {@code alternateHypothesis} is null.
101     * @see AlternativeHypothesis
102     */
103    public double binomialTest(int numberOfTrials, int numberOfSuccesses, double probability,
104                               AlternativeHypothesis alternativeHypothesis) {
105        if (numberOfTrials < 0) {
106            throw new NotPositiveException(numberOfTrials);
107        }
108        if (numberOfSuccesses < 0) {
109            throw new NotPositiveException(numberOfSuccesses);
110        }
111        if (probability < 0 || probability > 1) {
112            throw new OutOfRangeException(probability, 0, 1);
113        }
114        if (numberOfTrials < numberOfSuccesses) {
115            throw new MathIllegalArgumentException(
116                LocalizedFormats.BINOMIAL_INVALID_PARAMETERS_ORDER,
117                numberOfTrials, numberOfSuccesses);
118        }
119        if (alternativeHypothesis == null) {
120            throw new NullArgumentException();
121        }
122
123        final BinomialDistribution distribution = new BinomialDistribution(numberOfTrials, probability);
124        switch (alternativeHypothesis) {
125        case GREATER_THAN:
126            return 1 - distribution.cumulativeProbability(numberOfSuccesses - 1);
127        case LESS_THAN:
128            return distribution.cumulativeProbability(numberOfSuccesses);
129        case TWO_SIDED:
130            int criticalValueLow = 0;
131            int criticalValueHigh = numberOfTrials;
132            double pTotal = 0;
133
134            while (true) {
135                double pLow = distribution.probability(criticalValueLow);
136                double pHigh = distribution.probability(criticalValueHigh);
137
138                if (pLow == pHigh) {
139                    pTotal += 2 * pLow;
140                    criticalValueLow++;
141                    criticalValueHigh--;
142                } else if (pLow < pHigh) {
143                    pTotal += pLow;
144                    criticalValueLow++;
145                } else {
146                    pTotal += pHigh;
147                    criticalValueHigh--;
148                }
149
150                if (criticalValueLow > numberOfSuccesses || criticalValueHigh < numberOfSuccesses) {
151                    break;
152                }
153            }
154            return pTotal;
155        default:
156            throw new MathInternalError(LocalizedFormats. OUT_OF_RANGE_SIMPLE, alternativeHypothesis,
157                      AlternativeHypothesis.TWO_SIDED, AlternativeHypothesis.LESS_THAN);
158        }
159    }
160}