Library QBS.normal_algebra
From HB Require Import structures.
From mathcomp Require Import all_boot all_algebra reals ereal topology
normedtype numfun measure sequences lebesgue_integral lebesgue_integral_fubini
lebesgue_stieltjes_measure probability boolp trigo exp.
From mathcomp.algebra_tactics Require Import ring.
Import GRing.Theory Num.Def Num.Theory.
Set Implicit Arguments.
Unset Strict Implicit.
Unset Printing Implicit Defensive.
Import order.Order.POrderTheory.
Section normal_density_algebra.
Variable R : realType.
Local Open Scope ring_scope.
Implicit Types (m s x K : R) (n : nat).
Let half_neq0 : (2%:R^-1 : R) != 0.
Proof. by rewrite invr_neq0 // pnatr_eq0. Qed.
Let three_neq0 : (3%:R : R) != 0.
Proof. by rewrite pnatr_eq0. Qed.
Local Lemma sqrtr_pnat_div_neq0 (a b : nat) :
(0 < a)%N → (0 < b)%N → sqrtr (a%:R / b%:R : R) != 0.
Proof.
move⇒ ha hb; apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ltr0n.
Qed.
Section normal_pdf_general.
Lemma complete_the_square (a b c x : R) (ha : a != 0) :
a × x ^+ 2 + b × x + c =
a × (x + b / (a *+ 2)) ^+ 2 - (b ^+ 2 - a *+ 4 × c) / (a *+ 4).
Proof.
have ha4 : a *+ 4 != 0 by rewrite mulrn_eq0.
have ha2 : a *+ 2 != 0 by rewrite mulrn_eq0.
field. repeat split ⇒ //; rewrite ?mulrn_eq0 //.
Qed.
Lemma normal_pdf_times (m m' s s' x : R) :
s != 0 → s' != 0 →
normal_pdf m s x × normal_pdf m' s' x =
normal_peak (sqrtr (s ^+ 2 + s' ^+ 2)) ×
normal_fun m (sqrtr (s ^+ 2 + s' ^+ 2)) m' ×
normal_pdf ((m × s' ^+ 2 + m' × s ^+ 2) / (s ^+ 2 + s' ^+ 2))
(s × s' / sqrtr (s ^+ 2 + s' ^+ 2)) x.
Proof.
move⇒ hs hs'.
have hS : s ^+ 2 + s' ^+ 2 != 0.
by rewrite paddr_eq0 ?sqr_ge0 // !sqrf_eq0; apply/nandP; left.
have hSgt0 : 0 < s ^+ 2 + s' ^+ 2.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
by apply: lt0r_neq0; rewrite sqrtr_gt0.
have hsnew : s × s' / sqrtr (s ^+ 2 + s' ^+ 2) != 0.
by rewrite mulf_neq0 ?invr_neq0 // mulf_neq0.
have hpi : 0 ≤ pi := pi_ge0 R.
have hs2pi : 0 ≤ s ^+ 2 × pi := mulr_ge0 (sqr_ge0 _) hpi.
have hS2pi : 0 ≤ (s ^+ 2 + s' ^+ 2) × pi.
by apply: mulr_ge0 ⇒ //; exact: addr_ge0 (sqr_ge0 _) (sqr_ge0 _).
rewrite (normal_pdfE _ hs) (normal_pdfE _ hs') (normal_pdfE _ hsnew).
rewrite [LHS]mulrACA [RHS]mulrACA.
congr (_ × _).
-
rewrite /normal_peak -invfM -sqrtrM ?pmulrn_lge0 ?sqr_ge0 //.
rewrite -[RHS]invfM -sqrtrM ?pmulrn_lge0 ?sqr_ge0 ?sqr_sqrtr
?addr_ge0 ?sqr_ge0 //.
congr (_^-1); congr (sqrtr _).
rewrite exprMn exprVn sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
by field; exact: hS.
-
rewrite /normal_fun -expRD -[RHS]expRD.
congr (expR _).
rewrite sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
rewrite exprMn exprVn sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
field.
by rewrite hS hs hs'.
Qed.
Definition gaussian_prod_mu (m m' s s' : R) :=
(m × s' ^+ 2 + m' × s ^+ 2) / (s ^+ 2 + s' ^+ 2).
Definition gaussian_prod_sigma (s s' : R) :=
s × s' / sqrtr (s ^+ 2 + s' ^+ 2).
Definition gaussian_prod_scalar (m m' s s' : R) :=
normal_peak (sqrtr (s ^+ 2 + s' ^+ 2)) ×
normal_fun m (sqrtr (s ^+ 2 + s' ^+ 2)) m'.
Lemma normal_pdf_times' (m m' s s' x : R) :
s != 0 → s' != 0 →
normal_pdf m s x × normal_pdf m' s' x =
gaussian_prod_scalar m m' s s' ×
normal_pdf (gaussian_prod_mu m m' s s') (gaussian_prod_sigma s s') x.
Proof.
move⇒ hs hs'.
by rewrite /gaussian_prod_scalar /gaussian_prod_mu /gaussian_prod_sigma
normal_pdf_times.
Qed.
Lemma normal_pdf_times_chain (K m m' s s' x : R) :
s != 0 → s' != 0 →
K × normal_pdf m s x × normal_pdf m' s' x =
K × gaussian_prod_scalar m m' s s' ×
normal_pdf (gaussian_prod_mu m m' s s') (gaussian_prod_sigma s s') x.
Proof.
move⇒ hs hs'.
by rewrite -mulrA normal_pdf_times' // mulrA.
Qed.
Lemma gaussian_prod_sigma_sqr (s s' : R) :
s != 0 → s' != 0 →
(gaussian_prod_sigma s s') ^+ 2 =
s ^+ 2 × s' ^+ 2 / (s ^+ 2 + s' ^+ 2).
Proof.
move⇒ hs hs'.
have hS : s ^+ 2 + s' ^+ 2 != 0.
by rewrite paddr_eq0 ?sqr_ge0 // !sqrf_eq0; apply/nandP; left.
rewrite /gaussian_prod_sigma exprMn exprVn.
rewrite sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
by field.
Qed.
Lemma gaussian_prod_sigma_neq0 (s s' : R) :
s != 0 → s' != 0 → gaussian_prod_sigma s s' != 0.
Proof.
move⇒ hs hs'.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
by rewrite /gaussian_prod_sigma mulf_neq0 ?invr_neq0 // mulf_neq0.
Qed.
Lemma gaussian_prod_scalar_gt0 (m m' s s' : R) :
s != 0 → s' != 0 → 0 < gaussian_prod_scalar m m' s s'.
Proof.
move⇒ hs hs'.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
rewrite /gaussian_prod_scalar.
apply: mulr_gt0; first exact: (normal_peak_gt0 hsqS).
rewrite /normal_fun; exact: expR_gt0.
Qed.
End normal_pdf_general.
Section bayesian_phase1_intercept.
Lemma var_recurrence (V : R) :
0 < V →
sqrtr V ^+ 2 × (2%:R ^-1) ^+ 2 / (sqrtr V ^+ 2 + (2%:R ^-1) ^+ 2) =
V / (V *+ 4 + 1).
Proof.
move⇒ hV.
have hV0 : V != 0 by apply: lt0r_neq0.
have hVn4 : V *+ 4 + 1 != 0.
apply: lt0r_neq0.
by rewrite addr_gt0 // ?ltr01 // mulrn_wgt0.
rewrite sqr_sqrtr; last exact: (ltW hV).
by field; rewrite mulr_natr.
Qed.
Lemma phase1_step01_sigma2 :
3%:R ^+ 2 × (2%:R^-1 : R) ^+ 2 / (3%:R ^+ 2 + (2%:R^-1) ^+ 2) =
9%:R / 37%:R.
Proof. by field. Qed.
Lemma phase1_step01_mu (s : R) :
gaussian_prod_mu 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
(90%:R - 36%:R × s) / 37%:R.
Proof. rewrite /gaussian_prod_mu. by field. Qed.
Lemma phase1_step12_sigma2 :
(9%:R / 37%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 37%:R + (2%:R^-1) ^+ 2) = 9%:R / 73%:R.
Proof. by field. Qed.
Lemma phase1_step12_mu (s : R) :
gaussian_prod_mu ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R))
(2%:R^-1) =
(1134%:R - 540%:R × s) / 365%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_37 : (0 ≤ 9%:R / 37%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step23_sigma2 :
(9%:R / 73%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 73%:R + (2%:R^-1) ^+ 2) = 9%:R / 109%:R.
Proof. by field. Qed.
Lemma phase1_step23_mu (s : R) :
gaussian_prod_mu ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R))
(2%:R^-1) =
(1944%:R - 1080%:R × s) / 545%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_73 : (0 ≤ 9%:R / 73%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step34_sigma2 :
(9%:R / 109%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 109%:R + (2%:R^-1) ^+ 2) = 9%:R / 145%:R.
Proof. by field. Qed.
Lemma phase1_step34_mu (s : R) :
gaussian_prod_mu ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R))
(2%:R^-1) =
(612%:R - 360%:R × s) / 145%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_109 : (0 ≤ 9%:R / 109%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step45_sigma2 :
(9%:R / 145%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 145%:R + (2%:R^-1) ^+ 2) = 9%:R / 181%:R.
Proof. by field. Qed.
Lemma phase1_step45_mu (s : R) :
gaussian_prod_mu ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R))
(2%:R^-1) =
(900%:R - 540%:R × s) / 181%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_145 : (0 ≤ 9%:R / 145%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma nine_div_gt0 (n : nat) : (0 < n)%N → 0 < 9%:R / n%:R :> R.
Proof. by move⇒ hn; apply: divr_gt0; rewrite ?ltr0n. Qed.
Lemma sqrtr_nine_div_neq0 (n : nat) : (0 < n)%N →
sqrtr (9%:R / n%:R) != (0 : R).
Proof.
by move⇒ hn; apply: lt0r_neq0; rewrite sqrtr_gt0;
exact: nine_div_gt0.
Qed.
Lemma phase1_sqrtr_var1_neq0 : sqrtr (9%:R / 37%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var2_neq0 : sqrtr (9%:R / 73%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var3_neq0 : sqrtr (9%:R / 109%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var4_neq0 : sqrtr (9%:R / 145%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var5_neq0 : sqrtr (9%:R / 181%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_combine3 (s b : R) :
normal_pdf 0 3%:R b ×
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b =
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
normal_pdf ((1944%:R - 1080%:R × s) / 545%:R)
(sqrtr (9%:R / 109%:R)) b.
Proof.
rewrite (normal_pdf_times' _ _ _ three_neq0 half_neq0).
rewrite phase1_step01_mu.
have sigma01 : gaussian_prod_sigma 3%:R (2%:R^-1 : R) = sqrtr (9%:R / 37%:R).
have hge0 : (0 ≤ gaussian_prod_sigma 3%:R (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?pnatr_eq0 ?invr_neq0 ?pnatr_eq0 //.
by field.
rewrite sigma01.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var1_neq0 half_neq0).
rewrite phase1_step12_mu.
have sigma12 : gaussian_prod_sigma (sqrtr (9%:R / 37%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 73%:R).
have hge0 : (0 ≤ gaussian_prod_sigma (sqrtr (9%:R / 37%:R)) (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma12.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var2_neq0 half_neq0).
rewrite phase1_step23_mu.
have sigma23 : gaussian_prod_sigma (sqrtr (9%:R / 73%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 109%:R).
have hge0 : (0 ≤ gaussian_prod_sigma (sqrtr (9%:R / 73%:R)) (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma23.
by rewrite !mulrA.
Qed.
Lemma phase1_combine5 (s b : R) :
normal_pdf 0 3%:R b ×
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b ×
normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b ×
normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b =
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) ×
normal_pdf ((900%:R - 540%:R × s) / 181%:R)
(sqrtr (9%:R / 181%:R)) b.
Proof.
set N0 := normal_pdf 0 3%:R b.
set N1 := normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b.
set N2 := normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b.
set N3 := normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b.
set N4 := normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b.
set N5 := normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b.
rewrite /N4 /N5 /N0 /N1 /N2 /N3 (phase1_combine3 s b).
set K3 := gaussian_prod_scalar 0 _ _ _ ×
gaussian_prod_scalar _ _ _ _ ×
gaussian_prod_scalar _ (9%:R / 2%:R - 3%:R × s) _ _.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var3_neq0 half_neq0).
rewrite phase1_step34_mu.
have sigma34 : gaussian_prod_sigma (sqrtr (9%:R / 109%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 145%:R).
have hge0 : (0 ≤ gaussian_prod_sigma
(sqrtr (9%:R / 109%:R)) (2%:R^-1 : R)).
by apply: divr_ge0;
rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma34.
set K4 := K3 × gaussian_prod_scalar _ _ _ _.
rewrite (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var4_neq0 half_neq0).
rewrite phase1_step45_mu.
have sigma45 : gaussian_prod_sigma (sqrtr (9%:R / 145%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 181%:R).
have hge0 : (0 ≤ gaussian_prod_sigma
(sqrtr (9%:R / 145%:R)) (2%:R^-1 : R)).
by apply: divr_ge0;
rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma45.
rewrite /gaussian_prod_scalar.
by rewrite !mulrA.
Qed.
End bayesian_phase1_intercept.
Section bayesian_phase2_slope.
Lemma normal_fun_sym (m s x : R) :
normal_fun m s x = normal_fun x s m.
Proof.
rewrite /normal_fun; congr (expR _).
rewrite -sqrrN; congr (- _ / _); congr (_ ^+ 2).
ring.
Qed.
Lemma normal_fun_sub (m s a : R) :
normal_fun 0 s (a - m) = normal_fun a s m.
Proof.
rewrite /normal_fun; congr (expR _).
congr (- _ / _).
rewrite -sqrrN; congr (_ ^+ 2).
ring.
Qed.
Lemma normal_fun_linear (m m' sigma a B c s : R) :
B != 0 → c != 0 → sigma != 0 →
m' - m = (a - B × s) / c →
normal_fun m sigma m' = normal_fun (a / B) (sigma × c / B) s.
Proof.
move⇒ hB hc hsigma hdiff.
suff heq : (m' - m) ^+ 2 / (sigma ^+ 2 *+ 2) =
(s - a / B) ^+ 2 / ((sigma × c / B) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
rewrite !mulNr; congr (- _).
exact: heq.
rewrite hdiff.
have h : (a - B × s) / c = - (B / c) × (s - a / B).
field. by rewrite hB hc.
rewrite h -mulNr exprMn.
rewrite -(mulr_natr (sigma ^+ 2) 2) -(mulr_natr ((sigma × c / B) ^+ 2) 2).
rewrite !exprMn exprVn.
field.
by rewrite hsigma hB hc.
Qed.
Lemma phase2_scalar01_diff (s : R) :
(5%:R / 2%:R - s) - 0 = (5%:R / 2%:R - 1 × s) / 1.
Proof. by field. Qed.
Lemma phase2_scalar01_fun (s : R) :
normal_fun 0 (sqrtr (3%:R ^+ 2 + (2%:R^-1) ^+ 2))
(5%:R / 2%:R - s) =
normal_fun (5%:R / 2%:R) (sqrtr (3%:R ^+ 2 + (2%:R^-1) ^+ 2)) s.
Proof. by rewrite normal_fun_sub. Qed.
Lemma phase2_scalar12_diff (s : R) :
(19%:R / 5%:R - 2%:R × s) - (90%:R - 36%:R × s) / 37%:R =
(253%:R - 190%:R × s) / 185%:R.
Proof. by field. Qed.
Lemma phase2_scalar23_diff (s : R) :
(9%:R / 2%:R - 3%:R × s) - (1134%:R - 540%:R × s) / 365%:R =
(1017%:R - 1110%:R × s) / 730%:R.
Proof. by field. Qed.
Lemma phase2_scalar34_diff (s : R) :
(31%:R / 5%:R - 4%:R × s) - (1944%:R - 1080%:R × s) / 545%:R =
(287%:R - 220%:R × s) / 109%:R.
Proof. by field. Qed.
Lemma phase2_scalar45_diff (s : R) :
(8%:R - 5%:R × s) - (612%:R - 360%:R × s) / 145%:R =
(548%:R - 365%:R × s) / 145%:R.
Proof. by field. Qed.
Lemma phase2_scalar01E (s : R) :
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
normal_pdf (5%:R / 2%:R) (sqrtr (37%:R / 4%:R)) s.
Proof.
rewrite /gaussian_prod_scalar.
have hsq : 3%:R ^+ 2 + (2%:R^-1 : R) ^+ 2 = 37%:R / 4%:R by field.
rewrite hsq normal_fun_sub.
have hsqrt_neq0 : sqrtr (37%:R / 4%:R : R) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
by rewrite (normal_pdfE _ hsqrt_neq0).
Qed.
Lemma phase2_step0 (s : R) :
normal_pdf 0 3%:R s ×
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R (sqrtr (37%:R / 4%:R)) ×
normal_pdf (90%:R / 73%:R)
(gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R))) s.
Proof.
rewrite phase2_scalar01E.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
rewrite normal_pdf_times' //.
have → : gaussian_prod_mu 0 (5%:R / 2%:R) 3%:R
(sqrtr (37%:R / 4%:R : R)) = 90%:R / 73%:R.
rewrite /gaussian_prod_mu.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
done.
Qed.
Lemma phase2_step0_sigma2 :
(gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R : R))) ^+ 2 =
333%:R / 73%:R.
Proof.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
rewrite gaussian_prod_sigma_sqr ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_step0_sigma_neq0 :
gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R : R)) != 0.
Proof.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
exact: gaussian_prod_sigma_neq0.
Qed.
Lemma phase2_S12 :
sqrtr (9%:R / 37%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 73%:R / 148%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar12_fun (s : R) :
normal_fun ((90%:R - 36%:R × s) / 37%:R)
(sqrtr (73%:R / 148%:R))
(19%:R / 5%:R - 2%:R × s) =
normal_fun (253%:R / 190%:R)
(sqrtr (73%:R / 148%:R) × 37%:R / 38%:R) s.
Proof.
set S12 := sqrtr _.
have hS : S12 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 73) (b := 148) isT isT.
suff heq : ((19%:R / 5%:R - 2%:R × s) -
(90%:R - 36%:R × s) / 37%:R) ^+ 2 /
(S12 ^+ 2 *+ 2) =
(s - 253%:R / 190%:R) ^+ 2 /
((S12 × 37%:R / 38%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S12 ^+ 2) 2)
-(mulr_natr ((S12 × 37%:R / 38%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W12 : R := sqrtr (73%:R / 148%:R) × 37%:R / 38%:R.
Lemma phase2_W12_neq0 : W12 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W12.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar12E (s : R) :
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) =
normal_peak (sqrtr (73%:R / 148%:R)) /
normal_peak W12 ×
normal_pdf (253%:R / 190%:R) W12 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S12 /W12 phase2_scalar12_fun -/W12.
rewrite (normal_pdfE _ phase2_W12_neq0).
have hpeak : normal_peak W12 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W12_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma0 : R := gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R)).
Lemma phase2_step1 (s : R) :
let K0 := gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R
(sqrtr (37%:R / 4%:R)) in
K0 × normal_pdf (90%:R / 73%:R) sigma0 s ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) =
K0 × (normal_peak (sqrtr (73%:R / 148%:R)) / normal_peak W12) ×
gaussian_prod_scalar (90%:R / 73%:R) (253%:R / 190%:R) sigma0 W12 ×
normal_pdf (gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R) sigma0 W12)
(gaussian_prod_sigma sigma0 W12) s.
Proof.
rewrite /= phase2_scalar12E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step0_sigma_neq0 phase2_W12_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step1_sigma2 :
(gaussian_prod_sigma sigma0 W12) ^+ 2 =
sigma0 ^+ 2 × W12 ^+ 2 / (sigma0 ^+ 2 + W12 ^+ 2).
Proof.
exact: gaussian_prod_sigma_sqr phase2_step0_sigma_neq0 phase2_W12_neq0.
Qed.
Lemma phase2_step1_sigma_neq0 :
gaussian_prod_sigma sigma0 W12 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step0_sigma_neq0 phase2_W12_neq0.
Qed.
Lemma phase2_S23 :
sqrtr (9%:R / 73%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 109%:R / 292%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar23_fun (s : R) :
normal_fun ((1134%:R - 540%:R × s) / 365%:R)
(sqrtr (109%:R / 292%:R))
(9%:R / 2%:R - 3%:R × s) =
normal_fun (1017%:R / 1110%:R)
(sqrtr (109%:R / 292%:R) × 73%:R / 111%:R) s.
Proof.
set S23 := sqrtr _.
have hS : S23 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 109) (b := 292) isT isT.
suff heq : ((9%:R / 2%:R - 3%:R × s) -
(1134%:R - 540%:R × s) / 365%:R) ^+ 2 /
(S23 ^+ 2 *+ 2) =
(s - 1017%:R / 1110%:R) ^+ 2 /
((S23 × 73%:R / 111%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S23 ^+ 2) 2)
-(mulr_natr ((S23 × 73%:R / 111%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W23 : R := sqrtr (109%:R / 292%:R) × 73%:R / 111%:R.
Lemma phase2_W23_neq0 : W23 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W23.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar23E (s : R) :
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) =
normal_peak (sqrtr (109%:R / 292%:R)) /
normal_peak W23 ×
normal_pdf (1017%:R / 1110%:R) W23 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S23 /W23 phase2_scalar23_fun -/W23.
rewrite (normal_pdfE _ phase2_W23_neq0).
have hpeak : normal_peak W23 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W23_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma1 : R := gaussian_prod_sigma sigma0 W12.
Lemma phase2_step2 (s : R) :
∀ K1 : R,
let mu1 := gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12 in
K1 × normal_pdf mu1 sigma1 s ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) =
K1 × (normal_peak (sqrtr (109%:R / 292%:R)) / normal_peak W23) ×
gaussian_prod_scalar mu1 (1017%:R / 1110%:R) sigma1 W23 ×
normal_pdf (gaussian_prod_mu mu1 (1017%:R / 1110%:R) sigma1 W23)
(gaussian_prod_sigma sigma1 W23) s.
Proof.
move⇒ K1.
rewrite /= phase2_scalar23E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step1_sigma_neq0 phase2_W23_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step2_sigma_neq0 :
gaussian_prod_sigma sigma1 W23 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step1_sigma_neq0 phase2_W23_neq0.
Qed.
Lemma phase2_S34 :
sqrtr (9%:R / 109%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 145%:R / 436%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar34_fun (s : R) :
normal_fun ((1944%:R - 1080%:R × s) / 545%:R)
(sqrtr (145%:R / 436%:R))
(31%:R / 5%:R - 4%:R × s) =
normal_fun (287%:R / 220%:R)
(sqrtr (145%:R / 436%:R) × 109%:R / 220%:R) s.
Proof.
set S34 := sqrtr _.
have hS : S34 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 145) (b := 436) isT isT.
suff heq : ((31%:R / 5%:R - 4%:R × s) -
(1944%:R - 1080%:R × s) / 545%:R) ^+ 2 /
(S34 ^+ 2 *+ 2) =
(s - 287%:R / 220%:R) ^+ 2 /
((S34 × 109%:R / 220%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S34 ^+ 2) 2)
-(mulr_natr ((S34 × 109%:R / 220%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W34 : R := sqrtr (145%:R / 436%:R) × 109%:R / 220%:R.
Lemma phase2_W34_neq0 : W34 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W34.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar34E (s : R) :
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) =
normal_peak (sqrtr (145%:R / 436%:R)) /
normal_peak W34 ×
normal_pdf (287%:R / 220%:R) W34 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S34 /W34 phase2_scalar34_fun -/W34.
rewrite (normal_pdfE _ phase2_W34_neq0).
have hpeak : normal_peak W34 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W34_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma2 : R := gaussian_prod_sigma sigma1 W23.
Lemma phase2_step3 (s : R) :
∀ K2 : R,
let mu2 := gaussian_prod_mu
(gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12)
(1017%:R / 1110%:R) sigma1 W23 in
K2 × normal_pdf mu2 sigma2 s ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) =
K2 × (normal_peak (sqrtr (145%:R / 436%:R)) / normal_peak W34) ×
gaussian_prod_scalar mu2 (287%:R / 220%:R) sigma2 W34 ×
normal_pdf (gaussian_prod_mu mu2 (287%:R / 220%:R) sigma2 W34)
(gaussian_prod_sigma sigma2 W34) s.
Proof.
move⇒ K2.
rewrite /= phase2_scalar34E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step2_sigma_neq0 phase2_W34_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step3_sigma_neq0 :
gaussian_prod_sigma sigma2 W34 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step2_sigma_neq0 phase2_W34_neq0.
Qed.
Lemma phase2_S45 :
sqrtr (9%:R / 145%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 181%:R / 580%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar45_fun (s : R) :
normal_fun ((612%:R - 360%:R × s) / 145%:R)
(sqrtr (181%:R / 580%:R))
(8%:R - 5%:R × s) =
normal_fun (548%:R / 365%:R)
(sqrtr (181%:R / 580%:R) × 29%:R / 73%:R) s.
Proof.
set S45 := sqrtr _.
have hS : S45 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 181) (b := 580) isT isT.
suff heq : ((8%:R - 5%:R × s) -
(612%:R - 360%:R × s) / 145%:R) ^+ 2 /
(S45 ^+ 2 *+ 2) =
(s - 548%:R / 365%:R) ^+ 2 /
((S45 × 29%:R / 73%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S45 ^+ 2) 2)
-(mulr_natr ((S45 × 29%:R / 73%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W45 : R := sqrtr (181%:R / 580%:R) × 29%:R / 73%:R.
Lemma phase2_W45_neq0 : W45 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W45.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar45E (s : R) :
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
normal_peak (sqrtr (181%:R / 580%:R)) /
normal_peak W45 ×
normal_pdf (548%:R / 365%:R) W45 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S45 /W45 phase2_scalar45_fun -/W45.
rewrite (normal_pdfE _ phase2_W45_neq0).
have hpeak : normal_peak W45 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W45_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma3 : R := gaussian_prod_sigma sigma2 W34.
Lemma phase2_step4 (s : R) :
∀ K3 : R,
let mu3 := gaussian_prod_mu
(gaussian_prod_mu
(gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12)
(1017%:R / 1110%:R) sigma1 W23)
(287%:R / 220%:R) sigma2 W34 in
K3 × normal_pdf mu3 sigma3 s ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
K3 × (normal_peak (sqrtr (181%:R / 580%:R)) / normal_peak W45) ×
gaussian_prod_scalar mu3 (548%:R / 365%:R) sigma3 W45 ×
normal_pdf (gaussian_prod_mu mu3 (548%:R / 365%:R) sigma3 W45)
(gaussian_prod_sigma sigma3 W45) s.
Proof.
move⇒ K3.
rewrite /= phase2_scalar45E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step3_sigma_neq0 phase2_W45_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step4_sigma_neq0 :
gaussian_prod_sigma sigma3 W45 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step3_sigma_neq0 phase2_W45_neq0.
Qed.
Let mu0_s : R := 90%:R / 73%:R.
Let mu1_s : R := gaussian_prod_mu mu0_s (253%:R / 190%:R) sigma0 W12.
Let mu2_s : R := gaussian_prod_mu mu1_s (1017%:R / 1110%:R) sigma1 W23.
Let mu3_s : R := gaussian_prod_mu mu2_s (287%:R / 220%:R) sigma2 W34.
Definition phase2_final_mu : R :=
gaussian_prod_mu mu3_s (548%:R / 365%:R) sigma3 W45.
Definition phase2_final_sigma : R :=
gaussian_prod_sigma sigma3 W45.
Definition phase2_const : R :=
gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R (sqrtr (37%:R / 4%:R)) ×
(normal_peak (sqrtr (73%:R / 148%:R)) / normal_peak W12) ×
gaussian_prod_scalar mu0_s (253%:R / 190%:R) sigma0 W12 ×
(normal_peak (sqrtr (109%:R / 292%:R)) / normal_peak W23) ×
gaussian_prod_scalar mu1_s (1017%:R / 1110%:R) sigma1 W23 ×
(normal_peak (sqrtr (145%:R / 436%:R)) / normal_peak W34) ×
gaussian_prod_scalar mu2_s (287%:R / 220%:R) sigma2 W34 ×
(normal_peak (sqrtr (181%:R / 580%:R)) / normal_peak W45) ×
gaussian_prod_scalar mu3_s (548%:R / 365%:R) sigma3 W45.
Lemma phase2_final_sigma_neq0 : phase2_final_sigma != 0.
Proof. exact: phase2_step4_sigma_neq0. Qed.
Lemma phase2_const_gt0 : 0 < phase2_const.
Proof.
rewrite /phase2_const.
have sqrtr37_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
have sqrtr73_neq0 := sqrtr_pnat_div_neq0 (a := 73) (b := 148) isT isT.
have sqrtr109_neq0 := sqrtr_pnat_div_neq0 (a := 109) (b := 292) isT isT.
have sqrtr145_neq0 := sqrtr_pnat_div_neq0 (a := 145) (b := 436) isT isT.
have sqrtr181_neq0 := sqrtr_pnat_div_neq0 (a := 181) (b := 580) isT isT.
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0 phase2_step3_sigma_neq0 phase2_W45_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr181_neq0|
exact: normal_peak_gt0 phase2_W45_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0
phase2_step2_sigma_neq0 phase2_W34_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr145_neq0|
exact: normal_peak_gt0 phase2_W34_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0
phase2_step1_sigma_neq0 phase2_W23_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr109_neq0|
exact: normal_peak_gt0 phase2_W23_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0 phase2_step0_sigma_neq0 phase2_W12_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr73_neq0|
exact: normal_peak_gt0 phase2_W12_neq0].
exact: gaussian_prod_scalar_gt0 three_neq0 sqrtr37_neq0.
Qed.
Lemma phase2_combine5 (s : R) :
normal_pdf 0 3%:R s ×
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
phase2_const × normal_pdf phase2_final_mu phase2_final_sigma s.
Proof.
rewrite phase2_step0.
rewrite (phase2_step1 s).
rewrite (phase2_step2 s).
rewrite (phase2_step3 s).
rewrite (phase2_step4 s).
rewrite /phase2_const /phase2_final_mu /phase2_final_sigma.
by rewrite !mulrA.
Qed.
Lemma normal_pdf_recenter (y s k sigma b : R) (hs : sigma != 0) :
normal_pdf (s × k + b) sigma y = normal_pdf (y - k × s) sigma b.
Proof.
rewrite !(normal_pdfE _ hs); congr (_ × _).
rewrite /normal_fun; congr (expR _).
congr (- _ / _).
suff → : (y - (s × k + b)) = - (b - (y - k × s)) by rewrite sqrrN.
ring.
Qed.
Arguments normal_pdf_recenter : clear implicits.
Lemma obs_rewrite (s b : R) :
normal_pdf (s × 1 + b) (2%:R^-1) (5%:R / 2%:R) ×
normal_pdf (s × 2%:R + b) (2%:R^-1) (19%:R / 5%:R) ×
normal_pdf (s × 3%:R + b) (2%:R^-1) (9%:R / 2%:R) ×
normal_pdf (s × 4%:R + b) (2%:R^-1) (31%:R / 5%:R) ×
normal_pdf (s × 5%:R + b) (2%:R^-1) 8%:R =
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b ×
normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b ×
normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b.
Proof.
rewrite (normal_pdf_recenter (5%:R/2%:R) s 1 (2%:R^-1) b half_neq0)
(normal_pdf_recenter (19%:R/5%:R) s 2%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter (9%:R/2%:R) s 3%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter (31%:R/5%:R) s 4%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter 8%:R s 5%:R (2%:R^-1) b half_neq0).
congr (_ × _ × _ × _ × _); rewrite ?mul1r //.
Qed.
End bayesian_phase2_slope.
End normal_density_algebra.
From mathcomp Require Import all_boot all_algebra reals ereal topology
normedtype numfun measure sequences lebesgue_integral lebesgue_integral_fubini
lebesgue_stieltjes_measure probability boolp trigo exp.
From mathcomp.algebra_tactics Require Import ring.
Import GRing.Theory Num.Def Num.Theory.
Set Implicit Arguments.
Unset Strict Implicit.
Unset Printing Implicit Defensive.
Import order.Order.POrderTheory.
Section normal_density_algebra.
Variable R : realType.
Local Open Scope ring_scope.
Implicit Types (m s x K : R) (n : nat).
Let half_neq0 : (2%:R^-1 : R) != 0.
Proof. by rewrite invr_neq0 // pnatr_eq0. Qed.
Let three_neq0 : (3%:R : R) != 0.
Proof. by rewrite pnatr_eq0. Qed.
Local Lemma sqrtr_pnat_div_neq0 (a b : nat) :
(0 < a)%N → (0 < b)%N → sqrtr (a%:R / b%:R : R) != 0.
Proof.
move⇒ ha hb; apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ltr0n.
Qed.
Section normal_pdf_general.
Lemma complete_the_square (a b c x : R) (ha : a != 0) :
a × x ^+ 2 + b × x + c =
a × (x + b / (a *+ 2)) ^+ 2 - (b ^+ 2 - a *+ 4 × c) / (a *+ 4).
Proof.
have ha4 : a *+ 4 != 0 by rewrite mulrn_eq0.
have ha2 : a *+ 2 != 0 by rewrite mulrn_eq0.
field. repeat split ⇒ //; rewrite ?mulrn_eq0 //.
Qed.
Lemma normal_pdf_times (m m' s s' x : R) :
s != 0 → s' != 0 →
normal_pdf m s x × normal_pdf m' s' x =
normal_peak (sqrtr (s ^+ 2 + s' ^+ 2)) ×
normal_fun m (sqrtr (s ^+ 2 + s' ^+ 2)) m' ×
normal_pdf ((m × s' ^+ 2 + m' × s ^+ 2) / (s ^+ 2 + s' ^+ 2))
(s × s' / sqrtr (s ^+ 2 + s' ^+ 2)) x.
Proof.
move⇒ hs hs'.
have hS : s ^+ 2 + s' ^+ 2 != 0.
by rewrite paddr_eq0 ?sqr_ge0 // !sqrf_eq0; apply/nandP; left.
have hSgt0 : 0 < s ^+ 2 + s' ^+ 2.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
by apply: lt0r_neq0; rewrite sqrtr_gt0.
have hsnew : s × s' / sqrtr (s ^+ 2 + s' ^+ 2) != 0.
by rewrite mulf_neq0 ?invr_neq0 // mulf_neq0.
have hpi : 0 ≤ pi := pi_ge0 R.
have hs2pi : 0 ≤ s ^+ 2 × pi := mulr_ge0 (sqr_ge0 _) hpi.
have hS2pi : 0 ≤ (s ^+ 2 + s' ^+ 2) × pi.
by apply: mulr_ge0 ⇒ //; exact: addr_ge0 (sqr_ge0 _) (sqr_ge0 _).
rewrite (normal_pdfE _ hs) (normal_pdfE _ hs') (normal_pdfE _ hsnew).
rewrite [LHS]mulrACA [RHS]mulrACA.
congr (_ × _).
-
rewrite /normal_peak -invfM -sqrtrM ?pmulrn_lge0 ?sqr_ge0 //.
rewrite -[RHS]invfM -sqrtrM ?pmulrn_lge0 ?sqr_ge0 ?sqr_sqrtr
?addr_ge0 ?sqr_ge0 //.
congr (_^-1); congr (sqrtr _).
rewrite exprMn exprVn sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
by field; exact: hS.
-
rewrite /normal_fun -expRD -[RHS]expRD.
congr (expR _).
rewrite sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
rewrite exprMn exprVn sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
field.
by rewrite hS hs hs'.
Qed.
Definition gaussian_prod_mu (m m' s s' : R) :=
(m × s' ^+ 2 + m' × s ^+ 2) / (s ^+ 2 + s' ^+ 2).
Definition gaussian_prod_sigma (s s' : R) :=
s × s' / sqrtr (s ^+ 2 + s' ^+ 2).
Definition gaussian_prod_scalar (m m' s s' : R) :=
normal_peak (sqrtr (s ^+ 2 + s' ^+ 2)) ×
normal_fun m (sqrtr (s ^+ 2 + s' ^+ 2)) m'.
Lemma normal_pdf_times' (m m' s s' x : R) :
s != 0 → s' != 0 →
normal_pdf m s x × normal_pdf m' s' x =
gaussian_prod_scalar m m' s s' ×
normal_pdf (gaussian_prod_mu m m' s s') (gaussian_prod_sigma s s') x.
Proof.
move⇒ hs hs'.
by rewrite /gaussian_prod_scalar /gaussian_prod_mu /gaussian_prod_sigma
normal_pdf_times.
Qed.
Lemma normal_pdf_times_chain (K m m' s s' x : R) :
s != 0 → s' != 0 →
K × normal_pdf m s x × normal_pdf m' s' x =
K × gaussian_prod_scalar m m' s s' ×
normal_pdf (gaussian_prod_mu m m' s s') (gaussian_prod_sigma s s') x.
Proof.
move⇒ hs hs'.
by rewrite -mulrA normal_pdf_times' // mulrA.
Qed.
Lemma gaussian_prod_sigma_sqr (s s' : R) :
s != 0 → s' != 0 →
(gaussian_prod_sigma s s') ^+ 2 =
s ^+ 2 × s' ^+ 2 / (s ^+ 2 + s' ^+ 2).
Proof.
move⇒ hs hs'.
have hS : s ^+ 2 + s' ^+ 2 != 0.
by rewrite paddr_eq0 ?sqr_ge0 // !sqrf_eq0; apply/nandP; left.
rewrite /gaussian_prod_sigma exprMn exprVn.
rewrite sqr_sqrtr ?addr_ge0 ?sqr_ge0 //.
by field.
Qed.
Lemma gaussian_prod_sigma_neq0 (s s' : R) :
s != 0 → s' != 0 → gaussian_prod_sigma s s' != 0.
Proof.
move⇒ hs hs'.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
by rewrite /gaussian_prod_sigma mulf_neq0 ?invr_neq0 // mulf_neq0.
Qed.
Lemma gaussian_prod_scalar_gt0 (m m' s s' : R) :
s != 0 → s' != 0 → 0 < gaussian_prod_scalar m m' s s'.
Proof.
move⇒ hs hs'.
have hsqS : sqrtr (s ^+ 2 + s' ^+ 2) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: addr_gt0; rewrite exprn_even_gt0 //=.
rewrite /gaussian_prod_scalar.
apply: mulr_gt0; first exact: (normal_peak_gt0 hsqS).
rewrite /normal_fun; exact: expR_gt0.
Qed.
End normal_pdf_general.
Section bayesian_phase1_intercept.
Lemma var_recurrence (V : R) :
0 < V →
sqrtr V ^+ 2 × (2%:R ^-1) ^+ 2 / (sqrtr V ^+ 2 + (2%:R ^-1) ^+ 2) =
V / (V *+ 4 + 1).
Proof.
move⇒ hV.
have hV0 : V != 0 by apply: lt0r_neq0.
have hVn4 : V *+ 4 + 1 != 0.
apply: lt0r_neq0.
by rewrite addr_gt0 // ?ltr01 // mulrn_wgt0.
rewrite sqr_sqrtr; last exact: (ltW hV).
by field; rewrite mulr_natr.
Qed.
Lemma phase1_step01_sigma2 :
3%:R ^+ 2 × (2%:R^-1 : R) ^+ 2 / (3%:R ^+ 2 + (2%:R^-1) ^+ 2) =
9%:R / 37%:R.
Proof. by field. Qed.
Lemma phase1_step01_mu (s : R) :
gaussian_prod_mu 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
(90%:R - 36%:R × s) / 37%:R.
Proof. rewrite /gaussian_prod_mu. by field. Qed.
Lemma phase1_step12_sigma2 :
(9%:R / 37%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 37%:R + (2%:R^-1) ^+ 2) = 9%:R / 73%:R.
Proof. by field. Qed.
Lemma phase1_step12_mu (s : R) :
gaussian_prod_mu ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R))
(2%:R^-1) =
(1134%:R - 540%:R × s) / 365%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_37 : (0 ≤ 9%:R / 37%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step23_sigma2 :
(9%:R / 73%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 73%:R + (2%:R^-1) ^+ 2) = 9%:R / 109%:R.
Proof. by field. Qed.
Lemma phase1_step23_mu (s : R) :
gaussian_prod_mu ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R))
(2%:R^-1) =
(1944%:R - 1080%:R × s) / 545%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_73 : (0 ≤ 9%:R / 73%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step34_sigma2 :
(9%:R / 109%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 109%:R + (2%:R^-1) ^+ 2) = 9%:R / 145%:R.
Proof. by field. Qed.
Lemma phase1_step34_mu (s : R) :
gaussian_prod_mu ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R))
(2%:R^-1) =
(612%:R - 360%:R × s) / 145%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_109 : (0 ≤ 9%:R / 109%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma phase1_step45_sigma2 :
(9%:R / 145%:R) × (2%:R^-1 : R) ^+ 2 /
(9%:R / 145%:R + (2%:R^-1) ^+ 2) = 9%:R / 181%:R.
Proof. by field. Qed.
Lemma phase1_step45_mu (s : R) :
gaussian_prod_mu ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R))
(2%:R^-1) =
(900%:R - 540%:R × s) / 181%:R.
Proof.
rewrite /gaussian_prod_mu.
have h9_145 : (0 ≤ 9%:R / 145%:R :> R).
by apply: divr_ge0; rewrite ?ler0n.
rewrite sqr_sqrtr //.
by field.
Qed.
Lemma nine_div_gt0 (n : nat) : (0 < n)%N → 0 < 9%:R / n%:R :> R.
Proof. by move⇒ hn; apply: divr_gt0; rewrite ?ltr0n. Qed.
Lemma sqrtr_nine_div_neq0 (n : nat) : (0 < n)%N →
sqrtr (9%:R / n%:R) != (0 : R).
Proof.
by move⇒ hn; apply: lt0r_neq0; rewrite sqrtr_gt0;
exact: nine_div_gt0.
Qed.
Lemma phase1_sqrtr_var1_neq0 : sqrtr (9%:R / 37%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var2_neq0 : sqrtr (9%:R / 73%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var3_neq0 : sqrtr (9%:R / 109%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var4_neq0 : sqrtr (9%:R / 145%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_sqrtr_var5_neq0 : sqrtr (9%:R / 181%:R) != (0 : R).
Proof. exact: sqrtr_nine_div_neq0. Qed.
Lemma phase1_combine3 (s b : R) :
normal_pdf 0 3%:R b ×
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b =
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
normal_pdf ((1944%:R - 1080%:R × s) / 545%:R)
(sqrtr (9%:R / 109%:R)) b.
Proof.
rewrite (normal_pdf_times' _ _ _ three_neq0 half_neq0).
rewrite phase1_step01_mu.
have sigma01 : gaussian_prod_sigma 3%:R (2%:R^-1 : R) = sqrtr (9%:R / 37%:R).
have hge0 : (0 ≤ gaussian_prod_sigma 3%:R (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?pnatr_eq0 ?invr_neq0 ?pnatr_eq0 //.
by field.
rewrite sigma01.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var1_neq0 half_neq0).
rewrite phase1_step12_mu.
have sigma12 : gaussian_prod_sigma (sqrtr (9%:R / 37%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 73%:R).
have hge0 : (0 ≤ gaussian_prod_sigma (sqrtr (9%:R / 37%:R)) (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma12.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var2_neq0 half_neq0).
rewrite phase1_step23_mu.
have sigma23 : gaussian_prod_sigma (sqrtr (9%:R / 73%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 109%:R).
have hge0 : (0 ≤ gaussian_prod_sigma (sqrtr (9%:R / 73%:R)) (2%:R^-1 : R)).
by apply: divr_ge0; rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma23.
by rewrite !mulrA.
Qed.
Lemma phase1_combine5 (s b : R) :
normal_pdf 0 3%:R b ×
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b ×
normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b ×
normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b =
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) ×
normal_pdf ((900%:R - 540%:R × s) / 181%:R)
(sqrtr (9%:R / 181%:R)) b.
Proof.
set N0 := normal_pdf 0 3%:R b.
set N1 := normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b.
set N2 := normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b.
set N3 := normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b.
set N4 := normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b.
set N5 := normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b.
rewrite /N4 /N5 /N0 /N1 /N2 /N3 (phase1_combine3 s b).
set K3 := gaussian_prod_scalar 0 _ _ _ ×
gaussian_prod_scalar _ _ _ _ ×
gaussian_prod_scalar _ (9%:R / 2%:R - 3%:R × s) _ _.
rewrite mulrA (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var3_neq0 half_neq0).
rewrite phase1_step34_mu.
have sigma34 : gaussian_prod_sigma (sqrtr (9%:R / 109%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 145%:R).
have hge0 : (0 ≤ gaussian_prod_sigma
(sqrtr (9%:R / 109%:R)) (2%:R^-1 : R)).
by apply: divr_ge0;
rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma34.
set K4 := K3 × gaussian_prod_scalar _ _ _ _.
rewrite (normal_pdf_times_chain _ _ _ _ phase1_sqrtr_var4_neq0 half_neq0).
rewrite phase1_step45_mu.
have sigma45 : gaussian_prod_sigma (sqrtr (9%:R / 145%:R)) (2%:R^-1 : R) =
sqrtr (9%:R / 181%:R).
have hge0 : (0 ≤ gaussian_prod_sigma
(sqrtr (9%:R / 145%:R)) (2%:R^-1 : R)).
by apply: divr_ge0;
rewrite ?mulr_ge0 ?invr_ge0 ?ler0n ?sqrtr_ge0.
rewrite -[LHS](ger0_norm hge0) -sqrtr_sqr.
congr (sqrtr _).
rewrite gaussian_prod_sigma_sqr ?invr_neq0 ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
rewrite sigma45.
rewrite /gaussian_prod_scalar.
by rewrite !mulrA.
Qed.
End bayesian_phase1_intercept.
Section bayesian_phase2_slope.
Lemma normal_fun_sym (m s x : R) :
normal_fun m s x = normal_fun x s m.
Proof.
rewrite /normal_fun; congr (expR _).
rewrite -sqrrN; congr (- _ / _); congr (_ ^+ 2).
ring.
Qed.
Lemma normal_fun_sub (m s a : R) :
normal_fun 0 s (a - m) = normal_fun a s m.
Proof.
rewrite /normal_fun; congr (expR _).
congr (- _ / _).
rewrite -sqrrN; congr (_ ^+ 2).
ring.
Qed.
Lemma normal_fun_linear (m m' sigma a B c s : R) :
B != 0 → c != 0 → sigma != 0 →
m' - m = (a - B × s) / c →
normal_fun m sigma m' = normal_fun (a / B) (sigma × c / B) s.
Proof.
move⇒ hB hc hsigma hdiff.
suff heq : (m' - m) ^+ 2 / (sigma ^+ 2 *+ 2) =
(s - a / B) ^+ 2 / ((sigma × c / B) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
rewrite !mulNr; congr (- _).
exact: heq.
rewrite hdiff.
have h : (a - B × s) / c = - (B / c) × (s - a / B).
field. by rewrite hB hc.
rewrite h -mulNr exprMn.
rewrite -(mulr_natr (sigma ^+ 2) 2) -(mulr_natr ((sigma × c / B) ^+ 2) 2).
rewrite !exprMn exprVn.
field.
by rewrite hsigma hB hc.
Qed.
Lemma phase2_scalar01_diff (s : R) :
(5%:R / 2%:R - s) - 0 = (5%:R / 2%:R - 1 × s) / 1.
Proof. by field. Qed.
Lemma phase2_scalar01_fun (s : R) :
normal_fun 0 (sqrtr (3%:R ^+ 2 + (2%:R^-1) ^+ 2))
(5%:R / 2%:R - s) =
normal_fun (5%:R / 2%:R) (sqrtr (3%:R ^+ 2 + (2%:R^-1) ^+ 2)) s.
Proof. by rewrite normal_fun_sub. Qed.
Lemma phase2_scalar12_diff (s : R) :
(19%:R / 5%:R - 2%:R × s) - (90%:R - 36%:R × s) / 37%:R =
(253%:R - 190%:R × s) / 185%:R.
Proof. by field. Qed.
Lemma phase2_scalar23_diff (s : R) :
(9%:R / 2%:R - 3%:R × s) - (1134%:R - 540%:R × s) / 365%:R =
(1017%:R - 1110%:R × s) / 730%:R.
Proof. by field. Qed.
Lemma phase2_scalar34_diff (s : R) :
(31%:R / 5%:R - 4%:R × s) - (1944%:R - 1080%:R × s) / 545%:R =
(287%:R - 220%:R × s) / 109%:R.
Proof. by field. Qed.
Lemma phase2_scalar45_diff (s : R) :
(8%:R - 5%:R × s) - (612%:R - 360%:R × s) / 145%:R =
(548%:R - 365%:R × s) / 145%:R.
Proof. by field. Qed.
Lemma phase2_scalar01E (s : R) :
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
normal_pdf (5%:R / 2%:R) (sqrtr (37%:R / 4%:R)) s.
Proof.
rewrite /gaussian_prod_scalar.
have hsq : 3%:R ^+ 2 + (2%:R^-1 : R) ^+ 2 = 37%:R / 4%:R by field.
rewrite hsq normal_fun_sub.
have hsqrt_neq0 : sqrtr (37%:R / 4%:R : R) != 0.
apply: lt0r_neq0; rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
by rewrite (normal_pdfE _ hsqrt_neq0).
Qed.
Lemma phase2_step0 (s : R) :
normal_pdf 0 3%:R s ×
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) =
gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R (sqrtr (37%:R / 4%:R)) ×
normal_pdf (90%:R / 73%:R)
(gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R))) s.
Proof.
rewrite phase2_scalar01E.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
rewrite normal_pdf_times' //.
have → : gaussian_prod_mu 0 (5%:R / 2%:R) 3%:R
(sqrtr (37%:R / 4%:R : R)) = 90%:R / 73%:R.
rewrite /gaussian_prod_mu.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
done.
Qed.
Lemma phase2_step0_sigma2 :
(gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R : R))) ^+ 2 =
333%:R / 73%:R.
Proof.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
rewrite gaussian_prod_sigma_sqr ?pnatr_eq0 //.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_step0_sigma_neq0 :
gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R : R)) != 0.
Proof.
have hsqrt_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
exact: gaussian_prod_sigma_neq0.
Qed.
Lemma phase2_S12 :
sqrtr (9%:R / 37%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 73%:R / 148%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar12_fun (s : R) :
normal_fun ((90%:R - 36%:R × s) / 37%:R)
(sqrtr (73%:R / 148%:R))
(19%:R / 5%:R - 2%:R × s) =
normal_fun (253%:R / 190%:R)
(sqrtr (73%:R / 148%:R) × 37%:R / 38%:R) s.
Proof.
set S12 := sqrtr _.
have hS : S12 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 73) (b := 148) isT isT.
suff heq : ((19%:R / 5%:R - 2%:R × s) -
(90%:R - 36%:R × s) / 37%:R) ^+ 2 /
(S12 ^+ 2 *+ 2) =
(s - 253%:R / 190%:R) ^+ 2 /
((S12 × 37%:R / 38%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S12 ^+ 2) 2)
-(mulr_natr ((S12 × 37%:R / 38%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W12 : R := sqrtr (73%:R / 148%:R) × 37%:R / 38%:R.
Lemma phase2_W12_neq0 : W12 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W12.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar12E (s : R) :
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) =
normal_peak (sqrtr (73%:R / 148%:R)) /
normal_peak W12 ×
normal_pdf (253%:R / 190%:R) W12 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S12 /W12 phase2_scalar12_fun -/W12.
rewrite (normal_pdfE _ phase2_W12_neq0).
have hpeak : normal_peak W12 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W12_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma0 : R := gaussian_prod_sigma 3%:R (sqrtr (37%:R / 4%:R)).
Lemma phase2_step1 (s : R) :
let K0 := gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R
(sqrtr (37%:R / 4%:R)) in
K0 × normal_pdf (90%:R / 73%:R) sigma0 s ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) =
K0 × (normal_peak (sqrtr (73%:R / 148%:R)) / normal_peak W12) ×
gaussian_prod_scalar (90%:R / 73%:R) (253%:R / 190%:R) sigma0 W12 ×
normal_pdf (gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R) sigma0 W12)
(gaussian_prod_sigma sigma0 W12) s.
Proof.
rewrite /= phase2_scalar12E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step0_sigma_neq0 phase2_W12_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step1_sigma2 :
(gaussian_prod_sigma sigma0 W12) ^+ 2 =
sigma0 ^+ 2 × W12 ^+ 2 / (sigma0 ^+ 2 + W12 ^+ 2).
Proof.
exact: gaussian_prod_sigma_sqr phase2_step0_sigma_neq0 phase2_W12_neq0.
Qed.
Lemma phase2_step1_sigma_neq0 :
gaussian_prod_sigma sigma0 W12 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step0_sigma_neq0 phase2_W12_neq0.
Qed.
Lemma phase2_S23 :
sqrtr (9%:R / 73%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 109%:R / 292%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar23_fun (s : R) :
normal_fun ((1134%:R - 540%:R × s) / 365%:R)
(sqrtr (109%:R / 292%:R))
(9%:R / 2%:R - 3%:R × s) =
normal_fun (1017%:R / 1110%:R)
(sqrtr (109%:R / 292%:R) × 73%:R / 111%:R) s.
Proof.
set S23 := sqrtr _.
have hS : S23 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 109) (b := 292) isT isT.
suff heq : ((9%:R / 2%:R - 3%:R × s) -
(1134%:R - 540%:R × s) / 365%:R) ^+ 2 /
(S23 ^+ 2 *+ 2) =
(s - 1017%:R / 1110%:R) ^+ 2 /
((S23 × 73%:R / 111%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S23 ^+ 2) 2)
-(mulr_natr ((S23 × 73%:R / 111%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W23 : R := sqrtr (109%:R / 292%:R) × 73%:R / 111%:R.
Lemma phase2_W23_neq0 : W23 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W23.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar23E (s : R) :
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) =
normal_peak (sqrtr (109%:R / 292%:R)) /
normal_peak W23 ×
normal_pdf (1017%:R / 1110%:R) W23 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S23 /W23 phase2_scalar23_fun -/W23.
rewrite (normal_pdfE _ phase2_W23_neq0).
have hpeak : normal_peak W23 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W23_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma1 : R := gaussian_prod_sigma sigma0 W12.
Lemma phase2_step2 (s : R) :
∀ K1 : R,
let mu1 := gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12 in
K1 × normal_pdf mu1 sigma1 s ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) =
K1 × (normal_peak (sqrtr (109%:R / 292%:R)) / normal_peak W23) ×
gaussian_prod_scalar mu1 (1017%:R / 1110%:R) sigma1 W23 ×
normal_pdf (gaussian_prod_mu mu1 (1017%:R / 1110%:R) sigma1 W23)
(gaussian_prod_sigma sigma1 W23) s.
Proof.
move⇒ K1.
rewrite /= phase2_scalar23E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step1_sigma_neq0 phase2_W23_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step2_sigma_neq0 :
gaussian_prod_sigma sigma1 W23 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step1_sigma_neq0 phase2_W23_neq0.
Qed.
Lemma phase2_S34 :
sqrtr (9%:R / 109%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 145%:R / 436%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar34_fun (s : R) :
normal_fun ((1944%:R - 1080%:R × s) / 545%:R)
(sqrtr (145%:R / 436%:R))
(31%:R / 5%:R - 4%:R × s) =
normal_fun (287%:R / 220%:R)
(sqrtr (145%:R / 436%:R) × 109%:R / 220%:R) s.
Proof.
set S34 := sqrtr _.
have hS : S34 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 145) (b := 436) isT isT.
suff heq : ((31%:R / 5%:R - 4%:R × s) -
(1944%:R - 1080%:R × s) / 545%:R) ^+ 2 /
(S34 ^+ 2 *+ 2) =
(s - 287%:R / 220%:R) ^+ 2 /
((S34 × 109%:R / 220%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S34 ^+ 2) 2)
-(mulr_natr ((S34 × 109%:R / 220%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W34 : R := sqrtr (145%:R / 436%:R) × 109%:R / 220%:R.
Lemma phase2_W34_neq0 : W34 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W34.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar34E (s : R) :
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) =
normal_peak (sqrtr (145%:R / 436%:R)) /
normal_peak W34 ×
normal_pdf (287%:R / 220%:R) W34 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S34 /W34 phase2_scalar34_fun -/W34.
rewrite (normal_pdfE _ phase2_W34_neq0).
have hpeak : normal_peak W34 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W34_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma2 : R := gaussian_prod_sigma sigma1 W23.
Lemma phase2_step3 (s : R) :
∀ K2 : R,
let mu2 := gaussian_prod_mu
(gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12)
(1017%:R / 1110%:R) sigma1 W23 in
K2 × normal_pdf mu2 sigma2 s ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) =
K2 × (normal_peak (sqrtr (145%:R / 436%:R)) / normal_peak W34) ×
gaussian_prod_scalar mu2 (287%:R / 220%:R) sigma2 W34 ×
normal_pdf (gaussian_prod_mu mu2 (287%:R / 220%:R) sigma2 W34)
(gaussian_prod_sigma sigma2 W34) s.
Proof.
move⇒ K2.
rewrite /= phase2_scalar34E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step2_sigma_neq0 phase2_W34_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step3_sigma_neq0 :
gaussian_prod_sigma sigma2 W34 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step2_sigma_neq0 phase2_W34_neq0.
Qed.
Lemma phase2_S45 :
sqrtr (9%:R / 145%:R) ^+ 2 + (2%:R^-1 : R) ^+ 2 = 181%:R / 580%:R.
Proof.
rewrite sqr_sqrtr; last by apply: divr_ge0; rewrite ?ler0n.
by field.
Qed.
Lemma phase2_scalar45_fun (s : R) :
normal_fun ((612%:R - 360%:R × s) / 145%:R)
(sqrtr (181%:R / 580%:R))
(8%:R - 5%:R × s) =
normal_fun (548%:R / 365%:R)
(sqrtr (181%:R / 580%:R) × 29%:R / 73%:R) s.
Proof.
set S45 := sqrtr _.
have hS : S45 != (0 : R)
:= sqrtr_pnat_div_neq0 (a := 181) (b := 580) isT isT.
suff heq : ((8%:R - 5%:R × s) -
(612%:R - 360%:R × s) / 145%:R) ^+ 2 /
(S45 ^+ 2 *+ 2) =
(s - 548%:R / 365%:R) ^+ 2 /
((S45 × 29%:R / 73%:R) ^+ 2 *+ 2).
rewrite /normal_fun; congr (expR _).
by rewrite !mulNr; congr (- _).
rewrite -(mulr_natr (S45 ^+ 2) 2)
-(mulr_natr ((S45 × 29%:R / 73%:R) ^+ 2) 2).
rewrite !exprMn exprVn.
by field; rewrite hS.
Qed.
Let W45 : R := sqrtr (181%:R / 580%:R) × 29%:R / 73%:R.
Lemma phase2_W45_neq0 : W45 != 0.
Proof.
apply: lt0r_neq0.
rewrite /W45.
apply: divr_gt0; last by rewrite ltr0n.
apply: mulr_gt0; last by rewrite ltr0n.
rewrite sqrtr_gt0.
by apply: divr_gt0; rewrite ?ltr0n.
Qed.
Lemma phase2_scalar45E (s : R) :
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
normal_peak (sqrtr (181%:R / 580%:R)) /
normal_peak W45 ×
normal_pdf (548%:R / 365%:R) W45 s.
Proof.
rewrite /gaussian_prod_scalar phase2_S45 /W45 phase2_scalar45_fun -/W45.
rewrite (normal_pdfE _ phase2_W45_neq0).
have hpeak : normal_peak W45 != 0.
by apply: lt0r_neq0; exact: (normal_peak_gt0 phase2_W45_neq0).
by rewrite mulrA divfK.
Qed.
Let sigma3 : R := gaussian_prod_sigma sigma2 W34.
Lemma phase2_step4 (s : R) :
∀ K3 : R,
let mu3 := gaussian_prod_mu
(gaussian_prod_mu
(gaussian_prod_mu (90%:R / 73%:R) (253%:R / 190%:R)
sigma0 W12)
(1017%:R / 1110%:R) sigma1 W23)
(287%:R / 220%:R) sigma2 W34 in
K3 × normal_pdf mu3 sigma3 s ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
K3 × (normal_peak (sqrtr (181%:R / 580%:R)) / normal_peak W45) ×
gaussian_prod_scalar mu3 (548%:R / 365%:R) sigma3 W45 ×
normal_pdf (gaussian_prod_mu mu3 (548%:R / 365%:R) sigma3 W45)
(gaussian_prod_sigma sigma3 W45) s.
Proof.
move⇒ K3.
rewrite /= phase2_scalar45E.
rewrite -mulrA (mulrCA (normal_pdf _ _ _)) !mulrA.
rewrite (normal_pdf_times_chain _ _ _ _
phase2_step3_sigma_neq0 phase2_W45_neq0).
by rewrite !mulrA.
Qed.
Lemma phase2_step4_sigma_neq0 :
gaussian_prod_sigma sigma3 W45 != 0.
Proof.
exact: gaussian_prod_sigma_neq0 phase2_step3_sigma_neq0 phase2_W45_neq0.
Qed.
Let mu0_s : R := 90%:R / 73%:R.
Let mu1_s : R := gaussian_prod_mu mu0_s (253%:R / 190%:R) sigma0 W12.
Let mu2_s : R := gaussian_prod_mu mu1_s (1017%:R / 1110%:R) sigma1 W23.
Let mu3_s : R := gaussian_prod_mu mu2_s (287%:R / 220%:R) sigma2 W34.
Definition phase2_final_mu : R :=
gaussian_prod_mu mu3_s (548%:R / 365%:R) sigma3 W45.
Definition phase2_final_sigma : R :=
gaussian_prod_sigma sigma3 W45.
Definition phase2_const : R :=
gaussian_prod_scalar 0 (5%:R / 2%:R) 3%:R (sqrtr (37%:R / 4%:R)) ×
(normal_peak (sqrtr (73%:R / 148%:R)) / normal_peak W12) ×
gaussian_prod_scalar mu0_s (253%:R / 190%:R) sigma0 W12 ×
(normal_peak (sqrtr (109%:R / 292%:R)) / normal_peak W23) ×
gaussian_prod_scalar mu1_s (1017%:R / 1110%:R) sigma1 W23 ×
(normal_peak (sqrtr (145%:R / 436%:R)) / normal_peak W34) ×
gaussian_prod_scalar mu2_s (287%:R / 220%:R) sigma2 W34 ×
(normal_peak (sqrtr (181%:R / 580%:R)) / normal_peak W45) ×
gaussian_prod_scalar mu3_s (548%:R / 365%:R) sigma3 W45.
Lemma phase2_final_sigma_neq0 : phase2_final_sigma != 0.
Proof. exact: phase2_step4_sigma_neq0. Qed.
Lemma phase2_const_gt0 : 0 < phase2_const.
Proof.
rewrite /phase2_const.
have sqrtr37_neq0 := sqrtr_pnat_div_neq0 (a := 37) (b := 4) isT isT.
have sqrtr73_neq0 := sqrtr_pnat_div_neq0 (a := 73) (b := 148) isT isT.
have sqrtr109_neq0 := sqrtr_pnat_div_neq0 (a := 109) (b := 292) isT isT.
have sqrtr145_neq0 := sqrtr_pnat_div_neq0 (a := 145) (b := 436) isT isT.
have sqrtr181_neq0 := sqrtr_pnat_div_neq0 (a := 181) (b := 580) isT isT.
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0 phase2_step3_sigma_neq0 phase2_W45_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr181_neq0|
exact: normal_peak_gt0 phase2_W45_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0
phase2_step2_sigma_neq0 phase2_W34_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr145_neq0|
exact: normal_peak_gt0 phase2_W34_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0
phase2_step1_sigma_neq0 phase2_W23_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr109_neq0|
exact: normal_peak_gt0 phase2_W23_neq0].
apply: mulr_gt0; last first.
exact: gaussian_prod_scalar_gt0 phase2_step0_sigma_neq0 phase2_W12_neq0.
apply: mulr_gt0; last first.
by apply: divr_gt0;
[exact: normal_peak_gt0 sqrtr73_neq0|
exact: normal_peak_gt0 phase2_W12_neq0].
exact: gaussian_prod_scalar_gt0 three_neq0 sqrtr37_neq0.
Qed.
Lemma phase2_combine5 (s : R) :
normal_pdf 0 3%:R s ×
gaussian_prod_scalar 0 (5%:R / 2%:R - s) 3%:R (2%:R^-1) ×
gaussian_prod_scalar ((90%:R - 36%:R × s) / 37%:R)
(19%:R / 5%:R - 2%:R × s)
(sqrtr (9%:R / 37%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1134%:R - 540%:R × s) / 365%:R)
(9%:R / 2%:R - 3%:R × s)
(sqrtr (9%:R / 73%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((1944%:R - 1080%:R × s) / 545%:R)
(31%:R / 5%:R - 4%:R × s)
(sqrtr (9%:R / 109%:R)) (2%:R^-1) ×
gaussian_prod_scalar ((612%:R - 360%:R × s) / 145%:R)
(8%:R - 5%:R × s)
(sqrtr (9%:R / 145%:R)) (2%:R^-1) =
phase2_const × normal_pdf phase2_final_mu phase2_final_sigma s.
Proof.
rewrite phase2_step0.
rewrite (phase2_step1 s).
rewrite (phase2_step2 s).
rewrite (phase2_step3 s).
rewrite (phase2_step4 s).
rewrite /phase2_const /phase2_final_mu /phase2_final_sigma.
by rewrite !mulrA.
Qed.
Lemma normal_pdf_recenter (y s k sigma b : R) (hs : sigma != 0) :
normal_pdf (s × k + b) sigma y = normal_pdf (y - k × s) sigma b.
Proof.
rewrite !(normal_pdfE _ hs); congr (_ × _).
rewrite /normal_fun; congr (expR _).
congr (- _ / _).
suff → : (y - (s × k + b)) = - (b - (y - k × s)) by rewrite sqrrN.
ring.
Qed.
Arguments normal_pdf_recenter : clear implicits.
Lemma obs_rewrite (s b : R) :
normal_pdf (s × 1 + b) (2%:R^-1) (5%:R / 2%:R) ×
normal_pdf (s × 2%:R + b) (2%:R^-1) (19%:R / 5%:R) ×
normal_pdf (s × 3%:R + b) (2%:R^-1) (9%:R / 2%:R) ×
normal_pdf (s × 4%:R + b) (2%:R^-1) (31%:R / 5%:R) ×
normal_pdf (s × 5%:R + b) (2%:R^-1) 8%:R =
normal_pdf (5%:R / 2%:R - s) (2%:R^-1) b ×
normal_pdf (19%:R / 5%:R - 2%:R × s) (2%:R^-1) b ×
normal_pdf (9%:R / 2%:R - 3%:R × s) (2%:R^-1) b ×
normal_pdf (31%:R / 5%:R - 4%:R × s) (2%:R^-1) b ×
normal_pdf (8%:R - 5%:R × s) (2%:R^-1) b.
Proof.
rewrite (normal_pdf_recenter (5%:R/2%:R) s 1 (2%:R^-1) b half_neq0)
(normal_pdf_recenter (19%:R/5%:R) s 2%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter (9%:R/2%:R) s 3%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter (31%:R/5%:R) s 4%:R (2%:R^-1) b half_neq0)
(normal_pdf_recenter 8%:R s 5%:R (2%:R^-1) b half_neq0).
congr (_ × _ × _ × _ × _); rewrite ?mul1r //.
Qed.
End bayesian_phase2_slope.
End normal_density_algebra.