查看: 343|回复: 2

因子分析：数学建模竞赛中的"降维 + 评价"利器

字体大小: 正常放大

1198 主题	4 听众	2978 积分

该用户从未签到

电梯直达

1^#

发表于 2026-8-19 15:31 |只看该作者 |倒序浏览

|招呼Ta 关注Ta

因子分析：数学建模竞赛中的"降维 + 评价"利器

写给参加全国大学生数学建模竞赛（CUMCM）的同学 · 基础概念 / 模型公式 / 可视化 / 可运行代码全干货

一句话定位： 因子分析（Factor Analysis, FA）是一种从一堆高度相关的指标中，提取出少数几个"公共因子"的统计方法。它既能降维去冗余，又能把指标归类成构念，是评价类、聚类类赛题的高频建模工具——但用错一步，结论就会"假了"。本文把算法讲透，并附可直接套用的 Python 流程。

为什么竞赛要用因子分析（典型赛题）核心概念：公共因子、载荷、方差贡献
数学模型与关键公式（含协方差分解、Varimax、因子得分）
标准算法流程（八步，含 KMO / Bartlett 检验）
四张可视化图表（结构图 / 碎石图 / 载荷热图 / 得分散点）
Python 实战代码（factor_analyzer）
竞赛案例：城市创新能力评价（简化）
常见误区与提分技巧
赛场自查清单

! R7 l9 q* B2 H5 l

一、为什么竞赛要用因子分析

全国大学生数学建模竞赛中，评价类问题几乎年年出现：给一堆指标（如水质、城市竞争力、上市公司绩效、区域创新能力），让你给对象排名或分类。这类题有四个典型痛点，因子分析正好对症：

指标高度相关：比如"GDP、人均GDP、财政收入"彼此高度共线，直接加权会重复计数。
指标太多：30 个指标直接打分，权重主观且不可解释。
需要客观赋权：评委偏爱"数据驱动"的权重，而不是拍脑袋。
需要可解释的分组：把指标归成"经济因子""环境因子"比罗列 30 个变量更有说服力。

$ E" L+ P! Q& [: c! g# L: S

⚠️ 先区分两个易混方法： 主成分分析（PCA）只做"最大方差压缩"，是综合指标；因子分析额外建模了误差项，目标是解释变量间的相关结构、还原潜在构念。评价建模若想"命名因子"，优先 FA；若只需压缩且解释力足够，PCA 也可。

[td]

维度	主成分分析 PCA	因子分析 FA
目标	保留最大方差（线性组合）	解释变量间相关结构（潜在构念）
模型	无误差项，主成分 = 原变量线性组合	含特殊因子 ε，公共因子建模
处理方差	全部方差（含独特方差）	仅"共同方差"
结果命名	难解释（是混合指标）	易命名（是潜在维度）
竞赛适用	纯降维 / 输入压缩	指标归类 / 构念提取 / 评价赋权

二、核心概念

因子分析假设：我们观测到的 $$p$$ 个标准化指标，其实是由少数 $$m$$ 个不可直接观测的公共因子加上每个指标各自的特殊因子共同决定的。

公共因子 $$F_j$$：被多个原指标共享的潜在维度，例如"经济活力""创新投入"。
特殊因子 $$\varepsilon_i$$：只属于第 i 个指标自身的误差/特异波动，不被其他指标共享。
因子载荷 $$a_{ij}$$：第 i 个指标在第 j 个因子上的相关系数，绝对值越大说明该指标越"听"这个因子的话。
公因子方差（共同度） $$h_i^2$$：第 i 个指标被公共因子解释的比例，$h_i^2$ 越接近 1 越好。
方差贡献 $$V_j$$：第 j 个因子对所有指标的总解释力，用于排权重。
5 _+ g& p; ^9 J

直觉图： 下面这张结构图把"6 个观测指标"通过载荷连到了"2 个潜在因子"。因子分析的本质，就是反过来——从指标的相关关系推断这张隐藏的网络。

图1　因子分析结构示意图F₁公共因子1F₂公共因子2X₁X₂X₃X₄X₅X₆载荷 aᵢⱼ（相关系数，粗细≈大小）

图1　因子分析结构示意图：观测变量经载荷连向潜在公共因子

三、数学模型与关键公式3.1 基本模型

