Discovery.Li­NGAM

Linear non-Gaussian acyclic model identification: permutation uniqueness, generalized-permutation algebra, kurtosis-based column support, and the LiNGAM identifiability theorem.

Kurtosis 3 core · 0 supporting This is the Marcinkiewicz-free route to LiNGAM identifiability. ★ kurt★ cross_fourth_cumulant_eq_sum★ colSupport_of_kurtosis

Kurtosis-based column support for ICA / LiNGAM

This is the Marcinkiewicz-free route to LiNGAM identifiability. Instead of the general Darmois–Skitovich theorem (which needs Marcinkiewicz), we assume the sources have non-zero fourth cumulant of one common sign (a standard ICA assumption: all super-Gaussian or all sub-Gaussian). Then the key column-support fact Wᵢⱼ · Wₖⱼ = 0 — the input to genPerm_of_det_ne_zero_of_colSupport — follows from a single fourth-cumulant identity and a sum-of-same-sign argument, with no characteristic-function functional equation at all.

cross_fourth_cumulant_eq_sum is the multilinear identity cum(yᵢ,yᵢ,yₖ,yₖ) = Σⱼ Wᵢⱼ² Wₖⱼ² κ₄(eⱼ) for yᵢ = Σⱼ Wᵢⱼ eⱼ with independent centered sources; colSupport_of_kurtosis combines it with independence of yᵢ, yₖ (which makes the cross-cumulant vanish) and the same-sign assumption.

def kurt reviewed
Causalean.Discovery.LiNGAM

Excess kurtosis, the fourth cumulant of a real random variable. For a centered variable X it is E[X⁴] − 3 · (E[X²])²; Gaussian variables have zero excess kurtosis, while LiNGAM's kurtosis route assumes every source has nonzero excess kurtosis of one common sign.

Definition (Lean source)
noncomputable def kurt (X : Ω → ℝ) (P : Measure Ω) : ℝ := (∫ ω, (X ω) ^ 4 ∂P) - 3 * (∫ ω, (X ω) ^ 2 ∂P) ^ 2
Causalean.Discovery.LiNGAM.kurt · Causalean/Discovery/LiNGAM/Kurtosis.lean:35
theorem cross_fourth_cumulant_eq_sum reviewed
Causalean.Discovery.LiNGAM

Fourth cross-cumulant identity (Isserlis / cumulant multilinearity). Let e be a family of real sources on a probability space such that each coordinate eⱼ is measurable, the coordinates are mutually independent, each has finite fourth moment, and each is centered. Then for any coefficient vectors a, b and the linear forms yₐ = Σⱼ aⱼ eⱼ, y_b = Σⱼ bⱼ eⱼ, the joint fourth cumulant cum(yₐ,yₐ,y_b,y_b) equals Σⱼ aⱼ² bⱼ² κ₄(eⱼ).

