Library Bernoulli

Bernoulli.v: Simulating Bernoulli and Binomial distributions

Require Export Cover.
Require Export Misc.

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 : ∀ x : U, {x < ½ }+{½ ≤ 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)).

Lemma Fbern_simpl : ∀ 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).

Definition bernoulli : U → distr bool := Mfix Fbern.

Properties of the Bernoulli program

Proofs using fixpoint rules

Definition Mubern (q: bool → U) : MF U -m > MF U.
Defined.

Lemma Mubern_simpl : ∀ (q: bool → U) f p,
  Mubern q f p = if dec_demi p then ½ *(q false )+½ *(f (p+p))
                           else ½ *(f (p&p)) + ½ *(q true ).

Lemma Mubern_eq : ∀ (q: bool → U) (f:U → distr bool) (p:U),
             μ (Fbern f p) q ≈ Mubern q (fun y ⇒ μ (f y) q) p.
Hint Resolve Mubern_eq.

Lemma Bern_eq :
    ∀ q : bool → U, ∀ p, μ (bernoulli p) q ≈ mufix (Mubern q) p.
Hint Resolve Bern_eq.

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

Lemma Bern_term : ∀ p, μ (bernoulli p) (fone bool) ≈ 1.
Hint Resolve Bern_term.

p is an invariant of Mubern qtrue

Lemma MuBern_true : ∀ p, Mubern B2U (fun q ⇒ q) p ≈ p.
Hint Resolve MuBern_true.

Lemma MuBern_false : ∀ p, Mubern (finv B2U) (finv (fun q ⇒ q)) p ≈ [1-]p.
Hint Resolve MuBern_false.

Lemma Bern_true : ∀ p, μ (bernoulli p) B2U ≈ p.

Lemma Bern_false : ∀ p, μ (bernoulli p) NB2U ≈ [1-]p.

Lemma Mubern_inv : ∀ (q: bool → U) (f:U → U) (p:U),
Mubern (finv q) (finv f) p ≈ [1-] Mubern q f p.

Proofs using lubs

Invariant pmin p pmin p n = p - ½ ^ n

Property : ∀ p, ok p (bernoulli p) χ (.=true)

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

Lemma bernoulli_true : okfun (fun p ⇒ p) bernoulli qtrue.

Property : ∀ p, ok (1-p) (bernoulli p) (χ (.=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 : ∀ (b:bool) (x:U), qtrue x b ≈ [1-] (qfalse x b).

Lemma bernoulli_eq_true : ∀ p, μ (bernoulli p) (qtrue p) ≈ p.

Lemma bernoulli_eq_false : ∀ p, μ (bernoulli p) (qfalse p) ≈ [1-]p.

Lemma bernoulli_eq : ∀ p f, μ (bernoulli p) f ≈ p * f true + ([1-]p) * f false.

Lemma bernoulli_total : ∀ p , μ (bernoulli p) (fone bool) ≈ 1.

Binomial distribution

(binomial p n) gives k with probability C(k,n) p^k(1-p)^(n-k)

Definition and properties of binomial coefficients

Fixpoint comb (k n:nat) {struct n} : nat :=
         match n with O ⇒ match k with O ⇒ (1%nat) | (S l) ⇒ O end
                | (S m) ⇒ match k with O ⇒ (1%nat)
                                                    | (S l) ⇒ ((comb l m) + (comb k m))%nat end
         end.

Lemma comb_0_n : ∀ n, comb 0 n = 1%nat.

Lemma comb_not_le : ∀ n k, le (S n) k → comb k n=0%nat.

Lemma comb_Sn_n : ∀ n, comb (S n) n =0%nat.

Lemma comb_n_n : ∀ n, comb n n = (1%nat).

Lemma comb_1_Sn : ∀ n, comb 1 (S n) = (S n).

Lemma comb_inv : ∀ n k, (k ≤ n)%nat → comb k n = comb (n-k) n.

Lemma comb_n_Sn : ∀ n, comb n (S n) = (S n).

Definition fc (p:U)(n k:nat) := (comb k n) */ (p^k * ([1-]p)^(n-k)).

Lemma fcp_0 : ∀ p n, fc p n O ≈ ([1-]p)^n.

Lemma fcp_n : ∀ p n, fc p n n ≈ p^n.

Lemma fcp_not_le : ∀ p n k, (S n ≤ k)%nat → fc p n k ≈ 0.

Lemma fc0 : ∀ n k, fc 0 n (S k) ≈ 0.
Hint Resolve fc0.

Add Morphism fc with signature Oeq ⇾ eq ⇾ eq ⇾ Oeq
as fc_eq_compat.
Save.

Hint Resolve fc_eq_compat.

Lemma sigma_fc0 : ∀ n k, sigma (fc 0 n) (S k) ≈ 1.

Lemma retract_class : ∀ f n, class (retract f n).
Hint Resolve retract_class.

Lemma fc_retract :
     ∀ p n, ([1-]p^n ≈ sigma (fc p n) n) → retract (fc p n) (S n).
Hint Resolve fc_retract.

Lemma fc_Nmult_def :
     ∀ p n k, ([1-]p^n ≈ sigma (fc p n) n) → Nmult_def (comb k n) (p^k * ([1-]p) ^(n-k)).
Hint Resolve fc_Nmult_def.

Lemma fcp_S :
    ∀ p n k, ([1-]p^n ≈ sigma (fc p n) n) → fc p (S n) (S k) ≈ p * (fc p n k) + ([1-]p) * (fc p n (S k)).

Lemma sigma_fc_1 : ∀ p n, ([1-]p^n ≈ sigma (fc p n) n) ->1 ≈ sigma (fc p n) (S n).
Hint Resolve sigma_fc_1.

Lemma Uinv_exp : ∀ p n, [1-](p^n) ≈ sigma (fc p n) n.

Hint Resolve Uinv_exp.

Lemma Nmult_comb : ∀ p n k, Nmult_def (comb k n) (p ^ k * ([1-] p) ^ (n - k)).
Hint Resolve Nmult_comb.

Definition of binomial distribution

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 binomial distribution

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

Lemma binomial_le_n :
   ∀ p n, 1 ≤ μ (binomial p n) (carac_le n).

Lemma binomial_le_S : ∀ p n k,
   μ (binomial p (S n)) (carac_le (S k)) ≈
   p * (μ (binomial p n) (carac_le k)) + ([1-]p) * (μ (binomial p n) (carac_le (S k))).

Lemma binomial_lt_S : ∀ p n k,
   μ (binomial p (S n)) (carac_lt (S k)) ≈
   p * (μ (binomial p n) (carac_lt k)) + ([1-]p) * (μ (binomial p n) (carac_lt (S k))).