AutoSamplingTheory.ExampleCases.SmoothedPicardHMC.GaussianPowerMoment

- The actual Gaussian likelihood has finite qth moment with exact coefficient q(q-1)/(2 tau), and its normalized log-moment is q norm(x-y)^2/(2 tau).

theorem gaussian_power_moment {E : Type*} [NormedAddCommGroup E] [InnerProductSpace ℝ E]
    [FiniteDimensional ℝ E] [MeasurableSpace E] [BorelSpace E]
    (x y : E) (τ q : ℝ) (hτ : 0 < τ) (hq : 1 < q) :
    let G := fun a : E => (stdGaussian E).map (fun z => a + Real.sqrt τ • z)
    let L := fun z : E => Real.exp (inner ℝ (x-y) (z-y) / τ - ‖x-y‖ ^ 2 / (2*τ))
    Measurable L ∧ G x = (G y).withDensity (fun z => ENNReal.ofReal (L z)) ∧
      G x ≪ G y ∧ (G x).rnDeriv (G y) =ᵐ[G y] (fun z => ENNReal.ofReal (L z)) ∧
      Integrable (fun z => ((G x).rnDeriv (G y) z).toReal ^ q) (G y) ∧
      (∫ z, ((G x).rnDeriv (G y) z).toReal ^ q ∂G y) =
        Real.exp (q*(q-1)*‖x-y‖^2/(2*τ)) ∧
      Real.log (∫ z, ((G x).rnDeriv (G y) z).toReal ^ q ∂G y) / (q-1) =
        q*‖x-y‖^2/(2*τ) := by
  have linear_exp
      (h : E) (a : ℝ) :
      Integrable (fun z : E => Real.exp (a * inner ℝ h z)) (stdGaussian E) ∧
      (∫ z : E, Real.exp (a * inner ℝ h z) ∂stdGaussian E) =
        Real.exp (a^2 * ‖h‖^2 / 2) := by
    let l : StrongDual ℝ E := innerSL ℝ h
    have hmap : (stdGaussian E).map l = gaussianReal 0 (‖h‖^2).toNNReal := by
      rw [IsGaussian.map_eq_gaussianReal l,integral_strongDual_stdGaussian,
        variance_dual_stdGaussian]
      simp only [l,innerSL_apply_norm]
    have hI := integrable_exp_mul_gaussianReal (μ := 0) (v := (‖h‖^2).toNNReal) a
    rw [← hmap] at hI
    have hI' := (integrable_map_measure (by fun_prop) (by fun_prop : AEMeasurable l (stdGaussian E))).mp hI
    refine ⟨hI',?_⟩
    have hm := mgf_gaussianReal hmap a
    simpa [mgf,l,Real.toNNReal_of_nonneg (sq_nonneg ‖h‖),mul_comm] using hm
  have power_moment
      (x y : E) (τ q : ℝ) (hτ : 0 < τ) (hq : 1 < q) :
      let G := (stdGaussian E).map (fun z => y + Real.sqrt τ • z)
      let L := fun z : E => Real.exp (inner ℝ (x-y) (z-y) / τ - ‖x-y‖ ^ 2 / (2*τ))
      Integrable (fun z => L z ^ q) G ∧
        (∫ z, L z ^ q ∂G) = Real.exp (q*(q-1)*‖x-y‖^2/(2*τ)) := by
    let G := (stdGaussian E).map (fun z => y + Real.sqrt τ • z)
    let L := fun z : E => Real.exp (inner ℝ (x-y) (z-y) / τ - ‖x-y‖ ^ 2 / (2*τ))
    let s := Real.sqrt τ
    have hs : 0 < s := Real.sqrt_pos.mpr hτ
    have hs2 : s^2 = τ := Real.sq_sqrt hτ.le
    let a := q / s
    let b := -q * ‖x-y‖^2 / (2*τ)
    have hq0 : 0 ≤ q := by linarith
    have hfunc : (fun z : E => L (y + Real.sqrt τ • z)^q) =
        (fun z => Real.exp b * Real.exp (a * inner ℝ (x-y) z)) := by
-- Source excerpt truncated; follow the exact source link.

AutoSamplingTheory/ExampleCases/SmoothedPicardHMC/GaussianPowerMoment.lean:25published source at 0e31a3cda412

Excerpt truncated; the exact source link is authoritative.