Skip to content

Instantly share code, notes, and snippets.

@komiya-atsushi
Created January 18, 2015 15:55
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save komiya-atsushi/141eeb2272367a5c72d0 to your computer and use it in GitHub Desktop.
Save komiya-atsushi/141eeb2272367a5c72d0 to your computer and use it in GitHub Desktop.
割合が 0% or 100% でも信頼区間を計算できるようにした Clopper-Pearson の実装
package org.apache.commons.math3.stat.interval;
import org.apache.commons.math3.distribution.FDistribution;
import org.apache.commons.math3.exception.NotPositiveException;
import org.apache.commons.math3.exception.NotStrictlyPositiveException;
import org.apache.commons.math3.exception.NumberIsTooLargeException;
import org.apache.commons.math3.exception.OutOfRangeException;
/**
* numberOfSuccesses が 0 もしくは numberOfTrials と等しい場合でも、
* とりあえず区間の一方のエンドポイントを算出できるようにした
* Clopper-Pearson 'Exact' の実装です。
*/
public class ModifiedClopperPearsonInterval implements BinomialConfidenceInterval {
@Override
public ConfidenceInterval createInterval(int numberOfTrials, int numberOfSuccesses, double confidenceLevel) throws NotStrictlyPositiveException, NotPositiveException, NumberIsTooLargeException, OutOfRangeException {
double lowerBound = 0;
double upperBound = 1;
final double alpha = (1.0 - confidenceLevel) / 2.0;
if (numberOfSuccesses > 0) {
final FDistribution distributionLowerBound = new FDistribution(2 * (numberOfTrials - numberOfSuccesses + 1),
2 * numberOfSuccesses);
final double fValueLowerBound = distributionLowerBound.inverseCumulativeProbability(1 - alpha);
lowerBound = numberOfSuccesses /
(numberOfSuccesses + (numberOfTrials - numberOfSuccesses + 1) * fValueLowerBound);
}
if (numberOfSuccesses < numberOfTrials) {
final FDistribution distributionUpperBound = new FDistribution(2 * (numberOfSuccesses + 1),
2 * (numberOfTrials - numberOfSuccesses));
final double fValueUpperBound = distributionUpperBound.inverseCumulativeProbability(1 - alpha);
upperBound = (numberOfSuccesses + 1) * fValueUpperBound /
(numberOfTrials - numberOfSuccesses + (numberOfSuccesses + 1) * fValueUpperBound);
}
return new ConfidenceInterval(lowerBound, upperBound, confidenceLevel);
}
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment