From b9c46b71f61fa77fb2416937045b53b329f64748 Mon Sep 17 00:00:00 2001 From: Remy Degenne Date: Sun, 15 Nov 2026 15:38:04 +0100 Subject: [PATCH 1/2] feat: the numeraire of two exponential or geometric observations against a common parameter MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit For a statistic whose conditional law given it is the same kernel under every member of the null, the likelihood ratio of the alternative to the law with its law of the statistic and that conditional law is an e-variable, and the numeraire when that law is in the effective null. For two exponential observations and the null of the products Exp(l) ⊗ Exp(l), this law (the sum under the alternative, then the first observation uniform on [0, s]) is the reverse information projection, the mixture of the members of the null whose mean is uniform between the two means, and no finite combination of members of the null; the same construction gives the numeraire for two geometric observations. EValues roadmap, Layer 1, item 5. Co-Authored-By: Claude Opus 5.5 Claude-Session: https://claude.ai/code/session_0196ftMVbgjSJzPQEwdek3kM --- LMLExplore/EValues/Numeraire/Conditional.lean | 133 ++++ LMLExplore/EValues/TwoSample/Exponential.lean | 580 ++++++++++++++++++ LMLExplore/EValues/TwoSample/Geometric.lean | 392 ++++++++++++ 3 files changed, 1105 insertions(+) create mode 100644 LMLExplore/EValues/Numeraire/Conditional.lean create mode 100644 LMLExplore/EValues/TwoSample/Exponential.lean create mode 100644 LMLExplore/EValues/TwoSample/Geometric.lean diff --git a/LMLExplore/EValues/Numeraire/Conditional.lean b/LMLExplore/EValues/Numeraire/Conditional.lean new file mode 100644 index 00000000..77819749 --- /dev/null +++ b/LMLExplore/EValues/Numeraire/Conditional.lean @@ -0,0 +1,133 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import LMLExplore.EValues.Numeraire + +/-! +# Numeraires of nulls with a common conditional law given a statistic + +Let `T : Ω → 𝓨` be a statistic and `κ : Kernel 𝓨 Ω` a kernel such that `κ y` is carried by the +fiber `T ⁻¹' {y}`. Suppose that `κ` is the conditional law given `T` of every member `P` of a +null `S`, in the sense that `P = κ ∘ₘ P.map T`. For an alternative `Q`, the measure +`κ ∘ₘ Q.map T` has the law of `T` under `Q` and the conditional law `κ` of the null given `T`. + +The likelihood ratio of `Q` to `κ ∘ₘ Q.map T` integrates to one against `κ y` for `Q.map T`-almost +every `y`: its conditional expectation given `T` is one. So it is an e-variable for `S`, with +expectation exactly `1` under every member `P` of `S` such that the law of `T` under `P` is +absolutely continuous with respect to its law under `Q`. This is the conditional e-variable of +[HGLLA24]. When moreover `κ ∘ₘ Q.map T` belongs to the effective null of `S` (for example when it +is a mixture of members of `S`, by `ProbabilityTheory.comp_mem_effectiveSet`), this likelihood +ratio is the numeraire of `S` under `Q` by the criterion of [LRR25, Theorem 4.1] +(`ProbabilityTheory.isNumeraire_of_ae_eq_rnDeriv`), and `κ ∘ₘ Q.map T` is the reverse information +projection of `Q` on `S` when it is equivalent to `Q` +(`ProbabilityTheory.ripr_eq_of_isEVar_rnDeriv`). + +## Main statements + +* `ProbabilityTheory.map_comp_eq_of_ae_eq`: `κ ∘ₘ ν` has law `ν` for `T`. +* `ProbabilityTheory.ae_lintegral_rnDeriv_comp_map_eq_one`: the likelihood ratio of `Q` to + `κ ∘ₘ Q.map T` has conditional expectation one given `T`. +* `ProbabilityTheory.lintegral_rnDeriv_comp_map`, `ProbabilityTheory.isEVar_rnDeriv_comp_map`: it + is an exact e-variable for every null of which `κ` is the conditional law given `T`. + +## References + +* [HGLLA24] Y. Hao, P. Grünwald, T. Lardy, L. Long and R. Adams, *E-values for k-sample tests + with exponential families*, arXiv:2303.00471. +* [LRR25] M. Larsson, A. Ramdas and J. Ruf, *The numeraire e-variable and reverse information + projection*, Ann. Statist. 53 (2025), arXiv:2402.18810v4, Theorem 4.1. +-/ + +public section + +open MeasureTheory Set + +open scoped ENNReal + +namespace ProbabilityTheory + +variable {Ω 𝓨 : Type*} {mΩ : MeasurableSpace Ω} {m𝓨 : MeasurableSpace 𝓨} {T : Ω → 𝓨} + {κ : Kernel 𝓨 Ω} {S : Set (Measure Ω)} {Q P : Measure Ω} + +/-- If `κ y` is a probability measure carried by the fiber `T ⁻¹' {y}` for `ν`-almost every `y`, +the law of `T` under `κ ∘ₘ ν` is `ν`. -/ +theorem map_comp_eq_of_ae_eq [IsMarkovKernel κ] (hT : Measurable T) {ν : Measure 𝓨} + (hκ : ∀ᵐ y ∂ν, ∀ᵐ x ∂κ y, T x = y) : (κ ∘ₘ ν).map T = ν := by + ext A hA + rw [Measure.map_apply hT hA, Measure.bind_apply (hT hA) κ.aemeasurable, + ← lintegral_indicator_one hA] + refine lintegral_congr_ae ?_ + filter_upwards [hκ] with y hκ + by_cases hy : y ∈ A + · rw [indicator_of_mem hy, Pi.one_apply, ← measure_univ (μ := κ y)] + refine measure_congr ?_ + filter_upwards [hκ] with x hx + simp [hx, hy] + · rw [indicator_of_notMem hy] + refine measure_eq_zero_iff_ae_notMem.mpr ?_ + filter_upwards [hκ] with x hx + simpa [hx] using hy + +/-- **The conditional expectation of the likelihood ratio given the statistic is one.** If `κ y` +is carried by the fiber `T ⁻¹' {y}` for `Q.map T`-almost every `y` and the finite measure `Q` is +absolutely continuous with respect to `κ ∘ₘ Q.map T`, the likelihood ratio of `Q` to +`κ ∘ₘ Q.map T` has integral one against `κ y` for `Q.map T`-almost every `y`. -/ +theorem ae_lintegral_rnDeriv_comp_map_eq_one [IsFiniteMeasure Q] [IsFiniteKernel κ] + (hT : Measurable T) (hκ : ∀ᵐ y ∂Q.map T, ∀ᵐ x ∂κ y, T x = y) (hQ : Q ≪ κ ∘ₘ Q.map T) : + ∀ᵐ y ∂Q.map T, ∫⁻ x, Q.rnDeriv (κ ∘ₘ Q.map T) x ∂κ y = 1 := by + set ν := Q.map T + set X := Q.rnDeriv (κ ∘ₘ ν) + have hX : Measurable X := Measure.measurable_rnDeriv _ _ + have hφ : Measurable fun y ↦ ∫⁻ x, X x ∂κ y := hX.lintegral_kernel + refine (withDensity_eq_iff_of_sigmaFinite (g := 1) hφ.aemeasurable aemeasurable_one).mp ?_ + rw [withDensity_one] + ext A hA + rw [withDensity_apply _ hA, Measure.map_apply hT hA, ← lintegral_indicator hA] + calc ∫⁻ y, A.indicator (fun y ↦ ∫⁻ x, X x ∂κ y) y ∂ν + _ = ∫⁻ y, ∫⁻ x, (T ⁻¹' A).indicator X x ∂κ y ∂ν := by + refine lintegral_congr_ae ?_ + filter_upwards [hκ] with y hκ + by_cases hy : y ∈ A + · rw [indicator_of_mem hy] + refine lintegral_congr_ae ?_ + filter_upwards [hκ] with x hx + rw [indicator_of_mem (show x ∈ T ⁻¹' A by simpa [hx])] + · rw [indicator_of_notMem hy, lintegral_congr_ae (g := 0), Pi.zero_def, lintegral_zero] + filter_upwards [hκ] with x hx + rw [indicator_of_notMem (show x ∉ T ⁻¹' A by simpa [hx]), Pi.zero_apply] + _ = ∫⁻ x, (T ⁻¹' A).indicator X x ∂(κ ∘ₘ ν) := + (Measure.lintegral_bind κ.aemeasurable (hX.indicator (hT hA)).aemeasurable).symm + _ = Q (T ⁻¹' A) := by + rw [lintegral_indicator (hT hA), ← withDensity_apply _ (hT hA), + Measure.withDensity_rnDeriv_eq _ _ hQ] + +/-- **The likelihood ratio to `κ ∘ₘ Q.map T` has expectation one under the null.** Under a +measure `P` whose conditional law given `T` is `κ` and such that the law of `T` under `P` is +absolutely continuous with respect to its law under `Q`, the likelihood ratio of `Q` to +`κ ∘ₘ Q.map T` has expectation the mass of `P`. -/ +theorem lintegral_rnDeriv_comp_map [IsFiniteMeasure Q] [IsFiniteKernel κ] (hT : Measurable T) + (hκ : ∀ᵐ y ∂Q.map T, ∀ᵐ x ∂κ y, T x = y) (hQ : Q ≪ κ ∘ₘ Q.map T) (hP : κ ∘ₘ P.map T = P) + (hPQ : P.map T ≪ Q.map T) : + ∫⁻ x, Q.rnDeriv (κ ∘ₘ Q.map T) x ∂P = P univ := by + conv_lhs => rw [← hP] + rw [Measure.lintegral_bind κ.aemeasurable (Measure.measurable_rnDeriv _ _).aemeasurable, + lintegral_congr_ae (hPQ.ae_le (ae_lintegral_rnDeriv_comp_map_eq_one hT hκ hQ)), + lintegral_one, Measure.map_apply hT .univ, preimage_univ] + +/-- **The conditional e-variable** [HGLLA24]. If `κ` is the conditional law given `T` of every +member of a null `S` of finite measures, under each of which the law of `T` is absolutely +continuous with respect to its law under `Q`, the likelihood ratio of `Q` to `κ ∘ₘ Q.map T` is an +e-variable for `S`, exact by `lintegral_rnDeriv_comp_map`. -/ +theorem isEVar_rnDeriv_comp_map [IsFiniteMeasure Q] [IsFiniteKernel κ] (hT : Measurable T) + (hκ : ∀ᵐ y ∂Q.map T, ∀ᵐ x ∂κ y, T x = y) (hQ : Q ≪ κ ∘ₘ Q.map T) + (hS : ∀ P ∈ S, IsFiniteMeasure P) (hSκ : ∀ P ∈ S, κ ∘ₘ P.map T = P) + (hSQ : ∀ P ∈ S, P.map T ≪ Q.map T) : + IsEVar (Q.rnDeriv (κ ∘ₘ Q.map T)) S := + .of_lintegral_le_measure_univ hS (Measure.measurable_rnDeriv _ _) fun P hP ↦ + (lintegral_rnDeriv_comp_map hT hκ hQ (hSκ P hP) (hSQ P hP)).le + +end ProbabilityTheory diff --git a/LMLExplore/EValues/TwoSample/Exponential.lean b/LMLExplore/EValues/TwoSample/Exponential.lean new file mode 100644 index 00000000..f64b11d8 --- /dev/null +++ b/LMLExplore/EValues/TwoSample/Exponential.lean @@ -0,0 +1,580 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import LMLExplore.EValues.Numeraire.Conditional +public import Mathlib.Probability.Distributions.Exponential +public import TauCeti.Probability.Distributions.Uniform.Basic +import TauCeti.MeasureTheory.Measure.Mixture.Basic +import TauCeti.Probability.Distributions.Beta.Uniform +import TauCeti.Probability.Distributions.Gamma.Beta +import TauCeti.Probability.Distributions.Gamma.Measurability +import TauCeti.Probability.Distributions.Gamma.PDF +import TauCeti.Probability.Distributions.Hypoexponential.PDF + +/-! +# The numeraire of two exponential observations against a common rate + +Two independent observations `X₁, X₂` have the laws `Exp(l₁)` and `Exp(l₂)`, and the null +hypothesis is that they have a common, unknown rate: it is the set of the products +`Exp(l) ⊗ Exp(l)`, `l > 0`. This is the two-sample test with a common nuisance parameter +of [HGLLA24] for the exponential family. + +Under every member of the null, the sum `s = X₁ + X₂` has the law `Gamma(2, l)` and, given `s`, +the first observation is uniform on `[0, s]` and the second one is `s - X₁`: the conditional law +given the sum is `ProbabilityTheory.Kernel.segmentUniform s` for every member. For the alternative +`Q = Exp(l₁) ⊗ Exp(l₂)`, let `Qᶜ = Kernel.segmentUniform ∘ₘ Q.map (x.1 + x.2)` be the law under +which the sum has its law under `Q` and, given the sum `s`, the first observation is uniform on +`[0, s]`. Then: + +* the likelihood ratio of `Q` to `Qᶜ`, the conditional e-variable of [HGLLA24], is the numeraire + of the null under `Q`, and `Qᶜ` is the reverse information projection of `Q` on the null; +* for `l₁ ≠ l₂`, `Qᶜ` is the mixture of the members `Exp(1 / m) ⊗ Exp(1 / m)` of the null whose + mean `m` is uniform between `1 / l₁` and `1 / l₂`; +* no finite combination of members of the null is equal to it, against the conjecture of + [HGLLA24, §5.2]. + +## Main definitions + +* `ProbabilityTheory.Kernel.segmentUniform`: the law of `(u * s, (1 - u) * s)` for `u` uniform + on `[0, 1]`. + +## Main statements + +* `ProbabilityTheory.isNumeraire_rnDeriv_segmentUniform_comp`: the likelihood ratio of `Q` to + `Qᶜ` is the numeraire. +* `ProbabilityTheory.ripr_prod_expMeasure`: `Qᶜ` is the reverse information projection. +* `ProbabilityTheory.segmentUniform_comp_map_add_prod_expMeasure`: `Qᶜ` is the mixture of the + members of the null whose mean is uniform between the two means. +* `ProbabilityTheory.sum_smul_prod_expMeasure_ne_ripr`: no finite combination of members of the + null is the reverse information projection. + +## Implementation notes + +That the law of the sum and the conditional law given the sum are as described under the null is +the beta–gamma algebra of Tau Ceti (`TauCeti.Probability.map_div_add_prod_gammaMeasure`). The law +of the sum under `Q` is Tau Ceti's hypoexponential law, with density +`l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)`, which is also the density of the +mixture of the laws `Gamma(2, 1 / m)` over the uniform mean `m`. The numeraire is then identified +through the effective null: the likelihood ratio is an e-variable by +`ProbabilityTheory.isEVar_rnDeriv_comp_map`, and `Qᶜ` is a mixture of members of the null, so +that `ProbabilityTheory.isNumeraire_of_ae_eq_rnDeriv` and +`ProbabilityTheory.ripr_eq_of_isEVar_rnDeriv` apply. For the last statement, a finite combination +of the laws `Gamma(2, l)` has a density `∑ c l * l ^ 2 * x * exp (- l x)`, which decays either too +slowly or too fast at infinity compared to `exp (- min l₁ l₂ * x)`. + +## References + +* [HGLLA24] Y. Hao, P. Grünwald, T. Lardy, L. Long and R. Adams, *E-values for k-sample tests + with exponential families*, arXiv:2303.00471, §5.2. +* [LRR25] M. Larsson, A. Ramdas and J. Ruf, *The numeraire e-variable and reverse information + projection*, Ann. Statist. 53 (2025), arXiv:2402.18810v4, Theorem 4.1. +-/ + +public section + +open MeasureTheory Set Filter Real TauCeti.Probability + +open scoped ENNReal Topology + +namespace ProbabilityTheory + +/-! ### The uniform law on a segment -/ + +/-- **The uniform law on a segment**: for `s : ℝ`, the law of `(u * s, (1 - u) * s)` for `u` +uniform on `[0, 1]`. For `s ≥ 0`, its first coordinate is uniform on `[0, s]` and the two +coordinates add up to `s`. -/ +noncomputable def Kernel.segmentUniform : Kernel ℝ (ℝ × ℝ) := + (Kernel.id ×ₖ Kernel.const ℝ (uniformMeasure 0 1)).map fun p ↦ (p.2 * p.1, (1 - p.2) * p.1) + +/-- The uniform law on a segment is a Markov kernel. -/ +instance Kernel.isMarkovKernel_segmentUniform : IsMarkovKernel Kernel.segmentUniform := by + unfold Kernel.segmentUniform + exact Kernel.IsMarkovKernel.map _ (by fun_prop) + +/-- The uniform law on the segment of the pairs with sum `s`. -/ +theorem Kernel.segmentUniform_apply (s : ℝ) : + Kernel.segmentUniform s = (uniformMeasure 0 1).map fun u ↦ (u * s, (1 - u) * s) := by + rw [Kernel.segmentUniform, Kernel.map_apply _ (by fun_prop), Kernel.prod_apply, Kernel.id_apply, + Kernel.const_apply, Measure.dirac_prod, Measure.map_map (by fun_prop) (by fun_prop)] + rfl + +/-- The mixture of the uniform laws on the segments against a law `ν` of the sum is the image of +`uniformMeasure 0 1 ⊗ ν` under `(u, s) ↦ (u * s, (1 - u) * s)`. -/ +theorem Kernel.segmentUniform_comp (ν : Measure ℝ) [SFinite ν] : + Kernel.segmentUniform ∘ₘ ν = + ((uniformMeasure 0 1).prod ν).map fun p ↦ (p.1 * p.2, (1 - p.1) * p.2) := by + rw [Kernel.segmentUniform, ← Measure.map_comp _ _ (by fun_prop), ← Measure.compProd_eq_comp_prod, + Measure.compProd_const, ← Measure.prod_swap, Measure.map_map (by fun_prop) measurable_swap] + rfl + +/-- Under the uniform law on the segment of the pairs with sum `s`, the coordinates add up to +`s`. -/ +theorem Kernel.ae_add_segmentUniform (s : ℝ) : ∀ᵐ x ∂Kernel.segmentUniform s, x.1 + x.2 = s := by + rw [Kernel.segmentUniform_apply, + ae_map_iff (by fun_prop) (measurableSet_eq_fun (by fun_prop) measurable_const)] + exact ae_of_all _ fun u ↦ by ring + +/-! ### The null of two exponential observations with a common rate -/ + +section Null + +variable {l l₁ l₂ : ℝ} + +/-- Exponential laws are s-finite. -/ +private instance sFinite_expMeasure (l : ℝ) : SFinite (expMeasure l) := + inferInstanceAs (SFinite (gammaMeasure 1 l)) + +/-- Exponential laws with positive rates are equivalent. -/ +private lemma expMeasure_absolutelyContinuous (l : ℝ) {l' : ℝ} (hl' : 0 < l') : + expMeasure l ≪ expMeasure l' := by + rw [expMeasure_eq_withDensity_restrict_Ioi l, expMeasure_eq_withDensity_restrict_Ioi l'] + refine (withDensity_absolutelyContinuous _ _).trans (withDensity_absolutelyContinuous' + (measurable_exponentialPDFReal l').ennreal_ofReal.aemeasurable ?_) + filter_upwards [ae_restrict_mem measurableSet_Ioi] with x hx + exact (ENNReal.ofReal_pos.mpr (exponentialPDFReal_pos hl' hx)).ne' + +/-- Products of exponential laws with positive rates are equivalent. -/ +private lemma prod_expMeasure_absolutelyContinuous (l₁ l₂ : ℝ) {l₁' l₂' : ℝ} (h₁ : 0 < l₁') + (h₂ : 0 < l₂') : + (expMeasure l₁).prod (expMeasure l₂) ≪ (expMeasure l₁').prod (expMeasure l₂') := + (expMeasure_absolutelyContinuous l₁ h₁).prod (expMeasure_absolutelyContinuous l₂ h₂) + +/-- The sum of two independent `Exp(l)` variables has the law `Gamma(2, l)`. -/ +private lemma map_add_prod_expMeasure_self (hl : 0 < l) : + ((expMeasure l).prod (expMeasure l)).map (fun x ↦ x.1 + x.2) = gammaMeasure 2 l := by + have : IsProbabilityMeasure (betaMeasure 1 1) := by rw [betaMeasure_one_one]; infer_instance + have h := congrArg (Measure.map Prod.snd) (map_div_add_prod_gammaMeasure one_pos one_pos hl) + rw [Measure.map_map measurable_snd (by fun_prop), Measure.map_snd_prod, measure_univ, + one_smul, one_add_one_eq_two] at h + exact h + +/-- **The conditional law given the sum under the null**: two independent `Exp(l)` variables are +`(u * s, (1 - u) * s)` for a sum `s` of law `Gamma(2, l)` and an independent uniform `u`. -/ +private lemma segmentUniform_comp_gammaMeasure (hl : 0 < l) : + Kernel.segmentUniform ∘ₘ gammaMeasure 2 l = (expMeasure l).prod (expMeasure l) := by + have h := map_div_add_prod_gammaMeasure one_pos one_pos hl + rw [betaMeasure_one_one, one_add_one_eq_two] at h + rw [Kernel.segmentUniform_comp, ← h, Measure.map_map (by fun_prop) (by fun_prop)] + conv_rhs => rw [← Measure.map_id (μ := (expMeasure l).prod (expMeasure l))] + refine Measure.map_congr ?_ + filter_upwards [ae_mem_prod_Ioi_gammaMeasure (a := 1) (r := l) (b := 1) (s := l)] with z hz + have hz0 : z.1 + z.2 ≠ 0 := (add_pos hz.1 hz.2).ne' + refine Prod.ext ?_ ?_ + · simp [field] + · simp [field] + +/-- Every member of the null is its own `Qᶜ`. -/ +private lemma segmentUniform_comp_map_add_prod_expMeasure_self (hl : 0 < l) : + Kernel.segmentUniform ∘ₘ ((expMeasure l).prod (expMeasure l)).map (fun x ↦ x.1 + x.2) = + (expMeasure l).prod (expMeasure l) := by + rw [map_add_prod_expMeasure_self hl, segmentUniform_comp_gammaMeasure hl] + +/-- `Exp(l₁) ⊗ Exp(l₂)` and its `Qᶜ` are equivalent: both are equivalent to `Exp(1) ⊗ Exp(1)`, +which is its own `Qᶜ`. -/ +private lemma absolutelyContinuous_segmentUniform_comp (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + (expMeasure l₁).prod (expMeasure l₂) ≪ + Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) ∧ + Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) ≪ + (expMeasure l₁).prod (expMeasure l₂) := by + have h1 := segmentUniform_comp_map_add_prod_expMeasure_self (l := 1) one_pos + refine ⟨(prod_expMeasure_absolutelyContinuous l₁ l₂ one_pos one_pos).trans ?_, ?_⟩ + · rw [← h1] + exact ((prod_expMeasure_absolutelyContinuous 1 1 h₁ h₂).map (by fun_prop)).comp_right _ + · refine (((prod_expMeasure_absolutelyContinuous l₁ l₂ one_pos one_pos).map + (by fun_prop)).comp_right _).trans ?_ + rw [h1] + exact prod_expMeasure_absolutelyContinuous 1 1 h₁ h₂ + +/-- The members of the null are probability measures. -/ +private lemma isProbabilityMeasure_of_mem_null {P : Measure (ℝ × ℝ)} + (hP : P ∈ Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) : + IsProbabilityMeasure P := by + obtain ⟨⟨l, hl⟩, rfl⟩ := hP + have := isProbabilityMeasure_expMeasure hl + infer_instance + +/-- Every member of the null has the conditional law `Kernel.segmentUniform` given the sum. -/ +private lemma segmentUniform_comp_map_add_of_mem {P : Measure (ℝ × ℝ)} + (hP : P ∈ Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) : + Kernel.segmentUniform ∘ₘ P.map (fun x ↦ x.1 + x.2) = P := by + obtain ⟨⟨l, hl⟩, rfl⟩ := hP + exact segmentUniform_comp_map_add_prod_expMeasure_self hl + +/-- Under every member of the null, the law of the sum is absolutely continuous with respect to +its law under `Exp(l₁) ⊗ Exp(l₂)`. -/ +private lemma map_add_absolutelyContinuous_of_mem (h₁ : 0 < l₁) (h₂ : 0 < l₂) + {P : Measure (ℝ × ℝ)} + (hP : P ∈ Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) : + P.map (fun x ↦ x.1 + x.2) ≪ + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) := by + obtain ⟨⟨l, hl⟩, rfl⟩ := hP + exact (prod_expMeasure_absolutelyContinuous l l h₁ h₂).map (by fun_prop) + +/-! ### The law of the sum and the mixture -/ + +/-- The density on `(0, ∞)` of the sum of two independent exponential variables with distinct +rates `l₁, l₂`. -/ +private noncomputable def sumDensity (l₁ l₂ x : ℝ) : ℝ := + l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁) + +/-- The density of the sum is nonnegative on `[0, ∞)`. -/ +private lemma sumDensity_nonneg (h₁ : 0 < l₁) (h₂ : 0 < l₂) {x : ℝ} (hx : 0 ≤ x) : + 0 ≤ sumDensity l₁ l₂ x := by + unfold sumDensity + rcases le_total l₁ l₂ with h | h + · exact div_nonneg (mul_nonneg (by positivity) (sub_nonneg.mpr (exp_le_exp.mpr + (by nlinarith)))) (sub_nonneg.mpr h) + · exact div_nonneg_of_nonpos (mul_nonpos_of_nonneg_of_nonpos (by positivity) + (sub_nonpos.mpr (exp_le_exp.mpr (by nlinarith)))) (sub_nonpos.mpr h) + +/-- **The law of the sum of two exponential variables** with distinct rates `l₁, l₂` has density +`l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)`. -/ +private lemma map_add_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + (volume.restrict (Ioi 0)).withDensity fun x ↦ ENNReal.ofReal (sumDensity l₁ l₂ x) := by + have := isProbabilityMeasure_expMeasure h₁ + have := isProbabilityMeasure_expMeasure h₂ + have hr : ∀ i, 0 < ![l₁, l₂] i := by simp [Fin.forall_fin_two, h₁, h₂] + have hinj : Function.Injective ![l₁, l₂] := by + intro i j hij + fin_cases i <;> fin_cases j <;> simp_all [eq_comm] + have hpi : (fun i ↦ expMeasure (![l₁, l₂] i)) = ![expMeasure l₁, expMeasure l₂] := by + ext i : 1 + fin_cases i <;> rfl + have hmap : ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + hypoexponentialMeasure ![l₁, l₂] := by + rw [hypoexponentialMeasure_eq_map hr, hpi, + ← (measurePreserving_finTwoArrow_vec (expMeasure l₁) (expMeasure l₂)).map_eq, + Measure.map_map (by fun_prop) (MeasurableEquiv.measurable _)] + congr 1 + ext x + simp [Fin.sum_univ_two, MeasurableEquiv.finTwoArrow] + rw [hmap, ← Measure.restrict_eq_self_of_ae_mem (s := Ioi 0) + (μ := hypoexponentialMeasure ![l₁, l₂]) (ae_pos_hypoexponentialMeasure _), + hypoexponentialMeasure_eq_withDensity, restrict_withDensity measurableSet_Ioi] + refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr + (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) + rw [hypoexponentialPDF_eq_ofReal, hypoexponentialPDFReal_of_injective hr hinj hx] + have h01 : (Finset.univ.erase (0 : Fin 2)) = {1} := by decide + have h10 : (Finset.univ.erase (1 : Fin 2)) = {0} := by decide + have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm + have h21 : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h + simp only [Fin.sum_univ_two, h01, h10, Finset.prod_singleton, Matrix.cons_val_zero, + Matrix.cons_val_one, sumDensity] + congr 1 + field_simp + ring + +/-- The density of `Gamma(2, r)` on `[0, ∞)`. -/ +private lemma gammaPDF_two {r x : ℝ} (hx : 0 ≤ x) : + gammaPDF 2 r x = ENNReal.ofReal (r ^ 2 * x * exp (-(r * x))) := by + rw [gammaPDF_of_nonneg hx, Real.Gamma_two, show (2 : ℝ) - 1 = 1 by norm_num, Real.rpow_one, + Real.rpow_two, div_one] + +/-- The uniform law on `(a, b]` lives on `(a, b]`. -/ +private lemma ae_mem_Ioc_uniformMeasure (a b : ℝ) : ∀ᵐ m ∂uniformMeasure a b, m ∈ Ioc a b := by + rw [uniformMeasure_eq_smul] + exact Measure.ae_smul_measure (ae_restrict_mem measurableSet_Ioc) _ + +/-- The density of the mixture of the laws `Gamma(2, 1 / m)` over the uniform mean `m`. -/ +private lemma lintegral_gammaPDF_uniformMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) + {x : ℝ} (hx : 0 < x) : + ∫⁻ m, gammaPDF 2 m⁻¹ x ∂uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹) = + ENNReal.ofReal (sumDensity l₁ l₂ x) := by + set α := min l₁⁻¹ l₂⁻¹ + set β := max l₁⁻¹ l₂⁻¹ + have hα : 0 < α := lt_min (inv_pos.mpr h₁) (inv_pos.mpr h₂) + have hαβ : α < β := min_lt_max.mpr (inv_injective.ne h) + have hpos : ∀ m ∈ uIcc α β, 0 < m := fun m hm ↦ hα.trans_le (by + rw [uIcc_of_le hαβ.le] at hm; exact hm.1) + have hderiv : ∀ m ∈ uIcc α β, HasDerivAt (fun m ↦ exp (-(m⁻¹ * x))) + (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) m := + fun m hm ↦ (((hasDerivAt_inv (hpos m hm).ne').mul_const x).neg.exp).congr_deriv + (by simp only [Pi.neg_apply]; ring) + have hcont : ContinuousOn (fun m ↦ exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) (uIcc α β) := by + refine ContinuousOn.mul (by fun_prop (disch := exact fun m hm ↦ (hpos m hm).ne')) ?_ + fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (hpos m hm).ne') + have hpdf : ∀ m ∈ Ioc α β, gammaPDF 2 m⁻¹ x = + ENNReal.ofReal (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) := by + intro m hm + rw [gammaPDF_two hx.le] + congr 1 + ring + rw [uniformMeasure_eq_smul, lintegral_smul_measure, setLIntegral_congr_fun measurableSet_Ioc hpdf, + ← ofReal_integral_eq_lintegral_ofReal, ← intervalIntegral.integral_of_le hαβ.le, + intervalIntegral.integral_eq_sub_of_hasDerivAt hderiv hcont.intervalIntegrable, + smul_eq_mul, ← ENNReal.ofReal_inv_of_pos (sub_pos.mpr hαβ), + ← ENNReal.ofReal_mul (inv_nonneg.mpr (sub_pos.mpr hαβ).le)] + · congr 1 + have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm + rcases h.lt_or_gt with h' | h' + · have hinv : l₂⁻¹ ≤ l₁⁻¹ := (inv_strictAnti₀ h₁ h').le + simp only [α, β, min_eq_right hinv, max_eq_left hinv, inv_inv, sumDensity] + field_simp + · have hinv : l₁⁻¹ ≤ l₂⁻¹ := (inv_strictAnti₀ h₂ h').le + simp only [α, β, min_eq_left hinv, max_eq_right hinv, inv_inv, sumDensity] + field_simp + ring + · rw [uIcc_of_le hαβ.le] at hcont + exact hcont.integrableOn_Icc.mono_set Ioc_subset_Icc_self + · refine ae_restrict_of_forall_mem measurableSet_Ioc fun m hm ↦ ?_ + have := hpos m (by rw [uIcc_of_le hαβ.le]; exact Ioc_subset_Icc_self hm) + positivity + +/-- **The sum of two exponential variables is a mixture of gamma laws**: for distinct rates, the +law of the sum is the mixture of the laws `Gamma(2, 1 / m)` over the mean `m` uniform between +`1 / l₁` and `1 / l₂`. -/ +private lemma map_add_prod_expMeasure_eq_bind (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + (uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)).bind fun m ↦ gammaMeasure 2 m⁻¹ := by + have := isProbabilityMeasure_uniformMeasure + (min_lt_max.mpr (inv_injective.ne h) : min l₁⁻¹ l₂⁻¹ < max l₁⁻¹ l₂⁻¹) + simp_rw [gammaMeasure_eq_withDensity_restrict_Ioi] + have hf : Measurable (Function.uncurry fun m x ↦ gammaPDF 2 m⁻¹ x) := + measurable_uncurry_gammaPDF.comp (f := fun p : ℝ × ℝ ↦ (((2 : ℝ), p.1⁻¹), p.2)) (by fun_prop) + rw [map_add_prod_expMeasure h₁ h₂ h, Measure.bind_withDensity _ _ hf] + refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr + (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) + exact (lintegral_gammaPDF_uniformMeasure h₁ h₂ h hx).symm + +/-- The laws `Gamma(2, 1 / m)`, as a kernel. -/ +private noncomputable def gammaTwoInvKernel : Kernel ℝ ℝ where + toFun m := gammaMeasure 2 m⁻¹ + measurable' := measurable_gammaMeasure.comp (f := fun m : ℝ ↦ ((2 : ℝ), m⁻¹)) (by fun_prop) + +/-- `Qᶜ` is in the effective null. -/ +private lemma segmentUniform_comp_mem_effectiveSet (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) ∈ + effectiveSet (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) := by + have hS : ∀ P ∈ Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l), + P univ ≤ 1 := fun P hP ↦ (isProbabilityMeasure_of_mem_null hP).measure_univ.le + rcases eq_or_ne l₁ l₂ with rfl | h + · rw [segmentUniform_comp_map_add_prod_expMeasure_self h₁] + exact subset_effectiveSet hS ⟨⟨l₁, h₁⟩, rfl⟩ + have hα : 0 < min l₁⁻¹ l₂⁻¹ := lt_min (inv_pos.mpr h₁) (inv_pos.mpr h₂) + have := isProbabilityMeasure_uniformMeasure + (min_lt_max.mpr (inv_injective.ne h) : min l₁⁻¹ l₂⁻¹ < max l₁⁻¹ l₂⁻¹) + rw [map_add_prod_expMeasure_eq_bind h₁ h₂ h] + change Kernel.segmentUniform ∘ₘ (gammaTwoInvKernel ∘ₘ _) ∈ _ + rw [Measure.comp_assoc] + refine comp_mem_effectiveSet prob_le_one ?_ + filter_upwards [ae_mem_Ioc_uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)] with m hm + have hm0 : 0 < m⁻¹ := inv_pos.mpr (hα.trans hm.1) + rw [Kernel.comp_apply] + change Kernel.segmentUniform ∘ₘ gammaMeasure 2 m⁻¹ ∈ _ + rw [segmentUniform_comp_gammaMeasure hm0] + exact subset_effectiveSet hS ⟨⟨m⁻¹, hm0⟩, rfl⟩ + +/-- The likelihood ratio of `Exp(l₁) ⊗ Exp(l₂)` to its `Qᶜ` is an e-variable for the null. -/ +private lemma isEVar_rnDeriv_segmentUniform_comp (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + IsEVar (((expMeasure l₁).prod (expMeasure l₂)).rnDeriv + (Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map fun x ↦ x.1 + x.2)) + (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) := by + have := isProbabilityMeasure_expMeasure h₁ + have := isProbabilityMeasure_expMeasure h₂ + exact isEVar_rnDeriv_comp_map (by fun_prop) (ae_of_all _ Kernel.ae_add_segmentUniform) + (absolutelyContinuous_segmentUniform_comp h₁ h₂).1 + (fun P hP ↦ have := isProbabilityMeasure_of_mem_null hP; inferInstance) + (fun _ ↦ segmentUniform_comp_map_add_of_mem) + (fun _ ↦ map_add_absolutelyContinuous_of_mem h₁ h₂) + +end Null + +/-! ### The numeraire and the reverse information projection -/ + +section Main + +variable {l₁ l₂ : ℝ} + +/-- **The numeraire of two exponential observations against a common rate.** For positive rates +`l₁, l₂`, the likelihood ratio of `Exp(l₁) ⊗ Exp(l₂)` to the law under which the sum of the two +observations has its law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation +is uniform on `[0, s]`, is the numeraire of the null of the products `Exp(l) ⊗ Exp(l)`, `l > 0`: +it is an e-variable, and no e-variable has a larger expected logarithm under +`Exp(l₁) ⊗ Exp(l₂)` (`ProbabilityTheory.IsNumeraire.eintegral_log_le`). -/ +theorem isNumeraire_rnDeriv_segmentUniform_comp (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + IsNumeraire (((expMeasure l₁).prod (expMeasure l₂)).rnDeriv + (Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map fun x ↦ x.1 + x.2)) + (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) + ((expMeasure l₁).prod (expMeasure l₂)) := by + have := isProbabilityMeasure_expMeasure h₁ + have := isProbabilityMeasure_expMeasure h₂ + exact isNumeraire_of_ae_eq_rnDeriv (isEVar_rnDeriv_segmentUniform_comp h₁ h₂) + (segmentUniform_comp_mem_effectiveSet h₁ h₂) (absolutelyContinuous_segmentUniform_comp h₁ h₂).1 + .rfl + +/-- **The reverse information projection of two exponential observations** on the null of the +products `Exp(l) ⊗ Exp(l)`, `l > 0`, is the law under which the sum of the two observations has +its law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation is uniform on +`[0, s]`. -/ +theorem ripr_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + ripr ((expMeasure l₁).prod (expMeasure l₂)) + (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) = + Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map fun x ↦ x.1 + x.2 := by + have := isProbabilityMeasure_expMeasure h₁ + have := isProbabilityMeasure_expMeasure h₂ + exact ripr_eq_of_isEVar_rnDeriv + (fun P hP ↦ have := isProbabilityMeasure_of_mem_null hP; inferInstance) + (segmentUniform_comp_mem_effectiveSet h₁ h₂) (absolutelyContinuous_segmentUniform_comp h₁ h₂).1 + (absolutelyContinuous_segmentUniform_comp h₁ h₂).2 (isEVar_rnDeriv_segmentUniform_comp h₁ h₂) + +/-- **The reverse information projection is a continuous mixture of members of the null.** For +distinct positive rates `l₁, l₂`, the law under which the sum of the two observations has its +law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation is uniform on +`[0, s]`, is the mixture of the members `Exp(1 / m) ⊗ Exp(1 / m)` of the null whose mean `m` is +uniform between `1 / l₁` and `1 / l₂`. -/ +theorem segmentUniform_comp_map_add_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : + Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + (uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)).bind + fun m ↦ (expMeasure m⁻¹).prod (expMeasure m⁻¹) := by + have hα : 0 < min l₁⁻¹ l₂⁻¹ := lt_min (inv_pos.mpr h₁) (inv_pos.mpr h₂) + rw [map_add_prod_expMeasure_eq_bind h₁ h₂ h, Measure.bind_bind (f := fun m ↦ gammaMeasure 2 m⁻¹) + gammaTwoInvKernel.measurable.aemeasurable Kernel.segmentUniform.aemeasurable] + refine Measure.bind_congr_right ?_ + filter_upwards [ae_mem_Ioc_uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)] with m hm + exact segmentUniform_comp_gammaMeasure (inv_pos.mpr (hα.trans hm.1)) + +/-- A nonnegative finite combination of the functions `x ↦ l ^ 2 * x * exp (- l x)`, `l > 0`, is +not `C (exp (- a x) - exp (- b x))` on `(0, ∞)` for `a < b` and `C > 0`: at infinity, the terms +with `l ≤ a` decay too slowly and the other ones too fast. -/ +private lemma not_forall_sum_eq_mul_sub {a b C : ℝ} (hab : a < b) (hC : 0 < C) (s : Finset ℝ) + (c : ℝ → ℝ) (hc : ∀ l ∈ s, 0 ≤ c l) (hs : ∀ l ∈ s, 0 < l) : + ¬ ∀ x > 0, ∑ l ∈ s, c l * (l ^ 2 * x * exp (-(l * x))) = + C * (exp (-(a * x)) - exp (-(b * x))) := by + intro h + have hterm : ∀ x > 0, ∀ l ∈ s, 0 ≤ c l * (l ^ 2 * x * exp (-(l * x))) := fun x hx l hl ↦ by + have := hc l hl + positivity + -- the coefficients of the rates `l ≤ a` vanish + have h0 : ∀ l ∈ s, l ≤ a → c l = 0 := by + intro l hl hla + by_contra hcl + have hcl' : 0 < c l * l ^ 2 := mul_pos ((hc l hl).lt_of_ne' hcl) (pow_pos (hs l hl) 2) + set x := 2 * C / (c l * l ^ 2) with hx_def + have hx : 0 < x := by positivity + have h1 := Finset.single_le_sum (hterm x hx) hl + rw [h x hx] at h1 + have h2 : C * exp (-(a * x)) ≤ C * exp (-(l * x)) := + mul_le_mul_of_nonneg_left (exp_le_exp.mpr (by nlinarith)) hC.le + have h3 : 0 < C * exp (-(b * x)) := mul_pos hC (exp_pos _) + have h4 : 0 < C * exp (-(l * x)) := mul_pos hC (exp_pos _) + have h5 : c l * (l ^ 2 * x * exp (-(l * x))) = 2 * (C * exp (-(l * x))) := by + have := (hs l hl).ne' + rw [hx_def] + field_simp + linarith + -- multiply by `exp (a x)` and let `x` tend to infinity + have hlim : Tendsto (fun x ↦ ∑ l ∈ s, c l * l ^ 2 * (x * exp (-((l - a) * x)))) atTop + (𝓝 0) := by + rw [show (0 : ℝ) = ∑ _l ∈ s, 0 by simp] + refine tendsto_finsetSum _ fun l hl ↦ ?_ + by_cases hla : l ≤ a + · simp [h0 l hl hla] + have hε : 0 < l - a := sub_pos.mpr (not_le.mp hla) + have := ((tendsto_pow_mul_exp_neg_atTop_nhds_zero 1).comp + (Tendsto.const_mul_atTop hε tendsto_id)).const_mul (c l * l ^ 2 / (l - a)) + rw [mul_zero] at this + refine this.congr fun x ↦ ?_ + simp only [Function.comp_apply, id, pow_one] + field_simp + have hlim' : Tendsto (fun x ↦ C * (1 - exp (-((b - a) * x)))) atTop (𝓝 C) := by + have : Tendsto (fun x ↦ exp (-((b - a) * x))) atTop (𝓝 0) := + tendsto_exp_neg_atTop_nhds_zero.comp (Tendsto.const_mul_atTop (sub_pos.mpr hab) tendsto_id) + simpa using (this.const_sub 1).const_mul C + refine hC.ne (tendsto_nhds_unique (hlim.congr' ?_) hlim') + filter_upwards [eventually_gt_atTop 0] with x hx + have key : ∀ l, exp (-((l - a) * x)) = exp (-(l * x)) * exp (a * x) := fun l ↦ by + rw [← exp_add] + ring_nf + have ha : exp (-(a * x)) * exp (a * x) = 1 := by + rw [← exp_add] + simp + calc ∑ l ∈ s, c l * l ^ 2 * (x * exp (-((l - a) * x))) + _ = (∑ l ∈ s, c l * (l ^ 2 * x * exp (-(l * x)))) * exp (a * x) := by + rw [Finset.sum_mul] + refine Finset.sum_congr rfl fun l _ ↦ ?_ + rw [key] + ring + _ = C * (1 - exp (-((b - a) * x))) := by + rw [h x hx, key b] + linear_combination C * ha + +/-- **No finite combination of members of the null is the reverse information projection.** For +distinct positive rates `l₁, l₂`, no finite combination of the laws `Exp(l) ⊗ Exp(l)` is the +reverse information projection of `Exp(l₁) ⊗ Exp(l₂)` on the null of the products +`Exp(l) ⊗ Exp(l)`, `l > 0`. In particular it is not a mixture of two members of the null, against +the conjecture of [HGLLA24, §5.2]. -/ +theorem sum_smul_prod_expMeasure_ne_ripr (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) + (s : Finset ℝ) (w : ℝ → ℝ≥0∞) : + ∑ l ∈ s, w l • (expMeasure l).prod (expMeasure l) ≠ + ripr ((expMeasure l₁).prod (expMeasure l₂)) + (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) := by + have := isProbabilityMeasure_expMeasure h₁ + have := isProbabilityMeasure_expMeasure h₂ + rw [ripr_prod_expMeasure h₁ h₂] + intro heq + set s' := s.filter (0 < ·) with hs' + have hs'pos : ∀ l ∈ s', 0 < l := fun l hl ↦ (Finset.mem_filter.mp hl).2 + -- the laws of the sum agree + have hsum : ∀ A, MeasurableSet A → ∑ l ∈ s', w l * gammaMeasure 2 l A = + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) A := by + intro A hA + have h1 := congrArg (fun μ : Measure (ℝ × ℝ) ↦ μ ((fun x ↦ x.1 + x.2) ⁻¹' A)) heq + simp only [FunLike.coe_sum, Finset.sum_apply, smul_apply, smul_eq_mul] at h1 + rw [← Measure.map_apply (by fun_prop) hA, + map_comp_eq_of_ae_eq (by fun_prop) (ae_of_all _ Kernel.ae_add_segmentUniform)] at h1 + rw [← h1, hs', Finset.sum_filter] + refine Finset.sum_congr rfl fun l _ ↦ ?_ + split_ifs with hl + · rw [← map_add_prod_expMeasure_self hl, Measure.map_apply (by fun_prop) hA] + · simp [expMeasure_of_nonpos (not_lt.mp hl)] + -- so the weights are finite + have hw : ∀ l ∈ s', w l ≠ ∞ := by + have hmass := hsum univ MeasurableSet.univ + rw [measure_univ, Finset.sum_congr rfl fun l hl ↦ by + have := isProbabilityMeasure_gammaMeasure two_pos (hs'pos l hl) + rw [measure_univ, mul_one]] at hmass + exact fun l hl ↦ ne_top_of_le_ne_top ENNReal.one_ne_top + (hmass ▸ Finset.single_le_sum (fun _ _ ↦ zero_le) hl) + -- and the densities of the laws of the sum agree on `(0, ∞)` + have hae : (fun x ↦ ∑ l ∈ s', w l * gammaPDF 2 l x) =ᵐ[volume.restrict (Ioi 0)] + fun x ↦ ENNReal.ofReal (sumDensity l₁ l₂ x) := by + refine ae_eq_of_forall_setLIntegral_eq_of_sigmaFinite + (Finset.measurable_fun_sum _ fun l _ ↦ (measurable_gammaPDF 2 l).const_mul _) + (by unfold sumDensity; fun_prop) fun A hA _ ↦ ?_ + rw [lintegral_finsetSum _ fun l _ ↦ (measurable_gammaPDF 2 l).const_mul _, + ← withDensity_apply _ hA, ← map_add_prod_expMeasure h₁ h₂ h, ← hsum A hA] + refine Finset.sum_congr rfl fun l _ ↦ ?_ + rw [lintegral_const_mul _ (measurable_gammaPDF 2 l), ← withDensity_apply _ hA, + ← gammaMeasure_eq_withDensity_restrict_Ioi] + have hreal : ∀ x ∈ Ioi (0 : ℝ), + ∑ l ∈ s', (w l).toReal * (l ^ 2 * x * exp (-(l * x))) = sumDensity l₁ l₂ x := by + refine Measure.eqOn_open_of_ae_eq (μ := volume) ?_ isOpen_Ioi (by fun_prop) + (by unfold sumDensity; fun_prop) + filter_upwards [hae, ae_restrict_mem measurableSet_Ioi] with x hx (hx0 : 0 < x) + rw [← ENNReal.ofReal_eq_ofReal_iff (Finset.sum_nonneg fun l hl ↦ + have := hs'pos l hl; by positivity) (sumDensity_nonneg h₁ h₂ hx0.le), ← hx, + ENNReal.ofReal_sum_of_nonneg fun l hl ↦ have := hs'pos l hl; by positivity] + refine Finset.sum_congr rfl fun l hl ↦ ?_ + rw [ENNReal.ofReal_mul ENNReal.toReal_nonneg, ENNReal.ofReal_toReal (hw l hl), + gammaPDF_two hx0.le] + rcases h.lt_or_gt with h12 | h12 + · refine not_forall_sum_eq_mul_sub h12 (div_pos (mul_pos h₁ h₂) (sub_pos.mpr h12)) s' + (fun l ↦ (w l).toReal) (fun _ _ ↦ ENNReal.toReal_nonneg) hs'pos fun x hx ↦ ?_ + rw [hreal x hx, sumDensity] + ring + · refine not_forall_sum_eq_mul_sub h12 (div_pos (mul_pos h₁ h₂) (sub_pos.mpr h12)) s' + (fun l ↦ (w l).toReal) (fun _ _ ↦ ENNReal.toReal_nonneg) hs'pos fun x hx ↦ ?_ + have h21 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm + have h12' : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h + rw [hreal x hx, sumDensity] + field_simp + ring + +end Main + +end ProbabilityTheory diff --git a/LMLExplore/EValues/TwoSample/Geometric.lean b/LMLExplore/EValues/TwoSample/Geometric.lean new file mode 100644 index 00000000..cf611578 --- /dev/null +++ b/LMLExplore/EValues/TwoSample/Geometric.lean @@ -0,0 +1,392 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import LMLExplore.EValues.Numeraire.Conditional +public import Mathlib.Probability.Distributions.Geometric +public import Mathlib.Probability.UniformOn +import TauCeti.Probability.Distributions.Geometric.Measurability +import TauCeti.Probability.Distributions.Uniform.Basic + +/-! +# The numeraire of two geometric observations against a common success probability + +Two independent observations `X₁, X₂` in `ℕ` have the geometric laws `Geom(p₁)` and `Geom(p₂)` +of Mathlib's `ProbabilityTheory.geometricMeasure`, and the null hypothesis is that they have a +common, unknown success probability: it is the set of the products `Geom(p) ⊗ Geom(p)`. Under +every member of the null the mass of a pair `(j, k)` depends only on `j + k`, so, given the sum +`s`, the first observation is uniform on `{0, …, s}`: the conditional law given the sum is +`ProbabilityTheory.Kernel.antidiagonalUniform s` for every member. + +For the alternative `Q = Geom(p₁) ⊗ Geom(p₂)`, the likelihood ratio of `Q` to the law +`Kernel.antidiagonalUniform ∘ₘ Q.map (x.1 + x.2)`, under which the sum has its law under `Q` and, +given the sum `s`, the first observation is uniform on `{0, …, s}`, is the numeraire of the null +under `Q` (`ProbabilityTheory.isNumeraire_rnDeriv_antidiagonalUniform_comp`). As for two +exponential observations (`ProbabilityTheory.isNumeraire_rnDeriv_segmentUniform_comp`), that law +is the mixture of the members `Geom(1 / m) ⊗ Geom(1 / m)` of the null over `m` uniform between +`1 / p₁` and `1 / p₂`, that is, over the mean `m - 1` uniform between the two means, so that it +belongs to the effective null. + +## Main definitions + +* `ProbabilityTheory.Kernel.antidiagonalUniform`: the uniform law on the pairs `(j, k)` of natural + numbers with `j + k = s`. + +## Main statements + +* `ProbabilityTheory.isNumeraire_rnDeriv_antidiagonalUniform_comp`: the numeraire of two geometric + observations against a common success probability. + +## Implementation notes + +Mathlib's `geometricMeasure 0` is the point mass at `0`, as is `geometricMeasure 1`, so the +statement holds for every pair of success probabilities except those for which `Q` is the point +mass at `(0, 0)`. The mixture is computed on the masses of the points, through the antiderivative +`(1 - 1 / m) ^ (s + 1) / (s + 1)` of `(1 - 1 / m) ^ s / m ^ 2`. + +## References + +* [HGLLA24] Y. Hao, P. Grünwald, T. Lardy, L. Long and R. Adams, *E-values for k-sample tests + with exponential families*, arXiv:2303.00471. +* [LRR25] M. Larsson, A. Ramdas and J. Ruf, *The numeraire e-variable and reverse information + projection*, Ann. Statist. 53 (2025), arXiv:2402.18810v4, Theorem 4.1. +-/ + +public section + +open MeasureTheory Set Finset TauCeti.Probability + +open scoped ENNReal + +namespace ProbabilityTheory + +/-! ### The uniform law on an antidiagonal -/ + +/-- **The uniform law on an antidiagonal**: for `s : ℕ`, the uniform law on the `s + 1` pairs +`(j, k)` of natural numbers with `j + k = s`. Its first coordinate is uniform on `{0, …, s}`. -/ +noncomputable def Kernel.antidiagonalUniform : Kernel ℕ (ℕ × ℕ) := + Kernel.ofFunOfCountable fun s ↦ uniformOn (antidiagonal s : Set (ℕ × ℕ)) + +/-- The uniform law on the pairs with sum `s`. -/ +theorem Kernel.antidiagonalUniform_apply (s : ℕ) : + Kernel.antidiagonalUniform s = uniformOn (antidiagonal s : Set (ℕ × ℕ)) := by + unfold Kernel.antidiagonalUniform + rfl + +/-- The uniform law on an antidiagonal is a Markov kernel. -/ +instance Kernel.isMarkovKernel_antidiagonalUniform : IsMarkovKernel Kernel.antidiagonalUniform := + ⟨fun s ↦ by + rw [Kernel.antidiagonalUniform_apply] + exact isProbabilityMeasure_uniformOn (finite_toSet _) ⟨(0, s), by simp⟩⟩ + +/-- The mass of a point under the uniform law on the pairs with sum `s`. -/ +theorem Kernel.antidiagonalUniform_singleton (s : ℕ) (x : ℕ × ℕ) : + Kernel.antidiagonalUniform s {x} = if x.1 + x.2 = s then ((s : ℝ≥0∞) + 1)⁻¹ else 0 := by + classical + rw [Kernel.antidiagonalUniform_apply, ← coe_singleton, uniformOn_apply_finset, + Nat.card_antidiagonal] + split_ifs with hx + · rw [Finset.inter_singleton_of_mem (mem_antidiagonal.mpr hx)] + simp + · rw [Finset.inter_singleton_of_notMem (by simpa using hx)] + simp + +/-- Under the uniform law on the pairs with sum `s`, the coordinates add up to `s`. -/ +theorem Kernel.ae_add_antidiagonalUniform (s : ℕ) : + ∀ᵐ x ∂Kernel.antidiagonalUniform s, x.1 + x.2 = s := by + rw [ae_iff, Kernel.antidiagonalUniform_apply, ← uniformOn_inter_self (finite_toSet _)] + convert measure_empty (μ := uniformOn (antidiagonal s : Set (ℕ × ℕ))) + ext x + simp + +/-! ### Laws on pairs whose masses depend on the sum -/ + +section Pairs + +variable {μ : Measure (ℕ × ℕ)} + +/-- The law of the sum at a point is the sum of the masses of the pairs with that sum. -/ +private lemma map_add_singleton (μ : Measure (ℕ × ℕ)) (s : ℕ) : + μ.map (fun x ↦ x.1 + x.2) {s} = ∑ x ∈ antidiagonal s, μ {x} := by + rw [Measure.map_apply (measurable_of_countable _) (measurableSet_singleton s), + sum_measure_singleton] + congr 1 + ext x + simp + +/-- The mass of a point under `Kernel.antidiagonalUniform ∘ₘ ν`. -/ +private lemma antidiagonalUniform_comp_singleton (ν : Measure ℕ) (x : ℕ × ℕ) : + (Kernel.antidiagonalUniform ∘ₘ ν) {x} = + ν {x.1 + x.2} * (((x.1 + x.2 : ℕ) : ℝ≥0∞) + 1)⁻¹ := by + rw [Measure.bind_apply (measurableSet_singleton x) (Kernel.aemeasurable _), lintegral_countable'] + simp_rw [Kernel.antidiagonalUniform_singleton] + rw [tsum_eq_single (x.1 + x.2) fun s hs ↦ by simp [Ne.symm hs]] + simp [mul_comm] + +/-- A law on pairs whose mass at a point depends only on the sum has the conditional law +`Kernel.antidiagonalUniform` given the sum. -/ +private lemma antidiagonalUniform_comp_map_add_of_singleton {g : ℕ → ℝ≥0∞} + (h : ∀ x, μ {x} = g (x.1 + x.2)) : + Kernel.antidiagonalUniform ∘ₘ μ.map (fun x ↦ x.1 + x.2) = μ := by + refine Measure.ext_iff_singleton.mpr fun x ↦ ?_ + rw [antidiagonalUniform_comp_singleton, map_add_singleton, + sum_congr rfl fun y hy ↦ by rw [h y, mem_antidiagonal.mp hy], sum_const, + Nat.card_antidiagonal, nsmul_eq_mul, h x, mul_comm _ (g _), mul_assoc] + push_cast + rw [ENNReal.mul_inv_cancel (by positivity) (by simp), mul_one] + +/-- Every law on pairs is absolutely continuous with respect to the law with its law of the sum +and the conditional law `Kernel.antidiagonalUniform` given the sum. -/ +private lemma absolutelyContinuous_antidiagonalUniform_comp (μ : Measure (ℕ × ℕ)) : + μ ≪ Kernel.antidiagonalUniform ∘ₘ μ.map (fun x ↦ x.1 + x.2) := by + refine Measure.AbsolutelyContinuous.mk fun t _ ht ↦ ?_ + rw [← biUnion_of_singleton t, measure_biUnion_null_iff t.to_countable] at ht ⊢ + intro x hx + have hx' := ht x hx + rw [antidiagonalUniform_comp_singleton, mul_eq_zero, map_add_singleton, + Finset.sum_eq_zero_iff] at hx' + exact hx'.resolve_right (by simp) x (mem_antidiagonal.mpr rfl) + +end Pairs + +/-! ### Geometric laws -/ + +section Geometric + +variable {p₁ p₂ : unitInterval} + +/-- The geometric law with success probability `1` is the point mass at `0`, as is the junk value +`geometricMeasure 0`. -/ +private lemma geometricMeasure_one : geometricMeasure 1 = geometricMeasure 0 := by + have h0 : geometricMeasure 0 = Measure.dirac 0 := by simp [geometricMeasure] + refine Measure.ext_iff_singleton.mpr fun n ↦ ?_ + rw [geometricMeasure_singleton one_ne_zero n, h0] + rcases n with _ | n <;> simp + +/-- Every geometric law is one with a positive success probability. -/ +private lemma exists_ne_zero_geometricMeasure_eq (p : unitInterval) : + ∃ q : unitInterval, q ≠ 0 ∧ geometricMeasure q = geometricMeasure p := by + by_cases hp : p = 0 + · exact ⟨1, one_ne_zero, hp ▸ geometricMeasure_one⟩ + · exact ⟨p, hp, rfl⟩ + +/-- The mass of a point under a product of geometric laws with positive success +probabilities. -/ +private lemma prod_geometricMeasure_singleton (h₁ : p₁ ≠ 0) (h₂ : p₂ ≠ 0) (x : ℕ × ℕ) : + (geometricMeasure p₁).prod (geometricMeasure p₂) {x} = + ENNReal.ofReal ((1 - p₁) ^ x.1 * p₁ * ((1 - p₂) ^ x.2 * p₂)) := by + rw [← singleton_prod_singleton, Measure.prod_prod, geometricMeasure_singleton h₁, + geometricMeasure_singleton h₂, ← ENNReal.ofReal_mul (geometricMeasure_nonneg _ _)] + +/-- The mass of a point under `Geom(p) ⊗ Geom(p)` depends only on the sum. -/ +private lemma exists_prod_geometricMeasure_self_singleton (p : unitInterval) : + ∃ g : ℕ → ℝ≥0∞, ∀ x : ℕ × ℕ, + (geometricMeasure p).prod (geometricMeasure p) {x} = g (x.1 + x.2) := by + obtain ⟨q, hq, hqp⟩ := exists_ne_zero_geometricMeasure_eq p + refine ⟨fun s ↦ ENNReal.ofReal ((1 - q) ^ s * q ^ 2), fun x ↦ ?_⟩ + rw [← hqp, prod_geometricMeasure_singleton hq hq] + congr 1 + ring + +/-- The point `0` has positive mass under every geometric law. -/ +private lemma geometricMeasure_singleton_zero_pos (p : unitInterval) : + 0 < geometricMeasure p {0} := by + obtain ⟨q, hq, hqp⟩ := exists_ne_zero_geometricMeasure_eq p + rw [← hqp, geometricMeasure_singleton hq] + simpa using unitInterval.pos_iff_ne_zero.mpr hq + +/-- If one of the two geometric laws is not a point mass, the sum charges every point. -/ +private lemma map_add_prod_geometricMeasure_singleton_pos + (h : (p₁ ≠ 0 ∧ p₁ ≠ 1) ∨ (p₂ ≠ 0 ∧ p₂ ≠ 1)) (s : ℕ) : + 0 < ((geometricMeasure p₁).prod (geometricMeasure p₂)).map (fun x ↦ x.1 + x.2) {s} := by + rw [map_add_singleton] + rcases h with ⟨h0, h1⟩ | ⟨h0, h1⟩ + · refine (lt_of_lt_of_le ?_ (single_le_sum (fun _ _ ↦ zero_le) (mem_antidiagonal.mpr + (add_zero s) : (s, 0) ∈ antidiagonal s))) + rw [← singleton_prod_singleton, Measure.prod_prod] + refine ENNReal.mul_pos ?_ (geometricMeasure_singleton_zero_pos p₂).ne' + rw [geometricMeasure_singleton h0] + exact (ENNReal.ofReal_pos.mpr (geometricMeasure_pos h0 h1 s)).ne' + · refine (lt_of_lt_of_le ?_ (single_le_sum (fun _ _ ↦ zero_le) (mem_antidiagonal.mpr + (zero_add s) : (0, s) ∈ antidiagonal s))) + rw [← singleton_prod_singleton, Measure.prod_prod] + refine ENNReal.mul_pos (geometricMeasure_singleton_zero_pos p₁).ne' ?_ + rw [geometricMeasure_singleton h0] + exact (ENNReal.ofReal_pos.mpr (geometricMeasure_pos h0 h1 s)).ne' + +/-- The sum of the masses of `Geom(p₁) ⊗ Geom(p₂)` on the pairs with sum `s`. -/ +private lemma sum_antidiagonal_geom {p₁ p₂ : ℝ} (h : p₁ ≠ p₂) (s : ℕ) : + ∑ y ∈ antidiagonal s, (1 - p₁) ^ y.1 * p₁ * ((1 - p₂) ^ y.2 * p₂) = + p₁ * p₂ * ((1 - p₁) ^ (s + 1) - (1 - p₂) ^ (s + 1)) / (p₂ - p₁) := by + rw [eq_div_iff (sub_ne_zero.mpr h.symm), Nat.sum_antidiagonal_eq_sum_range_succ_mk] + have hg := geom_sum₂_mul (1 - p₁) (1 - p₂) (s + 1) + simp only [add_tsub_cancel_right] at hg + rw [show p₂ - p₁ = (1 - p₁) - (1 - p₂) by ring, mul_assoc, ← hg] + simp only [sum_mul, mul_sum] + refine sum_congr rfl fun i _ ↦ ?_ + ring + +/-- The geometric laws, as a kernel. -/ +private noncomputable def geomKernel : Kernel unitInterval ℕ := + ⟨geometricMeasure, measurable_geometricMeasure⟩ + +/-- The geometric laws are probability measures. -/ +private instance : IsMarkovKernel geomKernel := ⟨fun _ ↦ isProbabilityMeasure_geometricMeasure⟩ + +/-- The laws `Geom(1 / m) ⊗ Geom(1 / m)`, as a kernel (the success probability is clipped to +`[0, 1]`). -/ +private noncomputable def geomPairInvKernel : Kernel ℝ (ℕ × ℕ) := + (geomKernel ×ₖ geomKernel).comap (fun m : ℝ ↦ projIcc 0 1 zero_le_one m⁻¹) (by fun_prop) + +/-- The mass of a point under the mixture of the laws `Geom(1 / m) ⊗ Geom(1 / m)` over `m` +uniform on `(a, b]`, for `1 ≤ a < b`. -/ +private lemma geomPairInvKernel_comp_singleton {a b : ℝ} (ha : 1 ≤ a) (hab : a < b) + (x : ℕ × ℕ) : + (geomPairInvKernel ∘ₘ uniformMeasure a b) {x} = + ENNReal.ofReal ((b - a)⁻¹ * (((1 - b⁻¹) ^ (x.1 + x.2 + 1) - (1 - a⁻¹) ^ (x.1 + x.2 + 1)) / + (((x.1 + x.2 : ℕ) : ℝ) + 1))) := by + set s := x.1 + x.2 + have hpos : ∀ m ∈ uIcc a b, 0 < m := fun m hm ↦ by + rw [uIcc_of_le hab.le] at hm + linarith [hm.1] + have hderiv : ∀ m ∈ uIcc a b, HasDerivAt (fun m : ℝ ↦ (1 - m⁻¹) ^ (s + 1) / (s + 1)) + ((1 - m⁻¹) ^ s * (m ^ 2)⁻¹) m := fun m hm ↦ + ((((hasDerivAt_inv (hpos m hm).ne').const_sub 1).pow (s + 1)).div_const _).congr_deriv + (by field_simp; push_cast; ring) + have hcont : ContinuousOn (fun m : ℝ ↦ (1 - m⁻¹) ^ s * (m ^ 2)⁻¹) (uIcc a b) := by + refine ContinuousOn.mul (by fun_prop (disch := exact fun m hm ↦ (hpos m hm).ne')) ?_ + fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (hpos m hm).ne') + have hmass : ∀ m ∈ Ioc a b, geomPairInvKernel m {x} = + ENNReal.ofReal ((1 - m⁻¹) ^ s * (m ^ 2)⁻¹) := by + intro m hm + have hm0 : 0 < m := by linarith [hm.1] + have hinv : m⁻¹ ∈ Icc (0 : ℝ) 1 := ⟨by positivity, inv_le_one_of_one_le₀ (by linarith [hm.1])⟩ + have hne : projIcc (0 : ℝ) 1 zero_le_one m⁻¹ ≠ 0 := by + rw [ne_eq, ← Subtype.coe_inj, projIcc_of_mem _ hinv] + exact (inv_pos.mpr hm0).ne' + rw [geomPairInvKernel, Kernel.comap_apply, Kernel.prod_apply] + change (geometricMeasure _).prod (geometricMeasure _) {x} = _ + rw [prod_geometricMeasure_singleton hne hne, projIcc_of_mem _ hinv] + congr 1 + simp only [s] + ring + rw [Measure.bind_apply (measurableSet_singleton x) (Kernel.aemeasurable _), + uniformMeasure_eq_smul, lintegral_smul_measure, setLIntegral_congr_fun measurableSet_Ioc hmass, + ← ofReal_integral_eq_lintegral_ofReal, ← intervalIntegral.integral_of_le hab.le, + intervalIntegral.integral_eq_sub_of_hasDerivAt hderiv hcont.intervalIntegrable, + smul_eq_mul, ← ENNReal.ofReal_inv_of_pos (sub_pos.mpr hab), + ← ENNReal.ofReal_mul (inv_nonneg.mpr (sub_pos.mpr hab).le), sub_div] + · rw [uIcc_of_le hab.le] at hcont + exact hcont.integrableOn_Icc.mono_set Ioc_subset_Icc_self + · refine ae_restrict_of_forall_mem measurableSet_Ioc fun m hm ↦ ?_ + have : 0 ≤ 1 - m⁻¹ := sub_nonneg.mpr (inv_le_one_of_one_le₀ (by linarith [hm.1])) + positivity + +/-- **The conditional law given the sum is a mixture of members of the null**: for distinct +positive success probabilities, the law with the law of the sum under `Geom(p₁) ⊗ Geom(p₂)` and +the conditional law `Kernel.antidiagonalUniform` given the sum is the mixture of the laws +`Geom(1 / m) ⊗ Geom(1 / m)` over `m` uniform between `1 / p₁` and `1 / p₂`. -/ +private lemma antidiagonalUniform_comp_map_add_prod_geometricMeasure (h₁ : p₁ ≠ 0) + (h₂ : p₂ ≠ 0) (h : p₁ ≠ p₂) : + Kernel.antidiagonalUniform ∘ₘ + ((geometricMeasure p₁).prod (geometricMeasure p₂)).map (fun x ↦ x.1 + x.2) = + geomPairInvKernel ∘ₘ uniformMeasure (min (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹) + (max (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹) := by + have h₁' : (0 : ℝ) < p₁ := unitInterval.pos_iff_ne_zero.mpr h₁ + have h₂' : (0 : ℝ) < p₂ := unitInterval.pos_iff_ne_zero.mpr h₂ + have hne : (p₁ : ℝ) ≠ p₂ := Subtype.coe_injective.ne h + have ha : 1 ≤ min (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ := + le_min (one_le_inv₀ h₁' |>.mpr p₁.2.2) (one_le_inv₀ h₂' |>.mpr p₂.2.2) + have hab : min (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ < max (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ := + min_lt_max.mpr (inv_injective.ne hne) + refine Measure.ext_iff_singleton.mpr fun x ↦ ?_ + have hc : (((x.1 + x.2 : ℕ) : ℝ≥0∞) + 1)⁻¹ = ENNReal.ofReal (((x.1 + x.2 : ℕ) : ℝ) + 1)⁻¹ := by + rw [ENNReal.ofReal_inv_of_pos (by positivity), ENNReal.ofReal_add (by positivity) zero_le_one, + ENNReal.ofReal_natCast, ENNReal.ofReal_one] + have hnn : ∀ y ∈ antidiagonal (x.1 + x.2), + 0 ≤ (1 - (p₁ : ℝ)) ^ y.1 * p₁ * ((1 - p₂) ^ y.2 * p₂) := fun y _ ↦ + mul_nonneg (geometricMeasure_nonneg _ _) (geometricMeasure_nonneg _ _) + rw [antidiagonalUniform_comp_singleton, map_add_singleton, + sum_congr rfl fun y _ ↦ prod_geometricMeasure_singleton h₁ h₂ y, + ← ENNReal.ofReal_sum_of_nonneg hnn, hc, ← ENNReal.ofReal_mul (sum_nonneg hnn), + sum_antidiagonal_geom hne, geomPairInvKernel_comp_singleton ha hab] + congr 1 + have h21 : (p₂ : ℝ) - p₁ ≠ 0 := sub_ne_zero.mpr hne.symm + rcases hne.lt_or_gt with h' | h' + · have hinv : (p₂ : ℝ)⁻¹ ≤ (p₁ : ℝ)⁻¹ := (inv_strictAnti₀ h₁' h').le + rw [min_eq_right hinv, max_eq_left hinv, inv_inv, inv_inv] + field_simp + · have hinv : (p₁ : ℝ)⁻¹ ≤ (p₂ : ℝ)⁻¹ := (inv_strictAnti₀ h₂' h').le + rw [min_eq_left hinv, max_eq_right hinv, inv_inv, inv_inv] + field_simp + ring + +/-- A measure on `ℕ` is absolutely continuous with respect to a measure that charges every +point. -/ +private lemma absolutelyContinuous_of_singleton_pos {μ ν : Measure ℕ} (h : ∀ s, 0 < ν {s}) : + μ ≪ ν := + Measure.AbsolutelyContinuous.mk fun t _ ht ↦ by + rcases t.eq_empty_or_nonempty with rfl | ⟨s, hs⟩ + · exact measure_empty + · exact absurd ht ((h s).trans_le (measure_mono (singleton_subset_iff.mpr hs))).ne' + +/-- The law with the law of the sum under `Geom(p₁) ⊗ Geom(p₂)` and the conditional law +`Kernel.antidiagonalUniform` given the sum is in the effective null. -/ +private lemma antidiagonalUniform_comp_mem_effectiveSet (p₁ p₂ : unitInterval) : + Kernel.antidiagonalUniform ∘ₘ + ((geometricMeasure p₁).prod (geometricMeasure p₂)).map (fun x ↦ x.1 + x.2) ∈ + effectiveSet (Set.range fun p : unitInterval ↦ + (geometricMeasure p).prod (geometricMeasure p)) := by + have hS : ∀ P ∈ Set.range fun p : unitInterval ↦ (geometricMeasure p).prod (geometricMeasure p), + P univ ≤ 1 := by + rintro _ ⟨p, rfl⟩ + exact prob_le_one + obtain ⟨q₁, hq₁, hq₁p⟩ := exists_ne_zero_geometricMeasure_eq p₁ + obtain ⟨q₂, hq₂, hq₂p⟩ := exists_ne_zero_geometricMeasure_eq p₂ + rw [← hq₁p, ← hq₂p] + rcases eq_or_ne q₁ q₂ with rfl | h + · obtain ⟨g, hg⟩ := exists_prod_geometricMeasure_self_singleton q₁ + rw [antidiagonalUniform_comp_map_add_of_singleton hg] + exact subset_effectiveSet hS ⟨q₁, rfl⟩ + have := isProbabilityMeasure_uniformMeasure + (min_lt_max.mpr (inv_injective.ne (Subtype.coe_injective.ne h)) : + min (q₁ : ℝ)⁻¹ (q₂ : ℝ)⁻¹ < max (q₁ : ℝ)⁻¹ (q₂ : ℝ)⁻¹) + rw [antidiagonalUniform_comp_map_add_prod_geometricMeasure hq₁ hq₂ h] + refine comp_mem_effectiveSet prob_le_one (ae_of_all _ fun m ↦ subset_effectiveSet hS ?_) + refine ⟨projIcc 0 1 zero_le_one m⁻¹, ?_⟩ + rw [geomPairInvKernel, Kernel.comap_apply, Kernel.prod_apply] + rfl + +end Geometric + +/-! ### The numeraire -/ + +/-- **The numeraire of two geometric observations against a common success probability.** If one +of the laws `Geom(p₁)` and `Geom(p₂)` is not a point mass, the likelihood ratio of +`Geom(p₁) ⊗ Geom(p₂)` to the law under which the sum of the two observations has its law under +`Geom(p₁) ⊗ Geom(p₂)` and, given the sum `s`, the first observation is uniform on `{0, …, s}`, is +the numeraire of the null of the products `Geom(p) ⊗ Geom(p)`: it is an e-variable, and no +e-variable has a larger expected logarithm under `Geom(p₁) ⊗ Geom(p₂)` +(`ProbabilityTheory.IsNumeraire.eintegral_log_le`). -/ +theorem isNumeraire_rnDeriv_antidiagonalUniform_comp {p₁ p₂ : unitInterval} + (h : (p₁ ≠ 0 ∧ p₁ ≠ 1) ∨ (p₂ ≠ 0 ∧ p₂ ≠ 1)) : + IsNumeraire (((geometricMeasure p₁).prod (geometricMeasure p₂)).rnDeriv + (Kernel.antidiagonalUniform ∘ₘ + ((geometricMeasure p₁).prod (geometricMeasure p₂)).map fun x ↦ x.1 + x.2)) + (Set.range fun p : unitInterval ↦ (geometricMeasure p).prod (geometricMeasure p)) + ((geometricMeasure p₁).prod (geometricMeasure p₂)) := by + refine isNumeraire_of_ae_eq_rnDeriv (isEVar_rnDeriv_comp_map (measurable_of_countable _) + (ae_of_all _ Kernel.ae_add_antidiagonalUniform) + (absolutelyContinuous_antidiagonalUniform_comp _) (by rintro _ ⟨p, rfl⟩; infer_instance) ?_ ?_) + (antidiagonalUniform_comp_mem_effectiveSet p₁ p₂) + (absolutelyContinuous_antidiagonalUniform_comp _) .rfl + · rintro _ ⟨p, rfl⟩ + obtain ⟨g, hg⟩ := exists_prod_geometricMeasure_self_singleton p + exact antidiagonalUniform_comp_map_add_of_singleton hg + · rintro _ ⟨p, rfl⟩ + exact absolutelyContinuous_of_singleton_pos (map_add_prod_geometricMeasure_singleton_pos h) + +end ProbabilityTheory From be4e40434b8969f665176f811ff24ea1be68d244 Mon Sep 17 00:00:00 2001 From: Remy Degenne Date: Sun, 15 Nov 2026 16:59:01 +0100 Subject: [PATCH 2/2] refactor: relocate the general lemmas of the two-sample numeraire and cite [EVR26] MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Move the general facts used by the two-sample instances out of the EValues files, as public lemmas: the law of a statistic under a composition with a kernel carried by its fibers (Probability/Kernel/Composition/MeasureComp), exponential laws and the sum of two exponential variables (Probability/Distributions/Exponential), the density of Gamma(2, r) (Probability/Distributions/Gamma), geometric laws with a positive success probability (Probability/Distributions/Geometric), integrating a derivative against the uniform law (Probability/Distributions/Uniform), measures through their masses on points (MeasureTheory/Measure/Singletons) and the geometric sum over an antidiagonal (Algebra/GeomSum). Reuse Tau Ceti's `geometricMeasure_one`, `geometricMeasure_zero` and `ae_mem_Ioc_uniformMeasure` instead of local copies. Expose the `Kernel.antidiagonalUniform ∘ₘ ν` API, tag the kernel evaluation lemmas `@[simp]`, replace the `change`/`rfl` unfoldings of the private kernels by apply lemmas, factor the min/max case splits through the symmetric uniform-law lemma, cite [EVR26, Note 28] for the instance results, and correct two docstrings. Co-Authored-By: Claude Code --- LMLExplore/Algebra/GeomSum.lean | 34 +++ LMLExplore/EValues/Numeraire/Conditional.lean | 36 +-- LMLExplore/EValues/TwoSample/Exponential.lean | 228 ++++-------------- LMLExplore/EValues/TwoSample/Geometric.lean | 212 +++++++--------- .../MeasureTheory/Measure/Singletons.lean | 52 ++++ .../Distributions/Exponential.lean | 173 +++++++++++++ .../Probability/Distributions/Gamma.lean | 33 +++ .../Probability/Distributions/Geometric.lean | 57 +++++ .../Probability/Distributions/Uniform.lean | 73 ++++++ .../Kernel/Composition/MeasureComp.lean | 50 ++++ 10 files changed, 611 insertions(+), 337 deletions(-) create mode 100644 LMLExplore/Algebra/GeomSum.lean create mode 100644 LMLExplore/MeasureTheory/Measure/Singletons.lean create mode 100644 LMLExplore/Probability/Distributions/Exponential.lean create mode 100644 LMLExplore/Probability/Distributions/Gamma.lean create mode 100644 LMLExplore/Probability/Distributions/Geometric.lean create mode 100644 LMLExplore/Probability/Distributions/Uniform.lean create mode 100644 LMLExplore/Probability/Kernel/Composition/MeasureComp.lean diff --git a/LMLExplore/Algebra/GeomSum.lean b/LMLExplore/Algebra/GeomSum.lean new file mode 100644 index 00000000..94e718d8 --- /dev/null +++ b/LMLExplore/Algebra/GeomSum.lean @@ -0,0 +1,34 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.Algebra.BigOperators.NatAntidiagonal +public import Mathlib.Algebra.Ring.GeomSum + +/-! +# Geometric sums over an antidiagonal + +The sum of `a ^ i * b ^ j` over the pairs `(i, j)` with `i + j = n` is the two-variable geometric +sum of Mathlib's `geom_sum₂_mul`: multiplied by `a - b`, it is `a ^ (n + 1) - b ^ (n + 1)`. + +## Main statements + +* `Finset.Nat.sum_antidiagonal_pow_mul_pow_mul_sub`: + `(∑ x ∈ antidiagonal n, a ^ x.1 * b ^ x.2) * (a - b) = a ^ (n + 1) - b ^ (n + 1)`. +-/ + +public section + +namespace Finset.Nat + +/-- **The geometric sum over an antidiagonal**, a form of `geom_sum₂_mul`: +`(∑ i + j = n, a ^ i * b ^ j) * (a - b) = a ^ (n + 1) - b ^ (n + 1)`. -/ +theorem sum_antidiagonal_pow_mul_pow_mul_sub {R : Type*} [CommRing R] (a b : R) (n : ℕ) : + (∑ x ∈ antidiagonal n, a ^ x.1 * b ^ x.2) * (a - b) = a ^ (n + 1) - b ^ (n + 1) := by + rw [sum_antidiagonal_eq_sum_range_succ_mk, ← geom_sum₂_mul] + simp + +end Finset.Nat diff --git a/LMLExplore/EValues/Numeraire/Conditional.lean b/LMLExplore/EValues/Numeraire/Conditional.lean index 77819749..491c5b39 100644 --- a/LMLExplore/EValues/Numeraire/Conditional.lean +++ b/LMLExplore/EValues/Numeraire/Conditional.lean @@ -17,18 +17,17 @@ null `S`, in the sense that `P = κ ∘ₘ P.map T`. For an alternative `Q`, the The likelihood ratio of `Q` to `κ ∘ₘ Q.map T` integrates to one against `κ y` for `Q.map T`-almost every `y`: its conditional expectation given `T` is one. So it is an e-variable for `S`, with -expectation exactly `1` under every member `P` of `S` such that the law of `T` under `P` is -absolutely continuous with respect to its law under `Q`. This is the conditional e-variable of -[HGLLA24]. When moreover `κ ∘ₘ Q.map T` belongs to the effective null of `S` (for example when it -is a mixture of members of `S`, by `ProbabilityTheory.comp_mem_effectiveSet`), this likelihood -ratio is the numeraire of `S` under `Q` by the criterion of [LRR25, Theorem 4.1] -(`ProbabilityTheory.isNumeraire_of_ae_eq_rnDeriv`), and `κ ∘ₘ Q.map T` is the reverse information -projection of `Q` on `S` when it is equivalent to `Q` -(`ProbabilityTheory.ripr_eq_of_isEVar_rnDeriv`). +expectation exactly the mass `P univ` (one for probability measures) under every member `P` of `S` +such that the law of `T` under `P` is absolutely continuous with respect to its law under `Q`. This +is the conditional e-variable of [HGLLA24]. When moreover `κ ∘ₘ Q.map T` belongs to the effective +null of `S` (for example when it is a mixture of members of `S`, by +`ProbabilityTheory.comp_mem_effectiveSet`), this likelihood ratio is the numeraire of `S` under +`Q` by the criterion of [LRR25, Theorem 4.1] (`ProbabilityTheory.isNumeraire_of_ae_eq_rnDeriv`), +and `κ ∘ₘ Q.map T` is the reverse information projection of `Q` on `S` when it is equivalent to +`Q` (`ProbabilityTheory.ripr_eq_of_isEVar_rnDeriv`). ## Main statements -* `ProbabilityTheory.map_comp_eq_of_ae_eq`: `κ ∘ₘ ν` has law `ν` for `T`. * `ProbabilityTheory.ae_lintegral_rnDeriv_comp_map_eq_one`: the likelihood ratio of `Q` to `κ ∘ₘ Q.map T` has conditional expectation one given `T`. * `ProbabilityTheory.lintegral_rnDeriv_comp_map`, `ProbabilityTheory.isEVar_rnDeriv_comp_map`: it @@ -53,25 +52,6 @@ namespace ProbabilityTheory variable {Ω 𝓨 : Type*} {mΩ : MeasurableSpace Ω} {m𝓨 : MeasurableSpace 𝓨} {T : Ω → 𝓨} {κ : Kernel 𝓨 Ω} {S : Set (Measure Ω)} {Q P : Measure Ω} -/-- If `κ y` is a probability measure carried by the fiber `T ⁻¹' {y}` for `ν`-almost every `y`, -the law of `T` under `κ ∘ₘ ν` is `ν`. -/ -theorem map_comp_eq_of_ae_eq [IsMarkovKernel κ] (hT : Measurable T) {ν : Measure 𝓨} - (hκ : ∀ᵐ y ∂ν, ∀ᵐ x ∂κ y, T x = y) : (κ ∘ₘ ν).map T = ν := by - ext A hA - rw [Measure.map_apply hT hA, Measure.bind_apply (hT hA) κ.aemeasurable, - ← lintegral_indicator_one hA] - refine lintegral_congr_ae ?_ - filter_upwards [hκ] with y hκ - by_cases hy : y ∈ A - · rw [indicator_of_mem hy, Pi.one_apply, ← measure_univ (μ := κ y)] - refine measure_congr ?_ - filter_upwards [hκ] with x hx - simp [hx, hy] - · rw [indicator_of_notMem hy] - refine measure_eq_zero_iff_ae_notMem.mpr ?_ - filter_upwards [hκ] with x hx - simpa [hx] using hy - /-- **The conditional expectation of the likelihood ratio given the statistic is one.** If `κ y` is carried by the fiber `T ⁻¹' {y}` for `Q.map T`-almost every `y` and the finite measure `Q` is absolutely continuous with respect to `κ ∘ₘ Q.map T`, the likelihood ratio of `Q` to diff --git a/LMLExplore/EValues/TwoSample/Exponential.lean b/LMLExplore/EValues/TwoSample/Exponential.lean index f64b11d8..556318c7 100644 --- a/LMLExplore/EValues/TwoSample/Exponential.lean +++ b/LMLExplore/EValues/TwoSample/Exponential.lean @@ -6,14 +6,14 @@ Authors: LMLExplore contributors module public import LMLExplore.EValues.Numeraire.Conditional -public import Mathlib.Probability.Distributions.Exponential +public import LMLExplore.Probability.Distributions.Exponential public import TauCeti.Probability.Distributions.Uniform.Basic -import TauCeti.MeasureTheory.Measure.Mixture.Basic +import LMLExplore.Probability.Distributions.Gamma +import LMLExplore.Probability.Kernel.Composition.MeasureComp import TauCeti.Probability.Distributions.Beta.Uniform import TauCeti.Probability.Distributions.Gamma.Beta import TauCeti.Probability.Distributions.Gamma.Measurability import TauCeti.Probability.Distributions.Gamma.PDF -import TauCeti.Probability.Distributions.Hypoexponential.PDF /-! # The numeraire of two exponential observations against a common rate @@ -37,6 +37,8 @@ which the sum has its law under `Q` and, given the sum `s`, the first observatio * no finite combination of members of the null is equal to it, against the conjecture of [HGLLA24, §5.2]. +These identifications, the mixture and the failure of finite combinations are [EVR26, Note 28]. + ## Main definitions * `ProbabilityTheory.Kernel.segmentUniform`: the law of `(u * s, (1 - u) * s)` for `u` uniform @@ -56,9 +58,11 @@ which the sum has its law under `Q` and, given the sum `s`, the first observatio That the law of the sum and the conditional law given the sum are as described under the null is the beta–gamma algebra of Tau Ceti (`TauCeti.Probability.map_div_add_prod_gammaMeasure`). The law -of the sum under `Q` is Tau Ceti's hypoexponential law, with density +of the sum under `Q` is Tau Ceti's hypoexponential law +(`ProbabilityTheory.map_add_prod_expMeasure_eq_withDensity`), with density `l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)`, which is also the density of the -mixture of the laws `Gamma(2, 1 / m)` over the uniform mean `m`. The numeraire is then identified +mixture of the laws `Gamma(2, 1 / m)` over the uniform mean `m` +(`ProbabilityTheory.map_add_prod_expMeasure_eq_bind`). The numeraire is then identified through the effective null: the likelihood ratio is an e-variable by `ProbabilityTheory.isEVar_rnDeriv_comp_map`, and `Qᶜ` is a mixture of members of the null, so that `ProbabilityTheory.isNumeraire_of_ae_eq_rnDeriv` and @@ -72,6 +76,8 @@ slowly or too fast at infinity compared to `exp (- min l₁ l₂ * x)`. with exponential families*, arXiv:2303.00471, §5.2. * [LRR25] M. Larsson, A. Ramdas and J. Ruf, *The numeraire e-variable and reverse information projection*, Ann. Statist. 53 (2025), arXiv:2402.18810v4, Theorem 4.1. +* [EVR26] *Answers to open questions in e-values and anytime-valid inference*, paper project + `RemyDegenne/evalues_research` (2026), Note 28. -/ public section @@ -96,6 +102,7 @@ instance Kernel.isMarkovKernel_segmentUniform : IsMarkovKernel Kernel.segmentUni exact Kernel.IsMarkovKernel.map _ (by fun_prop) /-- The uniform law on the segment of the pairs with sum `s`. -/ +@[simp] theorem Kernel.segmentUniform_apply (s : ℝ) : Kernel.segmentUniform s = (uniformMeasure 0 1).map fun u ↦ (u * s, (1 - u) * s) := by rw [Kernel.segmentUniform, Kernel.map_apply _ (by fun_prop), Kernel.prod_apply, Kernel.id_apply, @@ -124,34 +131,6 @@ section Null variable {l l₁ l₂ : ℝ} -/-- Exponential laws are s-finite. -/ -private instance sFinite_expMeasure (l : ℝ) : SFinite (expMeasure l) := - inferInstanceAs (SFinite (gammaMeasure 1 l)) - -/-- Exponential laws with positive rates are equivalent. -/ -private lemma expMeasure_absolutelyContinuous (l : ℝ) {l' : ℝ} (hl' : 0 < l') : - expMeasure l ≪ expMeasure l' := by - rw [expMeasure_eq_withDensity_restrict_Ioi l, expMeasure_eq_withDensity_restrict_Ioi l'] - refine (withDensity_absolutelyContinuous _ _).trans (withDensity_absolutelyContinuous' - (measurable_exponentialPDFReal l').ennreal_ofReal.aemeasurable ?_) - filter_upwards [ae_restrict_mem measurableSet_Ioi] with x hx - exact (ENNReal.ofReal_pos.mpr (exponentialPDFReal_pos hl' hx)).ne' - -/-- Products of exponential laws with positive rates are equivalent. -/ -private lemma prod_expMeasure_absolutelyContinuous (l₁ l₂ : ℝ) {l₁' l₂' : ℝ} (h₁ : 0 < l₁') - (h₂ : 0 < l₂') : - (expMeasure l₁).prod (expMeasure l₂) ≪ (expMeasure l₁').prod (expMeasure l₂') := - (expMeasure_absolutelyContinuous l₁ h₁).prod (expMeasure_absolutelyContinuous l₂ h₂) - -/-- The sum of two independent `Exp(l)` variables has the law `Gamma(2, l)`. -/ -private lemma map_add_prod_expMeasure_self (hl : 0 < l) : - ((expMeasure l).prod (expMeasure l)).map (fun x ↦ x.1 + x.2) = gammaMeasure 2 l := by - have : IsProbabilityMeasure (betaMeasure 1 1) := by rw [betaMeasure_one_one]; infer_instance - have h := congrArg (Measure.map Prod.snd) (map_div_add_prod_gammaMeasure one_pos one_pos hl) - rw [Measure.map_map measurable_snd (by fun_prop), Measure.map_snd_prod, measure_univ, - one_smul, one_add_one_eq_two] at h - exact h - /-- **The conditional law given the sum under the null**: two independent `Exp(l)` variables are `(u * s, (1 - u) * s)` for a sum `s` of law `Gamma(2, l)` and an independent uniform `u`. -/ private lemma segmentUniform_comp_gammaMeasure (hl : 0 < l) : @@ -216,136 +195,29 @@ private lemma map_add_absolutelyContinuous_of_mem (h₁ : 0 < l₁) (h₂ : 0 < /-! ### The law of the sum and the mixture -/ -/-- The density on `(0, ∞)` of the sum of two independent exponential variables with distinct -rates `l₁, l₂`. -/ -private noncomputable def sumDensity (l₁ l₂ x : ℝ) : ℝ := - l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁) - /-- The density of the sum is nonnegative on `[0, ∞)`. -/ private lemma sumDensity_nonneg (h₁ : 0 < l₁) (h₂ : 0 < l₂) {x : ℝ} (hx : 0 ≤ x) : - 0 ≤ sumDensity l₁ l₂ x := by - unfold sumDensity + 0 ≤ l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁) := by rcases le_total l₁ l₂ with h | h · exact div_nonneg (mul_nonneg (by positivity) (sub_nonneg.mpr (exp_le_exp.mpr (by nlinarith)))) (sub_nonneg.mpr h) · exact div_nonneg_of_nonpos (mul_nonpos_of_nonneg_of_nonpos (by positivity) (sub_nonpos.mpr (exp_le_exp.mpr (by nlinarith)))) (sub_nonpos.mpr h) -/-- **The law of the sum of two exponential variables** with distinct rates `l₁, l₂` has density -`l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)`. -/ -private lemma map_add_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : - ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = - (volume.restrict (Ioi 0)).withDensity fun x ↦ ENNReal.ofReal (sumDensity l₁ l₂ x) := by - have := isProbabilityMeasure_expMeasure h₁ - have := isProbabilityMeasure_expMeasure h₂ - have hr : ∀ i, 0 < ![l₁, l₂] i := by simp [Fin.forall_fin_two, h₁, h₂] - have hinj : Function.Injective ![l₁, l₂] := by - intro i j hij - fin_cases i <;> fin_cases j <;> simp_all [eq_comm] - have hpi : (fun i ↦ expMeasure (![l₁, l₂] i)) = ![expMeasure l₁, expMeasure l₂] := by - ext i : 1 - fin_cases i <;> rfl - have hmap : ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = - hypoexponentialMeasure ![l₁, l₂] := by - rw [hypoexponentialMeasure_eq_map hr, hpi, - ← (measurePreserving_finTwoArrow_vec (expMeasure l₁) (expMeasure l₂)).map_eq, - Measure.map_map (by fun_prop) (MeasurableEquiv.measurable _)] - congr 1 - ext x - simp [Fin.sum_univ_two, MeasurableEquiv.finTwoArrow] - rw [hmap, ← Measure.restrict_eq_self_of_ae_mem (s := Ioi 0) - (μ := hypoexponentialMeasure ![l₁, l₂]) (ae_pos_hypoexponentialMeasure _), - hypoexponentialMeasure_eq_withDensity, restrict_withDensity measurableSet_Ioi] - refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr - (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) - rw [hypoexponentialPDF_eq_ofReal, hypoexponentialPDFReal_of_injective hr hinj hx] - have h01 : (Finset.univ.erase (0 : Fin 2)) = {1} := by decide - have h10 : (Finset.univ.erase (1 : Fin 2)) = {0} := by decide - have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm - have h21 : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h - simp only [Fin.sum_univ_two, h01, h10, Finset.prod_singleton, Matrix.cons_val_zero, - Matrix.cons_val_one, sumDensity] - congr 1 - field_simp - ring - -/-- The density of `Gamma(2, r)` on `[0, ∞)`. -/ -private lemma gammaPDF_two {r x : ℝ} (hx : 0 ≤ x) : - gammaPDF 2 r x = ENNReal.ofReal (r ^ 2 * x * exp (-(r * x))) := by - rw [gammaPDF_of_nonneg hx, Real.Gamma_two, show (2 : ℝ) - 1 = 1 by norm_num, Real.rpow_one, - Real.rpow_two, div_one] - -/-- The uniform law on `(a, b]` lives on `(a, b]`. -/ -private lemma ae_mem_Ioc_uniformMeasure (a b : ℝ) : ∀ᵐ m ∂uniformMeasure a b, m ∈ Ioc a b := by - rw [uniformMeasure_eq_smul] - exact Measure.ae_smul_measure (ae_restrict_mem measurableSet_Ioc) _ - -/-- The density of the mixture of the laws `Gamma(2, 1 / m)` over the uniform mean `m`. -/ -private lemma lintegral_gammaPDF_uniformMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) - {x : ℝ} (hx : 0 < x) : - ∫⁻ m, gammaPDF 2 m⁻¹ x ∂uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹) = - ENNReal.ofReal (sumDensity l₁ l₂ x) := by - set α := min l₁⁻¹ l₂⁻¹ - set β := max l₁⁻¹ l₂⁻¹ - have hα : 0 < α := lt_min (inv_pos.mpr h₁) (inv_pos.mpr h₂) - have hαβ : α < β := min_lt_max.mpr (inv_injective.ne h) - have hpos : ∀ m ∈ uIcc α β, 0 < m := fun m hm ↦ hα.trans_le (by - rw [uIcc_of_le hαβ.le] at hm; exact hm.1) - have hderiv : ∀ m ∈ uIcc α β, HasDerivAt (fun m ↦ exp (-(m⁻¹ * x))) - (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) m := - fun m hm ↦ (((hasDerivAt_inv (hpos m hm).ne').mul_const x).neg.exp).congr_deriv - (by simp only [Pi.neg_apply]; ring) - have hcont : ContinuousOn (fun m ↦ exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) (uIcc α β) := by - refine ContinuousOn.mul (by fun_prop (disch := exact fun m hm ↦ (hpos m hm).ne')) ?_ - fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (hpos m hm).ne') - have hpdf : ∀ m ∈ Ioc α β, gammaPDF 2 m⁻¹ x = - ENNReal.ofReal (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) := by - intro m hm - rw [gammaPDF_two hx.le] - congr 1 - ring - rw [uniformMeasure_eq_smul, lintegral_smul_measure, setLIntegral_congr_fun measurableSet_Ioc hpdf, - ← ofReal_integral_eq_lintegral_ofReal, ← intervalIntegral.integral_of_le hαβ.le, - intervalIntegral.integral_eq_sub_of_hasDerivAt hderiv hcont.intervalIntegrable, - smul_eq_mul, ← ENNReal.ofReal_inv_of_pos (sub_pos.mpr hαβ), - ← ENNReal.ofReal_mul (inv_nonneg.mpr (sub_pos.mpr hαβ).le)] - · congr 1 - have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm - rcases h.lt_or_gt with h' | h' - · have hinv : l₂⁻¹ ≤ l₁⁻¹ := (inv_strictAnti₀ h₁ h').le - simp only [α, β, min_eq_right hinv, max_eq_left hinv, inv_inv, sumDensity] - field_simp - · have hinv : l₁⁻¹ ≤ l₂⁻¹ := (inv_strictAnti₀ h₂ h').le - simp only [α, β, min_eq_left hinv, max_eq_right hinv, inv_inv, sumDensity] - field_simp - ring - · rw [uIcc_of_le hαβ.le] at hcont - exact hcont.integrableOn_Icc.mono_set Ioc_subset_Icc_self - · refine ae_restrict_of_forall_mem measurableSet_Ioc fun m hm ↦ ?_ - have := hpos m (by rw [uIcc_of_le hαβ.le]; exact Ioc_subset_Icc_self hm) - positivity - -/-- **The sum of two exponential variables is a mixture of gamma laws**: for distinct rates, the -law of the sum is the mixture of the laws `Gamma(2, 1 / m)` over the mean `m` uniform between -`1 / l₁` and `1 / l₂`. -/ -private lemma map_add_prod_expMeasure_eq_bind (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : - ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = - (uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)).bind fun m ↦ gammaMeasure 2 m⁻¹ := by - have := isProbabilityMeasure_uniformMeasure - (min_lt_max.mpr (inv_injective.ne h) : min l₁⁻¹ l₂⁻¹ < max l₁⁻¹ l₂⁻¹) - simp_rw [gammaMeasure_eq_withDensity_restrict_Ioi] - have hf : Measurable (Function.uncurry fun m x ↦ gammaPDF 2 m⁻¹ x) := - measurable_uncurry_gammaPDF.comp (f := fun p : ℝ × ℝ ↦ (((2 : ℝ), p.1⁻¹), p.2)) (by fun_prop) - rw [map_add_prod_expMeasure h₁ h₂ h, Measure.bind_withDensity _ _ hf] - refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr - (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) - exact (lintegral_gammaPDF_uniformMeasure h₁ h₂ h hx).symm - /-- The laws `Gamma(2, 1 / m)`, as a kernel. -/ private noncomputable def gammaTwoInvKernel : Kernel ℝ ℝ where toFun m := gammaMeasure 2 m⁻¹ measurable' := measurable_gammaMeasure.comp (f := fun m : ℝ ↦ ((2 : ℝ), m⁻¹)) (by fun_prop) +/-- The kernel `gammaTwoInvKernel` maps `m` to `Gamma(2, 1 / m)`. -/ +private lemma gammaTwoInvKernel_apply (m : ℝ) : gammaTwoInvKernel m = gammaMeasure 2 m⁻¹ := rfl + +/-- The mixture of the laws `Gamma(2, 1 / m)` over a law of `m`, as a composition with +`gammaTwoInvKernel`. -/ +private lemma bind_gammaMeasure_eq_comp (μ : Measure ℝ) : + (μ.bind fun m ↦ gammaMeasure 2 m⁻¹) = gammaTwoInvKernel ∘ₘ μ := by + simp_rw [← gammaTwoInvKernel_apply] + /-- `Qᶜ` is in the effective null. -/ private lemma segmentUniform_comp_mem_effectiveSet (h₁ : 0 < l₁) (h₂ : 0 < l₂) : Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) ∈ @@ -358,15 +230,11 @@ private lemma segmentUniform_comp_mem_effectiveSet (h₁ : 0 < l₁) (h₂ : 0 < have hα : 0 < min l₁⁻¹ l₂⁻¹ := lt_min (inv_pos.mpr h₁) (inv_pos.mpr h₂) have := isProbabilityMeasure_uniformMeasure (min_lt_max.mpr (inv_injective.ne h) : min l₁⁻¹ l₂⁻¹ < max l₁⁻¹ l₂⁻¹) - rw [map_add_prod_expMeasure_eq_bind h₁ h₂ h] - change Kernel.segmentUniform ∘ₘ (gammaTwoInvKernel ∘ₘ _) ∈ _ - rw [Measure.comp_assoc] + rw [map_add_prod_expMeasure_eq_bind h₁ h₂ h, bind_gammaMeasure_eq_comp, Measure.comp_assoc] refine comp_mem_effectiveSet prob_le_one ?_ - filter_upwards [ae_mem_Ioc_uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)] with m hm + filter_upwards [ae_mem_Ioc_uniformMeasure] with m hm have hm0 : 0 < m⁻¹ := inv_pos.mpr (hα.trans hm.1) - rw [Kernel.comp_apply] - change Kernel.segmentUniform ∘ₘ gammaMeasure 2 m⁻¹ ∈ _ - rw [segmentUniform_comp_gammaMeasure hm0] + rw [Kernel.comp_apply, gammaTwoInvKernel_apply, segmentUniform_comp_gammaMeasure hm0] exact subset_effectiveSet hS ⟨⟨m⁻¹, hm0⟩, rfl⟩ /-- The likelihood ratio of `Exp(l₁) ⊗ Exp(l₂)` to its `Qᶜ` is an e-variable for the null. -/ @@ -395,7 +263,7 @@ variable {l₁ l₂ : ℝ} observations has its law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation is uniform on `[0, s]`, is the numeraire of the null of the products `Exp(l) ⊗ Exp(l)`, `l > 0`: it is an e-variable, and no e-variable has a larger expected logarithm under -`Exp(l₁) ⊗ Exp(l₂)` (`ProbabilityTheory.IsNumeraire.eintegral_log_le`). -/ +`Exp(l₁) ⊗ Exp(l₂)` (`ProbabilityTheory.IsNumeraire.eintegral_log_le`) [EVR26, Note 28]. -/ theorem isNumeraire_rnDeriv_segmentUniform_comp (h₁ : 0 < l₁) (h₂ : 0 < l₂) : IsNumeraire (((expMeasure l₁).prod (expMeasure l₂)).rnDeriv (Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map fun x ↦ x.1 + x.2)) @@ -410,7 +278,7 @@ theorem isNumeraire_rnDeriv_segmentUniform_comp (h₁ : 0 < l₁) (h₂ : 0 < l /-- **The reverse information projection of two exponential observations** on the null of the products `Exp(l) ⊗ Exp(l)`, `l > 0`, is the law under which the sum of the two observations has its law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation is uniform on -`[0, s]`. -/ +`[0, s]` [EVR26, Note 28]. -/ theorem ripr_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) : ripr ((expMeasure l₁).prod (expMeasure l₂)) (Set.range fun l : Ioi (0 : ℝ) ↦ (expMeasure l).prod (expMeasure l)) = @@ -426,7 +294,7 @@ theorem ripr_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) : distinct positive rates `l₁, l₂`, the law under which the sum of the two observations has its law under `Exp(l₁) ⊗ Exp(l₂)` and, given the sum `s`, the first observation is uniform on `[0, s]`, is the mixture of the members `Exp(1 / m) ⊗ Exp(1 / m)` of the null whose mean `m` is -uniform between `1 / l₁` and `1 / l₂`. -/ +uniform between `1 / l₁` and `1 / l₂` [EVR26, Note 28]. -/ theorem segmentUniform_comp_map_add_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : Kernel.segmentUniform ∘ₘ ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = (uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)).bind @@ -435,7 +303,7 @@ theorem segmentUniform_comp_map_add_prod_expMeasure (h₁ : 0 < l₁) (h₂ : 0 rw [map_add_prod_expMeasure_eq_bind h₁ h₂ h, Measure.bind_bind (f := fun m ↦ gammaMeasure 2 m⁻¹) gammaTwoInvKernel.measurable.aemeasurable Kernel.segmentUniform.aemeasurable] refine Measure.bind_congr_right ?_ - filter_upwards [ae_mem_Ioc_uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)] with m hm + filter_upwards [ae_mem_Ioc_uniformMeasure] with m hm exact segmentUniform_comp_gammaMeasure (inv_pos.mpr (hα.trans hm.1)) /-- A nonnegative finite combination of the functions `x ↦ l ^ 2 * x * exp (- l x)`, `l > 0`, is @@ -507,7 +375,7 @@ private lemma not_forall_sum_eq_mul_sub {a b C : ℝ} (hab : a < b) (hC : 0 < C) distinct positive rates `l₁, l₂`, no finite combination of the laws `Exp(l) ⊗ Exp(l)` is the reverse information projection of `Exp(l₁) ⊗ Exp(l₂)` on the null of the products `Exp(l) ⊗ Exp(l)`, `l > 0`. In particular it is not a mixture of two members of the null, against -the conjecture of [HGLLA24, §5.2]. -/ +the conjecture of [HGLLA24, §5.2]; this refutation is [EVR26, Note 28]. -/ theorem sum_smul_prod_expMeasure_ne_ripr (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) (s : Finset ℝ) (w : ℝ → ℝ≥0∞) : ∑ l ∈ s, w l • (expMeasure l).prod (expMeasure l) ≠ @@ -542,19 +410,19 @@ theorem sum_smul_prod_expMeasure_ne_ripr (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h (hmass ▸ Finset.single_le_sum (fun _ _ ↦ zero_le) hl) -- and the densities of the laws of the sum agree on `(0, ∞)` have hae : (fun x ↦ ∑ l ∈ s', w l * gammaPDF 2 l x) =ᵐ[volume.restrict (Ioi 0)] - fun x ↦ ENNReal.ofReal (sumDensity l₁ l₂ x) := by + fun x ↦ ENNReal.ofReal (l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁)) := by refine ae_eq_of_forall_setLIntegral_eq_of_sigmaFinite (Finset.measurable_fun_sum _ fun l _ ↦ (measurable_gammaPDF 2 l).const_mul _) - (by unfold sumDensity; fun_prop) fun A hA _ ↦ ?_ + (by fun_prop) fun A hA _ ↦ ?_ rw [lintegral_finsetSum _ fun l _ ↦ (measurable_gammaPDF 2 l).const_mul _, - ← withDensity_apply _ hA, ← map_add_prod_expMeasure h₁ h₂ h, ← hsum A hA] + ← withDensity_apply _ hA, ← map_add_prod_expMeasure_eq_withDensity h₁ h₂ h, ← hsum A hA] refine Finset.sum_congr rfl fun l _ ↦ ?_ rw [lintegral_const_mul _ (measurable_gammaPDF 2 l), ← withDensity_apply _ hA, ← gammaMeasure_eq_withDensity_restrict_Ioi] have hreal : ∀ x ∈ Ioi (0 : ℝ), - ∑ l ∈ s', (w l).toReal * (l ^ 2 * x * exp (-(l * x))) = sumDensity l₁ l₂ x := by - refine Measure.eqOn_open_of_ae_eq (μ := volume) ?_ isOpen_Ioi (by fun_prop) - (by unfold sumDensity; fun_prop) + ∑ l ∈ s', (w l).toReal * (l ^ 2 * x * exp (-(l * x))) = + l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁) := by + refine Measure.eqOn_open_of_ae_eq (μ := volume) ?_ isOpen_Ioi (by fun_prop) (by fun_prop) filter_upwards [hae, ae_restrict_mem measurableSet_Ioi] with x hx (hx0 : 0 < x) rw [← ENNReal.ofReal_eq_ofReal_iff (Finset.sum_nonneg fun l hl ↦ have := hs'pos l hl; by positivity) (sumDensity_nonneg h₁ h₂ hx0.le), ← hx, @@ -562,18 +430,20 @@ theorem sum_smul_prod_expMeasure_ne_ripr (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h refine Finset.sum_congr rfl fun l hl ↦ ?_ rw [ENNReal.ofReal_mul ENNReal.toReal_nonneg, ENNReal.ofReal_toReal (hw l hl), gammaPDF_two hx0.le] - rcases h.lt_or_gt with h12 | h12 - · refine not_forall_sum_eq_mul_sub h12 (div_pos (mul_pos h₁ h₂) (sub_pos.mpr h12)) s' - (fun l ↦ (w l).toReal) (fun _ _ ↦ ENNReal.toReal_nonneg) hs'pos fun x hx ↦ ?_ - rw [hreal x hx, sumDensity] - ring - · refine not_forall_sum_eq_mul_sub h12 (div_pos (mul_pos h₁ h₂) (sub_pos.mpr h12)) s' - (fun l ↦ (w l).toReal) (fun _ _ ↦ ENNReal.toReal_nonneg) hs'pos fun x hx ↦ ?_ + -- the density of the law of the sum is symmetric in the two rates: order them + obtain ⟨a, b, hab, ha, hdens⟩ : ∃ a b, a < b ∧ 0 < a ∧ ∀ x, + l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁) = + a * b / (b - a) * (exp (-(a * x)) - exp (-(b * x))) := by + have h12 : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h have h21 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm - have h12' : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h - rw [hreal x hx, sumDensity] - field_simp - ring + rcases h.lt_or_gt with h' | h' + · exact ⟨l₁, l₂, h', h₁, fun x ↦ by ring⟩ + · refine ⟨l₂, l₁, h', h₂, fun x ↦ ?_⟩ + field_simp + ring + exact not_forall_sum_eq_mul_sub hab (div_pos (mul_pos ha (ha.trans hab)) (sub_pos.mpr hab)) s' + (fun l ↦ (w l).toReal) (fun _ _ ↦ ENNReal.toReal_nonneg) hs'pos fun x hx ↦ by + rw [hreal x hx, hdens] end Main diff --git a/LMLExplore/EValues/TwoSample/Geometric.lean b/LMLExplore/EValues/TwoSample/Geometric.lean index cf611578..ca077712 100644 --- a/LMLExplore/EValues/TwoSample/Geometric.lean +++ b/LMLExplore/EValues/TwoSample/Geometric.lean @@ -8,6 +8,10 @@ module public import LMLExplore.EValues.Numeraire.Conditional public import Mathlib.Probability.Distributions.Geometric public import Mathlib.Probability.UniformOn +import LMLExplore.Algebra.GeomSum +import LMLExplore.MeasureTheory.Measure.Singletons +import LMLExplore.Probability.Distributions.Geometric +import LMLExplore.Probability.Distributions.Uniform import TauCeti.Probability.Distributions.Geometric.Measurability import TauCeti.Probability.Distributions.Uniform.Basic @@ -24,11 +28,12 @@ every member of the null the mass of a pair `(j, k)` depends only on `j + k`, so For the alternative `Q = Geom(p₁) ⊗ Geom(p₂)`, the likelihood ratio of `Q` to the law `Kernel.antidiagonalUniform ∘ₘ Q.map (x.1 + x.2)`, under which the sum has its law under `Q` and, given the sum `s`, the first observation is uniform on `{0, …, s}`, is the numeraire of the null -under `Q` (`ProbabilityTheory.isNumeraire_rnDeriv_antidiagonalUniform_comp`). As for two -exponential observations (`ProbabilityTheory.isNumeraire_rnDeriv_segmentUniform_comp`), that law -is the mixture of the members `Geom(1 / m) ⊗ Geom(1 / m)` of the null over `m` uniform between -`1 / p₁` and `1 / p₂`, that is, over the mean `m - 1` uniform between the two means, so that it -belongs to the effective null. +under `Q` (`ProbabilityTheory.isNumeraire_rnDeriv_antidiagonalUniform_comp`), as computed in +[EVR26, Note 28]. As for two exponential observations +(`ProbabilityTheory.isNumeraire_rnDeriv_segmentUniform_comp`), for distinct positive success +probabilities `p₁ ≠ p₂` that law is the mixture of the members `Geom(1 / m) ⊗ Geom(1 / m)` of the +null over `m` uniform between `1 / p₁` and `1 / p₂`, that is, over the mean `m - 1` uniform between +the two means; when `p₁ = p₂` it is `Q` itself. In both cases it belongs to the effective null. ## Main definitions @@ -37,8 +42,12 @@ belongs to the effective null. ## Main statements +* `ProbabilityTheory.Kernel.antidiagonalUniform_comp_singleton`, + `ProbabilityTheory.Kernel.antidiagonalUniform_comp_map_add_of_singleton`: the masses of the + points under `Kernel.antidiagonalUniform ∘ₘ ν`, and a law on pairs whose masses depend only on + the sum has the conditional law `Kernel.antidiagonalUniform` given the sum. * `ProbabilityTheory.isNumeraire_rnDeriv_antidiagonalUniform_comp`: the numeraire of two geometric - observations against a common success probability. + observations against a common success probability [EVR26, Note 28]. ## Implementation notes @@ -53,6 +62,8 @@ mass at `(0, 0)`. The mixture is computed on the masses of the points, through t with exponential families*, arXiv:2303.00471. * [LRR25] M. Larsson, A. Ramdas and J. Ruf, *The numeraire e-variable and reverse information projection*, Ann. Statist. 53 (2025), arXiv:2402.18810v4, Theorem 4.1. +* [EVR26] *Answers to open questions in e-values and anytime-valid inference*, paper project + `RemyDegenne/evalues_research` (2026), Note 28. -/ public section @@ -71,6 +82,7 @@ noncomputable def Kernel.antidiagonalUniform : Kernel ℕ (ℕ × ℕ) := Kernel.ofFunOfCountable fun s ↦ uniformOn (antidiagonal s : Set (ℕ × ℕ)) /-- The uniform law on the pairs with sum `s`. -/ +@[simp] theorem Kernel.antidiagonalUniform_apply (s : ℕ) : Kernel.antidiagonalUniform s = uniformOn (antidiagonal s : Set (ℕ × ℕ)) := by unfold Kernel.antidiagonalUniform @@ -108,17 +120,8 @@ section Pairs variable {μ : Measure (ℕ × ℕ)} -/-- The law of the sum at a point is the sum of the masses of the pairs with that sum. -/ -private lemma map_add_singleton (μ : Measure (ℕ × ℕ)) (s : ℕ) : - μ.map (fun x ↦ x.1 + x.2) {s} = ∑ x ∈ antidiagonal s, μ {x} := by - rw [Measure.map_apply (measurable_of_countable _) (measurableSet_singleton s), - sum_measure_singleton] - congr 1 - ext x - simp - /-- The mass of a point under `Kernel.antidiagonalUniform ∘ₘ ν`. -/ -private lemma antidiagonalUniform_comp_singleton (ν : Measure ℕ) (x : ℕ × ℕ) : +theorem Kernel.antidiagonalUniform_comp_singleton (ν : Measure ℕ) (x : ℕ × ℕ) : (Kernel.antidiagonalUniform ∘ₘ ν) {x} = ν {x.1 + x.2} * (((x.1 + x.2 : ℕ) : ℝ≥0∞) + 1)⁻¹ := by rw [Measure.bind_apply (measurableSet_singleton x) (Kernel.aemeasurable _), lintegral_countable'] @@ -128,26 +131,26 @@ private lemma antidiagonalUniform_comp_singleton (ν : Measure ℕ) (x : ℕ × /-- A law on pairs whose mass at a point depends only on the sum has the conditional law `Kernel.antidiagonalUniform` given the sum. -/ -private lemma antidiagonalUniform_comp_map_add_of_singleton {g : ℕ → ℝ≥0∞} +theorem Kernel.antidiagonalUniform_comp_map_add_of_singleton {g : ℕ → ℝ≥0∞} (h : ∀ x, μ {x} = g (x.1 + x.2)) : Kernel.antidiagonalUniform ∘ₘ μ.map (fun x ↦ x.1 + x.2) = μ := by refine Measure.ext_iff_singleton.mpr fun x ↦ ?_ - rw [antidiagonalUniform_comp_singleton, map_add_singleton, - sum_congr rfl fun y hy ↦ by rw [h y, mem_antidiagonal.mp hy], sum_const, + rw [Kernel.antidiagonalUniform_comp_singleton, Measure.map_add_singleton_eq_sum_antidiagonal, + Finset.sum_congr rfl fun y hy ↦ by rw [h y, mem_antidiagonal.mp hy], Finset.sum_const, Nat.card_antidiagonal, nsmul_eq_mul, h x, mul_comm _ (g _), mul_assoc] push_cast rw [ENNReal.mul_inv_cancel (by positivity) (by simp), mul_one] /-- Every law on pairs is absolutely continuous with respect to the law with its law of the sum and the conditional law `Kernel.antidiagonalUniform` given the sum. -/ -private lemma absolutelyContinuous_antidiagonalUniform_comp (μ : Measure (ℕ × ℕ)) : +theorem Kernel.absolutelyContinuous_antidiagonalUniform_comp (μ : Measure (ℕ × ℕ)) : μ ≪ Kernel.antidiagonalUniform ∘ₘ μ.map (fun x ↦ x.1 + x.2) := by refine Measure.AbsolutelyContinuous.mk fun t _ ht ↦ ?_ rw [← biUnion_of_singleton t, measure_biUnion_null_iff t.to_countable] at ht ⊢ intro x hx have hx' := ht x hx - rw [antidiagonalUniform_comp_singleton, mul_eq_zero, map_add_singleton, - Finset.sum_eq_zero_iff] at hx' + rw [Kernel.antidiagonalUniform_comp_singleton, mul_eq_zero, + Measure.map_add_singleton_eq_sum_antidiagonal, Finset.sum_eq_zero_iff] at hx' exact hx'.resolve_right (by simp) x (mem_antidiagonal.mpr rfl) end Pairs @@ -158,29 +161,6 @@ section Geometric variable {p₁ p₂ : unitInterval} -/-- The geometric law with success probability `1` is the point mass at `0`, as is the junk value -`geometricMeasure 0`. -/ -private lemma geometricMeasure_one : geometricMeasure 1 = geometricMeasure 0 := by - have h0 : geometricMeasure 0 = Measure.dirac 0 := by simp [geometricMeasure] - refine Measure.ext_iff_singleton.mpr fun n ↦ ?_ - rw [geometricMeasure_singleton one_ne_zero n, h0] - rcases n with _ | n <;> simp - -/-- Every geometric law is one with a positive success probability. -/ -private lemma exists_ne_zero_geometricMeasure_eq (p : unitInterval) : - ∃ q : unitInterval, q ≠ 0 ∧ geometricMeasure q = geometricMeasure p := by - by_cases hp : p = 0 - · exact ⟨1, one_ne_zero, hp ▸ geometricMeasure_one⟩ - · exact ⟨p, hp, rfl⟩ - -/-- The mass of a point under a product of geometric laws with positive success -probabilities. -/ -private lemma prod_geometricMeasure_singleton (h₁ : p₁ ≠ 0) (h₂ : p₂ ≠ 0) (x : ℕ × ℕ) : - (geometricMeasure p₁).prod (geometricMeasure p₂) {x} = - ENNReal.ofReal ((1 - p₁) ^ x.1 * p₁ * ((1 - p₂) ^ x.2 * p₂)) := by - rw [← singleton_prod_singleton, Measure.prod_prod, geometricMeasure_singleton h₁, - geometricMeasure_singleton h₂, ← ENNReal.ofReal_mul (geometricMeasure_nonneg _ _)] - /-- The mass of a point under `Geom(p) ⊗ Geom(p)` depends only on the sum. -/ private lemma exists_prod_geometricMeasure_self_singleton (p : unitInterval) : ∃ g : ℕ → ℝ≥0∞, ∀ x : ℕ × ℕ, @@ -191,18 +171,11 @@ private lemma exists_prod_geometricMeasure_self_singleton (p : unitInterval) : congr 1 ring -/-- The point `0` has positive mass under every geometric law. -/ -private lemma geometricMeasure_singleton_zero_pos (p : unitInterval) : - 0 < geometricMeasure p {0} := by - obtain ⟨q, hq, hqp⟩ := exists_ne_zero_geometricMeasure_eq p - rw [← hqp, geometricMeasure_singleton hq] - simpa using unitInterval.pos_iff_ne_zero.mpr hq - /-- If one of the two geometric laws is not a point mass, the sum charges every point. -/ private lemma map_add_prod_geometricMeasure_singleton_pos (h : (p₁ ≠ 0 ∧ p₁ ≠ 1) ∨ (p₂ ≠ 0 ∧ p₂ ≠ 1)) (s : ℕ) : 0 < ((geometricMeasure p₁).prod (geometricMeasure p₂)).map (fun x ↦ x.1 + x.2) {s} := by - rw [map_add_singleton] + rw [Measure.map_add_singleton_eq_sum_antidiagonal] rcases h with ⟨h0, h1⟩ | ⟨h0, h1⟩ · refine (lt_of_lt_of_le ?_ (single_le_sum (fun _ _ ↦ zero_le) (mem_antidiagonal.mpr (add_zero s) : (s, 0) ∈ antidiagonal s))) @@ -217,18 +190,6 @@ private lemma map_add_prod_geometricMeasure_singleton_pos rw [geometricMeasure_singleton h0] exact (ENNReal.ofReal_pos.mpr (geometricMeasure_pos h0 h1 s)).ne' -/-- The sum of the masses of `Geom(p₁) ⊗ Geom(p₂)` on the pairs with sum `s`. -/ -private lemma sum_antidiagonal_geom {p₁ p₂ : ℝ} (h : p₁ ≠ p₂) (s : ℕ) : - ∑ y ∈ antidiagonal s, (1 - p₁) ^ y.1 * p₁ * ((1 - p₂) ^ y.2 * p₂) = - p₁ * p₂ * ((1 - p₁) ^ (s + 1) - (1 - p₂) ^ (s + 1)) / (p₂ - p₁) := by - rw [eq_div_iff (sub_ne_zero.mpr h.symm), Nat.sum_antidiagonal_eq_sum_range_succ_mk] - have hg := geom_sum₂_mul (1 - p₁) (1 - p₂) (s + 1) - simp only [add_tsub_cancel_right] at hg - rw [show p₂ - p₁ = (1 - p₁) - (1 - p₂) by ring, mul_assoc, ← hg] - simp only [sum_mul, mul_sum] - refine sum_congr rfl fun i _ ↦ ?_ - ring - /-- The geometric laws, as a kernel. -/ private noncomputable def geomKernel : Kernel unitInterval ℕ := ⟨geometricMeasure, measurable_geometricMeasure⟩ @@ -236,54 +197,55 @@ private noncomputable def geomKernel : Kernel unitInterval ℕ := /-- The geometric laws are probability measures. -/ private instance : IsMarkovKernel geomKernel := ⟨fun _ ↦ isProbabilityMeasure_geometricMeasure⟩ +/-- The geometric kernel maps `p` to `Geom(p)`. -/ +private lemma geomKernel_apply (p : unitInterval) : geomKernel p = geometricMeasure p := rfl + /-- The laws `Geom(1 / m) ⊗ Geom(1 / m)`, as a kernel (the success probability is clipped to `[0, 1]`). -/ private noncomputable def geomPairInvKernel : Kernel ℝ (ℕ × ℕ) := (geomKernel ×ₖ geomKernel).comap (fun m : ℝ ↦ projIcc 0 1 zero_le_one m⁻¹) (by fun_prop) +/-- The kernel `geomPairInvKernel` maps `m` to `Geom(1 / m) ⊗ Geom(1 / m)`, with the success +probability clipped to `[0, 1]`. -/ +private lemma geomPairInvKernel_apply (m : ℝ) : + geomPairInvKernel m = (geometricMeasure (projIcc 0 1 zero_le_one m⁻¹)).prod + (geometricMeasure (projIcc 0 1 zero_le_one m⁻¹)) := by + rw [geomPairInvKernel, Kernel.comap_apply, Kernel.prod_apply, geomKernel_apply] + /-- The mass of a point under the mixture of the laws `Geom(1 / m) ⊗ Geom(1 / m)` over `m` -uniform on `(a, b]`, for `1 ≤ a < b`. -/ -private lemma geomPairInvKernel_comp_singleton {a b : ℝ} (ha : 1 ≤ a) (hab : a < b) +uniform between two distinct `a, b ≥ 1`. -/ +private lemma geomPairInvKernel_comp_singleton {a b : ℝ} (ha : 1 ≤ a) (hb : 1 ≤ b) (hab : a ≠ b) (x : ℕ × ℕ) : - (geomPairInvKernel ∘ₘ uniformMeasure a b) {x} = + (geomPairInvKernel ∘ₘ uniformMeasure (min a b) (max a b)) {x} = ENNReal.ofReal ((b - a)⁻¹ * (((1 - b⁻¹) ^ (x.1 + x.2 + 1) - (1 - a⁻¹) ^ (x.1 + x.2 + 1)) / (((x.1 + x.2 : ℕ) : ℝ) + 1))) := by set s := x.1 + x.2 - have hpos : ∀ m ∈ uIcc a b, 0 < m := fun m hm ↦ by - rw [uIcc_of_le hab.le] at hm - linarith [hm.1] + have hone : ∀ m ∈ uIcc a b, 1 ≤ m := fun m hm ↦ (le_inf ha hb).trans hm.1 have hderiv : ∀ m ∈ uIcc a b, HasDerivAt (fun m : ℝ ↦ (1 - m⁻¹) ^ (s + 1) / (s + 1)) ((1 - m⁻¹) ^ s * (m ^ 2)⁻¹) m := fun m hm ↦ - ((((hasDerivAt_inv (hpos m hm).ne').const_sub 1).pow (s + 1)).div_const _).congr_deriv - (by field_simp; push_cast; ring) + ((((hasDerivAt_inv (zero_lt_one.trans_le (hone m hm)).ne').const_sub 1).pow (s + 1)).div_const + _).congr_deriv (by field_simp; push_cast; ring) have hcont : ContinuousOn (fun m : ℝ ↦ (1 - m⁻¹) ^ s * (m ^ 2)⁻¹) (uIcc a b) := by - refine ContinuousOn.mul (by fun_prop (disch := exact fun m hm ↦ (hpos m hm).ne')) ?_ - fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (hpos m hm).ne') - have hmass : ∀ m ∈ Ioc a b, geomPairInvKernel m {x} = + refine ContinuousOn.mul + (by fun_prop (disch := exact fun m hm ↦ (zero_lt_one.trans_le (hone m hm)).ne')) ?_ + fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (zero_lt_one.trans_le (hone m hm)).ne') + have hmass : ∀ᵐ m ∂uniformMeasure (min a b) (max a b), geomPairInvKernel m {x} = ENNReal.ofReal ((1 - m⁻¹) ^ s * (m ^ 2)⁻¹) := by - intro m hm - have hm0 : 0 < m := by linarith [hm.1] - have hinv : m⁻¹ ∈ Icc (0 : ℝ) 1 := ⟨by positivity, inv_le_one_of_one_le₀ (by linarith [hm.1])⟩ + filter_upwards [ae_mem_Ioc_uniformMeasure] with m hm + have hm1 : 1 ≤ m := (le_min ha hb).trans hm.1.le + have hinv : m⁻¹ ∈ Icc (0 : ℝ) 1 := ⟨by positivity, inv_le_one_of_one_le₀ hm1⟩ have hne : projIcc (0 : ℝ) 1 zero_le_one m⁻¹ ≠ 0 := by rw [ne_eq, ← Subtype.coe_inj, projIcc_of_mem _ hinv] - exact (inv_pos.mpr hm0).ne' - rw [geomPairInvKernel, Kernel.comap_apply, Kernel.prod_apply] - change (geometricMeasure _).prod (geometricMeasure _) {x} = _ - rw [prod_geometricMeasure_singleton hne hne, projIcc_of_mem _ hinv] + exact (inv_pos.mpr (zero_lt_one.trans_le hm1)).ne' + rw [geomPairInvKernel_apply, prod_geometricMeasure_singleton hne hne, projIcc_of_mem _ hinv] congr 1 simp only [s] ring rw [Measure.bind_apply (measurableSet_singleton x) (Kernel.aemeasurable _), - uniformMeasure_eq_smul, lintegral_smul_measure, setLIntegral_congr_fun measurableSet_Ioc hmass, - ← ofReal_integral_eq_lintegral_ofReal, ← intervalIntegral.integral_of_le hab.le, - intervalIntegral.integral_eq_sub_of_hasDerivAt hderiv hcont.intervalIntegrable, - smul_eq_mul, ← ENNReal.ofReal_inv_of_pos (sub_pos.mpr hab), - ← ENNReal.ofReal_mul (inv_nonneg.mpr (sub_pos.mpr hab).le), sub_div] - · rw [uIcc_of_le hab.le] at hcont - exact hcont.integrableOn_Icc.mono_set Ioc_subset_Icc_self - · refine ae_restrict_of_forall_mem measurableSet_Ioc fun m hm ↦ ?_ - have : 0 ≤ 1 - m⁻¹ := sub_nonneg.mpr (inv_le_one_of_one_le₀ (by linarith [hm.1])) - positivity + lintegral_congr_ae hmass, + lintegral_ofReal_uniformMeasure_min_max_of_hasDerivAt hab hderiv hcont fun m hm ↦ ?_, sub_div] + have : 0 ≤ 1 - m⁻¹ := sub_nonneg.mpr (inv_le_one_of_one_le₀ (hone m (uIoc_subset_uIcc hm))) + positivity /-- **The conditional law given the sum is a mixture of members of the null**: for distinct positive success probabilities, the law with the law of the sum under `Geom(p₁) ⊗ Geom(p₂)` and @@ -298,40 +260,30 @@ private lemma antidiagonalUniform_comp_map_add_prod_geometricMeasure (h₁ : p have h₁' : (0 : ℝ) < p₁ := unitInterval.pos_iff_ne_zero.mpr h₁ have h₂' : (0 : ℝ) < p₂ := unitInterval.pos_iff_ne_zero.mpr h₂ have hne : (p₁ : ℝ) ≠ p₂ := Subtype.coe_injective.ne h - have ha : 1 ≤ min (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ := - le_min (one_le_inv₀ h₁' |>.mpr p₁.2.2) (one_le_inv₀ h₂' |>.mpr p₂.2.2) - have hab : min (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ < max (p₁ : ℝ)⁻¹ (p₂ : ℝ)⁻¹ := - min_lt_max.mpr (inv_injective.ne hne) + have h12 : (p₁ : ℝ) - p₂ ≠ 0 := sub_ne_zero.mpr hne + have h21 : (p₂ : ℝ) - p₁ ≠ 0 := sub_ne_zero.mpr hne.symm refine Measure.ext_iff_singleton.mpr fun x ↦ ?_ - have hc : (((x.1 + x.2 : ℕ) : ℝ≥0∞) + 1)⁻¹ = ENNReal.ofReal (((x.1 + x.2 : ℕ) : ℝ) + 1)⁻¹ := by + set s := x.1 + x.2 + have hc : ((s : ℝ≥0∞) + 1)⁻¹ = ENNReal.ofReal ((s : ℝ) + 1)⁻¹ := by rw [ENNReal.ofReal_inv_of_pos (by positivity), ENNReal.ofReal_add (by positivity) zero_le_one, ENNReal.ofReal_natCast, ENNReal.ofReal_one] - have hnn : ∀ y ∈ antidiagonal (x.1 + x.2), + have hnn : ∀ y ∈ antidiagonal s, 0 ≤ (1 - (p₁ : ℝ)) ^ y.1 * p₁ * ((1 - p₂) ^ y.2 * p₂) := fun y _ ↦ mul_nonneg (geometricMeasure_nonneg _ _) (geometricMeasure_nonneg _ _) - rw [antidiagonalUniform_comp_singleton, map_add_singleton, + -- the law of the sum at `s`, a geometric sum over the antidiagonal + have hsum : ∑ y ∈ antidiagonal s, (1 - (p₁ : ℝ)) ^ y.1 * p₁ * ((1 - p₂) ^ y.2 * p₂) = + p₁ * p₂ * ((1 - p₁) ^ (s + 1) - (1 - p₂) ^ (s + 1)) / (p₂ - p₁) := by + rw [eq_div_iff h21, ← Nat.sum_antidiagonal_pow_mul_pow_mul_sub, sum_mul, sum_mul, mul_sum] + refine sum_congr rfl fun y _ ↦ ?_ + ring + rw [Kernel.antidiagonalUniform_comp_singleton, Measure.map_add_singleton_eq_sum_antidiagonal, sum_congr rfl fun y _ ↦ prod_geometricMeasure_singleton h₁ h₂ y, - ← ENNReal.ofReal_sum_of_nonneg hnn, hc, ← ENNReal.ofReal_mul (sum_nonneg hnn), - sum_antidiagonal_geom hne, geomPairInvKernel_comp_singleton ha hab] + ← ENNReal.ofReal_sum_of_nonneg hnn, hc, ← ENNReal.ofReal_mul (sum_nonneg hnn), hsum, + geomPairInvKernel_comp_singleton (one_le_inv₀ h₁' |>.mpr p₁.2.2) + (one_le_inv₀ h₂' |>.mpr p₂.2.2) (inv_injective.ne hne), inv_inv, inv_inv] congr 1 - have h21 : (p₂ : ℝ) - p₁ ≠ 0 := sub_ne_zero.mpr hne.symm - rcases hne.lt_or_gt with h' | h' - · have hinv : (p₂ : ℝ)⁻¹ ≤ (p₁ : ℝ)⁻¹ := (inv_strictAnti₀ h₁' h').le - rw [min_eq_right hinv, max_eq_left hinv, inv_inv, inv_inv] - field_simp - · have hinv : (p₁ : ℝ)⁻¹ ≤ (p₂ : ℝ)⁻¹ := (inv_strictAnti₀ h₂' h').le - rw [min_eq_left hinv, max_eq_right hinv, inv_inv, inv_inv] - field_simp - ring - -/-- A measure on `ℕ` is absolutely continuous with respect to a measure that charges every -point. -/ -private lemma absolutelyContinuous_of_singleton_pos {μ ν : Measure ℕ} (h : ∀ s, 0 < ν {s}) : - μ ≪ ν := - Measure.AbsolutelyContinuous.mk fun t _ ht ↦ by - rcases t.eq_empty_or_nonempty with rfl | ⟨s, hs⟩ - · exact measure_empty - · exact absurd ht ((h s).trans_le (measure_mono (singleton_subset_iff.mpr hs))).ne' + field_simp + ring /-- The law with the law of the sum under `Geom(p₁) ⊗ Geom(p₂)` and the conditional law `Kernel.antidiagonalUniform` given the sum is in the effective null. -/ @@ -349,16 +301,14 @@ private lemma antidiagonalUniform_comp_mem_effectiveSet (p₁ p₂ : unitInterva rw [← hq₁p, ← hq₂p] rcases eq_or_ne q₁ q₂ with rfl | h · obtain ⟨g, hg⟩ := exists_prod_geometricMeasure_self_singleton q₁ - rw [antidiagonalUniform_comp_map_add_of_singleton hg] + rw [Kernel.antidiagonalUniform_comp_map_add_of_singleton hg] exact subset_effectiveSet hS ⟨q₁, rfl⟩ have := isProbabilityMeasure_uniformMeasure (min_lt_max.mpr (inv_injective.ne (Subtype.coe_injective.ne h)) : min (q₁ : ℝ)⁻¹ (q₂ : ℝ)⁻¹ < max (q₁ : ℝ)⁻¹ (q₂ : ℝ)⁻¹) rw [antidiagonalUniform_comp_map_add_prod_geometricMeasure hq₁ hq₂ h] refine comp_mem_effectiveSet prob_le_one (ae_of_all _ fun m ↦ subset_effectiveSet hS ?_) - refine ⟨projIcc 0 1 zero_le_one m⁻¹, ?_⟩ - rw [geomPairInvKernel, Kernel.comap_apply, Kernel.prod_apply] - rfl + exact ⟨projIcc 0 1 zero_le_one m⁻¹, (geomPairInvKernel_apply m).symm⟩ end Geometric @@ -370,7 +320,8 @@ of the laws `Geom(p₁)` and `Geom(p₂)` is not a point mass, the likelihood ra `Geom(p₁) ⊗ Geom(p₂)` and, given the sum `s`, the first observation is uniform on `{0, …, s}`, is the numeraire of the null of the products `Geom(p) ⊗ Geom(p)`: it is an e-variable, and no e-variable has a larger expected logarithm under `Geom(p₁) ⊗ Geom(p₂)` -(`ProbabilityTheory.IsNumeraire.eintegral_log_le`). -/ +(`ProbabilityTheory.IsNumeraire.eintegral_log_le`). This is the geometric instance of +[EVR26, Note 28]. -/ theorem isNumeraire_rnDeriv_antidiagonalUniform_comp {p₁ p₂ : unitInterval} (h : (p₁ ≠ 0 ∧ p₁ ≠ 1) ∨ (p₂ ≠ 0 ∧ p₂ ≠ 1)) : IsNumeraire (((geometricMeasure p₁).prod (geometricMeasure p₂)).rnDeriv @@ -380,13 +331,14 @@ theorem isNumeraire_rnDeriv_antidiagonalUniform_comp {p₁ p₂ : unitInterval} ((geometricMeasure p₁).prod (geometricMeasure p₂)) := by refine isNumeraire_of_ae_eq_rnDeriv (isEVar_rnDeriv_comp_map (measurable_of_countable _) (ae_of_all _ Kernel.ae_add_antidiagonalUniform) - (absolutelyContinuous_antidiagonalUniform_comp _) (by rintro _ ⟨p, rfl⟩; infer_instance) ?_ ?_) - (antidiagonalUniform_comp_mem_effectiveSet p₁ p₂) - (absolutelyContinuous_antidiagonalUniform_comp _) .rfl + (Kernel.absolutelyContinuous_antidiagonalUniform_comp _) (by rintro _ ⟨p, rfl⟩; infer_instance) + ?_ ?_) (antidiagonalUniform_comp_mem_effectiveSet p₁ p₂) + (Kernel.absolutelyContinuous_antidiagonalUniform_comp _) .rfl · rintro _ ⟨p, rfl⟩ obtain ⟨g, hg⟩ := exists_prod_geometricMeasure_self_singleton p - exact antidiagonalUniform_comp_map_add_of_singleton hg + exact Kernel.antidiagonalUniform_comp_map_add_of_singleton hg · rintro _ ⟨p, rfl⟩ - exact absolutelyContinuous_of_singleton_pos (map_add_prod_geometricMeasure_singleton_pos h) + exact Measure.absolutelyContinuous_of_forall_singleton_pos + (map_add_prod_geometricMeasure_singleton_pos h) end ProbabilityTheory diff --git a/LMLExplore/MeasureTheory/Measure/Singletons.lean b/LMLExplore/MeasureTheory/Measure/Singletons.lean new file mode 100644 index 00000000..e04503a3 --- /dev/null +++ b/LMLExplore/MeasureTheory/Measure/Singletons.lean @@ -0,0 +1,52 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.MeasureTheory.Measure.AbsolutelyContinuous +public import Mathlib.MeasureTheory.Measure.Map +public import Mathlib.Algebra.Order.Antidiag.Prod + +/-! +# Measures through their masses on points + +Two facts about measures read off their masses on points: a measure that charges every point +dominates every measure, and on a countable additive monoid with an antidiagonal (such as `ℕ`), the +law of the sum of the coordinates of a measure on pairs gives a point `s` the total mass of the +pairs with sum `s`. + +## Main statements + +* `MeasureTheory.Measure.absolutelyContinuous_of_forall_singleton_pos`: every measure is absolutely + continuous with respect to a measure that charges every point. +* `MeasureTheory.Measure.map_add_singleton_eq_sum_antidiagonal`: the mass of `{s}` under the law of + `x.1 + x.2` is `∑ x ∈ antidiagonal s, μ {x}`. +-/ + +public section + +open Set Finset + +namespace MeasureTheory.Measure + +/-- Every measure is absolutely continuous with respect to a measure that charges every point. -/ +theorem absolutelyContinuous_of_forall_singleton_pos {α : Type*} {_ : MeasurableSpace α} + {μ ν : Measure α} (h : ∀ x, 0 < ν {x}) : μ ≪ ν := + AbsolutelyContinuous.mk fun t _ ht ↦ by + rcases t.eq_empty_or_nonempty with rfl | ⟨x, hx⟩ + · exact measure_empty + · exact absurd ht ((h x).trans_le (measure_mono (singleton_subset_iff.mpr hx))).ne' + +/-- On a countable additive monoid with an antidiagonal, the mass of a point `s` under the law of +the sum of the coordinates is the sum of the masses of the pairs with sum `s`. -/ +theorem map_add_singleton_eq_sum_antidiagonal {M : Type*} [AddMonoid M] [HasAntidiagonal M] + [MeasurableSpace M] [Countable M] [MeasurableSingletonClass M] (μ : Measure (M × M)) (s : M) : + μ.map (fun x ↦ x.1 + x.2) {s} = ∑ x ∈ antidiagonal s, μ {x} := by + rw [map_apply (measurable_of_countable _) (measurableSet_singleton s), sum_measure_singleton] + congr 1 + ext x + simp + +end MeasureTheory.Measure diff --git a/LMLExplore/Probability/Distributions/Exponential.lean b/LMLExplore/Probability/Distributions/Exponential.lean new file mode 100644 index 00000000..478b8b44 --- /dev/null +++ b/LMLExplore/Probability/Distributions/Exponential.lean @@ -0,0 +1,173 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.Probability.Distributions.Exponential +public import TauCeti.Probability.Distributions.Hypoexponential.Basic +public import TauCeti.Probability.Distributions.Uniform.Basic +import LMLExplore.Probability.Distributions.Gamma +import LMLExplore.Probability.Distributions.Uniform +import TauCeti.MeasureTheory.Measure.Mixture.Basic +import TauCeti.Probability.Distributions.Beta.Uniform +import TauCeti.Probability.Distributions.Gamma.Beta +import TauCeti.Probability.Distributions.Gamma.Measurability +import TauCeti.Probability.Distributions.Hypoexponential.PDF + +/-! +# Exponential laws and the sum of two independent exponential variables + +Facts about Mathlib's exponential laws `expMeasure l`: they are finite measures, and those with +positive rates are equivalent. The sum of two independent exponential variables of rates `l₁` and +`l₂` has Tau Ceti's hypoexponential law, which is `Gamma(2, l)` when `l₁ = l₂ = l` and has the +density `l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)` when `l₁ ≠ l₂`. In the +second case it is the mixture of the laws `Gamma(2, 1 / m)` over the mean `m` uniform between +`1 / l₁` and `1 / l₂`. + +## Main statements + +* `ProbabilityTheory.isFiniteMeasure_expMeasure`: exponential laws are finite measures. +* `ProbabilityTheory.expMeasure_absolutelyContinuous`, + `ProbabilityTheory.prod_expMeasure_absolutelyContinuous`: exponential laws, and their products, + with positive rates are equivalent. +* `ProbabilityTheory.map_add_prod_expMeasure_self`: `Exp(l) + Exp(l)` has the law `Gamma(2, l)`. +* `ProbabilityTheory.map_add_prod_expMeasure_eq_hypoexponentialMeasure`: `Exp(l₁) + Exp(l₂)` has + the hypoexponential law of rates `l₁, l₂`. +* `ProbabilityTheory.map_add_prod_expMeasure_eq_withDensity`: its density for `l₁ ≠ l₂`. +* `ProbabilityTheory.map_add_prod_expMeasure_eq_bind`: for `l₁ ≠ l₂`, it is the mixture of the + laws `Gamma(2, 1 / m)` over `m` uniform between `1 / l₁` and `1 / l₂`. +-/ + +public section + +open MeasureTheory Set Real TauCeti.Probability + +namespace ProbabilityTheory + +variable {l l₁ l₂ : ℝ} + +/-- Exponential laws are finite measures: probability measures for a positive rate, and the zero +measure otherwise. -/ +instance isFiniteMeasure_expMeasure (l : ℝ) : IsFiniteMeasure (expMeasure l) := by + rcases le_or_gt l 0 with hl | hl + · rw [expMeasure_of_nonpos hl] + infer_instance + · have := isProbabilityMeasure_expMeasure hl + infer_instance + +/-- Exponential laws with positive rates are equivalent. -/ +theorem expMeasure_absolutelyContinuous (l : ℝ) {l' : ℝ} (hl' : 0 < l') : + expMeasure l ≪ expMeasure l' := by + rw [expMeasure_eq_withDensity_restrict_Ioi l, expMeasure_eq_withDensity_restrict_Ioi l'] + refine (withDensity_absolutelyContinuous _ _).trans (withDensity_absolutelyContinuous' + (measurable_exponentialPDFReal l').ennreal_ofReal.aemeasurable ?_) + filter_upwards [ae_restrict_mem measurableSet_Ioi] with x hx + exact (ENNReal.ofReal_pos.mpr (exponentialPDFReal_pos hl' hx)).ne' + +/-- Products of exponential laws with positive rates are equivalent. -/ +theorem prod_expMeasure_absolutelyContinuous (l₁ l₂ : ℝ) {l₁' l₂' : ℝ} (h₁ : 0 < l₁') + (h₂ : 0 < l₂') : + (expMeasure l₁).prod (expMeasure l₂) ≪ (expMeasure l₁').prod (expMeasure l₂') := + (expMeasure_absolutelyContinuous l₁ h₁).prod (expMeasure_absolutelyContinuous l₂ h₂) + +/-- The sum of two independent `Exp(l)` variables has the law `Gamma(2, l)`. -/ +theorem map_add_prod_expMeasure_self (hl : 0 < l) : + ((expMeasure l).prod (expMeasure l)).map (fun x ↦ x.1 + x.2) = gammaMeasure 2 l := by + have : IsProbabilityMeasure (betaMeasure 1 1) := by rw [betaMeasure_one_one]; infer_instance + have h := congrArg (Measure.map Prod.snd) (map_div_add_prod_gammaMeasure one_pos one_pos hl) + rw [Measure.map_map measurable_snd (by fun_prop), Measure.map_snd_prod, measure_univ, + one_smul, one_add_one_eq_two] at h + exact h + +/-- **The sum of two independent exponential variables is hypoexponential**: for positive rates +`l₁, l₂`, the law of the sum under `Exp(l₁) ⊗ Exp(l₂)` is `hypoexponentialMeasure ![l₁, l₂]`. -/ +theorem map_add_prod_expMeasure_eq_hypoexponentialMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) : + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + hypoexponentialMeasure ![l₁, l₂] := by + have hr : ∀ i, 0 < ![l₁, l₂] i := by simp [Fin.forall_fin_two, h₁, h₂] + have hpi : (fun i ↦ expMeasure (![l₁, l₂] i)) = ![expMeasure l₁, expMeasure l₂] := by + ext i : 1 + fin_cases i <;> rfl + rw [hypoexponentialMeasure_eq_map hr, hpi, + ← (measurePreserving_finTwoArrow_vec (expMeasure l₁) (expMeasure l₂)).map_eq, + Measure.map_map (by fun_prop) (MeasurableEquiv.measurable _)] + congr 1 + ext x + simp [Fin.sum_univ_two, MeasurableEquiv.finTwoArrow] + +/-- **The density of the sum of two exponential variables** with distinct positive rates +`l₁, l₂` is `l₁ l₂ (exp (- l₁ x) - exp (- l₂ x)) / (l₂ - l₁)` on `(0, ∞)`. -/ +theorem map_add_prod_expMeasure_eq_withDensity (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + (volume.restrict (Ioi 0)).withDensity fun x ↦ + ENNReal.ofReal (l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁)) := by + have hr : ∀ i, 0 < ![l₁, l₂] i := by simp [Fin.forall_fin_two, h₁, h₂] + have hinj : Function.Injective ![l₁, l₂] := by + intro i j hij + fin_cases i <;> fin_cases j <;> simp_all [eq_comm] + rw [map_add_prod_expMeasure_eq_hypoexponentialMeasure h₁ h₂, + ← Measure.restrict_eq_self_of_ae_mem (s := Ioi 0) + (μ := hypoexponentialMeasure ![l₁, l₂]) (ae_pos_hypoexponentialMeasure _), + hypoexponentialMeasure_eq_withDensity, restrict_withDensity measurableSet_Ioi] + refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr + (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) + rw [hypoexponentialPDF_eq_ofReal, hypoexponentialPDFReal_of_injective hr hinj hx] + have h01 : (Finset.univ.erase (0 : Fin 2)) = {1} := by decide + have h10 : (Finset.univ.erase (1 : Fin 2)) = {0} := by decide + have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm + have h21 : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h + simp only [Fin.sum_univ_two, h01, h10, Finset.prod_singleton, Matrix.cons_val_zero, + Matrix.cons_val_one] + congr 1 + field_simp + ring + +/-- The density of the mixture of the laws `Gamma(2, 1 / m)` over the mean `m` uniform between +`1 / l₁` and `1 / l₂`, at `x > 0`. -/ +private lemma lintegral_gammaPDF_uniformMeasure (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) + {x : ℝ} (hx : 0 < x) : + ∫⁻ m, gammaPDF 2 m⁻¹ x ∂uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹) = + ENNReal.ofReal (l₁ * l₂ * (exp (-(l₁ * x)) - exp (-(l₂ * x))) / (l₂ - l₁)) := by + have hpos : ∀ m ∈ uIcc l₁⁻¹ l₂⁻¹, 0 < m := fun m hm ↦ + (lt_inf_iff.mpr ⟨inv_pos.mpr h₁, inv_pos.mpr h₂⟩).trans_le hm.1 + have hderiv : ∀ m ∈ uIcc l₁⁻¹ l₂⁻¹, HasDerivAt (fun m ↦ exp (-(m⁻¹ * x))) + (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) m := + fun m hm ↦ (((hasDerivAt_inv (hpos m hm).ne').mul_const x).neg.exp).congr_deriv + (by simp only [Pi.neg_apply]; ring) + have hcont : ContinuousOn (fun m ↦ exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) (uIcc l₁⁻¹ l₂⁻¹) := by + refine ContinuousOn.mul (by fun_prop (disch := exact fun m hm ↦ (hpos m hm).ne')) ?_ + fun_prop (disch := exact fun m hm ↦ pow_ne_zero 2 (hpos m hm).ne') + have hpdf : ∀ᵐ m ∂uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹), gammaPDF 2 m⁻¹ x = + ENNReal.ofReal (exp (-(m⁻¹ * x)) * ((m ^ 2)⁻¹ * x)) := by + filter_upwards [ae_mem_Ioc_uniformMeasure] with m hm + rw [gammaPDF_two hx.le] + congr 1 + ring + rw [lintegral_congr_ae hpdf, lintegral_ofReal_uniformMeasure_min_max_of_hasDerivAt + (inv_injective.ne h) hderiv hcont fun m hm ↦ ?_] + · have h12 : l₂ - l₁ ≠ 0 := sub_ne_zero.mpr h.symm + have h21 : l₁ - l₂ ≠ 0 := sub_ne_zero.mpr h + simp only [inv_inv] + congr 1 + field_simp + ring + · have := hpos m (uIoc_subset_uIcc hm) + positivity + +/-- **The sum of two exponential variables is a mixture of gamma laws**: for distinct positive +rates, the law of the sum is the mixture of the laws `Gamma(2, 1 / m)` over the mean `m` uniform +between `1 / l₁` and `1 / l₂`. -/ +theorem map_add_prod_expMeasure_eq_bind (h₁ : 0 < l₁) (h₂ : 0 < l₂) (h : l₁ ≠ l₂) : + ((expMeasure l₁).prod (expMeasure l₂)).map (fun x ↦ x.1 + x.2) = + (uniformMeasure (min l₁⁻¹ l₂⁻¹) (max l₁⁻¹ l₂⁻¹)).bind fun m ↦ gammaMeasure 2 m⁻¹ := by + simp_rw [gammaMeasure_eq_withDensity_restrict_Ioi] + have hf : Measurable (Function.uncurry fun m x ↦ gammaPDF 2 m⁻¹ x) := + measurable_uncurry_gammaPDF.comp (f := fun p : ℝ × ℝ ↦ (((2 : ℝ), p.1⁻¹), p.2)) (by fun_prop) + rw [map_add_prod_expMeasure_eq_withDensity h₁ h₂ h, Measure.bind_withDensity _ _ hf] + refine withDensity_congr_ae ((ae_restrict_iff' measurableSet_Ioi).mpr + (ae_of_all _ fun x (hx : 0 < x) ↦ ?_)) + exact (lintegral_gammaPDF_uniformMeasure h₁ h₂ h hx).symm + +end ProbabilityTheory diff --git a/LMLExplore/Probability/Distributions/Gamma.lean b/LMLExplore/Probability/Distributions/Gamma.lean new file mode 100644 index 00000000..c9d3868e --- /dev/null +++ b/LMLExplore/Probability/Distributions/Gamma.lean @@ -0,0 +1,33 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.Probability.Distributions.Gamma +import Mathlib.Analysis.SpecialFunctions.Gamma.BohrMollerup + +/-! +# The density of the gamma law of shape two + +The density of Mathlib's gamma law `gammaMeasure 2 r` at `x ≥ 0` is `r ^ 2 * x * exp (- r x)`. + +## Main statements + +* `ProbabilityTheory.gammaPDF_two`: the density of `Gamma(2, r)` on `[0, ∞)`. +-/ + +public section + +open Real + +namespace ProbabilityTheory + +/-- The density of `Gamma(2, r)` on `[0, ∞)` is `r ^ 2 * x * exp (- r x)`. -/ +theorem gammaPDF_two {r x : ℝ} (hx : 0 ≤ x) : + gammaPDF 2 r x = ENNReal.ofReal (r ^ 2 * x * exp (-(r * x))) := by + rw [gammaPDF_of_nonneg hx, Real.Gamma_two, show (2 : ℝ) - 1 = 1 by norm_num, Real.rpow_one, + Real.rpow_two, div_one] + +end ProbabilityTheory diff --git a/LMLExplore/Probability/Distributions/Geometric.lean b/LMLExplore/Probability/Distributions/Geometric.lean new file mode 100644 index 00000000..f80079e3 --- /dev/null +++ b/LMLExplore/Probability/Distributions/Geometric.lean @@ -0,0 +1,57 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.Probability.Distributions.Geometric +import TauCeti.Probability.Distributions.Geometric.Basic + +/-! +# Geometric laws with a positive success probability + +Mathlib's `geometricMeasure 0` is the point mass at `0`, as is `geometricMeasure 1` +(Tau Ceti's `TauCeti.Probability.geometricMeasure_zero` and +`TauCeti.Probability.geometricMeasure_one`), so every geometric law is one with a positive success +probability, for which Mathlib's `geometricMeasure_singleton` gives the masses of the points. + +## Main statements + +* `ProbabilityTheory.exists_ne_zero_geometricMeasure_eq`: every geometric law is one with a + positive success probability. +* `ProbabilityTheory.geometricMeasure_singleton_zero_pos`: every geometric law charges `0`. +* `ProbabilityTheory.prod_geometricMeasure_singleton`: the mass of a point under a product of two + geometric laws with positive success probabilities. +-/ + +public section + +open MeasureTheory Set + +namespace ProbabilityTheory + +/-- Every geometric law is one with a positive success probability. -/ +theorem exists_ne_zero_geometricMeasure_eq (p : unitInterval) : + ∃ q : unitInterval, q ≠ 0 ∧ geometricMeasure q = geometricMeasure p := by + by_cases hp : p = 0 + · exact ⟨1, one_ne_zero, by rw [hp, TauCeti.Probability.geometricMeasure_one, + TauCeti.Probability.geometricMeasure_zero]⟩ + · exact ⟨p, hp, rfl⟩ + +/-- The point `0` has positive mass under every geometric law. -/ +theorem geometricMeasure_singleton_zero_pos (p : unitInterval) : 0 < geometricMeasure p {0} := by + obtain ⟨q, hq, hqp⟩ := exists_ne_zero_geometricMeasure_eq p + rw [← hqp, geometricMeasure_singleton hq] + simpa using unitInterval.pos_iff_ne_zero.mpr hq + +/-- The mass of a point under a product of geometric laws with positive success +probabilities. -/ +theorem prod_geometricMeasure_singleton {p₁ p₂ : unitInterval} (h₁ : p₁ ≠ 0) (h₂ : p₂ ≠ 0) + (x : ℕ × ℕ) : + (geometricMeasure p₁).prod (geometricMeasure p₂) {x} = + ENNReal.ofReal ((1 - p₁) ^ x.1 * p₁ * ((1 - p₂) ^ x.2 * p₂)) := by + rw [← singleton_prod_singleton, Measure.prod_prod, geometricMeasure_singleton h₁, + geometricMeasure_singleton h₂, ← ENNReal.ofReal_mul (geometricMeasure_nonneg _ _)] + +end ProbabilityTheory diff --git a/LMLExplore/Probability/Distributions/Uniform.lean b/LMLExplore/Probability/Distributions/Uniform.lean new file mode 100644 index 00000000..d9edf34c --- /dev/null +++ b/LMLExplore/Probability/Distributions/Uniform.lean @@ -0,0 +1,73 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.MeasureTheory.Integral.IntervalIntegral.FundThmCalculus +public import TauCeti.Probability.Distributions.Uniform.Basic + +/-! +# Integrating a derivative against the uniform law + +Against Tau Ceti's uniform law `uniformMeasure a b` on `(a, b]`, the lower Lebesgue integral of a +nonnegative continuous derivative `f'` of `f` is the mean slope `(b - a)⁻¹ * (f b - f a)`, by the +fundamental theorem of calculus. The second form takes the uniform law between two distinct +points `a` and `b` in either order, `uniformMeasure (min a b) (max a b)`; the mean slope is +symmetric in `a` and `b`. + +## Main statements + +* `ProbabilityTheory.lintegral_ofReal_uniformMeasure_of_hasDerivAt`: for `a < b`, + `∫⁻ x, ENNReal.ofReal (f' x) ∂uniformMeasure a b = ENNReal.ofReal ((b - a)⁻¹ * (f b - f a))`. +* `ProbabilityTheory.lintegral_ofReal_uniformMeasure_min_max_of_hasDerivAt`: the same for the + uniform law between `a ≠ b`. +-/ + +public section + +open MeasureTheory Set TauCeti.Probability + +namespace ProbabilityTheory + +variable {a b : ℝ} {f f' : ℝ → ℝ} + +/-- **The integral of a derivative against the uniform law is the mean slope.** For `a < b` and a +function `f` with a continuous derivative `f'` on `[a, b]`, nonnegative on `(a, b]`, the integral +of `f'` against the uniform law on `(a, b]` is `(b - a)⁻¹ * (f b - f a)`. -/ +theorem lintegral_ofReal_uniformMeasure_of_hasDerivAt (hab : a < b) + (hf : ∀ x ∈ Icc a b, HasDerivAt f (f' x) x) (hf' : ContinuousOn f' (Icc a b)) + (hf'0 : ∀ x ∈ Ioc a b, 0 ≤ f' x) : + ∫⁻ x, ENNReal.ofReal (f' x) ∂uniformMeasure a b = + ENNReal.ofReal ((b - a)⁻¹ * (f b - f a)) := by + rw [← uIcc_of_le hab.le] at hf hf' + rw [uniformMeasure_eq_smul, lintegral_smul_measure, ← ofReal_integral_eq_lintegral_ofReal, + ← intervalIntegral.integral_of_le hab.le, + intervalIntegral.integral_eq_sub_of_hasDerivAt hf hf'.intervalIntegrable, smul_eq_mul, + ← ENNReal.ofReal_inv_of_pos (sub_pos.mpr hab), + ← ENNReal.ofReal_mul (inv_nonneg.mpr (sub_pos.mpr hab).le)] + · rw [uIcc_of_le hab.le] at hf' + exact hf'.integrableOn_Icc.mono_set Ioc_subset_Icc_self + · exact ae_restrict_of_forall_mem measurableSet_Ioc hf'0 + +/-- **The integral of a derivative against the uniform law between two points.** For `a ≠ b` and +a function `f` with a continuous derivative `f'` on `uIcc a b`, nonnegative on `uIoc a b`, the +integral of `f'` against the uniform law between `a` and `b` is `(b - a)⁻¹ * (f b - f a)`. -/ +theorem lintegral_ofReal_uniformMeasure_min_max_of_hasDerivAt (hab : a ≠ b) + (hf : ∀ x ∈ uIcc a b, HasDerivAt f (f' x) x) (hf' : ContinuousOn f' (uIcc a b)) + (hf'0 : ∀ x ∈ uIoc a b, 0 ≤ f' x) : + ∫⁻ x, ENNReal.ofReal (f' x) ∂uniformMeasure (min a b) (max a b) = + ENNReal.ofReal ((b - a)⁻¹ * (f b - f a)) := by + rcases hab.lt_or_gt with h | h + · rw [min_eq_left h.le, max_eq_right h.le] + rw [uIcc_of_le h.le] at hf hf' + rw [uIoc_of_le h.le] at hf'0 + exact lintegral_ofReal_uniformMeasure_of_hasDerivAt h hf hf' hf'0 + · rw [min_eq_right h.le, max_eq_left h.le, ← neg_sub a b, inv_neg, ← neg_sub (f a) (f b), + neg_mul_neg] + rw [uIcc_of_ge h.le] at hf hf' + rw [uIoc_of_ge h.le] at hf'0 + exact lintegral_ofReal_uniformMeasure_of_hasDerivAt h hf hf' hf'0 + +end ProbabilityTheory diff --git a/LMLExplore/Probability/Kernel/Composition/MeasureComp.lean b/LMLExplore/Probability/Kernel/Composition/MeasureComp.lean new file mode 100644 index 00000000..c901e496 --- /dev/null +++ b/LMLExplore/Probability/Kernel/Composition/MeasureComp.lean @@ -0,0 +1,50 @@ +/- +Copyright (c) 2026 LMLExplore contributors. All rights reserved. +Released under Apache 2.0 license as described in the file LICENSE. +Authors: LMLExplore contributors +-/ +module + +public import Mathlib.Probability.Kernel.Composition.MeasureComp + +/-! +# The law of a statistic under a composition with a kernel carried by its fibers + +If a Markov kernel `κ : Kernel 𝓨 Ω` is carried by the fibers of a statistic `T : Ω → 𝓨`, that is, +`κ y` is carried by `T ⁻¹' {y}`, then the law of `T` under `κ ∘ₘ ν` is `ν`: the composition with +`κ` draws a point of the fiber of each `y`. + +## Main statements + +* `ProbabilityTheory.map_comp_eq_of_ae_eq`: `(κ ∘ₘ ν).map T = ν`. +-/ + +public section + +open MeasureTheory Set + +namespace ProbabilityTheory + +variable {Ω 𝓨 : Type*} {mΩ : MeasurableSpace Ω} {m𝓨 : MeasurableSpace 𝓨} {T : Ω → 𝓨} + {κ : Kernel 𝓨 Ω} + +/-- If `κ y` is a probability measure carried by the fiber `T ⁻¹' {y}` for `ν`-almost every `y`, +the law of `T` under `κ ∘ₘ ν` is `ν`. -/ +theorem map_comp_eq_of_ae_eq [IsMarkovKernel κ] (hT : Measurable T) {ν : Measure 𝓨} + (hκ : ∀ᵐ y ∂ν, ∀ᵐ x ∂κ y, T x = y) : (κ ∘ₘ ν).map T = ν := by + ext A hA + rw [Measure.map_apply hT hA, Measure.bind_apply (hT hA) κ.aemeasurable, + ← lintegral_indicator_one hA] + refine lintegral_congr_ae ?_ + filter_upwards [hκ] with y hκ + by_cases hy : y ∈ A + · rw [indicator_of_mem hy, Pi.one_apply, ← measure_univ (μ := κ y)] + refine measure_congr ?_ + filter_upwards [hκ] with x hx + simp [hx, hy] + · rw [indicator_of_notMem hy] + refine measure_eq_zero_iff_ae_notMem.mpr ?_ + filter_upwards [hκ] with x hx + simpa [hx] using hy + +end ProbabilityTheory