lib/distribution/binomial/ruby.rb

ruby.rb

        # Returns the inverse-CDF or the quantile given the probability `prob`,
        # the total number of trials `n` and the number of successes `k`
        # Note: This is a binary search under the hood and is a candidate for
        #   updates once more stable techniques are found.
        #   This still needs to be optimized and benchmarked against 
        #   the vanilla quantile  
        #
        # @paran qn [Float] the cumulative function value to be inverted
        # @param n [Fixnum, Bignum] total number of trials
        # @param prob [Float] probabilty of success in a single independant trial
        #
        # @return [Fixnum, Bignum] the integer quantile `k` given cumulative value
        #
        # @raise RangeError if qn is from outside of the closed interval [0, 1]
        def quantile(qn, n, prob)
          fail RangeError, 'cdf value(qn) must be from [0, 1]. '\
            "Cannot find quantile for qn=#{qn}" if qn > 1 || qn < 0

          return n if qn == 1.0 || qn == 1
          return 0 if qn == 0.0 || qn == 0
          low, high = 0, n
          
          # Search predicate: qn <= CDF(inverse_val)
          predicate = lambda {|k| qn <= exact_cdf(k, n, prob) }
          mid_proportion = prob
          return __binary_search_inv(low, high, predicate, mid_proportion)
        end
        
        private
        # :nodoc:
        # Inverts a function predicate via a binary search
        # Accepts bounds and a predicate (lambda resulting in boolean values)
        # to look for the variable
        def __binary_search_inv(lower, upper, predicate, mid_proportion)
          # The smallest value in the given domain satisfies our needs
          return lower if predicate.call(lower)
          
          # Sucessively squeeze the bounds so that we find the first value 
          # in the domain where the predicate is satisfied
          while lower <= upper
            mid = lower + ((upper - lower) * mid_proportion).floor
            if predicate.call(mid)
              upper = mid - 1
            else
              lower = mid + 1
            end
          end
          
          return lower
        end