## Defining Bernoulli and Binomial distributions from flip

Definition and properties of Bernoulli and Binomial distributions using the ALEA library.

Authors: Christine Paulin-Mohring

Topics: Randomized Algorithms

Tools: Alea / Coq

see also the index (by topic, by tool, by reference, by year)

# Bernoulli.v: Simulating Bernoulli and Binomial distributions

Require Export Cover.
Require Export Misc.
Require Export BinCoeff.

## Program for computing a Bernoulli distribution

bernoulli p gives true with probability p and false with probability (1-p)
```let rec bernoulli p =
if flip
then (if p < 1/2 then false else bernoulli (2 p - 1))
else (if p < 1/2 then bernoulli (2 p) else true)
```

Hypothesis dec_demi : forall x : U, {x < [1/2]}+{[1/2] <= x}.

Instance Fbern_mon : monotonic
(fun (f:U -> distr bool) p =>
Mif Flip
(if dec_demi p then Munit false else f (p & p))
(if dec_demi p then f (p + p) else Munit true)).
Save.

Definition Fbern : (U -> distr bool) -m> (U -> distr bool)
:= mon (fun f p => Mif Flip
(if dec_demi p then Munit false else f (p & p))
(if dec_demi p then f (p + p) else Munit true)).

Definition bernoulli : U -> distr bool := Mfix Fbern.

Recursive definition of binomial distribution using bernoulli
(binomial p n) gives k with probability C(k,n) p^k (1-p)^(n-k)

Fixpoint binomial (p:U)(n:nat) {struct n}: distr nat :=
match n with O => Munit O
| S m => Mlet (binomial p m)
(fun x => Mif (bernoulli p) (Munit (S x)) (Munit x))
end.

## Properties of the Bernoulli program

Lemma Fbern_simpl : forall f p,
Fbern f p = Mif Flip
(if dec_demi p then Munit false else f (p & p))
(if dec_demi p then f (p + p) else Munit true).

### Proofs using fixpoint rules

Instance Mubern_mon : forall (q: bool -> U),
monotonic
(fun bern (p:U) => if dec_demi p then [1/2]*(q false)+[1/2]*(bern (p+p))
else [1/2]*(bern (p&p)) + [1/2]*(q true)).
Save.

Definition Mubern (q: bool -> U) : MF U -m> MF U
:= mon (fun bern (p:U) => if dec_demi p then [1/2]*(q false)+[1/2]*(bern (p+p))
else [1/2]*(bern (p&p)) + [1/2]*(q true)).

Lemma Mubern_simpl : forall (q: bool -> U) f p,
Mubern q f p = if dec_demi p then [1/2]*(q false)+[1/2]*(f (p+p))
else [1/2]*(f (p&p)) + [1/2]*(q true).

Mubern commutes with the measure of Fbern
Lemma Mubern_eq : forall (q: bool -> U) (f:U -> distr bool) (p:U),
mu (Fbern f p) q == Mubern q (fun y => mu (f y) q) p.
Hint Resolve Mubern_eq.

Lemma Bern_eq :
forall q : bool -> U, forall p, mu (bernoulli p) q == mufix (Mubern q) p.
Hint Resolve Bern_eq.

Lemma Bern_commute : forall q : bool -> U,
mu_muF_commute_le Fbern (fun (x:U) => q) (Mubern q).
Hint Resolve Bern_commute.

bernoulli terminates with probability 1
Lemma Bern_term : forall p, mu (bernoulli p) (fone bool) == 1.
Hint Resolve Bern_term.

### p is an invariant of Mubern qtrue

Lemma MuBern_true : forall p, Mubern B2U (fun q => q) p == p.
Hint Resolve MuBern_true.

Lemma MuBern_false : forall p, Mubern (finv B2U) (finv (fun q => q)) p == [1-]p.
Hint Resolve MuBern_false.

Lemma Mubern_inv : forall (q: bool -> U) (f:U -> U) (p:U),
Mubern (finv q) (finv f) p == [1-] Mubern q f p.

prob(bernoulli = true) = p
Lemma Bern_true : forall p, mu (bernoulli p) B2U == p.

prob(bernoulli = false) = 1-p
Lemma Bern_false : forall p, mu (bernoulli p) NB2U == [1-]p.

### Direct proofs using lubs

Invariant pmin p with pmin p n = p - [1/2] ^ n

Property : forall p, ok p (bernoulli p) chi (.=true)

Definition qtrue (p:U) := B2U.
Definition qfalse (p:U) := NB2U.

Lemma bernoulli_true : okfun (fun p => p) bernoulli qtrue.

Property : forall p, ok (1-p) (bernoulli p) (chi (.=false))

Lemma bernoulli_false : okfun (fun p => [1-] p) bernoulli qfalse.

Probability for the result of (bernoulli p) to be true is exactly p

Lemma qtrue_qfalse_inv : forall (b:bool) (x:U), qtrue x b == [1-] (qfalse x b).

Lemma bernoulli_eq_true : forall p, mu (bernoulli p) (qtrue p) == p.

Lemma bernoulli_eq_false : forall p, mu (bernoulli p) (qfalse p)== [1-]p.

Lemma bernoulli_eq : forall p f,
mu (bernoulli p) f == p * f true + ([1-]p) * f false.

Lemma bernoulli_total : forall p , mu (bernoulli p) (fone bool)==1.

## Properties of Binomial distribution

prob(binomial p n = k) = C(k,n) p ^ k (1-p)^(n-k)

Lemma binomial_eq_k :
p n k, μ (binomial p n) (carac_eq k) fc p n k.

prob(binomial p nn) = 1
Lemma binomial_le_n :
p n, 1 μ (binomial p n) (carac_le n).

prob(binomial p (S n) <= S k) = p prob(binomial p n <= k) + (1-p) prob(binomial p n <= S k)
Lemma binomial_le_S : forall p n k,
mu (binomial p (S n)) (carac_le (S k)) ==
p * (mu (binomial p n) (carac_le k)) + ([1-]p) * (mu (binomial p n) (carac_le (S k))).

prob(binomial p (S n) < S k) = p prob(binomial p n < k) + (1-p) prob(binomial p n < S k)
Lemma binomial_lt_S : forall p n k,
mu (binomial p (S n)) (carac_lt (S k)) ==
p * (mu (binomial p n) (carac_lt k)) + ([1-]p) * (mu (binomial p n) (carac_lt (S k))).