BinomialProportion.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.  *      https://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.  */

  18. /*
  19.  * This is not the original file distributed by the Apache Software Foundation
  20.  * It has been modified by the Hipparchus project
  21.  */
  22. package org.hipparchus.stat.interval;

  24. import org.hipparchus.distribution.continuous.FDistribution;
  25. import org.hipparchus.distribution.continuous.NormalDistribution;
  26. import org.hipparchus.exception.LocalizedCoreFormats;
  27. import org.hipparchus.exception.MathIllegalArgumentException;
  28. import org.hipparchus.stat.LocalizedStatFormats;
  29. import org.hipparchus.util.FastMath;
  30. import org.hipparchus.util.MathUtils;

  32. /**
  33.  * Utility methods to generate confidence intervals for a binomial proportion.
  34.  *
  35.  * @see
  36.  * <a href="http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval">
  37.  * Binomial proportion confidence interval (Wikipedia)</a>
  38.  */
  39. public class BinomialProportion {

  41.     /**
  42.      * The standard normal distribution to calculate the inverse cumulative probability.
  43.      * Accessed and used in a thread-safe way.
  44.      */
  45.     private static final NormalDistribution NORMAL_DISTRIBUTION = new NormalDistribution(0, 1);

  47.     /** Utility class, prevent instantiation. */
  48.     private BinomialProportion() {}

  50.     /**
  51.      * Create an Agresti-Coull binomial confidence interval for the true
  52.      * probability of success of an unknown binomial distribution with
  53.      * the given observed number of trials, probability of success and
  54.      * confidence level.
  55.      * <p>
  56.      * Preconditions:
  57.      * <ul>
  58.      * <li>{@code numberOfTrials} must be positive</li>
  59.      * <li>{@code probabilityOfSuccess} must be between 0 and 1 (inclusive)</li>
  60.      * <li>{@code confidenceLevel} must be strictly between 0 and 1 (exclusive)</li>
  61.      * </ul>
  62.      *
  63.      * @see
  64.      * <a href="http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Agresti-Coull_Interval">
  65.      * Agresti-Coull interval (Wikipedia)</a>
  66.      *
  67.      * @param numberOfTrials number of trials
  68.      * @param probabilityOfSuccess observed probability of success
  69.      * @param confidenceLevel desired probability that the true probability of
  70.      * success falls within the returned interval
  71.      * @return Confidence interval containing the probability of success with
  72.      * probability {@code confidenceLevel}
  73.      * @throws MathIllegalArgumentException if {@code numberOfTrials <= 0}.
  74.      * @throws MathIllegalArgumentException if {@code probabilityOfSuccess} is not in the interval [0, 1].
  75.      * @throws MathIllegalArgumentException if {@code confidenceLevel} is not in the interval (0, 1).
  76.      */
  77.     public static ConfidenceInterval getAgrestiCoullInterval(int numberOfTrials,
  78.                                                              double probabilityOfSuccess,
  79.                                                              double confidenceLevel)
  80.         throws MathIllegalArgumentException {

  82.         checkParameters(numberOfTrials, probabilityOfSuccess, confidenceLevel);

  84.         final int numberOfSuccesses = (int) (numberOfTrials * probabilityOfSuccess);

  86.         final double alpha = (1.0 - confidenceLevel) / 2;
  87.         final double z = NORMAL_DISTRIBUTION.inverseCumulativeProbability(1 - alpha);
  88.         final double zSquared = FastMath.pow(z, 2);
  89.         final double modifiedNumberOfTrials = numberOfTrials + zSquared;
  90.         final double modifiedSuccessesRatio = (1.0 / modifiedNumberOfTrials) *
  91.                                               (numberOfSuccesses + 0.5 * zSquared);
  92.         final double difference = z * FastMath.sqrt(1.0 / modifiedNumberOfTrials *
  93.                                                     modifiedSuccessesRatio *
  94.                                                     (1 - modifiedSuccessesRatio));

  96.         return new ConfidenceInterval(modifiedSuccessesRatio - difference,
  97.                                       modifiedSuccessesRatio + difference,
  98.                                       confidenceLevel);
  99.     }

  101.     /**
  102.      * Create a Clopper-Pearson binomial confidence interval for the true
  103.      * probability of success of an unknown binomial distribution with
  104.      * the given observed number of trials, probability of success and
  105.      * confidence level.
  106.      * <p>
  107.      * Preconditions:
  108.      * <ul>
  109.      * <li>{@code numberOfTrials} must be positive</li>
  110.      * <li>{@code probabilityOfSuccess} must be between 0 and 1 (inclusive)</li>
  111.      * <li>{@code confidenceLevel} must be strictly between 0 and 1 (exclusive)</li>
  112.      * </ul>
  113.      *
  114.      * @see
  115.      * <a href="http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Clopper-Pearson_interval">
  116.      * Clopper-Pearson interval (Wikipedia)</a>
  117.      *
  118.      * @param numberOfTrials number of trials
  119.      * @param probabilityOfSuccess observed probability of success
  120.      * @param confidenceLevel desired probability that the true probability of
  121.      * success falls within the returned interval
  122.      * @return Confidence interval containing the probability of success with
  123.      * probability {@code confidenceLevel}
  124.      * @throws MathIllegalArgumentException if {@code numberOfTrials <= 0}.
  125.      * @throws MathIllegalArgumentException if {@code probabilityOfSuccess} is not in the interval [0, 1].
  126.      * @throws MathIllegalArgumentException if {@code confidenceLevel} is not in the interval (0, 1).
  127.      */
  128.     public static ConfidenceInterval getClopperPearsonInterval(int numberOfTrials,
  129.                                                                double probabilityOfSuccess,
  130.                                                                double confidenceLevel)
  131.         throws MathIllegalArgumentException {

  133.         checkParameters(numberOfTrials, probabilityOfSuccess, confidenceLevel);

  135.         double lowerBound = 0;
  136.         double upperBound = 0;
  137.         final int numberOfSuccesses = (int) (numberOfTrials * probabilityOfSuccess);

  139.         if (numberOfSuccesses > 0) {
  140.             final double alpha = (1.0 - confidenceLevel) / 2.0;

  142.             final FDistribution distributionLowerBound =
  143.                     new FDistribution(2 * (numberOfTrials - numberOfSuccesses + 1),
  144.                                       2 * numberOfSuccesses);

  146.             final double fValueLowerBound =
  147.                     distributionLowerBound.inverseCumulativeProbability(1 - alpha);
  148.             lowerBound = numberOfSuccesses /
  149.                          (numberOfSuccesses + (numberOfTrials - numberOfSuccesses + 1) * fValueLowerBound);

  151.             final FDistribution distributionUpperBound =
  152.                     new FDistribution(2 * (numberOfSuccesses + 1),
  153.                                       2 * (numberOfTrials - numberOfSuccesses));

  155.             final double fValueUpperBound =
  156.                     distributionUpperBound.inverseCumulativeProbability(1 - alpha);
  157.             upperBound = (numberOfSuccesses + 1) * fValueUpperBound /
  158.                          (numberOfTrials - numberOfSuccesses + (numberOfSuccesses + 1) * fValueUpperBound);
  159.         }

  161.         return new ConfidenceInterval(lowerBound, upperBound, confidenceLevel);
  162.     }

  164.     /**
  165.      * Create a binomial confidence interval using normal approximation
  166.      * for the true probability of success of an unknown binomial distribution
  167.      * with the given observed number of trials, probability of success and
  168.      * confidence level.
  169.      * <p>
  170.      * Preconditions:
  171.      * <ul>
  172.      * <li>{@code numberOfTrials} must be positive</li>
  173.      * <li>{@code probabilityOfSuccess} must be between 0 and 1 (inclusive)</li>
  174.      * <li>{@code confidenceLevel} must be strictly between 0 and 1 (exclusive)</li>
  175.      * </ul>
  176.      *
  177.      * @see
  178.      * <a href="http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Normal_approximation_interval">
  179.      * Normal approximation interval (Wikipedia)</a>
  180.      *
  181.      * @param numberOfTrials number of trials
  182.      * @param probabilityOfSuccess observed probability of success
  183.      * @param confidenceLevel desired probability that the true probability of
  184.      * success falls within the returned interval
  185.      * @return Confidence interval containing the probability of success with
  186.      * probability {@code confidenceLevel}
  187.      * @throws MathIllegalArgumentException if {@code numberOfTrials <= 0}.
  188.      * @throws MathIllegalArgumentException if {@code probabilityOfSuccess} is not in the interval [0, 1].
  189.      * @throws MathIllegalArgumentException if {@code confidenceLevel} is not in the interval (0, 1).
  190.      */
  191.     public static ConfidenceInterval getNormalApproximationInterval(int numberOfTrials,
  192.                                                                     double probabilityOfSuccess,
  193.                                                                     double confidenceLevel)
  194.         throws MathIllegalArgumentException {

  196.         checkParameters(numberOfTrials, probabilityOfSuccess, confidenceLevel);

  198.         final double mean = probabilityOfSuccess;
  199.         final double alpha = (1.0 - confidenceLevel) / 2;

  201.         final double difference = NORMAL_DISTRIBUTION.inverseCumulativeProbability(1 - alpha) *
  202.                                   FastMath.sqrt(1.0 / numberOfTrials * mean * (1 - mean));
  203.         return new ConfidenceInterval(mean - difference, mean + difference, confidenceLevel);
  204.     }

  206.     /**
  207.      * Create an Wilson score binomial confidence interval for the true
  208.      * probability of success of an unknown binomial distribution with
  209.      * the given observed number of trials, probability of success and
  210.      * confidence level.
  211.      * <p>
  212.      * Preconditions:
  213.      * <ul>
  214.      * <li>{@code numberOfTrials} must be positive</li>
  215.      * <li>{@code probabilityOfSuccess} must be between 0 and 1 (inclusive)</li>
  216.      * <li>{@code confidenceLevel} must be strictly between 0 and 1 (exclusive)</li>
  217.      * </ul>
  218.      *
  219.      * @see
  220.      * <a href="http://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Wilson_score_interval">
  221.      * Wilson score interval (Wikipedia)</a>
  222.      *
  223.      * @param numberOfTrials number of trials
  224.      * @param probabilityOfSuccess observed probability of success
  225.      * @param confidenceLevel desired probability that the true probability of
  226.      * success falls within the returned interval
  227.      * @return Confidence interval containing the probability of success with
  228.      * probability {@code confidenceLevel}
  229.      * @throws MathIllegalArgumentException if {@code numberOfTrials <= 0}.
  230.      * @throws MathIllegalArgumentException if {@code probabilityOfSuccess} is not in the interval [0, 1].
  231.      * @throws MathIllegalArgumentException if {@code confidenceLevel} is not in the interval (0, 1).
  232.      */
  233.     public static ConfidenceInterval getWilsonScoreInterval(int numberOfTrials,
  234.                                                             double probabilityOfSuccess,
  235.                                                             double confidenceLevel)
  236.         throws MathIllegalArgumentException {

  238.         checkParameters(numberOfTrials, probabilityOfSuccess, confidenceLevel);

  240.         final double alpha = (1.0 - confidenceLevel) / 2;
  241.         final double z = NORMAL_DISTRIBUTION.inverseCumulativeProbability(1 - alpha);
  242.         final double zSquared = FastMath.pow(z, 2);
  243.         final double mean = probabilityOfSuccess;

  245.         final double factor = 1.0 / (1 + (1.0 / numberOfTrials) * zSquared);
  246.         final double modifiedSuccessRatio = mean + (1.0 / (2 * numberOfTrials)) * zSquared;
  247.         final double difference =
  248.                 z * FastMath.sqrt(1.0 / numberOfTrials * mean * (1 - mean) +
  249.                                   (1.0 / (4 * FastMath.pow(numberOfTrials, 2)) * zSquared));

  251.         final double lowerBound = factor * (modifiedSuccessRatio - difference);
  252.         final double upperBound = factor * (modifiedSuccessRatio + difference);
  253.         return new ConfidenceInterval(lowerBound, upperBound, confidenceLevel);
  254.     }

  256.     /**
  257.      * Verifies that parameters satisfy preconditions.
  258.      *
  259.      * @param numberOfTrials number of trials (must be positive)
  260.      * @param probabilityOfSuccess probability of successes (must be between 0 and 1)
  261.      * @param confidenceLevel confidence level (must be strictly between 0 and 1)
  262.      * @throws MathIllegalArgumentException if {@code numberOfTrials <= 0}.
  263.      * @throws MathIllegalArgumentException if {@code probabilityOfSuccess is not in the interval [0, 1]}.
  264.      * @throws MathIllegalArgumentException if {@code confidenceLevel} is not in the interval (0, 1)}.
  265.      */
  266.     private static void checkParameters(int numberOfTrials,
  267.                                         double probabilityOfSuccess,
  268.                                         double confidenceLevel) {
  269.         if (numberOfTrials <= 0) {
  270.             throw new MathIllegalArgumentException(LocalizedCoreFormats.NUMBER_OF_TRIALS,
  271.                                                    numberOfTrials);
  272.         }
  273.         MathUtils.checkRangeInclusive(probabilityOfSuccess, 0, 1);
  274.         if (confidenceLevel <= 0 || confidenceLevel >= 1) {
  275.             throw new MathIllegalArgumentException(LocalizedStatFormats.OUT_OF_BOUNDS_CONFIDENCE_LEVEL,
  276.                                                    confidenceLevel, 0, 1);
  277.         }
  278.     }

  280. }
Created with JaCoCo 0.8.11.202310140853