Stat.Minimax.Multinomial.Two­Sample­L1

lemma twoSampleLaw_singleton {k : ℕ} (n : ℕ) (R S : ProbabilitySimplex k) (x : (Fin n → Fin k) × (Fin n → Fin k)) : twoSampleLaw n (R, S) {x} = (∏ t : Fin n, (simplexPMF R) (x.1 t)) * (∏ t : Fin n, (simplexPMF S) (x.2 t)) := by have hsingle : ({x} : Set ((Fin n → Fin k) × (Fin n → Fin k))) = {x.1} ×ˢ {x.2} := by ext; simp [Prod.ext_iff] rw [twoSampleLaw, hsingle, Measure.prod_prod] simp only [simplexSampleLaw, Measure.pi_singleton] congr 1 <;> apply Finset.prod_congr rfl <;> intro t _ <;> exact PMF.toMeasure_apply_singleton _ _ (MeasurableSet.singleton _)

theorem exists_scalarMomentPriors (L : ℕ) (hL : 0 < L) : Nonempty (ScalarMomentPriors L) := by classical obtain ⟨D⟩ := exists_finiteMomentDual (f := fun x : ℝ => |x|) (r := -1) (s := 1) (by norm_num) (continuous_abs.continuousOn) L let T : ℝ := ∑ i, D.weights i * |D.nodes i| let a : Fin (L + 2) → ℝ := fun i => if 0 ≤ T then D.weights i else -D.weights i have ha_abs (i : Fin (L + 2)) : |a i| = |D.weights i| := by simp only [a] split_ifs <;> simp have ha_sum_abs : ∑ i, |a i| = 1 := by simpa only [ha_abs] using D.weights_normalized have ha_mom (j : ℕ) (hj : j ≤ L) : ∑ i, a i * D.nodes i ^ j = 0 := by by_cases hT : 0 ≤ T · simpa only [a, if_pos hT] using D.moments_zero j hj · simp only [a, if_neg hT, neg_mul, Finset.sum_neg_distrib, D.moments_zero j hj, neg_zero] have ha_target : ∑ i, a i * |D.nodes i| = bestUniformApproxErrorAbs L := by have hD : |T| = bestUniformApproxErrorAbs L := by simpa only [T, bestUniformApproxErrorAbs, uniformApproxErrorAbs, bestUniformApproxError, uniformApproxError, intervalSupNorm, symmUnitInterval] using D.target_abs_eq by_cases hT : 0 ≤ T · simp only [a, if_pos hT] change T = bestUniformApproxErrorAbs L simpa only [abs_of_nonneg hT] using hD · have hneg : T ≤ 0 := le_of_lt (lt_of_not_ge hT) simp only [a, if_neg hT, neg_mul, Finset.sum_neg_distrib] change -T = bestUniformApproxErrorAbs L simpa only [abs_of_nonpos hneg] using hD have ha_sum : ∑ i, a i = 0 := by simpa only [pow_zero, mul_one] using ha_mom 0 (zero_le L) have hparts_abs : (∑ i, (a i)⁺) + (∑ i, (a i)⁻) = 1 := by rw [← Finset.sum_add_distrib] simpa only [posPart_add_negPart] using ha_sum_abs have hparts_diff : (∑ i, (a i)⁺) - (∑ i, (a i)⁻) = 0 := by rw [← Finset.sum_sub_distrib] simpa only [posPart_sub_negPart] using ha_sum have hpos : ∑ i, 2 * (a i)⁺ = 1 := by rw [← Finset.mul_sum] linarith have hneg : ∑ i, 2 * (a i)⁻ = 1 := by rw [← Finset.mul_sum] linarith have hmoment (j : ℕ) (hj : j ≤ L) : ∑ i, (2 * (a i)⁻) * D.nodes i ^ j = ∑ i, (2 * (a i)⁺) * D.nodes i ^ j := by have h : ∑ i, ((2 * (a i)⁺) * D.nodes i ^ j - (2 * (a i)⁻) * D.nodes i ^ j) = 0 := by calc _ = 2 * ∑ i, a i * D.nodes i ^ j := by rw [Finset.mul_sum] apply Finset.sum_congr rfl intro i _ calc _ = 2 * ((a i)⁺ - (a i)⁻) * D.nodes i ^ j := by ring _ = 2 * (a i * D.nodes i ^ j) := by rw [posPart_sub_negPart]; ring _ = 0 := by rw [ha_mom j hj, mul_zero] rw [Finset.sum_sub_distrib] at h linarith have hgap : (1 / 50 : ℝ) / (L : ℝ) ≤ ∑ i, (2 * (a i)⁺) * |D.nodes i| - ∑ i, (2 * (a i)⁻) * |D.nodes i| := by have htarget : ∑ i, (2 * (a i)⁺) * |D.nodes i| - ∑ i, (2 * (a i)⁻) * |D.nodes i| = 2 * bestUniformApproxErrorAbs L := by calc _ = 2 * ∑ i, a i * |D.nodes i| := by rw [Finset.mul_sum, ← Finset.sum_sub_distrib] apply Finset.sum_congr rfl intro i _ calc _ = 2 * ((a i)⁺ - (a i)⁻) * |D.nodes i| := by ring _ = 2 * (a i * |D.nodes i|) := by rw [posPart_sub_negPart]; ring _ = 2 * bestUniformApproxErrorAbs L := by rw [ha_target] rw [htarget] calc (1 / 50 : ℝ) / (L : ℝ) = 2 * ((1 / 100 : ℝ) / (L : ℝ)) := by ring _ ≤ 2 * bestUniformApproxErrorAbs L := by gcongr exact bestUniformApproxErrorAbs_lower L hL exact ⟨{ m := L + 2 m_pos := by omega node := D.nodes node_mem := D.nodes_mem w₀ := fun i => 2 * (a i)⁻ w₁ := fun i => 2 * (a i)⁺ w₀_nonneg := fun i => mul_nonneg (by norm_num) (negPart_nonneg _) w₁_nonneg := fun i => mul_nonneg (by norm_num) (posPart_nonneg _) w₀_sum := hneg w₁_sum := hpos moments_eq := hmoment abs_gap := hgap }⟩

theorem pairedBase_tilt_l1 (b : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (u : Fin b → ℝ) (hu : ∀ j, |u j| ≤ 1) : simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 u hu) = t / (b : ℝ) * ∑ j : Fin b, |u j| := by have hb' : (b : ℝ) ≠ 0 := by exact_mod_cast (Nat.ne_of_gt hb) have hden : 0 < (b : ℝ) * 2 := by positivity unfold simplexL1 rw [← Equiv.sum_comp (finProdFinEquiv : Fin b × Fin 2 ≃ Fin (b * 2))] simp only [Fintype.sum_prod_type, Fin.sum_univ_two] simp only [pairedBaseVector, pairedTiltVector, Subtype.coe_mk, Equiv.symm_apply_apply, Fin.isValue, ↓reduceIte, one_ne_zero] rw [Finset.mul_sum] apply Finset.sum_congr rfl intro j hj have hfirst : 1 / ((b : ℝ) * 2) - (1 + t * u j) / ((b : ℝ) * 2) = -(t * u j) / ((b : ℝ) * 2) := by ring have hsecond : 1 / ((b : ℝ) * 2) - (1 + -(t * u j)) / ((b : ℝ) * 2) = (t * u j) / ((b : ℝ) * 2) := by ring rw [hfirst, hsecond, abs_div, abs_div, abs_neg, abs_mul, abs_of_nonneg ht, abs_of_pos hden] field_simp ring

theorem largeSample_fuzzyCertificate : ∃ a : ℝ, 0 < a ∧ ∀ k n : ℕ, 2 ≤ k → k ^ 2 < n → ∃ W : FuzzyCertificate n k, a * ((k : ℝ) / ((n : ℝ) * logAlphabet k)) ≤ W.delta ^ 2 := by obtain ⟨a₁, ha₁, K, hK, hlarge⟩ := largeAlphabet_fuzzyCertificate obtain ⟨a₂, ha₂, hbounded⟩ := boundedAlphabet_fuzzyCertificate K hK refine ⟨min a₁ a₂, lt_min ha₁ ha₂, ?_⟩ intro k n hk hn have hrate : 0 ≤ (k : ℝ) / ((n : ℝ) * logAlphabet k) := by have hnpos : (0 : ℝ) < n := by exact_mod_cast (by omega : 0 < n) have hlog : 0 < logAlphabet k := by have h := logAlphabet_one_le k (by omega) linarith positivity rcases lt_or_ge k K with hsmall | hlargeK · obtain ⟨W, hW⟩ := hbounded k n hk hsmall hn exact ⟨W, (mul_le_mul_of_nonneg_right (min_le_right a₁ a₂) hrate).trans hW⟩ · obtain ⟨W, hW⟩ := hlarge k n hlargeK hn exact ⟨W, (mul_le_mul_of_nonneg_right (min_le_left a₁ a₂) hrate).trans hW⟩