设标准化后的观测向量 $$\mathbf{X}=(X_1,\dots,X_p)^\top$$，公共因子 $$\mathbf{F}=(F_1,\dots,F_m)^\top$$，特殊因子 $$\boldsymbol{\varepsilon}=(\varepsilon_1,\dots,\varepsilon_p)^\top$$，则

$$
X_i = a_{i1}F_1 + a_{i2}F_2 + \cdots + a_{im}F_m + \varepsilon_i,\qquad i=1,\dots,p
$$

等价于矩阵形式：

$$
\mathbf{X} = \mathbf{\Lambda}\mathbf{F} + \boldsymbol{\varepsilon}
$$

其中 $$\mathbf{\Lambda}$$为 $$p\times m$$ 因子载荷矩阵，$$E[\mathbf{F}]=\mathbf{0},\ \mathrm{Cov}(\mathbf{F})=\mathbf{I}_m$$，$$E[\boldsymbol{\varepsilon}]=\mathbf{0}$$，且$$\mathbf{F}$$与 $$\boldsymbol{\varepsilon}$$ 不相关。

3.2 协方差分解（模型的核心）

由模型可得观测变量的协方差（标准化后即为相关矩阵 $$\mathbf{R}$$）：

$$
\mathbf{\Sigma} = \mathrm{Cov}(\mathbf{X}) = \mathbf{\Lambda}\mathbf{\Lambda}^\top + \mathbf{\Psi},\qquad \mathbf{\Psi}=\mathrm{diag}(\psi_1,\dots,\psi_p)
$$

$$\mathbf{\Lambda}\mathbf{\Lambda}^\top$$ 是公共因子贡献，$$\mathbf{\Psi}$ $是特殊方差（唯一性）。这正是 FA 与 PCA 的根本区别：FA 显式地把"误差方差"$$\Psi$$ 单独拆了出来。

3.3 公因子方差与方差贡献

公因子方差（共同度）：
共同度（变量 i 的共性方差）：\[h_i^2 = \sum_{j=1}^{m} a_{ij}^2 = 1 - \psi_i\]因子 $j$ 的方差贡献：\[V_j = \sum_{i=1}^{p} a_{ij}^2\]累计方差贡献率：\[\frac{\sum_{j=1}^{m} V_j}{\sum_{i=1}^{p} h_i^2} \approx \frac{\sum_{j=1}^{m} V_j}{p}\]

实际选因子时，常要求累计方差贡献率 ≥ 80% 且每个指标的共同度 $$h_i^2$$ 不宜过低（一般 ≥ 0.4）。

3.4 因子旋转：让因子"可命名"

初始载荷往往分散、难解释。通过正交旋转（保持因子不相关）让每个变量只在少数因子上有高载荷。最常用的是 Varimax（方差最大化旋转），其目标是最大化各列载荷平方的方差：

$$
\max \sum_{j=1}^{m} \left[ \frac{1}{p}\sum_{i=1}^{p} (\tilde a_{ij})^4 - \left(\frac{1}{p}\sum_{i=1}^{p} (\tilde a_{ij})^2\right)^2 \right]
$$

旋转不改变模型拟合优度与共同度，只改变载荷分布，方便命名（例如把 F1 命名为"规模因子"、F2 命名为"效率因子"）。

3.5 因子得分与综合得分

对每个样本估计其因子得分（回归法）：

$$
\hat{\mathbf{F}} = \big(\mathbf{\Lambda}^\top\mathbf{\Psi}^{-1}\mathbf{\Lambda}\big)^{-1}\mathbf{\Lambda}^\top\mathbf{\Psi}^{-1}\mathbf{X}
$$

用各因子方差贡献作权重，得到综合得分（排名依据）：

$$
S_k = \sum_{j=1}^{m} w_j\,\hat F_{kj}, \qquad w_j = \frac{V_j}{\sum_{j=1}^{m} V_j}
$$

3.6 做之前必须先做的两个检验

不是任何数据都适合因子分析。先过两关：

Bartlett 球形检验： 原假设"相关矩阵 = 单位阵（变量独立）"。检验统计量

$$
\chi^2 = -\left(n-1-\frac{2p+5}{6}\right)\ln|\mathbf{R}|
$$

要求 $p$ 值 < 0.05，才拒绝原假设、适合做 FA。

KMO 取样适切性量数：

$$
\mathrm{KMO} = \dfrac{\sum_{i\neq j} r_{ij}^2}{\sum_{i\neq j} r_{ij}^2 + \sum_{i\neq j} p_{ij}^2}
$$