Formal statement
n :
Ω → Fin n → ℝ
a b :
Fin n → ℝ
hmeas :
∀ j, Measurable (fun ω => e ω j)
hindep :
iIndepFun (fun j ω => e ω j) P
hL4 :
∀ j, MemLp (fun ω => e ω j) 4 P
hcent :
∀ j, ∫ ω, e ω j ∂P = 0
(∫ ω, (∑ j, a j * e ω j) ^ 2 * (∑ j, b j * e ω j) ^ 2 ∂P)
- (∫ ω, (∑ j, a j * e ω j) ^ 2 ∂P) * (∫ ω, (∑ j, b j * e ω j) ^ 2 ∂P)
- 2 * (∫ ω, (∑ j, a j * e ω j) * (∑ j, b j * e ω j) ∂P) ^ 2
= ∑ j, (a j) ^ 2 * (b j) ^ 2 * kurt (fun ω => e ω j) P
Proof (Lean source)
theorem cross_fourth_cumulant_eq_sum {n : ℕ} {e : Ω → Fin n → ℝ} (a b : Fin n → ℝ) (hmeas : ∀ j, Measurable (fun ω => e ω j)) (hindep : iIndepFun (fun j ω => e ω j) P) (hL4 : ∀ j, MemLp (fun ω => e ω j) 4 P) (hcent : ∀ j, ∫ ω, e ω j ∂P = 0) : (∫ ω, (∑ j, a j * e ω j) ^ 2 * (∑ j, b j * e ω j) ^ 2 ∂P) - (∫ ω, (∑ j, a j * e ω j) ^ 2 ∂P) * (∫ ω, (∑ j, b j * e ω j) ^ 2 ∂P) - 2 * (∫ ω, (∑ j, a j * e ω j) * (∑ j, b j * e ω j) ∂P) ^ 2 = ∑ j, (a j) ^ 2 * (b j) ^ 2 * kurt (fun ω => e ω j) P := by classical let m : Fin n → ℝ := fun j => ∫ ω, (e ω j) ^ 2 ∂P let q : Fin n → ℝ := fun j => ∫ ω, (e ω j) ^ 4 ∂P have hcross : ∀ c d : Fin n → ℝ, ∫ ω, (∑ j, c j * e ω j) * (∑ j, d j * e ω j) ∂P = ∑ j, c j * d j * m j := by intro c d have h2 : ∀ j, MemLp (fun ω => e ω j) 2 P := fun j => (hL4 j).mono_exponent (by norm_num) have hint_pair : ∀ i j, Integrable (fun ω => (c i * e ω i) * (d j * e ω j)) P := by intro i j have hmul : MemLp (fun ω => e ω i * e ω j) 1 P := by simpa only [Pi.mul_apply] using (h2 j).mul' (h2 i) convert (memLp_one_iff_integrable.mp hmul).const_mul (c i * d j) using 1 ext ω ring have hpair : ∀ i j, ∫ ω, e ω i * e ω j ∂P = if i = j then m i else 0 := by intro i j by_cases hij : i = j · subst j simp [m, pow_two] · have h_ind : IndepFun (fun ω => e ω i) (fun ω => e ω j) P := by simpa using hindep.indepFun (i := i) (j := j) hij have h := h_ind.integral_fun_mul_eq_mul_integral (h2 i).aestronglyMeasurable (h2 j).aestronglyMeasurable simp [hij, hcent i, hcent j, h] have hterm : ∀ i j, ∫ ω, (c i * e ω i) * (d j * e ω j) ∂P = c i * d j * ∫ ω, e ω i * e ω j ∂P := by intro i j rw [← MeasureTheory.integral_const_mul] congr 1 ext ω ring calc ∫ ω, (∑ j, c j * e ω j) * (∑ j, d j * e ω j) ∂P = ∫ ω, ∑ i, ∑ j, (c i * e ω i) * (d j * e ω j) ∂P := by congr 1 ext ω rw [Finset.sum_mul_sum] _ = ∑ i, ∑ j, ∫ ω, (c i * e ω i) * (d j * e ω j) ∂P := by rw [MeasureTheory.integral_finset_sum univ] · congr with i rw [MeasureTheory.integral_finset_sum univ] intro j _ exact hint_pair i j · intro i _ exact integrable_finset_sum univ (fun j _ => hint_pair i j) _ = ∑ i, ∑ j, c i * d j * (∫ ω, e ω i * e ω j ∂P) := by apply Finset.sum_congr rfl intro i _ apply Finset.sum_congr rfl intro j _ exact hterm i j _ = ∑ i, ∑ j, c i * d j * (if i = j then m i else 0) := by apply Finset.sum_congr rfl intro i _ apply Finset.sum_congr rfl intro j _ rw [hpair i j] _ = ∑ i, c i * d i * m i := by simp [m, mul_assoc] have hAA : ∫ ω, (∑ j, a j * e ω j) ^ 2 ∂P = ∑ j, (a j) ^ 2 * m j := by simpa [pow_two, mul_assoc] using hcross a a have hBB : ∫ ω, (∑ j, b j * e ω j) ^ 2 ∂P = ∑ j, (b j) ^ 2 * m j := by simpa [pow_two, mul_assoc] using hcross b b have hAB : ∫ ω, (∑ j, a j * e ω j) * (∑ j, b j * e ω j) ∂P = ∑ j, a j * b j * m j := hcross a b have h4 : ∫ ω, (∑ j, a j * e ω j) ^ 2 * (∑ j, b j * e ω j) ^ 2 ∂P = (∑ j, (a j) ^ 2 * (b j) ^ 2 * q j) + ((∑ j, (a j) ^ 2 * m j) * (∑ j, (b j) ^ 2 * m j) - ∑ j, (a j) ^ 2 * (b j) ^ 2 * (m j) ^ 2) + 2 * ((∑ j, a j * b j * m j) ^ 2 - ∑ j, (a j * b j * m j) ^ 2) := by -- Fourth-moment expansion under independence and centering: -- only the all-equal index and the three pair-partition patterns survive. have hI4 : ∀ i i' k k' : Fin n, ∫ ω, e ω i * e ω i' * e ω k * e ω k' ∂P = ((if i = i' ∧ i = k ∧ i = k' then q i else 0) + (if i = i' ∧ k = k' ∧ i ≠ k then m i * m k else 0) + (if i = k ∧ i' = k' ∧ i ≠ i' then m i * m i' else 0) + (if i = k' ∧ i' = k ∧ i ≠ i' then m i * m i' else 0)) := by have h2 : ∀ j, MemLp (fun ω => e ω j) 2 P := fun j => (hL4 j).mono_exponent (by norm_num) have hpair : ∀ x y, ∫ ω, e ω x * e ω y ∂P = if x = y then m x else 0 := by intro x y by_cases hxy : x = y · subst y simp [m, pow_two] · have h_ind : IndepFun (fun ω => e ω x) (fun ω => e ω y) P := by simpa using hindep.indepFun (i := x) (j := y) hxy have h := h_ind.integral_fun_mul_eq_mul_integral (h2 x).aestronglyMeasurable (h2 y).aestronglyMeasurable simp [hxy, hcent x, hcent y, h] have hpairpair_eval : ∀ x y z w, x ≠ z → x ≠ w → y ≠ z → y ≠ w → ∫ ω, e ω x * e ω y * e ω z * e ω w ∂P = (if x = y then m x else 0) * (if z = w then m z else 0) := by intro x y z w hxz hxw hyz hyw have h_ind : IndepFun (fun ω => e ω x * e ω y) (fun ω => e ω z * e ω w) P := by exact hindep.indepFun_mul_mul hmeas x y z w hxz hxw hyz hyw haveI h442 : ENNReal.HolderTriple 4 4 2 := by change ENNReal.HolderTriple ((4 : NNReal) : ENNReal) ((4 : NNReal) : ENNReal) ((2 : NNReal) : ENNReal) exact NNReal.HolderTriple.coe_ennreal (by norm_num) (by constructor <;> norm_num : NNReal.HolderTriple 4 4 2) have hxy_mem : MemLp (fun ω => e ω x * e ω y) 2 P := by simpa only [Pi.mul_apply] using (hL4 y).mul' (hL4 x) have hzw_mem : MemLp (fun ω => e ω z * e ω w) 2 P := by simpa only [Pi.mul_apply] using (hL4 w).mul' (hL4 z) calc ∫ ω, e ω x * e ω y * e ω z * e ω w ∂P = ∫ ω, (e ω x * e ω y) * (e ω z * e ω w) ∂P := by congr 1 ext ω ring _ = (∫ ω, e ω x * e ω y ∂P) * (∫ ω, e ω z * e ω w ∂P) := by exact h_ind.integral_fun_mul_eq_mul_integral hxy_mem.aestronglyMeasurable hzw_mem.aestronglyMeasurable _ = (if x = y then m x else 0) * (if z = w then m z else 0) := by rw [hpair x y, hpair z w] have htwoPair : ∀ x y, x ≠ y → ∫ ω, e ω x * e ω x * e ω y * e ω y ∂P = m x * m y := by intro x y hxy simpa using hpairpair_eval x x y y hxy hxy hxy hxy have hcube : ∀ x y, x ≠ y → ∫ ω, e ω x * e ω y * e ω y * e ω y ∂P = 0 := by intro x y hxy have h_ind : IndepFun (fun ω => e ω x) (fun ω => (e ω y) ^ 3) P := by have h0 : IndepFun (fun ω => e ω x) (fun ω => e ω y) P := by simpa using hindep.indepFun (i := x) (j := y) hxy simpa [Function.comp_def] using h0.comp measurable_id (measurable_id.pow_const 3) have h := h_ind.integral_fun_mul_eq_mul_integral (h2 x).aestronglyMeasurable ((hmeas y).pow_const 3).aestronglyMeasurable calc ∫ ω, e ω x * e ω y * e ω y * e ω y ∂P = ∫ ω, e ω x * (e ω y) ^ 3 ∂P := by congr 1 ext ω ring _ = 0 := by simpa [hcent x] using h have hall : ∀ x, ∫ ω, e ω x * e ω x * e ω x * e ω x ∂P = q x := by intro x change ∫ ω, e ω x * e ω x * e ω x * e ω x ∂P = ∫ ω, (e ω x) ^ 4 ∂P congr 1 ext ω ring intro i i' k k' by_cases hii' : i = i' · subst i' by_cases hik : i = k · subst k by_cases hik' : i = k' · subst k' calc ∫ ω, e ω i * e ω i * e ω i * e ω i ∂P = q i := hall i _ = ((if i = i ∧ i = i ∧ i = i then q i else 0) + (if i = i ∧ i = i ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i = i ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i = i ∧ i ≠ i then m i * m i else 0)) := by simp · calc ∫ ω, e ω i * e ω i * e ω i * e ω k' ∂P = 0 := by simpa [mul_comm, mul_left_comm, mul_assoc] using hcube k' i (Ne.symm hik') _ = ((if i = i ∧ i = i ∧ i = k' then q i else 0) + (if i = i ∧ i = k' ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i = k' ∧ i ≠ i then m i * m i else 0) + (if i = k' ∧ i = i ∧ i ≠ i then m i * m i else 0)) := by simp [hik'] · by_cases hkk' : k = k' · subst k' calc ∫ ω, e ω i * e ω i * e ω k * e ω k ∂P = m i * m k := htwoPair i k hik _ = ((if i = i ∧ i = k ∧ i = k then q i else 0) + (if i = i ∧ k = k ∧ i ≠ k then m i * m k else 0) + (if i = k ∧ i = k ∧ i ≠ i then m i * m i else 0) + (if i = k ∧ i = k ∧ i ≠ i then m i * m i else 0)) := by simp [hik] · by_cases hik' : i = k' · subst k' calc ∫ ω, e ω i * e ω i * e ω k * e ω i ∂P = 0 := by simpa [mul_comm, mul_left_comm, mul_assoc] using hcube k i (Ne.symm hik) _ = ((if i = i ∧ i = k ∧ i = i then q i else 0) + (if i = i ∧ k = i ∧ i ≠ k then m i * m k else 0) + (if i = k ∧ i = i ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i = k ∧ i ≠ i then m i * m i else 0)) := by simp [hik, hkk'] · calc ∫ ω, e ω i * e ω i * e ω k * e ω k' ∂P = 0 := by have h := hpairpair_eval k k' i i (Ne.symm hik) (Ne.symm hik) (Ne.symm hik') (Ne.symm hik') calc ∫ ω, e ω i * e ω i * e ω k * e ω k' ∂P = ∫ ω, e ω k * e ω k' * e ω i * e ω i ∂P := by congr 1 ext ω ring _ = 0 := by simpa [hkk'] using h _ = ((if i = i ∧ i = k ∧ i = k' then q i else 0) + (if i = i ∧ k = k' ∧ i ≠ k then m i * m k else 0) + (if i = k ∧ i = k' ∧ i ≠ i then m i * m i else 0) + (if i = k' ∧ i = k ∧ i ≠ i then m i * m i else 0)) := by simp [hik, hik', hkk'] · by_cases hik : i = k · subst k by_cases hi'k' : i' = k' · subst k' calc ∫ ω, e ω i * e ω i' * e ω i * e ω i' ∂P = m i * m i' := by simpa [mul_comm, mul_left_comm, mul_assoc] using htwoPair i i' hii' _ = ((if i = i' ∧ i = i ∧ i = i' then q i else 0) + (if i = i' ∧ i = i' ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i' = i' ∧ i ≠ i' then m i * m i' else 0) + (if i = i' ∧ i' = i ∧ i ≠ i' then m i * m i' else 0)) := by simp [hii'] · by_cases hik' : i = k' · subst k' calc ∫ ω, e ω i * e ω i' * e ω i * e ω i ∂P = 0 := by simpa [mul_comm, mul_left_comm, mul_assoc] using hcube i' i (Ne.symm hii') _ = ((if i = i' ∧ i = i ∧ i = i then q i else 0) + (if i = i' ∧ i = i ∧ i ≠ i then m i * m i else 0) + (if i = i ∧ i' = i ∧ i ≠ i' then m i * m i' else 0) + (if i = i ∧ i' = i ∧ i ≠ i' then m i * m i' else 0)) := by simp [hii', eq_comm] · calc ∫ ω, e ω i * e ω i' * e ω i * e ω k' ∂P = 0 := by have h := hpairpair_eval i' k' i i (Ne.symm hii') (Ne.symm hii') (Ne.symm hik') (Ne.symm hik') -- … truncated; follow the source link for the rest …
Causalean.Discovery.LiNGAM.cross_fourth_cumulant_eq_sum · Causalean/Discovery/LiNGAM/Kurtosis.lean:42 · uses kurt
theorem colSupport_of_kurtosis reviewed
Causalean.Discovery.LiNGAM

Kurtosis-based column support. Let e be a family of real sources on a probability space such that each coordinate eⱼ is measurable, the coordinates are mutually independent, each has finite fourth moment, each is centered, and the fourth cumulant (excess kurtosis) of every coordinate is nonzero and of one common sign, all positive or all negative. For two distinct row indices i ≠ k of a mixing matrix W, if the linear forms Σⱼ Wᵢⱼ eⱼ and Σⱼ Wₖⱼ eⱼ are independent, then every column j satisfies Wᵢⱼ · Wₖⱼ = 0. This is the input required by genPerm_of_det_ne_zero_of_colSupport.

Formal statement
n :
Ω → Fin n → ℝ
Matrix (Fin n) (Fin n) ℝ
hmeas :
∀ j, Measurable (fun ω => e ω j)
hindep :
iIndepFun (fun j ω => e ω j) P
hL4 :
∀ j, MemLp (fun ω => e ω j) 4 P
hcent :
∀ j, ∫ ω, e ω j ∂P = 0
hsign :
(∀ j, 0 < kurt (fun ω => e ω j) P) ∨ (∀ j, kurt (fun ω => e ω j) P < 0)
i k :
Fin n
hik :
i ≠ k
hyindep :
IndepFun (fun ω => ∑ j, W i j * e ω j) (fun ω => ∑ j, W k j * e ω j) P
∀ j, W i j * W k j = 0
Proof (Lean source)
theorem colSupport_of_kurtosis {n : ℕ} {e : Ω → Fin n → ℝ} {W : Matrix (Fin n) (Fin n) ℝ} (hmeas : ∀ j, Measurable (fun ω => e ω j)) (hindep : iIndepFun (fun j ω => e ω j) P) (hL4 : ∀ j, MemLp (fun ω => e ω j) 4 P) (hcent : ∀ j, ∫ ω, e ω j ∂P = 0) (hsign : (∀ j, 0 < kurt (fun ω => e ω j) P) ∨ (∀ j, kurt (fun ω => e ω j) P < 0)) {i k : Fin n} (hik : i ≠ k) (hyindep : IndepFun (fun ω => ∑ j, W i j * e ω j) (fun ω => ∑ j, W k j * e ω j) P) : ∀ j, W i j * W k j = 0 := by classical let y : Ω → ℝ := fun ω => ∑ j, W i j * e ω j let z : Ω → ℝ := fun ω => ∑ j, W k j * e ω j have hy_meas : Measurable y := by simp only [y] exact measurable_sum Finset.univ fun j _ => (hmeas j).const_mul _ have hz_meas : Measurable z := by simp only [z] exact measurable_sum Finset.univ fun j _ => (hmeas j).const_mul _ have hy_int : ∫ ω, y ω ∂P = 0 := by have hterm : ∀ j ∈ Finset.univ, Integrable (fun ω => W i j * e ω j) P := by intro j _ have h14 : (1 : ENNReal) ≤ (4 : ENNReal) := by norm_num exact (memLp_one_iff_integrable.mp ((hL4 j).mono_exponent h14)).const_mul _ calc ∫ ω, y ω ∂P = ∑ j, ∫ ω, W i j * e ω j ∂P := by simp only [y] rw [MeasureTheory.integral_finset_sum Finset.univ hterm] _ = ∑ j, W i j * ∫ ω, e ω j ∂P := by congr 1 ext j rw [MeasureTheory.integral_const_mul] _ = 0 := by simp [hcent] have hz_int : ∫ ω, z ω ∂P = 0 := by have hterm : ∀ j ∈ Finset.univ, Integrable (fun ω => W k j * e ω j) P := by intro j _ have h14 : (1 : ENNReal) ≤ (4 : ENNReal) := by norm_num exact (memLp_one_iff_integrable.mp ((hL4 j).mono_exponent h14)).const_mul _ calc ∫ ω, z ω ∂P = ∑ j, ∫ ω, W k j * e ω j ∂P := by simp only [z] rw [MeasureTheory.integral_finset_sum Finset.univ hterm] _ = ∑ j, W k j * ∫ ω, e ω j ∂P := by congr 1 ext j rw [MeasureTheory.integral_const_mul] _ = 0 := by simp [hcent] have hyz_int : ∫ ω, y ω * z ω ∂P = 0 := by have h := hyindep.integral_fun_mul_eq_mul_integral hy_meas.aestronglyMeasurable hz_meas.aestronglyMeasurable simpa [y, z, hy_int, hz_int] using h have hy2z2_int : (∫ ω, y ω ^ 2 * z ω ^ 2 ∂P) = (∫ ω, y ω ^ 2 ∂P) * (∫ ω, z ω ^ 2 ∂P) := by have hsq_indep : IndepFun (fun ω => y ω ^ 2) (fun ω => z ω ^ 2) P := by exact hyindep.comp (measurable_id.pow_const 2) (measurable_id.pow_const 2) exact hsq_indep.integral_fun_mul_eq_mul_integral (hy_meas.pow_const 2).aestronglyMeasurable (hz_meas.pow_const 2).aestronglyMeasurable have hsum0 : (∑ j, (W i j) ^ 2 * (W k j) ^ 2 * kurt (fun ω => e ω j) P) = 0 := by have hcum := cross_fourth_cumulant_eq_sum (P := P) (e := e) (fun j => W i j) (fun j => W k j) hmeas hindep hL4 hcent simpa [y, z, hy2z2_int, hyz_int] using hcum.symm intro j rcases hsign with hpos | hneg · have hterm0 : (W i j) ^ 2 * (W k j) ^ 2 * kurt (fun ω => e ω j) P = 0 := by have hnonneg : ∀ x ∈ (Finset.univ : Finset (Fin n)), 0 ≤ (W i x) ^ 2 * (W k x) ^ 2 * kurt (fun ω => e ω x) P := by intro x _ have hsq : 0 ≤ (W i x) ^ 2 * (W k x) ^ 2 := mul_nonneg (sq_nonneg _) (sq_nonneg _) exact mul_nonneg hsq (le_of_lt (hpos x)) exact (Finset.sum_eq_zero_iff_of_nonneg hnonneg).mp hsum0 j (Finset.mem_univ j) have hkurt_ne : kurt (fun ω => e ω j) P ≠ 0 := ne_of_gt (hpos j) have hsqsq : (W i j) ^ 2 * (W k j) ^ 2 = 0 := by exact (mul_eq_zero.mp hterm0).resolve_right hkurt_ne have hpowsq : (W i j * W k j) ^ 2 = 0 := by simpa [mul_pow] using hsqsq exact eq_zero_of_pow_eq_zero hpowsq · have hsum0_neg : (∑ x, -((W i x) ^ 2 * (W k x) ^ 2 * kurt (fun ω => e ω x) P)) = 0 := by rw [Finset.sum_neg_distrib, hsum0, neg_zero] have hterm0_neg : -((W i j) ^ 2 * (W k j) ^ 2 * kurt (fun ω => e ω j) P) = 0 := by have hnonneg : ∀ x ∈ (Finset.univ : Finset (Fin n)), 0 ≤ -((W i x) ^ 2 * (W k x) ^ 2 * kurt (fun ω => e ω x) P) := by intro x _ have hsq : 0 ≤ (W i x) ^ 2 * (W k x) ^ 2 := mul_nonneg (sq_nonneg _) (sq_nonneg _) have hk : kurt (fun ω => e ω x) P ≤ 0 := le_of_lt (hneg x) exact neg_nonneg.mpr (mul_nonpos_of_nonneg_of_nonpos hsq hk) exact (Finset.sum_eq_zero_iff_of_nonneg hnonneg).mp hsum0_neg j (Finset.mem_univ j) have hterm0 : (W i j) ^ 2 * (W k j) ^ 2 * kurt (fun ω => e ω j) P = 0 := by linarith have hkurt_ne : kurt (fun ω => e ω j) P ≠ 0 := ne_of_lt (hneg j) have hsqsq : (W i j) ^ 2 * (W k j) ^ 2 = 0 := by exact (mul_eq_zero.mp hterm0).resolve_right hkurt_ne have hpowsq : (W i j * W k j) ^ 2 = 0 := by simpa [mul_pow] using hsqsq exact eq_zero_of_pow_eq_zero hpowsq
Causalean.Discovery.LiNGAM.colSupport_of_kurtosis · Causalean/Discovery/LiNGAM/Kurtosis.lean:611 · uses kurt
Identifiability 1 core · 0 supporting This is the structural linear-algebra step in the LiNGAM identification argument (Shimizu et al. ★ lingam_identifiable

LiNGAM structural identifiability

This is the structural linear-algebra step in the LiNGAM identification argument (Shimizu et al. 2006): once two coefficient matrices are related by a generalized permutation, acyclicity and unit-diagonal normalization pin the coefficient matrix uniquely.

Write C = I − B for the coefficient matrix of a LiNGAM model x = B x + e (so x = C⁻¹ e, C has unit diagonal because Bᵢᵢ = 0, and C is lower triangular in the causal order because B is acyclic). ICA identifiability says that two observationally-equivalent models have mixing matrices related by a generalized permutation, i.e. C' = W C with W a permutation composed with a non-zero diagonal scaling: C' i j = d i · C (τ i) j with d i ≠ 0.

The public theorem lingam_identifiable specializes Causalean.Mathlib.LinearAlgebra.eq_of_genPerm_triangular_unitDiag to C = I − B. The kurtosis route in LiNGAMKurtosis.lean supplies the generalized-permutation relation through column support, so this file isolates the deterministic acyclicity-and-unit-diagonal pinning step.

theorem lingam_identifiable reviewed
Causalean.Discovery.LiNGAM

LiNGAM coefficient identifiability (modulo ICA). Let B and B' be two n × n real coefficient matrices such that B has zero diagonal, B' has zero diagonal, B is acyclic with respect to a causal order σ: B i j = 0 whenever σ i < σ j, and the matrices I − B and I − B' are related by a generalized permutation — a permutation τ composed with a nonzero diagonal scale d, so (I − B') i j = d i · (I − B) (τ i) j (as delivered by ICA identifiability). Then B equals B'. This packages eq_of_genPerm_triangular_unitDiag for C = I − B.

Formal statement
n :
Matrix (Fin n) (Fin n) ℝ
hBdiag :
∀ i, B i i = 0
hB'diag :
∀ i, B' i i = 0
σ :
Perm (Fin n)
hBacyc :
∀ i j
if
σ i < σ j
then
B i j = 0
τ :
Perm (Fin n)
d :
Fin n → ℝ
hICA :
∀ i j, (1 - B') i j = d i * (1 - B) (τ i) j
B = B'
Proof (Lean source)
theorem lingam_identifiable {n : ℕ} {B B' : Matrix (Fin n) (Fin n) ℝ} (hBdiag : ∀ i, B i i = 0) (hB'diag : ∀ i, B' i i = 0) {σ : Perm (Fin n)} (hBacyc : ∀ i j, σ i < σ j → B i j = 0) {τ : Perm (Fin n)} {d : Fin n → ℝ} (hICA : ∀ i j, (1 - B') i j = d i * (1 - B) (τ i) j) : B = B' := by have hC : (1 - B) = (1 - B') := by refine eq_of_genPerm_triangular_unitDiag (σ := σ) (τ := τ) (d := d) ?_ ?_ ?_ hICA · intro i; rw [Matrix.sub_apply, Matrix.one_apply_eq, hBdiag i, sub_zero] · intro i; rw [Matrix.sub_apply, Matrix.one_apply_eq, hB'diag i, sub_zero] · intro i j hlt have hij : i ≠ j := fun h => (ne_of_lt hlt) (congrArg σ h) rw [Matrix.sub_apply, Matrix.one_apply_ne hij, hBacyc i j hlt, sub_zero] exact sub_right_injective hC
Causalean.Discovery.LiNGAM.lingam_identifiable · Causalean/Discovery/LiNGAM/Identifiability.lean:36
Li­NGAMKurtosis 2 core · 0 supporting This file assembles a LiNGAM identification theorem under the standard ICA assumption that the disturbances have non-zero fourth cumulant of one common sign (all super-Gaussian or all sub-Gaussian). ★ ica_genPerm_relation★ lingam_identifiability_kurtosis

LiNGAM identification, kurtosis route

This file assembles a LiNGAM identification theorem under the standard ICA assumption that the disturbances have non-zero fourth cumulant of one common sign (all super-Gaussian or all sub-Gaussian). This kurtosis route bypasses the general Darmois–Skitovich theorem and Marcinkiewicz's theorem on entire functions:

* colSupport_of_kurtosis (Kurtosis.lean) supplies the column support Wᵢⱼ·Wₖⱼ = 0 from fourth cumulants instead of from Darmois–Skitovich; * ica_genPerm_relation (here) is the purely linear-algebraic core: column support plus invertibility give the generalized permutation relation between A⁻¹ and A'⁻¹; * eq_of_genPerm_triangular_unitDiag (Mathlib/LinearAlgebra/MonomialMatrix.lean) then forces A⁻¹ = A'⁻¹.

theorem ica_genPerm_relation reviewed
Causalean.Discovery.LiNGAM

Let A, A' be n × n real matrices such that A is invertible and A' is invertible. If the matrix A'⁻¹ A has at most one non-zero entry per column — the column-support fact — then there exist a permutation τ and a nowhere-zero scale vector d such that A'⁻¹ and A⁻¹ agree up to that generalized permutation: A'⁻¹ i j = d i · A⁻¹ (τ i) j for every i, j.

Formal statement
n :
Matrix (Fin n) (Fin n) ℝ
hAu :
IsUnit A.det
hA'u :
IsUnit A'.det
hcol :
∀ j i k
if
i ≠ k
then
(A'⁻¹ * A) i j = 0 ∨ (A'⁻¹ * A) k j = 0
∃ (τ : Perm (Fin n)) (d : Fin n → ℝ), ∀ i j, A'⁻¹ i j = d i * A⁻¹ (τ i) j
Proof (Lean source)
theorem ica_genPerm_relation {n : ℕ} {A A' : Matrix (Fin n) (Fin n) ℝ} (hAu : IsUnit A.det) (hA'u : IsUnit A'.det) (hcol : ∀ j i k, i ≠ k → (A'⁻¹ * A) i j = 0 ∨ (A'⁻¹ * A) k j = 0) : ∃ (τ : Perm (Fin n)) (d : Fin n → ℝ), ∀ i j, A'⁻¹ i j = d i * A⁻¹ (τ i) j := by classical let W : Matrix (Fin n) (Fin n) ℝ := A'⁻¹ * A have hWu : IsUnit W.det := by dsimp [W] rw [Matrix.det_mul] exact IsUnit.mul (Matrix.isUnit_nonsing_inv_det A' hA'u) hAu obtain ⟨τ, d, _hd, hWform⟩ := genPerm_of_det_ne_zero_of_colSupport hWu.ne_zero hcol refine ⟨τ, d, ?_⟩ intro i j have hAinv : A'⁻¹ = W * A⁻¹ := by calc A'⁻¹ = A'⁻¹ * 1 := by rw [Matrix.mul_one] _ = A'⁻¹ * (A * A⁻¹) := by rw [Matrix.mul_nonsing_inv _ hAu] _ = W * A⁻¹ := by change A'⁻¹ * (A * A⁻¹) = (A'⁻¹ * A) * A⁻¹ rw [Matrix.mul_assoc] calc A'⁻¹ i j = (W * A⁻¹) i j := by rw [hAinv] _ = ∑ c, W i c * A⁻¹ c j := by rw [Matrix.mul_apply] _ = d i * A⁻¹ (τ i) j := by rw [Finset.sum_eq_single (τ i)] · simp [hWform] · intro c _hc hc simp [hWform, hc] · simp
Causalean.Discovery.LiNGAM.ica_genPerm_relation · Causalean/Discovery/LiNGAM/LiNGAMKurtosis.lean:35
theorem lingam_identifiability_kurtosis reviewed
Causalean.Discovery.LiNGAM

LiNGAM identification theorem (kurtosis route, Marcinkiewicz-free). Let A, A' be n × n real mixing matrices, both invertible, with coefficient matrices A⁻¹ = I − B, A'⁻¹ = I − B' both having unit diagonal, and with A⁻¹ acyclic with respect to a causal order σ: (A⁻¹) i j = 0 whenever σ i < σ j. Let e, e' be families of real disturbances on a probability space, each coordinate measurable, each family's coordinates mutually independent, each coordinate of e with finite fourth moment, each coordinate of e centered, and the fourth cumulant of every coordinate of e nonzero and of one common sign. If the structural equations A·e and A'·e' produce the same observed law, then the coefficient matrices coincide, A⁻¹ = A'⁻¹ (i.e. B = B').

Formal statement
n :
Matrix (Fin n) (Fin n) ℝ
hAu :
IsUnit A.det
hA'u :
IsUnit A'.det
hCdiag :
∀ i, A⁻¹ i i = 1
hC'diag :
∀ i, A'⁻¹ i i = 1
σ :
Perm (Fin n)
hacyc :
∀ i j
if
σ i < σ j
then
A⁻¹ i j = 0
e e' :
Ω → Fin n → ℝ
hem :
∀ i, Measurable (fun ω => e ω i)
he'm :
∀ i, Measurable (fun ω => e' ω i)
heI :
iIndepFun (fun i ω => e ω i) P
he'I :
iIndepFun (fun i ω => e' ω i) P
heL4 :
∀ i, MemLp (fun ω => e ω i) 4 P
hcent :
∀ i, ∫ ω, e ω i ∂P = 0
hkurt :
(∀ j, 0 < kurt (fun ω => e ω j) P) ∨ (∀ j, kurt (fun ω => e ω j) P < 0)
hobs :
P.map (fun ω => A *ᵥ e ω) = P.map (fun ω => A' *ᵥ e' ω)
A⁻¹ = A'⁻¹
Proof (Lean source)
theorem lingam_identifiability_kurtosis {n : ℕ} {A A' : Matrix (Fin n) (Fin n) ℝ} (hAu : IsUnit A.det) (hA'u : IsUnit A'.det) (hCdiag : ∀ i, A⁻¹ i i = 1) (hC'diag : ∀ i, A'⁻¹ i i = 1) {σ : Perm (Fin n)} (hacyc : ∀ i j, σ i < σ j → A⁻¹ i j = 0) {e e' : Ω → Fin n → ℝ} (hem : ∀ i, Measurable (fun ω => e ω i)) (he'm : ∀ i, Measurable (fun ω => e' ω i)) (heI : iIndepFun (fun i ω => e ω i) P) (he'I : iIndepFun (fun i ω => e' ω i) P) (heL4 : ∀ i, MemLp (fun ω => e ω i) 4 P) (hcent : ∀ i, ∫ ω, e ω i ∂P = 0) (hkurt : (∀ j, 0 < kurt (fun ω => e ω j) P) ∨ (∀ j, kurt (fun ω => e ω j) P < 0)) (hobs : P.map (fun ω => A *ᵥ e ω) = P.map (fun ω => A' *ᵥ e' ω)) : A⁻¹ = A'⁻¹ := by classical let W : Matrix (Fin n) (Fin n) ℝ := A'⁻¹ * A have hyindep : ∀ i k, i ≠ k → IndepFun (fun ω => ∑ j, W i j * e ω j) (fun ω => ∑ j, W k j * e ω j) P := by let y : Ω → Fin n → ℝ := fun ω => W *ᵥ e ω have hy_meas : Measurable y := by fun_prop have he'_meas : Measurable e' := by rw [measurable_pi_iff] intro i exact he'm i let g : (Fin n → ℝ) → Fin n → ℝ := fun v => A'⁻¹ *ᵥ v have hg_meas : Measurable g := by fun_prop have hAe_meas : Measurable (fun ω => A *ᵥ e ω) := by fun_prop have hA'e'_meas : Measurable (fun ω => A' *ᵥ e' ω) := by fun_prop have hy : P.map y = P.map e' := by have hmap : Measure.map g (P.map (fun ω => A *ᵥ e ω)) = Measure.map g (P.map (fun ω => A' *ᵥ e' ω)) := by rw [hobs] rw [Measure.map_map hg_meas hAe_meas, Measure.map_map hg_meas hA'e'_meas] at hmap simpa [y, g, W, Function.comp_def, Matrix.mulVec_mulVec, Matrix.nonsing_inv_mul _ hA'u, Matrix.one_mulVec] using hmap have hy_coord : ∀ r, P.map (fun ω => y ω r) = P.map (fun ω => e' ω r) := by intro r let π : (Fin n → ℝ) → ℝ := fun v => v r have h := congrArg (Measure.map π) hy rw [Measure.map_map (measurable_pi_apply r) hy_meas, Measure.map_map (measurable_pi_apply r) he'_meas] at h simpa [π, Function.comp_def] using h have hy_pair : ∀ i k, P.map (fun ω => (y ω i, y ω k)) = P.map (fun ω => (e' ω i, e' ω k)) := by intro i k let proj : (Fin n → ℝ) → ℝ × ℝ := fun v => (v i, v k) have hproj : Measurable proj := (measurable_pi_apply i).prodMk (measurable_pi_apply k) have h := congrArg (Measure.map proj) hy rw [Measure.map_map hproj hy_meas, Measure.map_map hproj he'_meas] at h simpa [proj, Function.comp_def] using h have hy_indep : ∀ {i k : Fin n}, i ≠ k → IndepFun (fun ω => y ω i) (fun ω => y ω k) P := by intro i k hik refine (indepFun_iff_map_prod_eq_prod_map_map (hy_meas.eval (a := i)).aemeasurable (hy_meas.eval (a := k)).aemeasurable).mpr ?_ have he'_prod := (indepFun_iff_map_prod_eq_prod_map_map (he'm i).aemeasurable (he'm k).aemeasurable).mp (he'I.indepFun hik) calc P.map (fun ω => (y ω i, y ω k)) = P.map (fun ω => (e' ω i, e' ω k)) := hy_pair i k _ = (P.map (fun ω => e' ω i)).prod (P.map (fun ω => e' ω k)) := he'_prod _ = (P.map (fun ω => y ω i)).prod (P.map (fun ω => y ω k)) := by rw [← hy_coord i, ← hy_coord k] intro i k hik have h := hy_indep hik simpa [y, mulVec, dotProduct, Finset.mul_sum] using h have hcol : ∀ j i k, i ≠ k → W i j = 0 ∨ W k j = 0 := by intro j i k hik exact mul_eq_zero.mp (colSupport_of_kurtosis (W := W) (hmeas := hem) (hindep := heI) (hL4 := heL4) (hcent := hcent) (hsign := hkurt) hik (hyindep i k hik) j) have hcol' : ∀ j i k, i ≠ k → (A'⁻¹ * A) i j = 0 ∨ (A'⁻¹ * A) k j = 0 := by simpa [W] using hcol obtain ⟨τ, d, hgp⟩ := ica_genPerm_relation hAu hA'u hcol' exact eq_of_genPerm_triangular_unitDiag hCdiag hC'diag hacyc (τ := τ) (d := d) hgp
Causalean.Discovery.LiNGAM.lingam_identifiability_kurtosis · Causalean/Discovery/LiNGAM/LiNGAMKurtosis.lean:73 · uses kurt