guile-math-examples/probability/probability.scm

(library (probability (0 0 1))
  (export factorial
          n-choose-k
          probability-choose-k-out-of-n
          probability-min-choose-k-out-of-n
          probability-max-choose-k-out-of-n)
  (import
    (except (rnrs base) let-values)
    (only (guile)
          ;; lambda forms
          lambda* λ)
    (prefix (srfi srfi-1) srfi-1:)))


(define factorial
  (λ (n)
    "Calculate the factorial of the given number n."
    (define iter
      (λ (n product)
        (cond
         [(< n 2) product]
         [else
          (iter (- n 1) (* n product))])))
    (iter n 1)))


(define n-choose-k
  (λ (n k)
    "Calculate the number of possible ways to choose k out
of n objects."
    (/ (factorial n)
       (* (factorial (- n k))
          (factorial k)))))


(define probability-choose-k-out-of-n
  (lambda* (n k #:key (k-probability (/ 1 2)))
    "Calculate the probability to choose exactly K objects
of one specific type T out of N objects. The probability to
choose one object of type T is by default one half. It can
optionally be given as keyword argument K-PROBABILITY. The
calculation is done assuming, that the probability of
choosing an object of type T does not change, for example by
assuming, that chosen objects are put back into the pool of
objects, from which objects are chosen.

Another way of looking at this function is the following:
'Given a scenario, in which one chooses N items from an
inexhaustable source, where the probability of choosing an
item of type T is K-PROBABILITY, what is the chance of
choosing K items of type T?'

Example scenario: Roll a die N times, having a K-PROBABILITY
of 1/2 to roll an even number, what is the probability of
rolling K even numbers? -- The die is an inexhaustable
source of numbers from 1 to 6, the probability of rolling an
even number never changes."
    (* (n-choose-k n k)
       (expt k-probability k)
       (expt (- 1 k-probability) (- n k)))))


(define probability-min-choose-k-out-of-n
  (lambda* (n k #:key (k-probability (/ 1 2)))
    "Calculate the probability to choose at least k objects
out of n, optionally with given chance to pick a k object
out of the n objects."
    (let next-k ([k k] [sum 0])
      (cond
       [(< k n)
        (next-k (+ k 1)
                (+ sum
                   (probability-choose-k-out-of-n n k #:k-probability k-probability)))]
       [else
        (+ sum
           (expt k-probability n))]))))


(define probability-max-choose-k-out-of-n
  (lambda* (n k #:key (k-probability (/ 1 2)))
    "Calculate the probability to choose at most k objects
out of n, optionally with given chance to pick a k object
out of the n objects."
    (let next-k ([k k] [sum 0])
      (cond
       [(> k 0)
        (next-k (- k 1)
                (+ sum
                   (probability-choose-k-out-of-n n k #:k-probability k-probability)))]
       [else
        (+ sum
           (expt (- 1 k-probability) n))]))))