r{ij}为相关系数，p{ij}为偏相关系数。KMO ≥ 0.6 可用，≥ 0.8 很好。

四、标准算法流程（八步）

数据收集与清洗：缺失值处理、异常值甄别；指标方向统一（越大越优或越小越优）。
标准化：用 Z-score 消除量纲（FA 默认基于相关矩阵，强烈建议先标准化）。
适用性检验：KMO ≥ 0.6 且 Bartlett $p$ < 0.05 才继续。
确定因子数 $m$：特征值 > 1 准则 + 碎石图拐点 + 累计方差贡献率 ≥ 80% 三者结合。
提取因子：主成分法或极大似然法估计载荷矩阵。
因子旋转：Varimax 正交旋转，得到易解释结构。
命名与解释：根据高载荷指标给每个因子起名，写清"这个因子代表什么"。
因子得分 + 综合得分排名：得到每个样本的最终评分与排序。
; L1 g' [% j) V* s1 X

五、四张可视化图表图2　碎石图（Scree Plot）—— 决定保留几个因子

特征值按从大到小排列，落在"陡坡转平台"的拐点之后、且低于红线（特征值=1）的因子不再保留。本例保留 F1、F2。

因子序号特征值0123特征值 = 1保留F1(2.61)F2(1.74)F3(0.92)F4(0.41)F5(0.19)F6(0.13)

图2　碎石图：前两个因子特征值 > 1，构成拐点，予以保留

图3　因子载荷热图（旋转后）—— 看指标归到哪个因子

颜色越深表示载荷绝对值越大。理想旋转结果：每列出现少数高载荷（深色），其余接近 0，从而清晰命名。

F₁（规模因子）F₂（效率因子）指标X₁X₂X₃X₄X₅X₆0.860.810.150.080.790.120.100.220.880.83−0.180.77强正强正弱负近0

图3　旋转后载荷热图：X₁/X₂/X₅ 归 F₁，X₃/X₄/X₆ 归 F₂，命名清晰

图4　因子得分散点图（四象限）—— 样本分类

把每个样本的两个因子得分画成散点，可直观做聚类/定位：右上"双高"、左上"效率型"、右下"规模型"、左下"双低"。

F₁ 得分 →F₂ 得分 ↑双高型效率型规模型双低型ABCDEFGH

图4　因子得分散点图：按双因子得分对样本做四象限定位与聚类

六、Python 实战代码

推荐使用 factor_analyzer 库（专为因子分析设计，自动给出 KMO、Bartlett、载荷、方差贡献）。pip