theorem pairedProductTarget_mean {L : ℕ} (P : ScalarMomentPriors L) (b : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (side : Bool) : ∑ u : Fin b → Fin P.m, (pairedProductWeight P b side u).toReal * simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u)) = t * ∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i| := by classical let w : Fin P.m → ℝ := scalarPriorWeight P side have hw (i : Fin P.m) : 0 ≤ w i := by dsimp [w] cases side <;> simp [scalarPriorWeight, P.w₀_nonneg, P.w₁_nonneg] have hwsum : ∑ i : Fin P.m, w i = 1 := by dsimp [w] cases side <;> simp [scalarPriorWeight, P.w₀_sum, P.w₁_sum] have hweight (u : Fin b → Fin P.m) : (pairedProductWeight P b side u).toReal = ∏ j : Fin b, w (u j) := by simp [pairedProductWeight, ENNReal.toReal_prod, ENNReal.toReal_ofReal (hw _), w] have hmarg (j : Fin b) : ∑ u : Fin b → Fin P.m, (∏ k : Fin b, w (u k)) * |P.node (u j)| = ∑ i : Fin P.m, w i * |P.node i| := by have hfactor (u : Fin b → Fin P.m) : (∏ k : Fin b, w (u k)) * |P.node (u j)| = ∏ k : Fin b, (if k = j then w (u k) * |P.node (u k)| else w (u k)) := by calc _ = (∏ k : Fin b, w (u k)) * ∏ k : Fin b, (if k = j then |P.node (u k)| else 1) := by simp _ = ∏ k : Fin b, w (u k) * (if k = j then |P.node (u k)| else 1) := by rw [Finset.prod_mul_distrib] _ = _ := by simp calc _ = ∑ u : Fin b → Fin P.m, ∏ k : Fin b, (if k = j then w (u k) * |P.node (u k)| else w (u k)) := by apply Finset.sum_congr rfl intro u _ exact hfactor u _ = ∏ k : Fin b, ∑ i : Fin P.m, (if k = j then w i * |P.node i| else w i) := by exact (Fintype.prod_sum (fun k : Fin b => fun i : Fin P.m => if k = j then w i * |P.node i| else w i)).symm _ = ∑ i : Fin P.m, w i * |P.node i| := by simp [hwsum, Finset.prod_ite_eq'] simp_rw [pairedBase_tilt_l1, hweight] calc _ = (t / (b : ℝ)) * ∑ j : Fin b, ∑ u : Fin b → Fin P.m, (∏ k : Fin b, w (u k)) * |P.node (u j)| := by rw [Finset.mul_sum] calc _ = ∑ u : Fin b → Fin P.m, ∑ j : Fin b, t / (b : ℝ) * ((∏ k : Fin b, w (u k)) * |P.node (u j)|) := by apply Finset.sum_congr rfl intro u _ rw [Finset.mul_sum] rw [Finset.mul_sum] apply Finset.sum_congr rfl intro j _ simp only [pairedNodeVector] ring _ = ∑ j : Fin b, ∑ u : Fin b → Fin P.m, t / (b : ℝ) * ((∏ k : Fin b, w (u k)) * |P.node (u j)|) := Finset.sum_comm _ = _ := by apply Finset.sum_congr rfl intro j _ rw [← Finset.mul_sum] _ = t * ∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i| := by simp [hmarg, w] have hb' : (b : ℝ) ≠ 0 := by exact_mod_cast (Nat.ne_of_gt hb) field_simp

theorem pairedPoissonPredictive_tv_le {L : ℕ} (P : ScalarMomentPriors L) (b : ℕ) (lambda t : ℝ) (hlambda : 0 ≤ lambda) (ht : 0 ≤ t) (ht1 : t ≤ 1) (hscale : 100 * lambda * t ^ 2 ≤ (L : ℝ)) : tvDist (pairedPoissonPredictive P b lambda t false) (pairedPoissonPredictive P b lambda t true) ≤ (b : ℝ) * (2 : ℝ) ^ (-(L : ℝ) / 4) := by classical have hpair (u : ℝ) : IsProbabilityMeasure (scalarPoissonPairLaw lambda t u) := by unfold scalarPoissonPairLaw infer_instance have hscalar (side : Bool) : IsProbabilityMeasure (scalarPoissonPredictive P lambda t side) := by haveI (i : Fin P.m) : IsProbabilityMeasure (scalarPoissonPairLaw lambda t (P.node i)) := hpair _ unfold scalarPoissonPredictive apply mixture_isProbabilityMeasure cases side · simp only [Bool.false_eq_true, ↓reduceIte] rw [← ENNReal.ofReal_sum_of_nonneg (fun i _ => P.w₀_nonneg i), P.w₀_sum] simp · simp only [↓reduceIte] rw [← ENNReal.ofReal_sum_of_nonneg (fun i _ => P.w₁_nonneg i), P.w₁_sum] simp have hfactor (side : Bool) : pairedPoissonPredictive P b lambda t side = Measure.pi (fun _ : Fin b => scalarPoissonPredictive P lambda t side) := by apply Measure.ext_of_singleton intro x haveI (u : ℝ) : IsProbabilityMeasure (scalarPoissonPairLaw lambda t u) := hpair u haveI : IsProbabilityMeasure (scalarPoissonPredictive P lambda t side) := hscalar side rw [Measure.pi_singleton] simp only [pairedPoissonPredictive, mixture_apply, Measure.pi_singleton, pairedProductWeight, scalarPoissonPredictive, scalarPriorWeight] rw [Fintype.prod_sum] congr 1 funext u rw [← Finset.prod_mul_distrib] haveI : IsProbabilityMeasure (scalarPoissonPredictive P lambda t false) := hscalar false haveI : IsProbabilityMeasure (scalarPoissonPredictive P lambda t true) := hscalar true rw [hfactor false, hfactor true] exact (tvDist_pi_iid_le b _ _).trans (mul_le_mul_of_nonneg_left (scalarPoissonPredictive_tv_le P lambda t hlambda ht ht1 hscale) (Nat.cast_nonneg b))

theorem scalarPoissonPairLikelihood_inner (lambda t u v : ℝ) (hlambda : 0 ≤ lambda) : (∫ z, scalarPoissonPairLikelihood t u z * scalarPoissonPairLikelihood t v z ∂scalarPoissonPairLaw lambda t 0) = exp (2 * lambda * t ^ 2 * u * v) := by /- Factor the product-Poisson integral into two Poisson power MGFs. The one-count identity `poisson_power_mgf` is private in Causalean's MarkedPoisson/AggregatePoisson module, so prove or expose a local version. Algebra reduces the two exponents to `2*lambda*t^2*u*v`. -/ have hrate : Real.toNNReal (lambda * (1 + t * 0)) = Real.toNNReal lambda := by ring have hrate' : Real.toNNReal (lambda * (1 - t * 0)) = Real.toNNReal lambda := by ring have hcast : ((Real.toNNReal lambda : NNReal) : ℝ) = lambda := Real.coe_toNNReal lambda hlambda have hfactor (z : ℕ × ℕ) : scalarPoissonPairLikelihood t u z * scalarPoissonPairLikelihood t v z = ((1 + t * u) * (1 + t * v)) ^ z.1 * ((1 - t * u) * (1 - t * v)) ^ z.2 := by simp only [scalarPoissonPairLikelihood, mul_pow] ring simp_rw [hfactor] simp only [scalarPoissonPairLaw, hrate, hrate'] rw [integral_prod_mul, poisson_power_mgf, poisson_power_mgf, ← Real.exp_add] rw [hcast] congr 1 ring

theorem finitePoissonSampleLaw_map_prefix {X : Type*} [MeasurableSpace X] (P : Measure X) [IsProbabilityMeasure P] (lam : ℝ≥0) (x₀ : X) (n : ℕ) : Measure.map (finitePoissonPrefix x₀ n) (finitePoissonSampleLaw P lam) = (poissonMeasure lam) (Set.Ici n) • Measure.pi (fun _ : Fin n => P) + (poissonMeasure lam) (Set.Iio n) • Measure.dirac (fun _ : Fin n => x₀) := by let μ := finitePoissonSampleLaw P lam let good : Set (FiniteSample X) := FiniteSample.count ⁻¹' Set.Ici n let bad : Set (FiniteSample X) := FiniteSample.count ⁻¹' Set.Iio n have hgood_meas : MeasurableSet good := (by measurability : MeasurableSet (Set.Ici n)).preimage measurable_finiteSample_count have hbad_meas : MeasurableSet bad := (by measurability : MeasurableSet (Set.Iio n)).preimage measurable_finiteSample_count have hgood : Measure.map (finitePoissonPrefix x₀ n) (μ.restrict good) = (poissonMeasure lam) (Set.Ici n) • Measure.pi (fun _ : Fin n => P) := by have hdecomp : μ.restrict good = Measure.sum (fun m : (Set.Ici n : Set ℕ) => μ.restrict (FiniteSample.count ⁻¹' ({m.1} : Set ℕ))) := by dsimp [good] rw [← Set.biUnion_preimage_singleton] exact Measure.restrict_biUnion (Set.to_countable (Set.Ici n)) (pairwiseDisjoint_fiber FiniteSample.count (Set.Ici n)) (fun m => measurable_finiteSample_count (measurableSet_singleton m)) rw [hdecomp, Measure.map_sum (measurable_finitePoissonPrefix x₀ n).aemeasurable] simp_rw [show ∀ m : (Set.Ici n : Set ℕ), Measure.map (finitePoissonPrefix x₀ n) (μ.restrict (FiniteSample.count ⁻¹' ({m.1} : Set ℕ))) = (poissonMeasure lam) ({m.1} : Set ℕ) • Measure.pi (fun _ : Fin n => P) by intro m rw [show μ = finitePoissonSampleLaw P lam by rfl, finitePoissonSampleLaw_restrict_count_eq, Measure.map_smul, Measure.map_map (measurable_finitePoissonPrefix x₀ n) (measurable_fixedSizeEmbed m.1)] have hfun : finitePoissonPrefix x₀ n ∘ fixedSizeEmbed m.1 = fun x : Fin m.1 → X => fun i : Fin n => x (castLE m.2 i) := by have hm : n ≤ m.1 := m.2 funext x funext i simp [Function.comp_apply, finitePoissonPrefix, fixedSizeEmbed, FiniteSample.count, FiniteSample.points, hm] rw [hfun, map_prefixCoordinates_pi P m.2]] ext A hA simp only [Measure.sum_apply _ hA, Measure.smul_apply, smul_eq_mul] rw [ENNReal.tsum_mul_right] congr 1 simpa using (tsum_measure_preimage_singleton (μ := poissonMeasure lam) (s := Set.Ici n) (f := id) (Set.to_countable (Set.Ici n)) (fun m _ => measurableSet_singleton m)) have hbad : Measure.map (finitePoissonPrefix x₀ n) (μ.restrict bad) = (poissonMeasure lam) (Set.Iio n) • Measure.dirac (fun _ : Fin n => x₀) := by have heq : finitePoissonPrefix x₀ n =ᵐ[μ.restrict bad] fun _ => (fun _ : Fin n => x₀) := by apply ae_restrict_of_forall_mem hbad_meas intro s hs have hlt : s.count < n := hs simp [finitePoissonPrefix, not_le.mpr hlt] rw [Measure.map_congr heq, Measure.map_const] congr 1 rw [Measure.restrict_apply_univ] change μ (FiniteSample.count ⁻¹' Set.Iio n) = _ rw [show μ = finitePoissonSampleLaw P lam by rfl, ← finitePoissonSampleLaw_map_count P lam, Measure.map_apply measurable_finiteSample_count (by measurability)] have hcompl : bad = goodᶜ := by ext s simp [good, bad] change Measure.map (finitePoissonPrefix x₀ n) μ = _ rw [← (Measure.restrict_add_restrict_compl hgood_meas : μ.restrict good + μ.restrict goodᶜ = μ), ← hcompl, Measure.map_add _ _ (measurable_finitePoissonPrefix x₀ n), hgood, hbad]

theorem fixedSampleMixture_tv_le_finitePoisson (n k m : ℕ) (hk : 0 < k) (lam : ℝ≥0) (w₀ w₁ : Fin m → ℝ≥0∞) (hw₀ : ∑ i, w₀ i = 1) (hw₁ : ∑ i, w₁ i = 1) (R₀ R₁ : Fin m → ProbabilitySimplex k) : tvDist (mixture w₀ (fun i => simplexSampleLaw (R₀ i) n)) (mixture w₁ (fun i => simplexSampleLaw (R₁ i) n)) ≤ tvDist (mixture w₀ (fun i => finitePoissonSampleLaw (simplexPMF (R₀ i)).toMeasure lam)) (mixture w₁ (fun i => finitePoissonSampleLaw (simplexPMF (R₁ i)).toMeasure lam)) + 2 * (poissonMeasure lam (Iio n)).toReal := by let x₀ : Fin k := ⟨0, hk⟩ let f := finitePoissonPrefix x₀ n let q := poissonMeasure lam (Ici n) let p := poissonMeasure lam (Iio n) let δ : Measure (Fin n → Fin k) := Measure.dirac (fun _ => x₀) let Q (w : Fin m → ℝ≥0∞) (R : Fin m → ProbabilitySimplex k) := mixture w (fun i => finitePoissonSampleLaw (simplexPMF (R i)).toMeasure lam) have hqp : q + p = 1 := by change (poissonMeasure lam) (Ici n) + (poissonMeasure lam) (Iio n) = 1 rw [show Iio n = (Ici n)ᶜ by ext i; simp] rw [measure_add_measure_compl (by measurability), measure_univ] have hmap (w : Fin m → ℝ≥0∞) (hw : ∑ i, w i = 1) (R : Fin m → ProbabilitySimplex k) : Measure.map f (Q w R) = q • mixture w (fun i => simplexSampleLaw (R i) n) + p • δ := by change Measure.map f (∑ i, w i • finitePoissonSampleLaw (simplexPMF (R i)).toMeasure lam) = _ rw [Measure.map_finset_sum' (measurable_finitePoissonPrefix x₀ n).aemeasurable] simp_rw [Measure.map_smul, finitePoissonSampleLaw_map_prefix] simp_rw [smul_add, smul_smul] rw [Finset.sum_add_distrib] congr 1 · rw [mixture, Finset.smul_sum] congr 1 ext i simp [q, simplexSampleLaw, smul_smul, mul_comm] · rw [← Finset.sum_smul, ← Finset.sum_mul, hw, one_mul] have hprobP (w : Fin m → ℝ≥0∞) (hw : ∑ i, w i = 1) (R : Fin m → ProbabilitySimplex k) : IsProbabilityMeasure (mixture w (fun i => simplexSampleLaw (R i) n)) := mixture_isProbabilityMeasure w hw _ have hprobQ (w : Fin m → ℝ≥0∞) (hw : ∑ i, w i = 1) (R : Fin m → ProbabilitySimplex k) : IsProbabilityMeasure (Q w R) := by dsimp [Q] exact mixture_isProbabilityMeasure w hw _ let μ₀ := mixture w₀ (fun i => simplexSampleLaw (R₀ i) n) let μ₁ := mixture w₁ (fun i => simplexSampleLaw (R₁ i) n) let ρ₀ := Q w₀ R₀ let ρ₁ := Q w₁ R₁ letI : IsProbabilityMeasure μ₀ := hprobP w₀ hw₀ R₀ letI : IsProbabilityMeasure μ₁ := hprobP w₁ hw₁ R₁ letI : IsProbabilityMeasure ρ₀ := hprobQ w₀ hw₀ R₀ letI : IsProbabilityMeasure ρ₁ := hprobQ w₁ hw₁ R₁ letI : IsProbabilityMeasure (Measure.map f ρ₀) := Measure.isProbabilityMeasure_map (measurable_finitePoissonPrefix x₀ n).aemeasurable letI : IsProbabilityMeasure (Measure.map f ρ₁) := Measure.isProbabilityMeasure_map (measurable_finitePoissonPrefix x₀ n).aemeasurable have hmap₀ : Measure.map f ρ₀ = q • μ₀ + p • δ := hmap w₀ hw₀ R₀ have hmap₁ : Measure.map f ρ₁ = q • μ₁ + p • δ := hmap w₁ hw₁ R₁ have hqfin : q ≠ ⊤ := ne_top_of_le_ne_top ENNReal.one_ne_top (le_of_le_of_eq (le_add_right le_rfl) hqp) have hpfin : p ≠ ⊤ := ne_top_of_le_ne_top ENNReal.one_ne_top (le_of_le_of_eq (le_add_left le_rfl) hqp) have hqpReal : q.toReal + p.toReal = 1 := by rw [← ENNReal.toReal_add hqfin hpfin, hqp] rfl have hcontract : tvDist (Measure.map f ρ₀) (Measure.map f ρ₁) ≤ tvDist ρ₀ ρ₁ := by unfold tvDist apply ciSup_le rintro ⟨A, hA⟩ rw [map_measureReal_apply_of_aemeasurable (measurable_finitePoissonPrefix x₀ n).aemeasurable hA, map_measureReal_apply_of_aemeasurable (measurable_finitePoissonPrefix x₀ n).aemeasurable hA] exact le_ciSup bddAbove_tvDist_range ⟨f ⁻¹' A, hA.preimage (measurable_finitePoissonPrefix x₀ n)⟩ calc tvDist μ₀ μ₁ ≤ tvDist (Measure.map f ρ₀) (Measure.map f ρ₁) + 2 * p.toReal := by unfold tvDist apply ciSup_le rintro ⟨A, hA⟩ have hmass₀ : (Measure.map f ρ₀).real A = q.toReal * μ₀.real A + p.toReal * δ.real A := by rw [hmap₀, measureReal_add_apply (by rw [Measure.smul_apply, smul_eq_mul]; exact ENNReal.mul_ne_top hqfin (measure_ne_top μ₀ A)) (by rw [Measure.smul_apply, smul_eq_mul]; exact ENNReal.mul_ne_top hpfin (measure_ne_top δ A))] simp have hmass₁ : (Measure.map f ρ₁).real A = q.toReal * μ₁.real A + p.toReal * δ.real A := by rw [hmap₁, measureReal_add_apply (by rw [Measure.smul_apply, smul_eq_mul]; exact ENNReal.mul_ne_top hqfin (measure_ne_top μ₁ A)) (by rw [Measure.smul_apply, smul_eq_mul]; exact ENNReal.mul_ne_top hpfin (measure_ne_top δ A))] simp have hgap : |μ₀.real A - μ₁.real A| ≤ 1 := Causalean.Stat.abs_measureReal_sub_le_one A have hqnonneg : 0 ≤ q.toReal := ENNReal.toReal_nonneg have hpnonneg : 0 ≤ p.toReal := ENNReal.toReal_nonneg have habs : |(Measure.map f ρ₀).real A - (Measure.map f ρ₁).real A| = q.toReal * |μ₀.real A - μ₁.real A| := by rw [hmass₀, hmass₁] have heq : q.toReal * μ₀.real A + p.toReal * δ.real A - (q.toReal * μ₁.real A + p.toReal * δ.real A) = q.toReal * (μ₀.real A - μ₁.real A) := by ring rw [heq, abs_mul, abs_of_nonneg hqnonneg] have htv : q.toReal * |μ₀.real A - μ₁.real A| ≤ ⨆ A : {A : Set (Fin n → Fin k) // MeasurableSet A}, |(Measure.map f ρ₀).real A.1 - (Measure.map f ρ₁).real A.1| := by rw [← habs] exact le_ciSup bddAbove_tvDist_range ⟨A, hA⟩ have hpbound := mul_le_mul_of_nonneg_left hgap hpnonneg nlinarith _ ≤ tvDist ρ₀ ρ₁ + 2 * p.toReal := by gcongr

theorem FuzzyCertificate.minimax_lower {n k : ℕ} (W : FuzzyCertificate n k) : 11 * W.delta ^ 2 / 512 ≤ twoSampleL1MinimaxRisk n k := by classical let Θ := Fin W.m × Bool let X := (Fin n → Fin k) × (Fin n → Fin k) let model : Θ → ProbabilitySimplex k × ProbabilitySimplex k := fun θ => if θ.2 then W.θ₁ θ.1 else W.θ₀ θ.1 let K : Kernel Θ X := Kernel.ofFunOfCountable (fun θ => twoSampleLaw n (model θ)) let target : Θ → ℝ := fun θ => simplexL1 (model θ).1 (model θ).2 let π₀ : Measure Θ := ∑ i : Fin W.m, W.w₀ i • Measure.dirac (i, false) let π₁ : Measure Θ := ∑ i : Fin W.m, W.w₁ i • Measure.dirac (i, true) have hp₀ : IsProbabilityMeasure π₀ := by apply isProbabilityMeasure_iff.mpr simp [π₀, Measure.coe_finsetSum, W.w₀_sum] have hp₁ : IsProbabilityMeasure π₁ := by apply isProbabilityMeasure_iff.mpr simp [π₁, Measure.coe_finsetSum, W.w₁_sum] letI : IsProbabilityMeasure π₀ := hp₀ letI : IsProbabilityMeasure π₁ := hp₁ have hK : ∀ θ, IsProbabilityMeasure (K θ) := by intro θ change IsProbabilityMeasure (twoSampleLaw n (model θ)) infer_instance have htarget : Measurable target := measurable_of_finite target have hmass₀ (A : Set Θ) : π₀ A = ∑ i, if (i, false) ∈ A then W.w₀ i else 0 := by simp [π₀, Measure.coe_finsetSum, Measure.dirac_apply, indicator, smul_eq_mul] have hmass₁ (A : Set Θ) : π₁ A = ∑ i, if (i, true) ∈ A then W.w₁ i else 0 := by simp [π₁, Measure.coe_finsetSum, Measure.dirac_apply, indicator, smul_eq_mul] have hbad₀ : π₀.real {θ | W.delta / 4 < |target θ - W.center₀|} ≤ 1 / 8 := by rw [Measure.real, hmass₀] simp only [target, model] exact (ENNReal.toReal_mono (by norm_num) W.bad₀).trans_eq (by norm_num) have hbad₁ : π₁.real {θ | W.delta / 4 < |target θ - W.center₁|} ≤ 1 / 8 := by rw [Measure.real, hmass₁] simp only [target, model] exact (ENNReal.toReal_mono (by norm_num) W.bad₁).trans_eq (by norm_num) have hpredict₀ : priorPredictive π₀ K = mixture W.w₀ (fun i => twoSampleLaw n (W.θ₀ i)) := by ext A hA rw [priorPredictive_apply π₀ K hA] simp only [π₀, lintegral_finsetSum_measure, lintegral_smul_measure, lintegral_dirac, mixture_apply, smul_eq_mul] rfl have hpredict₁ : priorPredictive π₁ K = mixture W.w₁ (fun i => twoSampleLaw n (W.θ₁ i)) := by ext A hA rw [priorPredictive_apply π₁ K hA] simp only [π₁, lintegral_finsetSum_measure, lintegral_smul_measure, lintegral_dirac, mixture_apply, smul_eq_mul] rfl have htv : tvDist (priorPredictive π₀ K) (priorPredictive π₁ K) ≤ 1 / 16 := by rw [hpredict₀, hpredict₁] exact W.predictive_close let c : ℝ := 11 * W.delta ^ 2 / 512 have hbdd (est : TwoSampleEstimator n k) : BddAbove (range (twoSampleL1Risk n est)) := by let C : ℝ := ∑ z : X, (|est.1 z| + 2)^2 refine ⟨C, ?_⟩ rintro _ ⟨RS, rfl⟩ have hpoint (z : X) : (est.1 z - simplexL1 RS.1 RS.2)^2 ≤ C := by have hval := simplexL1_nonneg_le_two RS.1 RS.2 have habs : |est.1 z - simplexL1 RS.1 RS.2| ≤ |est.1 z| + 2 := by calc _ ≤ |est.1 z| + |simplexL1 RS.1 RS.2| := abs_sub _ _ _ ≤ |est.1 z| + 2 := by rw [abs_of_nonneg hval.1]; linarith have hsq : (est.1 z - simplexL1 RS.1 RS.2)^2 ≤ (|est.1 z| + 2)^2 := by rw [← sq_abs] exact (sq_le_sq₀ (abs_nonneg _) (by positivity)).2 habs exact hsq.trans (Finset.single_le_sum (fun x _ => sq_nonneg (|est.1 x| + 2)) (Finset.mem_univ z)) change (∫ z, (est.1 z - simplexL1 RS.1 RS.2)^2 ∂twoSampleLaw n RS) ≤ C calc _ ≤ ∫ _z, C ∂twoSampleLaw n RS := integral_mono_ae Integrable.of_finite (integrable_const C) (ae_of_all _ hpoint) _ = C := by simp have hrisk (est : TwoSampleEstimator n k) (θ : Θ) : squaredRisk K target est.1 θ = ofReal (twoSampleL1Risk n est (model θ)) := by unfold squaredRisk twoSampleL1Risk change (∫⁻ z, ofReal ((est.1 z - target θ)^2) ∂twoSampleLaw n (model θ)) = ofReal (∫ z, (est.1 z - target θ)^2 ∂twoSampleLaw n (model θ)) exact (ofReal_integral_eq_lintegral_ofReal Integrable.of_finite (ae_of_all _ fun z => sq_nonneg _)).symm have hall (est : TwoSampleEstimator n k) : c ≤ worstCaseRiskReal (twoSampleL1Risk n) est := by have hbound := Causalean.Stat.Minimax.FuzzyHypotheses.twoFuzzyHypotheses_worstCase_lower_standard π₀ π₁ K target hK est.1 est.2 htarget W.center₀ W.center₁ W.delta W.delta_pos W.center_sep hbad₀ hbad₁ htv have hsup : worstCaseSquaredRisk K target est.1 ≤ ofReal (worstCaseRiskReal (twoSampleL1Risk n) est) := by unfold worstCaseSquaredRisk apply iSup_le intro θ rw [hrisk] exact ENNReal.ofReal_le_ofReal (le_worstCaseRisk (hbdd est) (model θ)) have hnonneg : 0 ≤ worstCaseRiskReal (twoSampleL1Risk n) est := by have hθ := le_worstCaseRisk (hbdd est) (W.θ₀ ⟨0, W.m_pos⟩) have hr : 0 ≤ twoSampleL1Risk n est (W.θ₀ ⟨0, W.m_pos⟩) := by unfold twoSampleL1Risk exact integral_nonneg fun z => sq_nonneg _ exact hr.trans hθ exact (ENNReal.ofReal_le_ofReal_iff hnonneg).mp (hbound.trans hsup) haveI : Nonempty (TwoSampleEstimator n k) := ⟨⟨fun _ => 0, measurable_const⟩⟩ exact le_minimaxValue hall

theorem twoSampleL1MinimaxRisk_largeSample_lower : ∃ c : ℝ, 0 < c ∧ ∀ k n : ℕ, 2 ≤ k → k ^ 2 < n → c * ((k : ℝ) / ((n : ℝ) * logAlphabet k)) ≤ twoSampleL1MinimaxRisk n k := by obtain ⟨a, ha, hwitness⟩ := largeSample_fuzzyCertificate refine ⟨11 * a / 512, by positivity, ?_⟩ intro k n hk hn obtain ⟨W, hsep⟩ := hwitness k n hk hn calc 11 * a / 512 * ((k : ℝ) / ((n : ℝ) * logAlphabet k)) = (11 / 512 : ℝ) * (a * ((k : ℝ) / ((n : ℝ) * logAlphabet k))) := by ring _ ≤ (11 / 512 : ℝ) * W.delta ^ 2 := mul_le_mul_of_nonneg_left hsep (by norm_num) _ = 11 * W.delta ^ 2 / 512 := by ring _ ≤ twoSampleL1MinimaxRisk n k := W.minimax_lower

theorem twoSampleLaw_pad {r k n : ℕ} (h : r ≤ k) (RS : ProbabilitySimplex r × ProbabilitySimplex r) : Measure.map (fun z : (Fin n → Fin r) × (Fin n → Fin r) => ((fun i => castLE h (z.1 i)), (fun i => castLE h (z.2 i)))) (twoSampleLaw n RS) = twoSampleLaw n (padSimplex h RS.1, padSimplex h RS.2) := by let f : (Fin n → Fin r) → (Fin n → Fin k) := fun z i => castLE h (z i) change Measure.map (map f f) (twoSampleLaw n RS) = twoSampleLaw n (padSimplex h RS.1, padSimplex h RS.2) unfold twoSampleLaw rw [← Measure.map_prod_map (simplexSampleLaw RS.1 n) (simplexSampleLaw RS.2 n) (measurable_of_finite f) (measurable_of_finite f)] rw [simplexSampleLaw_pad h RS.1, simplexSampleLaw_pad h RS.2]

theorem pairedFixedPredictive_tv_le {L : ℕ} (P : ScalarMomentPriors L) (b n : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (hscale : 100 * ((n : ℝ) / (b : ℝ)) * t ^ 2 ≤ (L : ℝ)) : tvDist (pairedFixedPredictive P b n hb t ht ht1 false) (pairedFixedPredictive P b n hb t ht ht1 true) ≤ (b : ℝ) * (2 : ℝ) ^ (-(L : ℝ) / 4) + 2 * (poissonMeasure (toNNReal (2 * (n : ℝ))) (Iio n)).toReal := by classical let μ := simplexSampleLaw (pairedBaseVector b hb) n let R (u : Fin b → Fin P.m) := pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u) let ν (side : Bool) := mixture (pairedProductWeight P b side) (fun u => simplexSampleLaw (R u) n) have hfactor (side : Bool) : pairedFixedPredictive P b n hb t ht ht1 side = μ.prod (ν side) := by unfold pairedFixedPredictive mixture change (∑ u, pairedProductWeight P b side u • (μ.prod (simplexSampleLaw (R u) n))) = μ.prod (∑ u, pairedProductWeight P b side u • simplexSampleLaw (R u) n) rw [← Measure.sum_fintype] rw [← Measure.sum_fintype] rw [Measure.prod_sum_right] simp_rw [Measure.prod_smul_right] have hνprob (side : Bool) : IsProbabilityMeasure (ν side) := by change IsProbabilityMeasure (mixture (pairedProductWeight P b side) (fun u => simplexSampleLaw (R u) n)) exact mixture_isProbabilityMeasure (pairedProductWeight P b side) (pairedProductWeight_sum P b side) _ haveI : IsProbabilityMeasure (ν false) := hνprob false haveI : IsProbabilityMeasure (ν true) := hνprob true have hfixed : tvDist (ν false) (ν true) ≤ tvDist (pairedFinitePoissonPredictive P b n hb t ht ht1 false) (pairedFinitePoissonPredictive P b n hb t ht ht1 true) + 2 * (poissonMeasure (toNNReal (2 * (n : ℝ))) (Iio n)).toReal := by let e : Fin (card (Fin b → Fin P.m)) ≃ (Fin b → Fin P.m) := (Fintype.equivFin _).symm let w (side : Bool) (i : Fin (card (Fin b → Fin P.m))) := pairedProductWeight P b side (e i) let S (i : Fin (card (Fin b → Fin P.m))) := R (e i) have hw (side : Bool) : ∑ i, w side i = 1 := by exact (e.sum_comp (pairedProductWeight P b side)).trans (pairedProductWeight_sum P b side) have h := fixedSampleMixture_tv_le_finitePoisson n (b * 2) (card (Fin b → Fin P.m)) (by omega) (toNNReal (2 * (n : ℝ))) (w false) (w true) (hw false) (hw true) S S let F (side : Bool) (u : Fin b → Fin P.m) : Measure (Fin n → Fin (b * 2)) := pairedProductWeight P b side u • simplexSampleLaw (R u) n let G (side : Bool) (u : Fin b → Fin P.m) : Measure (FiniteSample (Fin (b * 2))) := pairedProductWeight P b side u • finitePoissonSampleLaw (simplexPMF (R u)).toMeasure (toNNReal (2 * (n : ℝ))) change tvDist (∑ i, F false (e i)) (∑ i, F true (e i)) ≤ tvDist (∑ i, G false (e i)) (∑ i, G true (e i)) + 2 * (poissonMeasure (toNNReal (2 * (n : ℝ))) (Iio n)).toReal at h rw [e.sum_comp, e.sum_comp, e.sum_comp, e.sum_comp] at h simpa only [ν, pairedFinitePoissonPredictive, mixture, F, G, R] using h have hcount := pairedFinitePoissonPredictive_tv_le_counts P b n hb t ht ht1 have hpoisson := pairedPoissonPredictive_tv_le P b ((n : ℝ) / (b : ℝ)) t (div_nonneg (Nat.cast_nonneg _) (Nat.cast_nonneg _)) ht ht1 hscale calc tvDist (pairedFixedPredictive P b n hb t ht ht1 false) (pairedFixedPredictive P b n hb t ht ht1 true) = tvDist (μ.prod (ν false)) (μ.prod (ν true)) := by rw [hfactor false, hfactor true] _ ≤ tvDist (ν false) (ν true) := by have h := tvDist_prod_le_add μ μ (ν false) (ν true) simpa [tvDist] using h _ ≤ _ := by linarith

theorem exists_pairedLargeParameters : ∃ a : ℝ, 0 < a ∧ ∃ K : ℕ, 2 ≤ K ∧ ∀ k n : ℕ, K ≤ k → k ^ 2 < n → ∃ b L : ℕ, ∃ t δ : ℝ, 0 < b ∧ 1 ≤ L ∧ b * 2 ≤ k ∧ 0 ≤ t ∧ t ≤ 1 ∧ 0 < δ ∧ 100 * ((n : ℝ) / (b : ℝ)) * t ^ 2 ≤ (L : ℝ) ∧ 128 * t ^ 2 ≤ (b : ℝ) * δ ^ 2 ∧ (b : ℝ) * (2 : ℝ) ^ (-(L : ℝ) / 4) + 2 * (poissonMeasure (toNNReal (2 * (n : ℝ))) (Iio n)).toReal ≤ 1 / 16 ∧ δ ≤ t * ((1 / 50 : ℝ) / (L : ℝ)) ∧ a * ((k : ℝ) / ((n : ℝ) * log (exp 1 * (k : ℝ)))) ≤ δ ^ 2 := by obtain ⟨K₀, hK₀⟩ := eventually_paired_degree_scale refine ⟨1 / 100000000, by norm_num, max K₀ 256, by omega, ?_⟩ intro k n hk hn have hkK : K₀ ≤ k := le_trans (Nat.le_max_left _ _) hk have hk256 : 256 ≤ k := le_trans (Nat.le_max_right _ _) hk have hk2 : 2 ≤ k := by omega have hk4 : 4 ≤ k := by omega have hn256 : 256 ≤ n := by have : k ^ 2 ≤ n := Nat.le_of_lt hn nlinarith let b : ℕ := k / 2 let L : ℕ := 16 * (Nat.log2 k + 1) have hbposN : 0 < b := by dsimp [b]; omega have hLposN : 0 < L := by dsimp [L]; omega have hbpos : (0 : ℝ) < b := by exact_mod_cast hbposN have hLpos : (0 : ℝ) < L := by exact_mod_cast hLposN have hnpos : (0 : ℝ) < n := by exact_mod_cast (show 0 < n by omega) have hkpos : (0 : ℝ) < k := by exact_mod_cast (show 0 < k by omega) have hscale : 320000 * (L : ℝ) ^ 2 ≤ (b : ℝ) := by simpa [L, b] using hK₀ k hkK have hLlog : (L : ℝ) ≤ 48 * log (exp 1 * (k : ℝ)) := by calc (L : ℝ) = 16 * ((Nat.log2 k + 1 : ℕ) : ℝ) := by simp [L] _ ≤ 16 * (3 * log (exp 1 * (k : ℝ))) := by gcongr exact paired_degree_log_bound k hk2 _ = _ := by ring have hlogpos : 0 < log (exp 1 * (k : ℝ)) := by apply Real.log_pos have he := Real.one_lt_exp_iff.mpr (show (0 : ℝ) < 1 by norm_num) have hk1 : (1 : ℝ) ≤ k := by exact_mod_cast (show 1 ≤ k by omega) nlinarith [Real.exp_pos 1] let t : ℝ := sqrt ((b : ℝ) * (L : ℝ) / (200 * (n : ℝ))) let δ : ℝ := t / (50 * (L : ℝ)) have ht0 : 0 ≤ t := Real.sqrt_nonneg _ have ht_sq : t ^ 2 = (b : ℝ) * (L : ℝ) / (200 * (n : ℝ)) := by exact Real.sq_sqrt (by positivity) have hδ_sq : δ ^ 2 = (b : ℝ) / (500000 * (n : ℝ) * (L : ℝ)) := by dsimp [δ] rw [div_pow, ht_sq] field_simp ring refine ⟨b, L, t, δ, hbposN, by omega, ?_, ht0, ?_, ?_, ?_, ?_, ?_, ?_, ?_⟩ · dsimp [b] omega · have hL1 : (1 : ℝ) ≤ L := by exact_mod_cast hLposN have hLb : (L : ℝ) ≤ b := by nlinarith [mul_nonneg hLpos.le (sub_nonneg.mpr hL1)] have hbK : (b : ℝ) ≤ k := by dsimp [b] exact_mod_cast (Nat.div_le_self k 2) have hnK : (k : ℝ) ^ 2 < n := by exact_mod_cast hn have hprod : (b : ℝ) * L ≤ (k : ℝ) ^ 2 := by calc (b : ℝ) * L ≤ b * b := by gcongr _ ≤ k * k := by gcongr _ = k ^ 2 := by ring have hs : t ^ 2 < 1 := by rw [ht_sq] apply (div_lt_iff₀ (by positivity : 0 < 200 * (n : ℝ))).2 nlinarith nlinarith · dsimp [δ, t] positivity · rw [ht_sq] have hh : 100 * ((n : ℝ) / (b : ℝ)) * ((b : ℝ) * (L : ℝ) / (200 * (n : ℝ))) = (L : ℝ) / 2 := by field_simp ring rw [hh] linarith · rw [ht_sq, hδ_sq] have hmult := mul_le_mul_of_nonneg_right hscale (mul_nonneg hbpos.le hnpos.le) have hh : 128 * ((b : ℝ) * (L : ℝ)) / (200 * (n : ℝ)) ≤ ((b : ℝ) * (b : ℝ)) / (500000 * (n : ℝ) * (L : ℝ)) := by apply (div_le_div_iff₀ (by positivity : 0 < 200 * (n : ℝ)) (by positivity : 0 < 500000 * (n : ℝ) * (L : ℝ))).2 nlinarith [hmult] convert hh using 1 <;> ring · have htail₁ := paired_degree_tail k hk4 have htail₂ := paired_poisson_tail n hn256 dsimp [b, L] linarith · dsimp [δ] ring_nf exact le_refl _ · rw [hδ_sq] have hbK : (k : ℝ) / 3 ≤ b := by dsimp [b] have hh : k ≤ 3 * (k / 2) := by omega have hhR : (k : ℝ) ≤ 3 * ((k / 2 : ℕ) : ℝ) := by exact_mod_cast hh linarith have hKL : (k : ℝ) * L ≤ 200 * (b : ℝ) * log (exp 1 * (k : ℝ)) := by have h1 := mul_le_mul_of_nonneg_left hLlog hkpos.le have h2 := mul_le_mul_of_nonneg_right hbK hlogpos.le nlinarith have hcross := mul_le_mul_of_nonneg_right hKL hnpos.le have he : (1 / 100000000 : ℝ) * ((k : ℝ) / ((n : ℝ) * log (exp 1 * (k : ℝ)))) = (k : ℝ) / (100000000 * (n : ℝ) * log (exp 1 * (k : ℝ))) := by ring rw [he] apply (div_le_div_iff₀ (by positivity : 0 < 100000000 * (n : ℝ) * log (exp 1 * (k : ℝ))) (by positivity : 0 < 500000 * (n : ℝ) * (L : ℝ))).2 nlinarith

theorem pairedCountVector_map_finitePoissonSampleLaw (b n : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (u : Fin b → ℝ) (hu : ∀ j, |u j| ≤ 1) : Measure.map (pairedCountVector b) (finitePoissonSampleLaw (simplexPMF (pairedTiltVector b hb t ht ht1 u hu)).toMeasure (toNNReal (2 * (n : ℝ)))) = Measure.pi (fun j : Fin b => scalarPoissonPairLaw ((n : ℝ) / (b : ℝ)) t (u j)) := by classical let P := (simplexPMF (pairedTiltVector b hb t ht ht1 u hu)).toMeasure let lam : ℝ≥0 := toNNReal (2 * (n : ℝ)) let H : FiniteSample (Fin (b * 2)) → (Fin (b * 2) → ℕ) := fun s => Causalean.Stat.FiniteRaoBlackwell.PairedPoissonHistogram.finiteSampleHistogram s.points let R : (Fin (b * 2) → ℕ) → (Fin b → ℕ × ℕ) := fun c j => (c (finProdFinEquiv (j, 0)), c (finProdFinEquiv (j, 1))) have hR : Measurable R := measurable_pi_lambda _ fun j => ((measurable_pi_apply _).prodMk (measurable_pi_apply _)) have hH : Measurable H := by apply measurable_to_countable' intro c rw [MeasurableSpace.measurableSet_iInf] intro m change MeasurableSet ((Sigma.mk m) ⁻¹' {s : FiniteSample (Fin (b * 2)) | H s = c}) exact (Set.to_countable _).measurableSet have hfactor : pairedCountVector b = R ∘ H := by funext s j unfold pairedCountVector R H change ((Finset.univ.filter fun i : Fin s.count => s.points i = finProdFinEquiv (j, 0)).card, (Finset.univ.filter fun i : Fin s.count => s.points i = finProdFinEquiv (j, 1)).card) = (Causalean.Stat.FiniteRaoBlackwell.PairedPoissonHistogram.finiteSampleHistogram s.points (finProdFinEquiv (j, 0)), Causalean.Stat.FiniteRaoBlackwell.PairedPoissonHistogram.finiteSampleHistogram s.points (finProdFinEquiv (j, 1))) simp [Causalean.Stat.FiniteRaoBlackwell.PairedPoissonHistogram.finiteSampleHistogram, Fintype.card_subtype] rw [hfactor, ← Measure.map_map hR hH] rw [finitePoissonSampleLaw_map_histogram P lam] unfold Causalean.Stat.FiniteRaoBlackwell.PairedPoissonHistogram.independentPoissonCountLaw apply Measure.ext_of_singleton intro z let c : Fin (b * 2) → ℕ := fun i => let p := finProdFinEquiv.symm i if p.2 = 0 then (z p.1).1 else (z p.1).2 have hc : R ⁻¹' {z} = {c} := by ext d simp only [Set.mem_preimage, Set.mem_singleton_iff] constructor · intro hd funext i have h := congrFun hd (finProdFinEquiv.symm i).1 change (d (finProdFinEquiv ((finProdFinEquiv.symm i).1, 0)), d (finProdFinEquiv ((finProdFinEquiv.symm i).1, 1))) = z (finProdFinEquiv.symm i).1 at h unfold c have hfin : (finProdFinEquiv.symm i).2 = 0 ∨ (finProdFinEquiv.symm i).2 = 1 := by omega rcases hfin with h0 | h1 · have hd0 := congrArg Prod.fst h rw [← h0, finProdFinEquiv.apply_symm_apply] at hd0 change i.modNat = 0 at h0 simpa [h0] using hd0 · have hd1 := congrArg snd h rw [← h1, finProdFinEquiv.apply_symm_apply] at hd1 change i.modNat = 1 at h1 simpa [h1] using hd1 · intro hd subst d funext j change (c (finProdFinEquiv (j, 0)), c (finProdFinEquiv (j, 1))) = z j simp only [c, Equiv.symm_apply_apply, ↓reduceIte, one_ne_zero] rw [Measure.map_apply hR (MeasurableSet.singleton z), hc, Measure.pi_singleton] simp only [scalarPoissonPairLaw, Measure.pi_singleton] have hrate0 (j : Fin b) : lam * (P {finProdFinEquiv (j, 0)}).toNNReal = toNNReal (((n : ℝ) / (b : ℝ)) * (1 + t * u j)) := by have hp : (P {finProdFinEquiv (j, 0)}).toNNReal = toNNReal ((1 + t * u j) / ((b : ℝ) * 2)) := by have he := finProdFinEquiv.symm_apply_apply (j, (0 : Fin 2)) have hdiv : (finProdFinEquiv (j, 0)).divNat = j := congrArg Prod.fst he have hmod : (finProdFinEquiv (j, 0)).modNat = 0 := congrArg snd he simp [P, simplexPMF, pairedTiltVector, ofReal, hdiv, hmod] rw [hp] change toNNReal (2 * (n : ℝ)) * toNNReal ((1 + t * u j) / ((b : ℝ) * 2)) = _ rw [← Real.toNNReal_mul (by positivity)] congr 1 have hb' : (b : ℝ) ≠ 0 := by exact_mod_cast (Nat.ne_of_gt hb) field_simp have hrate1 (j : Fin b) : lam * (P {finProdFinEquiv (j, 1)}).toNNReal = toNNReal (((n : ℝ) / (b : ℝ)) * (1 - t * u j)) := by have hp : (P {finProdFinEquiv (j, 1)}).toNNReal = toNNReal ((1 - t * u j) / ((b : ℝ) * 2)) := by have he := finProdFinEquiv.symm_apply_apply (j, (1 : Fin 2)) have hdiv : (finProdFinEquiv (j, 1)).divNat = j := congrArg Prod.fst he have hmod : (finProdFinEquiv (j, 1)).modNat = 1 := congrArg snd he simp [P, simplexPMF, pairedTiltVector, ofReal, hdiv, hmod, sub_eq_add_neg] rw [hp] change toNNReal (2 * (n : ℝ)) * toNNReal ((1 - t * u j) / ((b : ℝ) * 2)) = _ rw [← Real.toNNReal_mul (by positivity)] congr 1 have hb' : (b : ℝ) ≠ 0 := by exact_mod_cast (Nat.ne_of_gt hb) field_simp calc (∏ i : Fin (b * 2), poissonMeasure (lam * (P {i}).toNNReal) {c i}) = ∏ j : Fin b, poissonMeasure (lam * (P {finProdFinEquiv (j, 0)}).toNNReal) {(z j).1} * poissonMeasure (lam * (P {finProdFinEquiv (j, 1)}).toNNReal) {(z j).2} := by rw [← Equiv.prod_comp finProdFinEquiv] simp only [Fintype.prod_prod_type, Fin.prod_univ_two] apply Finset.prod_congr rfl intro j _ have hc0 : c (finProdFinEquiv (j, 0)) = (z j).1 := by have he := finProdFinEquiv.symm_apply_apply (j, (0 : Fin 2)) have hdiv : (finProdFinEquiv (j, 0)).divNat = j := congrArg Prod.fst he have hmod : (finProdFinEquiv (j, 0)).modNat = 0 := congrArg snd he simp [c, hdiv, hmod] have hc1 : c (finProdFinEquiv (j, 1)) = (z j).2 := by have he := finProdFinEquiv.symm_apply_apply (j, (1 : Fin 2)) have hdiv : (finProdFinEquiv (j, 1)).divNat = j := congrArg Prod.fst he have hmod : (finProdFinEquiv (j, 1)).modNat = 1 := congrArg snd he simp [c, hdiv, hmod] rw [hc0, hc1] _ = ∏ j : Fin b, (poissonMeasure (toNNReal (((n : ℝ) / (b : ℝ)) * (1 + t * u j)))).prod (poissonMeasure (toNNReal (((n : ℝ) / (b : ℝ)) * (1 - t * u j)))) {z j} := by apply Finset.prod_congr rfl intro j _ rw [hrate0 j, hrate1 j, ← Set.singleton_prod_singleton, Measure.prod_prod]

theorem pairedFinitePoissonPredictive_tv_le_counts {L : ℕ} (P : ScalarMomentPriors L) (b n : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) : tvDist (pairedFinitePoissonPredictive P b n hb t ht ht1 false) (pairedFinitePoissonPredictive P b n hb t ht ht1 true) ≤ tvDist (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false) (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true) := by classical let R : (Fin b → ℕ × ℕ) → (Fin (b * 2) → ℕ) := fun z i => if (finProdFinEquiv.symm i).2 = 0 then (z (finProdFinEquiv.symm i).1).1 else (z (finProdFinEquiv.symm i).1).2 let H : FiniteSample (Fin (b * 2)) → (Fin (b * 2) → ℕ) := fun s => finiteSampleHistogram s.points let K := histogramReconstructionKernel (Fin (b * 2)) have hR : Measurable R := measurable_of_countable _ have hH : Measurable H := by apply measurable_to_countable' intro c rw [MeasurableSpace.measurableSet_iInf] intro m change MeasurableSet ((Sigma.mk m) ⁻¹' {s : FiniteSample (Fin (b * 2)) | H s = c}) exact (Set.to_countable _).measurableSet have hC : Measurable (pairedCountVector b) := by apply measurable_to_countable' intro c rw [MeasurableSpace.measurableSet_iInf] intro m change MeasurableSet ((Sigma.mk m) ⁻¹' {s : FiniteSample (Fin (b * 2)) | pairedCountVector b s = c}) exact (Set.to_countable _).measurableSet have hfactor : H = R ∘ pairedCountVector b := by funext s i have he := finProdFinEquiv.apply_symm_apply i have hfin : (finProdFinEquiv.symm i).2 = 0 ∨ (finProdFinEquiv.symm i).2 = 1 := by omega rcases hfin with h0 | h1 · have hi : finProdFinEquiv ((finProdFinEquiv.symm i).1, 0) = i := by rw [← h0] exact he change i.modNat = 0 at h0 simp [H, R, pairedCountVector, finiteSampleHistogram, Fintype.card_subtype, h0] change finProdFinEquiv (i.divNat, 0) = i at hi rw [hi] · have hi : finProdFinEquiv ((finProdFinEquiv.symm i).1, 1) = i := by rw [← h1] exact he change i.modNat = 1 at h1 simp [H, R, pairedCountVector, finiteSampleHistogram, Fintype.card_subtype, h1] change finProdFinEquiv (i.divNat, 1) = i at hi rw [hi] have hcomponent (u : Fin b → Fin P.m) : K ∘ₘ Measure.map R (Measure.pi (fun j : Fin b => scalarPoissonPairLaw ((n : ℝ) / (b : ℝ)) t (P.node (u j)))) = finitePoissonSampleLaw (simplexPMF (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u))).toMeasure (toNNReal (2 * (n : ℝ))) := by let Q := (simplexPMF (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u))).toMeasure have hcount := pairedCountVector_map_finitePoissonSampleLaw b n hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u) change Measure.map (pairedCountVector b) (finitePoissonSampleLaw Q _) = _ at hcount change K ∘ₘ Measure.map R (Measure.pi (fun j : Fin b => scalarPoissonPairLaw ((n : ℝ) / (b : ℝ)) t (pairedNodeVector P u j))) = _ rw [← hcount, Measure.map_map hR hC] rw [← hfactor] rw [finitePoissonSampleLaw_map_histogram Q] exact independentPoissonCountLaw_comp_reconstruction Q _ have hmixture (side : Bool) : K ∘ₘ Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t side) = pairedFinitePoissonPredictive P b n hb t ht ht1 side := by unfold pairedPoissonPredictive pairedFinitePoissonPredictive mixture rw [Measure.map_finset_sum' hR.aemeasurable] simp_rw [Measure.map_smul] rw [← Measure.sum_fintype] rw [Measure.bind_sum _ _ K.aemeasurable] rw [Measure.sum_fintype] simp_rw [Measure.bind_smul] apply Finset.sum_congr rfl intro u _ rw [hcomponent u] have hprob (side : Bool) : IsProbabilityMeasure (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t side) := by haveI (u : Fin b → Fin P.m) : IsProbabilityMeasure (Measure.pi (fun j : Fin b => scalarPoissonPairLaw ((n : ℝ) / (b : ℝ)) t (P.node (u j)))) := by unfold scalarPoissonPairLaw infer_instance exact mixture_isProbabilityMeasure (pairedProductWeight P b side) (pairedProductWeight_sum P b side) _ letI : IsProbabilityMeasure (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false) := hprob false letI : IsProbabilityMeasure (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true) := hprob true letI : IsProbabilityMeasure (Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false)) := Measure.isProbabilityMeasure_map hR.aemeasurable letI : IsProbabilityMeasure (Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true)) := Measure.isProbabilityMeasure_map hR.aemeasurable let D := Kernel.deterministic R hR calc tvDist (pairedFinitePoissonPredictive P b n hb t ht ht1 false) (pairedFinitePoissonPredictive P b n hb t ht ht1 true) = tvDist (K ∘ₘ Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false)) (K ∘ₘ Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true)) := by rw [hmixture false, hmixture true] _ ≤ tvDist (Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false)) (Measure.map R (pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true)) := tvDist_bind_le _ _ K _ = tvDist (D ∘ₘ pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t false) (D ∘ₘ pairedPoissonPredictive P b ((n : ℝ) / (b : ℝ)) t true) := by rw [Measure.deterministic_comp_eq_map hR, Measure.deterministic_comp_eq_map hR] _ ≤ _ := tvDist_bind_le _ _ D

theorem pairedProductTarget_variance_le {L : ℕ} (P : ScalarMomentPriors L) (b : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (side : Bool) : ∑ u : Fin b → Fin P.m, (pairedProductWeight P b side u).toReal * (simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u)) - t * ∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i|) ^ 2 ≤ t ^ 2 / (b : ℝ) := by classical let w : Fin P.m → ℝ := scalarPriorWeight P side let μ : ℝ := ∑ i : Fin P.m, w i * |P.node i| let c : Fin P.m → ℝ := fun i => |P.node i| - μ let W : (Fin b → Fin P.m) → ℝ := fun u => ∏ j : Fin b, w (u j) have hw (i : Fin P.m) : 0 ≤ w i := by dsimp [w] cases side <;> simp [scalarPriorWeight, P.w₀_nonneg, P.w₁_nonneg] have hwsum : ∑ i : Fin P.m, w i = 1 := by dsimp [w] cases side <;> simp [scalarPriorWeight, P.w₀_sum, P.w₁_sum] have hcenter : ∑ i : Fin P.m, w i * c i = 0 := by simp only [c, mul_sub, Finset.sum_sub_distrib, ← Finset.sum_mul] simp [μ, hwsum] have hW (u : Fin b → Fin P.m) : 0 ≤ W u := by exact Finset.prod_nonneg (by intro i hi; exact hw _) have hweight (u : Fin b → Fin P.m) : (pairedProductWeight P b side u).toReal = W u := by simp [pairedProductWeight, ENNReal.toReal_prod, ENNReal.toReal_ofReal (hw _), W, w] have hpair (j k : Fin b) (hjk : j ≠ k) : ∑ u : Fin b → Fin P.m, W u * c (u j) * c (u k) = 0 := by let f : Fin b → Fin P.m → ℝ := fun l i => w i * (if l = j then c i else 1) * (if l = k then c i else 1) have hf (u : Fin b → Fin P.m) : W u * c (u j) * c (u k) = ∏ l : Fin b, f l (u l) := by simp only [f, Finset.prod_mul_distrib] simp [W, Finset.prod_ite_eq', mul_assoc] calc _ = ∑ u : Fin b → Fin P.m, ∏ l : Fin b, f l (u l) := by apply Finset.sum_congr rfl intro u _ exact hf u _ = ∏ l : Fin b, ∑ i : Fin P.m, f l i := by exact (Fintype.prod_sum f).symm _ = 0 := by apply Finset.prod_eq_zero (Finset.mem_univ j) simp [f, hjk, hcenter] have hdiag (j : Fin b) : ∑ u : Fin b → Fin P.m, W u * c (u j) ^ 2 ≤ 1 := by have hμlo : 0 ≤ μ := Finset.sum_nonneg (by intro i hi exact mul_nonneg (hw i) (abs_nonneg _)) have hμhi : μ ≤ 1 := by calc μ ≤ ∑ i : Fin P.m, w i := by apply Finset.sum_le_sum intro i hi exact mul_le_of_le_one_right (hw i) (pairedNodeVector_abs_le_one P (fun _ => i) j) _ = 1 := hwsum have hc (i : Fin P.m) : c i ^ 2 ≤ 1 := by have hai : |P.node i| ≤ 1 := by exact abs_le.mpr (P.node_mem i) have hlow : 0 ≤ |P.node i| := abs_nonneg _ nlinarith [sq_nonneg (c i - 1), sq_nonneg (c i + 1)] calc _ ≤ ∑ u : Fin b → Fin P.m, W u := by apply Finset.sum_le_sum intro u hu simpa [mul_comm] using mul_le_of_le_one_left (hW u) (hc (u j)) _ = 1 := by calc _ = ∏ _j : Fin b, ∑ i : Fin P.m, w i := by simpa [W] using (Fintype.prod_sum (fun _j : Fin b => w)).symm _ = 1 := by simp [hwsum] have hb' : (b : ℝ) ≠ 0 := by exact_mod_cast (Nat.ne_of_gt hb) have hscale (u : Fin b → Fin P.m) : simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u)) - t * ∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i| = (t / (b : ℝ)) * ∑ j : Fin b, c (u j) := by rw [pairedBase_tilt_l1] simp only [pairedNodeVector] have hs : ∑ j : Fin b, c (u j) = (∑ j : Fin b, |P.node (u j)|) - (b : ℝ) * μ := by simp [c, Finset.sum_sub_distrib] rw [hs] dsimp [μ, w] field_simp have hvar : ∑ u : Fin b → Fin P.m, W u * (∑ j : Fin b, c (u j)) ^ 2 ≤ (b : ℝ) := by have hsplit : ∑ u : Fin b → Fin P.m, W u * (∑ j : Fin b, c (u j)) ^ 2 = ∑ j : Fin b, ∑ u : Fin b → Fin P.m, W u * c (u j) ^ 2 := by calc _ = ∑ j : Fin b, ∑ k : Fin b, ∑ u : Fin b → Fin P.m, W u * c (u j) * c (u k) := by simp_rw [pow_two, Finset.sum_mul_sum, Finset.mul_sum] rw [Finset.sum_comm] apply Finset.sum_congr rfl intro j _ rw [Finset.sum_comm] apply Finset.sum_congr rfl intro k _ apply Finset.sum_congr rfl intro u _ ring _ = _ := by apply Finset.sum_congr rfl intro j _ calc _ = ∑ k : Fin b, if k = j then ∑ u : Fin b → Fin P.m, W u * c (u j) ^ 2 else 0 := by apply Finset.sum_congr rfl intro k _ split_ifs with h · subst k congr 1 ext u ring · exact hpair j k (Ne.symm h) _ = _ := by simp rw [hsplit] calc _ ≤ ∑ _j : Fin b, (1 : ℝ) := Finset.sum_le_sum (by intro j hj exact hdiag j) _ = b := by simp simp_rw [hweight, hscale] calc _ = (t / (b : ℝ)) ^ 2 * ∑ u : Fin b → Fin P.m, W u * (∑ j : Fin b, c (u j)) ^ 2 := by rw [Finset.mul_sum] apply Finset.sum_congr rfl intro u _ ring _ ≤ (t / (b : ℝ)) ^ 2 * (b : ℝ) := by exact mul_le_mul_of_nonneg_left hvar (sq_nonneg _) _ = t ^ 2 / (b : ℝ) := by field_simp

