private int inverseProbability()

in commons-statistics-distribution/src/main/java/org/apache/commons/statistics/distribution/AbstractDiscreteDistribution.java [148:201]

35 lines of code
10 McCabe index (conditional complexity)


    private int inverseProbability(double p, double q, boolean complement) {

        int lower = getSupportLowerBound();
        if (p == 0) {
            return lower;
        }
        int upper = getSupportUpperBound();
        if (q == 0) {
            return upper;
        }

        // The binary search sets the upper value to the mid-point
        // based on fun(x) >= 0. The upper value is returned.
        //
        // Create a function to search for x where the upper bound can be
        // lowered if:
        // cdf(x) >= p
        // sf(x)  <= q
        final IntUnaryOperator fun = complement ?
            x -> Double.compare(q, survivalProbability(x)) :
            x -> Double.compare(cumulativeProbability(x), p);

        if (lower == Integer.MIN_VALUE) {
            if (fun.applyAsInt(lower) >= 0) {
                return lower;
            }
        } else {
            // this ensures:
            // cumulativeProbability(lower) < p
            // survivalProbability(lower) > q
            // which is important for the solving step
            lower -= 1;
        }

        // use the one-sided Chebyshev inequality to narrow the bracket
        // cf. AbstractContinuousDistribution.inverseCumulativeProbability(double)
        final double mu = getMean();
        final double sig = Math.sqrt(getVariance());
        final boolean chebyshevApplies = Double.isFinite(mu) &&
                                         ArgumentUtils.isFiniteStrictlyPositive(sig);

        if (chebyshevApplies) {
            double tmp = mu - sig * Math.sqrt(q / p);
            if (tmp > lower) {
                lower = ((int) Math.ceil(tmp)) - 1;
            }
            tmp = mu + sig * Math.sqrt(p / q);
            if (tmp < upper) {
                upper = ((int) Math.ceil(tmp)) - 1;
            }
        }

        return solveInverseProbability(fun, lower, upper);
    }