import pandas as pd# o# |) w3 r( L, q3 @/ [; U
import numpy as np
! B# _5 M1 G# o- }, V) ^
from factor_analyzer import FactorAnalyzer& ~1 m- A& ?$ ]8 N( u6 }
from factor_analyzer.factor_analyzer import calculate_kmo, calculate_bartlett_sphericity( x, v0 j6 n9 m& f
2 c7 E; W( S# \2 r0 L
# 1) 读取数据（行=样本，列=指标），并 Z-score 标准化
M( H! U4 J7 [: R6 m) E
df = pd.read_excel('indicators.xlsx')
h. H( B8 h8 d$ Q; _\" r- Q
X = (df - df.mean()) / df.std()
3 o( I# P) u- B. k* `
; ?% l; w+ [$ ^0 h, H$ Y. J
# 2) 适用性检验4 u/ B3 C, |4 g' J% T3 a2 h8 @4 K
chi2, p_value = calculate_bartlett_sphericity(X)/ Z1 K4 Q, T G* g3 @. L1 o6 N. n& |# l
kmo_all, kmo_model = calculate_kmo(X)
9 s! ]4 F- w8 E5 ]) X- T
print(f'Bartlett p={p_value:.3g} KMO={kmo_model:.3f}') # 需 p<0.05 且 KMO≥0.6
' G5 f9 p$ N' n- q, B r# o: c
5 v2 w7 a/ S- L9 b. f1 v
# 3) 先看所有特征值，确定因子数 m（特征值>1 准则）) `0 q. f( G' \: z2 j2 ~
fa0 = FactorAnalyzer(n_factors=df.shape[1], rotation=None, method='principal')0 o. ^: N) } d' d8 n
fa0.fit(X)) a+ {! h+ |4 @( I P t9 T
ev, _ = fa0.get_eigenvalues()
* k4 j6 N\" l/ I2 `, X1 m( t
print(pd.Series(ev, name='eigenvalue').round(3))
4 J! L) S' F' |8 F
\" K% `- q6 T+ M8 i0 G\" x$ g D0 T' N
# 4) 提取 m 个因子 + Varimax 旋转
* M% @- i\" l# u
m = 29 a\" S$ p3 H) K: ^
fa = FactorAnalyzer(n_factors=m, rotation='varimax', method='principal')5 f6 ?3 E. J. Q4 F1 D( {
fa.fit(X)
7 D2 w0 @, O `3 Y& U3 C9 F
: F. ^! E! u9 T# {4 f0 s
loadings = pd.DataFrame(fa.loadings_, index=df.columns,( y: y/ q4 g* B
columns=[f'F{j+1}' for j in range(m)])7 c5 d! r: F3 H7 C$ p/ t1 l- W
print(loadings.round(3))
: T, x) Q$ u6 J
$ n7 e; z; A7 c$ K$ [
# 5) 因子得分 → 方差贡献加权得到综合得分
7 D* \, H) Z. A1 n) S6 H/ i8 Q
scores = fa.transform(X) # 形状 (n, m)
+ `) C5 g: |5 y' |0 ]- @- [
variance = fa.get_factor_variance()[0] # 各因子方差贡献
1 N$ Z3 t; ~) F# }' I# n0 J
weights = variance / variance.sum()/ h3 i& f\" y# P( ^& u; S. M
composite = (scores * weights).sum(axis=1)3 f; f/ j# U* v% l2 ]7 R
% f+ }: B: d1 Q, }6 Q0 L4 \# e
result = df.copy()\" Y3 D3 q! [0 `6 n( G( g
result[[f'F{j+1}' for j in range(m)]] = scores4 M+ u0 K% ^- {
result['综合得分'] = composite* H& d* p; K# |2 A! u' F$ Z7 B
result['排名'] = composite.rank(ascending=False).astype(int)- p7 t( o\" x# a3 A; ^2 K6 E' l. @3 I- v1 q
print(result.sort_values('排名').head(10))

复制代码

⚠️ 小技巧： 若数据量较小（n < 50 或 n < 5p），优先用 method='principal'（主成分法）而非极大似然法，后者对小样本不稳定。竞赛数据常见样本数有限，主成分法更稳健。

七、竞赛案例：城市创新能力评价（简化）

假设有 6 个指标、8 个城市，指标为：R&D 经费投入(X₁)、每万人专利数(X₂)、高新技术企业数(X₃)、技术市场成交额(X₄)、产业结构高级化(X₅)、创新环境指数(X₆)。

第一步检验： KMO=0.79、Bartlett $p$<0.001 → 适合 FA。
第二步定数： 碎石图显示前 2 个特征值 >1（2.61、1.74），累计贡献率 72.5%，加第三个仅到 87.8% 但 F3 难以命名；经权衡取 2 个因子（或可写"保守取 3 个"，这是评委加分点——要说明取舍理由）。
第三步旋转后解读： F₁ 在 X₁/X₂/X₅ 上高载荷 → 命名为"创新投入与产出规模因子"；F₂ 在 X₃/X₄/X₆ 上高载荷 → 命名为"创新转化与环境因子"。
第四步排名： 用方差贡献加权得综合得分，输出 8 城市排名，并对"双高型"城市写政策建议。

, S* n: Y8 a* c4 S; }' y

这一套路在评阅中很"稳"：每一步都有统计依据（KMO/Bartlett 证明能做、特征值证明取几个、载荷证明怎么命名、方差贡献证明权重怎么来），比"我拍脑袋定了 5 个一级指标"强太多。

八、常见误区与提分技巧[td]

误区	正确做法
跳过 KMO / Bartlett 直接做	先检验，不适配就换方法（如 PCA 或 TOPSIS）
指标方向不统一就标准化	成本类、污染类指标先"正向化"，否则因子含义反了
只看特征值>1 定因子数	结合碎石图拐点 + 累计贡献率 + 能否命名，三者互相印证
不旋转，载荷难解释	默认 Varimax；如需因子相关可用 Promax（斜交）
用因子得分直接排名但不写权重来源	明确写出 $w_j=V_j/\sum V_j$，说明客观赋权
共同度过低仍保留该指标	$h_i^2$<0.4 的指标考虑剔除或说明信息损失