theorem pairedProductTarget_bad_mass_le {L : ℕ} (P : ScalarMomentPriors L) (b : ℕ) (hb : 0 < b) (t : ℝ) (ht : 0 ≤ t) (ht1 : t ≤ 1) (δ : ℝ) (hδ : 0 < δ) (hscale : 128 * t ^ 2 ≤ (b : ℝ) * δ ^ 2) (side : Bool) : (∑ u : Fin b → Fin P.m, if δ / 4 < |simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u)) - t * (∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i|)| then pairedProductWeight P b side u else 0) ≤ (1 / 8 : ℝ≥0∞) := by classical let X : (Fin b → Fin P.m) → ℝ := fun u => simplexL1 (pairedBaseVector b hb) (pairedTiltVector b hb t ht ht1 (pairedNodeVector P u) (pairedNodeVector_abs_le_one P u)) - t * (∑ i : Fin P.m, scalarPriorWeight P side i * |P.node i|) let w : (Fin b → Fin P.m) → ℝ≥0∞ := pairedProductWeight P b side have hwfinite (u : Fin b → Fin P.m) : w u ≠ ⊤ := by have hu : w u ≤ ∑ v : Fin b → Fin P.m, w v := Finset.single_le_sum (fun _ _ => bot_le) (Finset.mem_univ u) rw [show (∑ v : Fin b → Fin P.m, w v) = 1 from pairedProductWeight_sum P b side] at hu exact ne_top_of_le_ne_top (by norm_num : (1 : ℝ≥0∞) ≠ ⊤) hu have hδsq : 0 < δ ^ 2 := sq_pos_of_pos hδ have hpoint (u : Fin b → Fin P.m) : (if δ / 4 < |X u| then (w u).toReal else 0) * δ ^ 2 ≤ 16 * ((w u).toReal * X u ^ 2) := by split_ifs with hbad · have hsquare : δ ^ 2 < 16 * X u ^ 2 := by nlinarith [sq_nonneg (|X u| - δ / 4), sq_abs (X u)] exact (mul_le_mul_of_nonneg_left hsquare.le (show 0 ≤ (w u).toReal from ENNReal.toReal_nonneg)).trans_eq (by ring) · nlinarith [show 0 ≤ (w u).toReal from ENNReal.toReal_nonneg, sq_nonneg (X u)] have hsum : (∑ u : Fin b → Fin P.m, if δ / 4 < |X u| then (w u).toReal else 0) * δ ^ 2 ≤ 16 * (t ^ 2 / (b : ℝ)) := by calc _ = ∑ u : Fin b → Fin P.m, (if δ / 4 < |X u| then (w u).toReal else 0) * δ ^ 2 := by rw [Finset.sum_mul] _ ≤ ∑ u : Fin b → Fin P.m, 16 * ((w u).toReal * X u ^ 2) := Finset.sum_le_sum (by intro u _; exact hpoint u) _ = 16 * ∑ u : Fin b → Fin P.m, (w u).toReal * X u ^ 2 := by rw [Finset.mul_sum] _ ≤ _ := by exact mul_le_mul_of_nonneg_left (pairedProductTarget_variance_le P b hb t ht ht1 side) (by norm_num) have hbpos : (0 : ℝ) < b := by exact_mod_cast hb have hratio : t ^ 2 / (b : ℝ) ≤ δ ^ 2 / 128 := by apply (div_le_iff₀ hbpos).2 nlinarith [hscale] have hreal : (∑ u : Fin b → Fin P.m, if δ / 4 < |X u| then (w u).toReal else 0) ≤ (1 / 8 : ℝ) := by have h := hsum.trans (mul_le_mul_of_nonneg_left hratio (by norm_num)) nlinarith have hmass : (∑ u : Fin b → Fin P.m, if δ / 4 < |X u| then w u else 0) = ofReal (∑ u : Fin b → Fin P.m, if δ / 4 < |X u| then (w u).toReal else 0) := by rw [ENNReal.ofReal_sum_of_nonneg (by intro u _ split_ifs <;> simp [ENNReal.toReal_nonneg])] apply Finset.sum_congr rfl intro u _ split_ifs <;> simp [ENNReal.ofReal_toReal (hwfinite u)] change (∑ u : Fin b → Fin P.m, if δ / 4 < |X u| then w u else 0) ≤ 1 / 8 rw [hmass] exact (ENNReal.ofReal_le_ofReal hreal).trans_eq (by rw [ENNReal.ofReal_div_of_pos (by norm_num : (0 : ℝ) < 8)] norm_num)

