Stein's lemma and James–Stein estimation

TheoremStein's lemma#✓fully proved

For X ~ N(μ, v) and an everywhere differentiable g with integrable derivative, (X - μ) g(X) is integrable and E[(X - μ) g(X)] = v E[g'(X)]. For v = 0 both sides vanish.

theorem integrable_and_integral_sub_mul_gaussianReal (μ : ℝ) (v : NNReal) {g g' : ℝ → ℝ}
    (hg : ∀ x, HasDerivAt g (g' x) x) (hg' : Integrable g' (gaussianReal μ v)) :
    Integrable (fun x => (x - μ) * g x) (gaussianReal μ v) ∧
      ∫ x, (x - μ) * g x ∂gaussianReal μ v = v * ∫ x, g' x ∂gaussianReal μ v
Used by 1

TheoremStein's lemma#✓fully proved

Stein's lemma, coordinate form: for X ~ N(θ, I_p) and a function h whose slices along the i-th axis are differentiable everywhere (for almost every base point) with integrable partial derivative dh, E[(Xᵢ - θᵢ) h(X)] = E[dh(X)] whenever (Xᵢ - θᵢ) h(X) is integrable. Unlike the one-dimensional case, this integrability does not follow from that of dh: h(x) = exp(x₁²) has ∂₀ h = 0 but x₀ h(x) is not integrable.

theorem integral_sub_mul_pi_gaussianReal {p : ℕ} (θ : Fin p → ℝ) (i : Fin p)
    {h dh : (Fin p → ℝ) → ℝ}
    (hd : ∀ᵐ x ∂(Measure.pi fun j => gaussianReal (θ j) 1),
      ∀ t, HasDerivAt (fun s => h (Function.update x i s)) (dh (Function.update x i t)) t)
    (hdi : Integrable dh (Measure.pi fun j => gaussianReal (θ j) 1))
    (hhi : Integrable (fun x => (x i - θ i) * h x) (Measure.pi fun j => gaussianReal (θ j) 1)) :
    ∫ x, (x i - θ i) * h x ∂(Measure.pi fun j => gaussianReal (θ j) 1) =
      ∫ x, dh x ∂(Measure.pi fun j => gaussianReal (θ j) 1)
Uses 1 · Used by 1

TheoremQuadratic risk of the usual Gaussian mean estimator#✓fully proved

The quadratic risk of the usual estimator X of the mean of N(θ, I_p) is p.

theorem integral_sum_sq_sub_pi_gaussianReal {p : ℕ} (θ : Fin p → ℝ) :
    ∫ x, ∑ i, (x i - θ i) ^ 2 ∂(Measure.pi fun j => gaussianReal (θ j) 1) = p
Used by 2

TheoremStein's unbiased risk estimate#✓fully proved

For X ~ N(θ, I_p) and an estimator X + g(X) whose coordinates gᵢ are square integrable and satisfy the hypotheses of Stein's lemma in the i-th coordinate with partial derivative dgᵢ, the quadratic risk is p + E[∑ gᵢ(X)² + 2 ∑ dgᵢ(X)].

theorem integral_sum_sq_add_sub_pi_gaussianReal {p : ℕ} (θ : Fin p → ℝ)
    {g dg : Fin p → (Fin p → ℝ) → ℝ}
    (hd : ∀ i, ∀ᵐ x ∂(Measure.pi fun j => gaussianReal (θ j) 1),
      ∀ t, HasDerivAt (fun s => g i (Function.update x i s)) (dg i (Function.update x i t)) t)
    (hdi : ∀ i, Integrable (dg i) (Measure.pi fun j => gaussianReal (θ j) 1))
    (hg2 : ∀ i, MemLp (g i) 2 (Measure.pi fun j => gaussianReal (θ j) 1)) :
    ∫ x, ∑ i, (x i + g i x - θ i) ^ 2 ∂(Measure.pi fun j => gaussianReal (θ j) 1) =
      p + ∫ x, ∑ i, (g i x ^ 2 + 2 * dg i x) ∂(Measure.pi fun j => gaussianReal (θ j) 1)
Uses 2 · Used by 1

DefinitionJames–Stein estimator#✓fully proved

The James–Stein estimator δ(x) = (1 - (p - 2) / ‖x‖²) x of the mean of N(θ, I_p). At x = 0 the factor is 1 by Lean's a / 0 = 0 convention, so δ(0) = 0.

noncomputable def jamesStein {p : ℕ} (x : Fin p → ℝ) : Fin p → ℝ :=
  fun i => (1 - ((p : ℝ) - 2) / ∑ j, x j ^ 2) * x i
Used by 2

TheoremJames–Stein#✓fully proved

James–Stein (1961): for p ≥ 3 and X ~ N(θ, I_p), ‖X‖⁻² is integrable and the risk of the James–Stein estimator is p - (p - 2)² E[‖X‖⁻²].

theorem integrable_and_integral_jamesStein_risk {p : ℕ} (hp : 3 ≤ p) (θ : Fin p → ℝ) :
    Integrable (fun x : Fin p → ℝ => (∑ j, x j ^ 2)⁻¹)
        (Measure.pi fun j => gaussianReal (θ j) 1) ∧
      ∫ x, ∑ i, (jamesStein x i - θ i) ^ 2 ∂(Measure.pi fun j => gaussianReal (θ j) 1) =
        p - ((p : ℝ) - 2) ^ 2 *
          ∫ x, (∑ j, x j ^ 2)⁻¹ ∂(Measure.pi fun j => gaussianReal (θ j) 1)
Uses 2 · Used by 1

TheoremJames–Stein dominance#✓fully proved

James–Stein dominance (Stein 1956, James–Stein 1961): for p ≥ 3, the James–Stein estimator has strictly smaller quadratic risk than X at every mean θ.

theorem integral_jamesStein_risk_lt {p : ℕ} (hp : 3 ≤ p) (θ : Fin p → ℝ) :
    ∫ x, ∑ i, (jamesStein x i - θ i) ^ 2 ∂(Measure.pi fun j => gaussianReal (θ j) 1) <
      ∫ x, ∑ i, (x i - θ i) ^ 2 ∂(Measure.pi fun j => gaussianReal (θ j) 1)
Uses 3
7 fully proved

The dependency explorer needs JavaScript and its graph data.

Download the complete graph data.