ObservableMomentAlgebra.lean

1 import KomorebiLSM.ObservableMomentConditional
2 import KomorebiLSM.FiniteMomentIdentification
3 import KomorebiLSM.RegressionGeometry
4 
5 /-! Actual basis moments, exact coefficient inversion and directional
6 confidence propagation for P82. -/
7 set_option autoImplicit false
8 namespace Komorebi.LSM.ObservableMoments
9 open MeasureTheory Matrix FiniteMomentIdentification
10 open scoped BigOperators Matrix
11 
12 section Moments
13 variable {E A Y : Type*} [NormedAddCommGroup E] [NormedSpace ℝ E]
14 variable [FiniteDimensional ℝ E] [MeasurableSpace E] [BorelSpace E]
15 variable [MeasurableSpace Y] {C : BaseNorm.ProbabilityCone E} {ν : Measure Y}
16 variable (D : ContinuousInstrument.DensityInstrument (A := A) C ν) {p : ℕ}
17 
18 /-- The known integration matrix of the specified test and observation bases. -/
19 noncomputable def basisMomentMatrix (g ψ : Fin p → Y → ℝ) : Matrix (Fin p) (Fin p) ℝ :=
20   fun j k => ∫ y,g j y*ψ k y∂ν
21 
22 omit [FiniteDimensional ℝ E] [MeasurableSpace E] [BorelSpace E] in
23 /-- Finite basis-product integrability supplies the actual weighted operator
24 integrals rather than presuming the moment equation. -/
25 theorem weighted_operator_integrable (B : A → Fin p → E →L[ℝ] E) (ψ : Fin p → Y → ℝ)
26     (hB : ∀ a y,D.op a y=∑ k,ψ k y • B a k)
27     (g : Y → ℝ) (hgψ : ∀ k,Integrable (fun y => g y*ψ k y) ν) (a : A) :
28     Integrable (fun y => g y • D.op a y) ν := by
29   have hi (k : Fin p) : Integrable (fun y => (g y*ψ k y) • B a k) ν :=
30     (hgψ k).smul_const (B a k)
31   have hsum : Integrable (fun y => ∑ k : Fin p,(g y*ψ k y) • B a k) ν :=
32     integrable_finsetSum (f := fun k : Fin p => fun y : Y => (g y*ψ k y) • B a k)
33       (Finset.univ) (fun k _ => hi k)
34   have he : (fun y => g y • D.op a y)=
35       (fun y => ∑ k : Fin p,(g y*ψ k y) • B a k) := by
36     funext y
37     rw [hB a y,Finset.smul_sum]
38     simp only [smul_smul]
39   rw [he]
40   exact hsum
41 
42 omit [MeasurableSpace E] [BorelSpace E] in
43 /-- Integrating the submitted finite observation expansion gives the literal
44 moment mixing equation coefficient by coefficient. -/
45 theorem momentOperator_mixing (B : A → Fin p → E →L[ℝ] E)
46     (g ψ : Fin p → Y → ℝ) (hB : ∀ a y,D.op a y=∑ k,ψ k y • B a k)
47     (hgψ : ∀ j k,Integrable (fun y => g j y*ψ k y) ν) (a : A) :
48     (fun j => momentOperator D (g j) a)=mix (basisMomentMatrix (ν := ν) g ψ) (B a) := by
49   funext j
50   unfold momentOperator mix basisMomentMatrix
51   have he : (fun y => g j y • D.op a y)=
52       (fun y => ∑ k : Fin p,(g j y*ψ k y) • B a k) := by
53     funext y
54     rw [hB a y,Finset.smul_sum]
55     simp only [smul_smul]
56   rw [he,integral_finsetSum _ (fun k _ => (hgψ j k).smul_const (B a k))]
57   apply Finset.sum_congr rfl
58   intro k _
59   exact integral_smul_const _ _
60 
61 /-- The actual inverse matrix reconstructs any finite family in a real vector
62 space; invertibility is an explicit condition on the known integration matrix. -/
63 theorem inverse_mix {F : Type*} [AddCommGroup F] [Module ℝ F]
64     (W : Matrix (Fin p) (Fin p) ℝ) (hW : IsUnit W) (B : Fin p → F) :
65     mix W⁻¹ (mix W B)=B := by
66   rw [mix_mul,Matrix.nonsing_inv_mul W ((Matrix.isUnit_iff_isUnit_det W).mp hW),mix_one]
67 
68 omit [MeasurableSpace E] [BorelSpace E] in
69 /-- P82's coefficient formula follows from actual moment integration and the
70 known invertible basis matrix, not from an assumed reconstruction relation. -/
71 theorem coefficient_from_moments (B : A → Fin p → E →L[ℝ] E)
72     (g ψ : Fin p → Y → ℝ) (hB : ∀ a y,D.op a y=∑ k,ψ k y • B a k)
73     (hgψ : ∀ j k,Integrable (fun y => g j y*ψ k y) ν)
74     (hW : IsUnit (basisMomentMatrix (ν := ν) g ψ)) (a : A) :
75     B a=mix (basisMomentMatrix (ν := ν) g ψ)⁻¹ (fun j => momentOperator D (g j) a) := by
76   rw [momentOperator_mixing D B g ψ hB hgψ a,inverse_mix _ hW]
77 
78 end Moments
79 
80 section Directional
81 variable {p r : ℕ}
82 
83 theorem mixed_row_apply (W : Matrix (Fin p) (Fin p) ℝ)
84     (M : Fin p → Matrix (Fin r) (Fin r) ℝ) (k : Fin p) (i : Fin r) (x : Fin r → ℝ) :
85     (mix W M k *ᵥ x) i=∑ j,W k j*(M j *ᵥ x) i := by
86   simp only [mix,Matrix.mulVec,dotProduct,Matrix.sum_apply,Matrix.smul_apply,
87     smul_eq_mul,Finset.sum_mul,Finset.mul_sum]
88   rw [Finset.sum_comm]
89   apply Finset.sum_congr rfl
90   intro j _
91   apply Finset.sum_congr rfl
92   intro l _
93   ring
94 
95 /-- Genuine weighted Cauchy--Schwarz propagates row ellipsoids through the
96 known inverse linear map in every submitted direction. -/
97 theorem mixed_directional_error (W : Matrix (Fin p) (Fin p) ℝ)
98     (Mhat M : Fin p → Matrix (Fin r) (Fin r) ℝ) (V : Matrix (Fin r) (Fin r) ℝ)
99     (hV : V.PosDef) (β : Fin p → Fin r → ℝ)
100     (hb : ∀ j i,RegressionConfidence.weightedNorm V (fun l => Mhat j i l-M j i l)≤β j i)
101     (k : Fin p) (i : Fin r) (x : Fin r → ℝ) :
102     |(mix W Mhat k *ᵥ x) i-(mix W M k *ᵥ x) i|≤
103       ∑ j,|W k j| * β j i*RegressionConfidence.weightedNorm V⁻¹ x := by
104   have he : (mix W Mhat k *ᵥ x) i-(mix W M k *ᵥ x) i=
105       ∑ j,W k j*((fun l => Mhat j i l-M j i l) ⬝ᵥ x) := by
106     rw [mixed_row_apply,mixed_row_apply,← Finset.sum_sub_distrib]
107     apply Finset.sum_congr rfl
108     intro j _
109     rw [← mul_sub]
110     congr 1
111     simp only [Matrix.mulVec,dotProduct,sub_mul,Finset.sum_sub_distrib]
112   rw [he]
113   calc
114     _ ≤ ∑ j,|W k j*((fun l => Mhat j i l-M j i l) ⬝ᵥ x)| :=
115       Finset.abs_sum_le_sum_abs _ _
116     _ ≤ _ := by
117       apply Finset.sum_le_sum
118       intro j _
119       rw [abs_mul]
120       have hcs := RegressionConfidence.weighted_cauchy_schwarz hV
121         (fun l => Mhat j i l-M j i l) x
122       have hh := hcs.trans (mul_le_mul_of_nonneg_right (hb j i) (Real.sqrt_nonneg _))
123       exact (mul_le_mul_of_nonneg_left hh (abs_nonneg _)).trans_eq (mul_assoc _ _ _).symm
124 
125 end Directional
126 end Komorebi.LSM.ObservableMoments

Back to source controls