theorem scalarPoissonPredictive_tv_le {L : ℕ} (P : ScalarMomentPriors L) (lambda t : ℝ) (hlambda : 0 ≤ lambda) (ht : 0 ≤ t) (ht1 : t ≤ 1) (hscale : 100 * lambda * t ^ 2 ≤ (L : ℝ)) : tvDist (scalarPoissonPredictive P lambda t false) (scalarPoissonPredictive P lambda t true) ≤ (2 : ℝ) ^ (-(L : ℝ) / 4) := by classical by_cases hzero : lambda = 0 · have hlaw (u : ℝ) : scalarPoissonPairLaw lambda t u = scalarPoissonPairLaw lambda t 0 := by simp [scalarPoissonPairLaw, hzero] have hpredict : scalarPoissonPredictive P lambda t false = scalarPoissonPredictive P lambda t true := by ext A hA simp only [scalarPoissonPredictive, mixture_apply, hlaw] rw [← Finset.sum_mul, ← Finset.sum_mul] have hw0 : ∑ i : Fin P.m, ofReal (P.w₀ i) = 1 := by rw [← ENNReal.ofReal_sum_of_nonneg (fun i _ => P.w₀_nonneg i), P.w₀_sum] norm_num have hw1 : ∑ i : Fin P.m, ofReal (P.w₁ i) = 1 := by rw [← ENNReal.ofReal_sum_of_nonneg (fun i _ => P.w₁_nonneg i), P.w₁_sum] norm_num simpa [hw0, hw1] rw [hpredict] simpa [tvDist] using (Real.rpow_nonneg (by norm_num : (0 : ℝ) ≤ 2) (-(L : ℝ) / 4)) have hlambda' : 0 < lambda := lt_of_le_of_ne hlambda (Ne.symm hzero) letI : IsSFiniteKernel (aggregatePoissonKernelOfRealRate (fun u => lambda * (1 + t * u))) := by unfold aggregatePoissonKernelOfRealRate apply Kernel.IsSFiniteKernel.withDensity simp letI : IsSFiniteKernel (aggregatePoissonKernelOfRealRate (fun u => lambda * (1 - t * u))) := by unfold aggregatePoissonKernelOfRealRate apply Kernel.IsSFiniteKernel.withDensity simp let K : Kernel ℝ (ℕ × ℕ) := (aggregatePoissonKernelOfRealRate (fun u => lambda * (1 + t * u))).prod (aggregatePoissonKernelOfRealRate (fun u => lambda * (1 - t * u))) have hkernel (r : ℝ → ℝ) (hr : Measurable r) (u : ℝ) : aggregatePoissonKernelOfRealRate r u = poissonMeasure (toNNReal (r u)) := by rw [aggregatePoissonKernelOfRealRate, Kernel.withDensity_apply] · apply Measure.ext_of_singleton intro n rw [withDensity_apply _ (measurableSet_singleton n), poissonMeasure_singleton] simp · apply measurable_from_prod_countable_left intro n fun_prop have hKfiber (u : ℝ) : K u = scalarPoissonPairLaw lambda t u := by dsimp [K, scalarPoissonPairLaw] rw [Kernel.prod_apply, hkernel _ (by fun_prop) u, hkernel _ (by fun_prop) u] let Q : Measure (ℕ × ℕ) := scalarPoissonPairLaw lambda t 0 let π₀ : Measure ℝ := scalarMomentPrior P false let π₁ : Measure ℝ := scalarMomentPrior P true letI : IsProbabilityMeasure π₀ := scalarMomentPrior_isProbability P false letI : IsProbabilityMeasure π₁ := scalarMomentPrior_isProbability P true have hQ : IsProbabilityMeasure Q := by dsimp [Q, scalarPoissonPairLaw] infer_instance letI : IsProbabilityMeasure Q := hQ have hK : ∀ u, IsProbabilityMeasure (K u) := by intro u rw [hKfiber] unfold scalarPoissonPairLaw infer_instance have hmeas : Measurable (fun p : ℝ × (ℕ × ℕ) => scalarPoissonPairLikelihood t p.1 p.2) := by unfold scalarPoissonPairLikelihood fun_prop have hnonneg : ∀ u : ℝ, |u| ≤ 1 → ∀ z, 0 ≤ scalarPoissonPairLikelihood t u z := by intro u hu z have htu : -1 ≤ t * u ∧ t * u ≤ 1 := by constructor <;> nlinarith [(abs_le.mp hu).1, (abs_le.mp hu).2] unfold scalarPoissonPairLikelihood exact mul_nonneg (pow_nonneg (by linarith) _) (pow_nonneg (by linarith) _) have hdensity : ∀ u : ℝ, |u| ≤ 1 → K u = Q.withDensity (fun z => ofReal (scalarPoissonPairLikelihood t u z)) := by intro u hu have htu : -1 ≤ t * u ∧ t * u ≤ 1 := by constructor <;> nlinarith [(abs_le.mp hu).1, (abs_le.mp hu).2] rw [hKfiber] exact scalarPoissonPairLaw_eq_withDensity lambda t u hlambda' (by linarith) (by linarith) have hinner : ∀ u : ℝ, |u| ≤ 1 → ∀ v : ℝ, |v| ≤ 1 → (∫ z, scalarPoissonPairLikelihood t u z * scalarPoissonPairLikelihood t v z ∂Q) = exp ((2 * lambda * t ^ 2) * u * v) := by intro u _ v _ exact scalarPoissonPairLikelihood_inner lambda t u v hlambda have hsupp0 : π₀ {u : ℝ | |u| ≤ 1} = 1 := scalarMomentPrior_supported P false have hsupp1 : π₁ {u : ℝ | |u| ≤ 1} = 1 := scalarMomentPrior_supported P true have hmom : ∀ n ≤ L, (∫ u, u ^ n ∂π₀) = ∫ u, u ^ n ∂π₁ := by intro n hn exact scalarMomentPrior_moments_eq P n hn have htv := momentMatchedMixture_tv_le_sqrt_tail_of_supported π₀ π₁ K Q (scalarPoissonPairLikelihood t) hK (2 * lambda * t ^ 2) 1 L (by positivity) (by norm_num) hmeas hnonneg hdensity hinner hsupp0 hsupp1 hmom have hpred (side : Bool) : priorPredictive (scalarMomentPrior P side) K = scalarPoissonPredictive P lambda t side := by ext A hA rw [priorPredictive_apply _ K hA] simp only [scalarMomentPrior, lintegral_finsetSum_measure, lintegral_smul_measure, lintegral_dirac, scalarPoissonPredictive, mixture_apply, smul_eq_mul] simp only [hKfiber] rw [hpred false, hpred true] at htv calc _ ≤ sqrt (exponentialSeriesTail L ((2 * lambda * t ^ 2) * 1 ^ 2)) := htv _ ≤ (2 : ℝ) ^ (-(L : ℝ) / 4) := by let z : ℝ := 2 * lambda * t ^ 2 have hz : 0 ≤ z := by dsimp [z]; positivity have hzL : z ≤ (L : ℝ) / 50 := by dsimp [z] linarith have hterm (n : ℕ) (hLn : L < n) : z ^ n / (n.factorial : ℝ) ≤ (1 / 4 : ℝ) ^ n := by have hn : 0 < n := lt_of_le_of_lt (zero_le L) hLn have hnR : 0 < (n : ℝ) := by exact_mod_cast hn have he : 0 < exp 1 := Real.exp_pos 1 have hsqrt : 1 ≤ sqrt (2 * pi * (n : ℝ)) := by rw [Real.one_le_sqrt] have hn1 : (1 : ℝ) ≤ n := by exact_mod_cast hn nlinarith [Real.pi_gt_three] have hbase : 0 ≤ ((n : ℝ) / exp 1) ^ n := by positivity have hfac : ((n : ℝ) / exp 1) ^ n ≤ (n.factorial : ℝ) := by calc _ ≤ sqrt (2 * pi * (n : ℝ)) * ((n : ℝ) / exp 1) ^ n := by nlinarith _ ≤ (n.factorial : ℝ) := Stirling.le_factorial_stirling n have hbasepos : 0 < ((n : ℝ) / exp 1) ^ n := by positivity have hratio_nonneg : 0 ≤ z * exp 1 / (n : ℝ) := by positivity have hratio : z * exp 1 / (n : ℝ) ≤ (1 / 4 : ℝ) := by have hLnR : (L : ℝ) ≤ n := by exact_mod_cast hLn.le have he3 : z * exp 1 ≤ 3 * z := by simpa [mul_comm] using (mul_le_mul_of_nonneg_left Real.exp_one_lt_three.le hz) have haux : z * exp 1 ≤ (1 / 4 : ℝ) * n := by linarith exact (div_le_iff₀ hnR).2 haux calc z ^ n / (n.factorial : ℝ) ≤ z ^ n / (((n : ℝ) / exp 1) ^ n) := div_le_div_of_nonneg_left (pow_nonneg hz n) hbasepos hfac _ = (z * exp 1 / (n : ℝ)) ^ n := by rw [← div_pow] congr 1 field_simp _ ≤ (1 / 4 : ℝ) ^ n := pow_le_pow_left₀ hratio_nonneg hratio n let f : ℕ → ℝ := fun n => if L < n then z ^ n / (n.factorial : ℝ) else 0 let g : ℕ → ℝ := fun n => if L < n then (1 / 4 : ℝ) ^ n else 0 have hq : Summable (fun n : ℕ => (1 / 4 : ℝ) ^ n) := summable_geometric_of_norm_lt_one (by norm_num) have hg : Summable g := by apply hq.of_nonneg_of_le · intro n dsimp [g] positivity · intro n dsimp [g] split_ifs · exact le_refl _ · positivity have hf_nonneg : ∀ n, 0 ≤ f n := by intro n dsimp [f] positivity have hfg : ∀ n, f n ≤ g n := by intro n dsimp [f, g] split_ifs with hn · exact hterm n hn · rfl have hf : Summable f := hg.of_nonneg_of_le hf_nonneg hfg have hgeo : (∑' n : ℕ, g n) ≤ (1 / 4 : ℝ) ^ L := by have hzero_sum : ((range (L + 1)).sum (fun n : ℕ => g n)) = 0 := by apply Finset.sum_eq_zero intro n hn have hnle : n ≤ L := le_of_lt_succ (Finset.mem_range.mp hn) simp [g, Nat.not_lt_of_ge hnle] rw [← hg.sum_add_tsum_nat_add (L + 1), hzero_sum, zero_add] simp only [show ∀ n : ℕ, L < n + (L + 1) by omega, if_true, g, pow_add] rw [tsum_mul_right, tsum_geometric_of_norm_lt_one (by norm_num)] norm_num [pow_succ] have hp : 0 ≤ (1 / 4 : ℝ) ^ L := pow_nonneg (by norm_num) L nlinarith have htail : exponentialSeriesTail L z ≤ (1 / 4 : ℝ) ^ L := by unfold exponentialSeriesTail change (∑' n, f n) ≤ _ exact (hf.tsum_le_tsum hfg hg).trans hgeo have hsqrt : sqrt (exponentialSeriesTail L z) ≤ (1 / 2 : ℝ) ^ L := by rw [Real.sqrt_le_iff] refine ⟨by positivity, htail.trans_eq ?_⟩ rw [show (1 / 4 : ℝ) = (1 / 2) ^ 2 by norm_num, ← pow_mul, ← pow_mul] congr 1 omega have hpow : (1 / 2 : ℝ) ^ L = (2 : ℝ) ^ (-(L : ℝ)) := by rw [show (1 / 2 : ℝ) = (2 : ℝ) ^ (-1 : ℝ) by norm_num, ← Real.rpow_natCast, ← Real.rpow_mul (by norm_num : (0 : ℝ) ≤ 2)] congr 1 ring simpa only [one_pow, mul_one, z] using hsqrt.trans (hpow ▸ Real.rpow_le_rpow_of_exponent_le (by norm_num : (1 : ℝ) ≤ 2) (by have := Nat.cast_nonneg (α := ℝ) L; linarith))

theorem scalarMomentPrior_moments_eq {L : ℕ} (P : ScalarMomentPriors L) (j : ℕ) (hj : j ≤ L) : (∫ u : ℝ, u ^ j ∂scalarMomentPrior P false) = ∫ u : ℝ, u ^ j ∂scalarMomentPrior P true := by classical have hIntegral (side : Bool) : (∫ u : ℝ, u ^ j ∂scalarMomentPrior P side) = ∑ i : Fin P.m, (if side then P.w₁ i else P.w₀ i) * P.node i ^ j := by cases side with | false => rw [scalarMomentPrior, integral_finsetSum_measure] · simp [integral_smul_measure, integral_dirac, P.w₀_nonneg] · intro i hi exact (integrable_dirac (by simp)).smul_measure ENNReal.ofReal_ne_top | true => rw [scalarMomentPrior, integral_finsetSum_measure] · simp [integral_smul_measure, integral_dirac, P.w₁_nonneg] · intro i hi exact (integrable_dirac (by simp)).smul_measure ENNReal.ofReal_ne_top rw [hIntegral false, hIntegral true] exact P.moments_eq j hj