AveragedCongruence.lean

1 import Mathlib.MeasureTheory.Function.L2Space
2 import Mathlib.MeasureTheory.SpecificCodomains.Pi
3 import Mathlib.LinearAlgebra.Matrix.PosDef
4 import Mathlib.Analysis.SpecialFunctions.Integrals.Basic
5 import Mathlib.Probability.ConditionalProbability
6 import Mathlib.Analysis.Matrix.Normed
7 
8 /-! The covariance correction for averaging a quadratic matrix action. -/
9 set_option autoImplicit false
10 
11 open MeasureTheory Matrix
12 open scoped BigOperators ENNReal MeasureTheory
13 
14 namespace Komorebi.LSM.AveragedCongruence
15 
16 variable {Ω : Type*} {mΩ : MeasurableSpace Ω} {ν : Measure Ω}
17   [IsProbabilityMeasure ν]
18 
19 theorem product_covariance {f g : Ω → ℝ} (hf : MemLp f 2 ν) (hg : MemLp g 2 ν) :
20     (∫ ω, f ω * g ω ∂ν) = (∫ ω, f ω ∂ν) * (∫ ω, g ω ∂ν) +
21       ∫ ω, (f ω - ∫ y, f y ∂ν) * (g ω - ∫ y, g y ∂ν) ∂ν := by
22   have hfi : Integrable f ν := hf.integrable (by norm_num)
23   have hgi : Integrable g ν := hg.integrable (by norm_num)
24   have hfg : Integrable (fun ω => f ω * g ω) ν := hf.integrable_mul hg
25   have heq : (fun ω => (f ω - ∫ y, f y ∂ν) * (g ω - ∫ y, g y ∂ν)) =
26       (fun ω => f ω * g ω - f ω * (∫ y, g y ∂ν) -
27         (∫ y, f y ∂ν) * g ω + (∫ y, f y ∂ν) * (∫ y, g y ∂ν)) := by
28     funext ω
29     ring
30   rw [heq, integral_add
31     (f := fun ω => f ω * g ω - f ω * (∫ y, g y ∂ν) - (∫ y, f y ∂ν) * g ω)
32     (g := fun _ => (∫ y, f y ∂ν) * (∫ y, g y ∂ν))
33     ((hfg.sub (hfi.mul_const _)).sub (hgi.const_mul _)) (integrable_const _),
34     integral_sub
35       (f := fun ω => f ω * g ω - f ω * (∫ y, g y ∂ν))
36       (g := fun ω => (∫ y, f y ∂ν) * g ω)
37       (hfg.sub (hfi.mul_const _)) (hgi.const_mul _),
38     integral_sub (f := fun ω => f ω * g ω)
39       (g := fun ω => f ω * (∫ y, g y ∂ν)) hfg (hfi.mul_const _),
40     integral_mul_const, integral_const_mul]
41   simp only [integral_const, probReal_univ, one_smul]
42   ring
43 
44 variable {M N : Type*} [Fintype M] [Fintype N]
45 
46 /-- Matrix integration is taken in its finite coordinate space. -/
47 noncomputable def matrixMean (K : Ω → Matrix M N ℝ) : Matrix M N ℝ :=
48   fun i j => ∫ ω, K ω i j ∂ν
49 
50 omit [Fintype M] in
51 theorem congruence_entry (K : Matrix M N ℝ) (X : Matrix N N ℝ) (i j : M) :
52     (K * X * K.transpose) i j = ∑ a, ∑ b, X a b * (K i a * K j b) := by
53   simp only [Matrix.mul_apply, Matrix.transpose_apply, Finset.sum_mul]
54   rw [Finset.sum_comm]
55   apply Finset.sum_congr rfl
56   intro a _
57   apply Finset.sum_congr rfl
58   intro b _
59   ring
60 
61 omit [IsProbabilityMeasure ν] [Fintype M] in
62 theorem congruence_entry_integrable {K : Ω → Matrix M N ℝ}
63     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν) (X : Matrix N N ℝ) (i j : M) :
64     Integrable (fun ω => (K ω * X * (K ω).transpose) i j) ν := by
65   simp_rw [congruence_entry]
66   apply integrable_finsetSum
67   intro a _
68   apply integrable_finsetSum
69   intro b _
70   exact ((hK i a).integrable_mul (hK j b)).const_mul (X a b)
71 
72 omit [IsProbabilityMeasure ν] [Fintype M] in
73 theorem mean_congruence_entry {K : Ω → Matrix M N ℝ}
74     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν) (X : Matrix N N ℝ) (i j : M) :
75     matrixMean (ν := ν) (fun ω => K ω * X * (K ω).transpose) i j =
76       ∑ a, ∑ b, X a b * ∫ ω, K ω i a * K ω j b ∂ν := by
77   unfold matrixMean
78   simp_rw [congruence_entry]
79   rw [integral_finsetSum]
80   · apply Finset.sum_congr rfl
81     intro a _
82     rw [integral_finsetSum]
83     · apply Finset.sum_congr rfl
84       intro b _
85       exact integral_const_mul _ _
86     · intro b _
87       exact ((hK i a).integrable_mul (hK j b)).const_mul (X a b)
88   · intro a _
89     exact integrable_finsetSum _ fun b _ =>
90       ((hK i a).integrable_mul (hK j b)).const_mul (X a b)
91 
92 omit [Fintype M] in
93 /-- P57 covariance identity. Matrix integrals are evaluated entry by entry,
94 the canonical finite-dimensional matrix integral. Symmetry of `X` is not
95 needed for the identity itself. -/
96 theorem averaged_congruence {K : Ω → Matrix M N ℝ}
97     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν) (X : Matrix N N ℝ) :
98     matrixMean (ν := ν) (fun ω => K ω * X * (K ω).transpose) =
99       matrixMean (ν := ν) K * X * (matrixMean (ν := ν) K).transpose +
100       matrixMean (ν := ν) (fun ω => (K ω - matrixMean (ν := ν) K) * X *
101         (K ω - matrixMean (ν := ν) K).transpose) := by
102   have hc : ∀ i j, MemLp (fun ω => (K ω - matrixMean (ν := ν) K) i j) 2 ν := by
103     intro i j
104     exact (hK i j).sub (memLp_const _)
105   ext i j
106   rw [Matrix.add_apply, mean_congruence_entry hK X i j,
107     congruence_entry, mean_congruence_entry hc X i j]
108   change (∑ a, ∑ b, X a b * ∫ ω, K ω i a * K ω j b ∂ν) =
109     (∑ a, ∑ b, X a b * ((∫ ω, K ω i a ∂ν) * (∫ ω, K ω j b ∂ν))) +
110       ∑ a, ∑ b, X a b *
111         ∫ ω, (K ω i a - ∫ y, K y i a ∂ν) * (K ω j b - ∫ y, K y j b ∂ν) ∂ν
112   calc
113     _ = ∑ a, ∑ b,
114         (X a b * ((∫ ω, K ω i a ∂ν) * (∫ ω, K ω j b ∂ν)) +
115           X a b * ∫ ω, (K ω i a - ∫ y, K y i a ∂ν) *
116             (K ω j b - ∫ y, K y j b ∂ν) ∂ν) := by
117       apply Finset.sum_congr rfl
118       intro a _
119       apply Finset.sum_congr rfl
120       intro b _
121       rw [product_covariance (hK i a) (hK j b), mul_add]
122     _ = _ := by simp only [Finset.sum_add_distrib]
123 
124 omit [IsProbabilityMeasure ν] [Fintype N] in
125 theorem mean_quadratic {A : Ω → Matrix M M ℝ}
126     (hA : ∀ i j, Integrable (fun ω => A ω i j) ν) (x : M → ℝ) :
127     star x ⬝ᵥ (matrixMean (ν := ν) A *ᵥ x) =
128       ∫ ω, star x ⬝ᵥ (A ω *ᵥ x) ∂ν := by
129   simp only [dotProduct, Matrix.mulVec, matrixMean, Finset.mul_sum]
130   rw [integral_finsetSum]
131   · apply Finset.sum_congr rfl
132     intro i _
133     rw [integral_finsetSum]
134     · apply Finset.sum_congr rfl
135       intro j _
136       rw [integral_const_mul, integral_mul_const]
137     · intro j _
138       exact ((hA i j).mul_const (x j)).const_mul (star (x i))
139   · intro i _
140     exact integrable_finsetSum _ fun j _ =>
141       ((hA i j).mul_const (x j)).const_mul (star (x i))
142 
143 omit [IsProbabilityMeasure ν] [Fintype N] in
144 /-- Integrating a finite real positive-semidefinite matrix family preserves
145 positive semidefiniteness. -/
146 theorem matrixMean_posSemidef {A : Ω → Matrix M M ℝ}
147     (hAi : ∀ i j, Integrable (fun ω => A ω i j) ν)
148     (hA : ∀ ω, (A ω).PosSemidef) : (matrixMean (ν := ν) A).PosSemidef := by
149   apply Matrix.PosSemidef.of_dotProduct_mulVec_nonneg
150   · ext i j
151     change (∫ ω, A ω j i ∂ν) = ∫ ω, A ω i j ∂ν
152     apply integral_congr_ae
153     filter_upwards with ω
154     have hs := congrArg (fun B : Matrix M M ℝ => B i j) (hA ω).isHermitian
155     simpa only [Matrix.conjTranspose_apply, star_trivial] using hs
156   · intro x
157     rw [mean_quadratic hAi x]
158     exact integral_nonneg fun ω => (hA ω).dotProduct_mulVec_nonneg x
159 
160 /-- The P57 covariance correction is positive semidefinite whenever `X` is. -/
161 theorem covariance_correction_posSemidef {K : Ω → Matrix M N ℝ}
162     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν)
163     {X : Matrix N N ℝ} (hX : X.PosSemidef) :
164     (matrixMean (ν := ν) (fun ω => (K ω - matrixMean (ν := ν) K) * X *
165       (K ω - matrixMean (ν := ν) K).transpose)).PosSemidef := by
166   have hc : ∀ i j, MemLp (fun ω => (K ω - matrixMean (ν := ν) K) i j) 2 ν := by
167     intro i j
168     exact (hK i j).sub (memLp_const _)
169   refine matrixMean_posSemidef (fun i j => congruence_entry_integrable hc X i j) ?_
170   intro ω
171   simpa only [Matrix.conjTranspose_eq_transpose_of_trivial] using
172     hX.mul_mul_conjTranspose_same (K ω - matrixMean (ν := ν) K)
173 
174 omit [Fintype M] [Fintype N] in
175 theorem mean_affine_matrix {u : Ω → ℝ} (hu : MemLp u 2 ν)
176     (K₀ K₁ : Matrix M N ℝ) :
177     matrixMean (ν := ν) (fun ω => K₀ + u ω • K₁) =
178       K₀ + (∫ ω, u ω ∂ν) • K₁ := by
179   ext i j
180   change (∫ ω, K₀ i j + u ω * K₁ i j ∂ν) = K₀ i j + (∫ ω, u ω ∂ν) * K₁ i j
181   rw [integral_add (integrable_const _) ((hu.integrable (by norm_num)).mul_const _),
182     integral_mul_const]
183   simp
184 
185 omit [Fintype M] [Fintype N] in
186 theorem affine_entries_memLp {u : Ω → ℝ} (hu : MemLp u 2 ν)
187     (K₀ K₁ : Matrix M N ℝ) (i : M) (j : N) :
188     MemLp (fun ω => (K₀ + u ω • K₁) i j) 2 ν :=
189   (memLp_const (K₀ i j)).add (hu.mul_const (K₁ i j))
190 
191 omit [Fintype M] in
192 /-- For an affine transformation family the correction is its scalar
193 variance times one fixed positive congruence. -/
194 theorem affine_covariance_correction {u : Ω → ℝ} (hu : MemLp u 2 ν)
195     (K₀ K₁ : Matrix M N ℝ) (X : Matrix N N ℝ) :
196     matrixMean (ν := ν) (fun ω =>
197       ((K₀ + u ω • K₁) - matrixMean (ν := ν) (fun y => K₀ + u y • K₁)) * X *
198       ((K₀ + u ω • K₁) - matrixMean (ν := ν) (fun y => K₀ + u y • K₁)).transpose) =
199       (∫ ω, (u ω - ∫ y, u y ∂ν) ^ 2 ∂ν) • (K₁ * X * K₁.transpose) := by
200   have hc (ω : Ω) : (K₀ + u ω • K₁) -
201       matrixMean (ν := ν) (fun y => K₀ + u y • K₁) =
202       (u ω - ∫ y, u y ∂ν) • K₁ := by
203     rw [mean_affine_matrix hu K₀ K₁]
204     ext i j
205     simp only [Matrix.sub_apply, Matrix.add_apply, Matrix.smul_apply, smul_eq_mul]
206     ring
207   simp_rw [hc, Matrix.transpose_smul, Matrix.smul_mul, Matrix.mul_smul, smul_smul,
208     ← pow_two]
209   ext i j
210   exact integral_mul_const _ _
211 
212 /-- The uniform probability measure on a nondegenerate real interval. -/
213 noncomputable def uniformInterval (a b : ℝ) : Measure ℝ :=
214   ProbabilityTheory.cond volume (Set.Ioc a b)
215 
216 theorem uniformInterval_isProbability {a b : ℝ} (hab : a < b) :
217     IsProbabilityMeasure (uniformInterval a b) := by
218   apply ProbabilityTheory.cond_isProbabilityMeasure_of_finite
219   · rw [Real.volume_Ioc]
220     exact ne_of_gt (ENNReal.ofReal_pos.mpr (sub_pos.mpr hab))
221   · rw [Real.volume_Ioc]
222     exact ENNReal.ofReal_ne_top
223 
224 theorem integral_uniformInterval {a b : ℝ} (hab : a < b) (f : ℝ → ℝ) :
225     (∫ x, f x ∂uniformInterval a b) = (b - a)⁻¹ * ∫ x in a..b, f x := by
226   rw [uniformInterval, ProbabilityTheory.cond, integral_smul_measure, Real.volume_Ioc,
227     ENNReal.toReal_inv, ENNReal.toReal_ofReal (sub_nonneg.mpr hab.le),
228     ← intervalIntegral.integral_of_le hab.le]
229   rfl
230 
231 theorem uniformInterval_mean {a b : ℝ} (hab : a < b) :
232     (∫ x, x ∂uniformInterval a b) = (a + b) / 2 := by
233   rw [integral_uniformInterval hab, integral_id]
234   field_simp [ne_of_gt (sub_pos.mpr hab)]
235   ring
236 
237 /-- The exact variance of an interval-uniform random variable is width
238 squared divided by twelve. -/
239 theorem uniformInterval_variance {a b : ℝ} (hab : a < b) :
240     (∫ x, (x - ∫ y, y ∂uniformInterval a b) ^ 2 ∂uniformInterval a b) =
241       (b - a) ^ 2 / 12 := by
242   rw [uniformInterval_mean hab, integral_uniformInterval hab]
243   have hpoly : (fun x : ℝ => (x - (a + b) / 2) ^ 2) =
244       (fun x => x ^ 2 - (a + b) * x + ((a + b) / 2) ^ 2) := by
245     funext x
246     ring
247   have hsq : IntervalIntegrable (fun x : ℝ => x ^ 2) volume a b :=
248     (continuous_id.pow 2).intervalIntegrable a b
249   have hlin : IntervalIntegrable (fun x : ℝ => (a + b) * x) volume a b :=
250     (continuous_const.mul continuous_id).intervalIntegrable a b
251   rw [hpoly, intervalIntegral.integral_add (f := fun x : ℝ => x ^ 2 - (a + b) * x)
252     (g := fun _ : ℝ => ((a + b) / 2) ^ 2) (hsq.sub hlin) intervalIntegrable_const,
253     intervalIntegral.integral_sub hsq hlin, integral_pow,
254     intervalIntegral.integral_const_mul, integral_id, intervalIntegral.integral_const]
255   simp only [smul_eq_mul]
256   field_simp [ne_of_gt (sub_pos.mpr hab)]
257   ring
258 
259 theorem uniformInterval_id_memLp {a b : ℝ} (hab : a < b) :
260     MemLp (fun x : ℝ => x) 2 (uniformInterval a b) := by
261   have hs : Integrable (fun x : ℝ => x ^ 2) (volume.restrict (Set.Ioc a b)) :=
262     ((continuous_id.pow 2).intervalIntegrable a b).1
263   have hi : MemLp (fun x : ℝ => x) 2 (volume.restrict (Set.Ioc a b)) :=
264     (memLp_two_iff_integrable_sq continuous_id.aestronglyMeasurable).mpr hs
265   apply hi.smul_measure
266   rw [ENNReal.inv_ne_top, Real.volume_Ioc]
267   exact ne_of_gt (ENNReal.ofReal_pos.mpr (sub_pos.mpr hab))
268 
269 omit [IsProbabilityMeasure ν] [Fintype M] in
270 /-- The interval-uniform specialization includes the derived `h²/12`
271 coefficient, with `h = b-a`, rather than assuming a variance oracle. -/
272 theorem uniformInterval_affine_correction {a b : ℝ} (hab : a < b)
273     (K₀ K₁ : Matrix M N ℝ) (X : Matrix N N ℝ) :
274     matrixMean (ν := uniformInterval a b) (fun u =>
275       ((K₀ + u • K₁) - matrixMean (ν := uniformInterval a b) (fun v => K₀ + v • K₁)) * X *
276       ((K₀ + u • K₁) - matrixMean (ν := uniformInterval a b) (fun v => K₀ + v • K₁)).transpose) =
277       ((b - a) ^ 2 / 12) • (K₁ * X * K₁.transpose) := by
278   let := uniformInterval_isProbability hab
279   rw [affine_covariance_correction (uniformInterval_id_memLp hab), uniformInterval_variance hab]
280 
281 section BochnerBridge
282 
283 open scoped Matrix.Norms.Elementwise
284 
285 set_option backward.isDefEq.respectTransparency false in
286 omit [IsProbabilityMeasure ν] in
287 /-- Entrywise integrability gives integrability in the actual finite matrix
288 space with its standard product topology and elementwise norm. -/
289 theorem matrix_integrable_of_entries {K : Ω → Matrix M N ℝ}
290     (hK : ∀ i j, Integrable (fun ω => K ω i j) ν) : Integrable K ν := by
291   change Integrable (fun ω i j => K ω i j) ν
292   exact Integrable.of_eval (fun i => Integrable.of_eval (hK i))
293 
294 set_option backward.isDefEq.respectTransparency false in
295 omit [IsProbabilityMeasure ν] in
296 /-- The coordinate matrix integral is exactly the Bochner matrix integral. -/
297 theorem matrixMean_eq_integral {K : Ω → Matrix M N ℝ}
298     (hK : ∀ i j, Integrable (fun ω => K ω i j) ν) :
299     matrixMean (ν := ν) K = ∫ ω, K ω ∂ν := by
300   ext i j
301   change (∫ ω, K ω i j ∂ν) = (∫ ω, K ω ∂ν) i j
302   rw [eval_integral (fun i => Integrable.of_eval (hK i)) i]
303   exact (eval_integral (hK i) j).symm
304 
305 set_option backward.isDefEq.respectTransparency false in
306 /-- The manuscript's L2-entry hypothesis supplies matrix integrability; it
307 is not an extra input to the covariance theorem. -/
308 theorem matrix_integrable_of_memLp_entries {K : Ω → Matrix M N ℝ}
309     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν) : Integrable K ν :=
310   matrix_integrable_of_entries (fun i j => (hK i j).integrable (by norm_num))
311 
312 set_option backward.isDefEq.respectTransparency false in
313 /-- P57 written directly with Bochner matrix integrals. -/
314 theorem averaged_congruence_integral {K : Ω → Matrix M N ℝ}
315     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν) (X : Matrix N N ℝ) :
316     (∫ ω, K ω * X * (K ω).transpose ∂ν) =
317       (∫ ω, K ω ∂ν) * X * (∫ ω, K ω ∂ν).transpose +
318       ∫ ω, (K ω - ∫ y, K y ∂ν) * X * (K ω - ∫ y, K y ∂ν).transpose ∂ν := by
319   have hKi : ∀ i j, Integrable (fun ω => K ω i j) ν :=
320     fun i j => (hK i j).integrable (by norm_num)
321   have hc : ∀ i j, MemLp (fun ω => (K ω - matrixMean (ν := ν) K) i j) 2 ν := by
322     intro i j
323     exact (hK i j).sub (memLp_const _)
324   have h := averaged_congruence hK X
325   rw [matrixMean_eq_integral (fun i j => congruence_entry_integrable hK X i j),
326     matrixMean_eq_integral (fun i j => congruence_entry_integrable hc X i j),
327     matrixMean_eq_integral hKi] at h
328   exact h
329 
330 set_option backward.isDefEq.respectTransparency false in
331 /-- The positive covariance correction in the same Bochner-integral
332 dictionary as the manuscript. -/
333 theorem covariance_integral_posSemidef {K : Ω → Matrix M N ℝ}
334     (hK : ∀ i j, MemLp (fun ω => K ω i j) 2 ν)
335     {X : Matrix N N ℝ} (hX : X.PosSemidef) :
336     (∫ ω, (K ω - ∫ y, K y ∂ν) * X *
337       (K ω - ∫ y, K y ∂ν).transpose ∂ν).PosSemidef := by
338   have hc : ∀ i j, MemLp (fun ω => (K ω - matrixMean (ν := ν) K) i j) 2 ν := by
339     intro i j
340     exact (hK i j).sub (memLp_const _)
341   have h := covariance_correction_posSemidef hK hX
342   rw [matrixMean_eq_integral (fun i j => congruence_entry_integrable hc X i j),
343     matrixMean_eq_integral (fun i j => (hK i j).integrable (by norm_num))] at h
344   exact h
345 
346 end BochnerBridge
347 
348 end Komorebi.LSM.AveragedCongruence

Back to source controls