企业官网建设流程全解析

简介：这份资源面向空间计量经济学的研究者与高年级学生，系统整理了截面数据下的主流空间回归估计方法，帮助解决空间依赖建模中模型选择与参数估计的实操难题。包内覆盖空间滞后模型SLM、空间误差模型SEM、空间杜宾模型SDM及其误差形式，并延伸至自变量空间滞后、Kelejian-Prucha模型、一般嵌套空间模型、空间扩展模型与地理加权回归等十余种估计程序，同时提供拉格朗日乘子检验LM用于模型诊断。资源附带原始数据，并集成jplv7与Elhorst_codes两套常用工具包，便于直接复现与二次开发。压缩包约6.15MB，以MATLAB脚本与数据文件为主，脚本对应各模型的估计流程，数据文件支撑实证演练。目前已有2073人学习下载，适合希望从经典模型过渡到扩展设定、并借助现成代码快速完成截面空间计量实证的读者参考。

1. 空间计量模型估计方法：截面数据里那些“看着像回归，其实差很远”的坑

如果你手头有一份带地理坐标的截面数据，比如各区县的人均产出、房价、污染浓度，直接跑 OLS 大概率会翻车。原因不复杂：相邻地区的变量往往互相影响，误差项也可能在空间上抱团。空间计量模型估计方法就是来解决这个问题的，核心思路是把“空间关系”显式写进模型——空间滞后模型（SAR）刻画因变量的溢出效应，空间误差模型（SEM）处理误差项的空间自相关，杜宾模型（SDM）则同时包含自变量和因变量的空间滞后项。截面数据是这套方法最经典的用武之地，因为每个观测单元只有一个时间点，空间权重矩阵直接决定了模型能不能用。这篇文章会从权重矩阵构造讲到三种模型的估计代码，再到参数解释和踩坑记录，适合手里有截面数据、想跑通空间回归的从业者。

2. 空间权重矩阵与三种模型：先搞清楚你的数据该用哪个

2.1 空间权重矩阵不是随便填个邻接就行

空间计量的第一步永远是构造空间权重矩阵 W。截面数据没有时间维度，W 就是 n×n 的矩阵，元素 w_ij 表示地区 i 和 j 的空间关系。最常见的三种构造方式：邻接矩阵（共享边界为 1，否则为 0）、距离阈值矩阵（距离小于 d 为 1）、反距离矩阵（w_ij = 1/d_ij）。我一般会先做邻接矩阵，因为它对边界敏感但解释直观；如果研究区域边界碎片化严重，再换反距离矩阵。

构造完必须做行标准化，让每行和为 1。这一步不做，后面估计出来的空间滞后系数 ρ 会偏得离谱。行标准化的含义是“每个地区的邻居影响取平均”，而不是简单加总。

import numpy as np import pandas as pd from scipy.spatial.distance import cdist # 假设 df 里有 lon, lat 两列 coords = df[['lon', 'lat']].values n = len(coords) # 反距离矩阵，对角线设为 0 D = cdist(coords, coords, metric='euclidean') W = 1.0 / (D + 1e-12) np.fill_diagonal(W, 0) # 距离阈值：超过 200 公里的邻居直接切断，避免远距离伪相关 W[D > 200] = 0 # 行标准化 row_sum = W.sum(axis=1, keepdims=True) W_std = W / row_sum

这段代码的关键参数是距离阈值 200 公里，你需要根据研究区域的实际尺度调整。阈值太小会导致孤立点（某行全为 0），太大则失去空间区分度。跑完检查np.any(row_sum == 0)，如果有孤立点，要么放宽阈值，要么手动指定最近邻。

2.2 SAR、SEM、SDM 的公式差异决定了你的假设

三个模型写出来差别不大，但假设完全不同。SAR 是 y = ρWy + Xβ + ε，它假设因变量之间存在溢出效应，比如相邻地区的房价会互相拉动。SEM 是 y = Xβ + u，u = λWu + ε，它认为溢出发生在误差项里，可能是遗漏了某些空间相关的变量。SDM 是 y = ρWy + Xβ + WXθ + ε，它把自变量也做了空间滞后，适合你想同时检验“本地自变量”和“邻居自变量”的影响。

选哪个模型不能靠拍脑袋。常见做法是先用 LM 检验（Lagrange Multiplier）判断残差是否存在空间自相关，再用 LR 或 Wald 检验判断 SDM 能否退化成 SAR 或 SEM。如果 SDM 的 θ 和 ρ 都显著，就别退化成 SAR，否则会遗漏自变量的空间溢出。

# 用 spreg 包做模型估计（PySAL 生态） from spreg import ML_Lag, ML_Error, ML_Lag_Regimes # 假设 y 是因变量，X 是自变量矩阵，W_std 是行标准化权重 # SAR 模型 sar = ML_Lag(y, X, w=W_std, name_y='y', name_x=['x1','x2']) print(sar.summary) # SEM 模型 sem = ML_Error(y, X, w=W_std, name_y='y', name_x=['x1','x2']) print(sem.summary) # SDM 模型：把 WX 也加入自变量 WX = W_std @ X X_sdm = np.hstack([X, WX]) sdm = ML_Lag(y, X_sdm, w=W_std, name_y='y', name_x=['x1','x2','Wx1','Wx2']) print(sdm.summary)

ML_Lag和ML_Error分别对应 SAR 和 SEM，ML_Lag加上 WX 就是 SDM。注意name_x的顺序必须和矩阵列顺序一致，否则输出表里的变量名会错位。估计方法默认是最大似然，截面数据样本量小于 500 时建议用method='full'，避免迭代不收敛。

2.3 参数解释：ρ、λ、θ 分别意味着什么

ρ 是空间滞后系数，衡量因变量的溢出强度。ρ 显著为正，说明邻居的 y 越高，本地的 y 也越高。λ 是误差项的空间自相关系数，显著说明遗漏变量在空间上聚集。θ 是自变量的空间滞后系数，显著意味着邻居的 X 会影响本地的 y。

解释时要注意直接效应和间接效应的分解。SDM 里一个自变量对 y 的总效应 = 直接效应 + 间接效应，直接效应不等于 β，间接效应也不等于 θ。LeSage 和 Pace 提出的分解方法需要用 (I - ρW)^{-1} 计算，spreg 包不直接输出，需要自己写。

# 计算 SDM 的直接效应和间接效应 I = np.eye(n) rho = sdm.betas[0] # 第一个系数是 rho beta = sdm.betas[1:3] # x1, x2 的系数 theta = sdm.betas[3:5] # Wx1, Wx2 的系数 # 总效应矩阵 S_total = np.linalg.inv(I - rho * W_std) @ (beta.reshape(-1,1) + theta.reshape(-1,1)) direct = np.mean(np.diag(S_total)) total = np.mean(S_total) indirect = total - direct print(f"直接效应: {direct:.4f}, 间接效应: {indirect:.4f}")

这段代码算的是平均效应，更严谨的做法是逐观测计算再取平均。ρ 接近 1 时 (I - ρW) 接近奇异，求逆会不稳定，这时候要检查模型是否过度拟合。

3. 截面数据实操：从原始数据到模型输出的完整链路

3.1 数据准备与缺失值处理

截面数据的空间计量对缺失值极其敏感。如果某个地区缺了自变量，整个 W 矩阵的那一行和那一列都得删掉，否则估计结果会有偏。我一般会先做三件事：检查 y 和 X 的缺失比例，缺失超过 10% 的变量直接换掉；检查空间权重矩阵的孤立点；检查 y 的分布，严重右偏就取对数。

# 数据清洗 df = df.dropna(subset=['y', 'x1', 'x2', 'lon', 'lat']) # 检查 y 的偏度 from scipy.stats import skew if abs(skew(df['y'])) > 2: df['y'] = np.log1p(df['y']) # 重新构造 W，确保和 df 的行顺序一致 coords = df[['lon', 'lat']].values D = cdist(coords, coords) W = 1.0 / (D + 1e-12) np.fill_diagonal(W, 0) W[D > 200] = 0 W_std = W / W.sum(axis=1, keepdims=True)

注意df的行顺序必须和 W 的行顺序严格对应。如果中间做了排序或筛选，W 必须重新构造。这个坑我踩过不止一次，模型跑出来 ρ 是负的，查了半天才发现是行错位。

3.2 模型选择：LM 检验与 LR 检验的代码实现

LM 检验用来判断要不要做空间模型。原假设是残差不存在空间自相关。如果 LM-Lag 显著而 LM-Error 不显著，选 SAR；反之选 SEM；两个都显著，选 SDM。

from spreg import OLS # 先跑 OLS 拿残差 ols = OLS(y, X, name_y='y', name_x=['x1','x2']) print(ols.summary) # LM 检验需要手动计算，或者用 GeoDa 等工具 # 这里给出简化版：用残差和 W 计算 Moran's I resid = ols.u moran_I = (n / W.sum()) * (resid.T @ W @ resid) / (resid.T @ resid) print(f"Moran's I: {moran_I:.4f}")

Moran's I 显著为正，说明残差有空间聚集，需要上空间模型。更严格的 LM 检验需要计算统计量，建议用 GeoDa 或 R 的 spdep 包做初步筛选，再用 Python 估计。

3.3 估计结果输出与稳健性检查

跑完模型别急着写结论，先做三件事：检查 ρ 或 λ 的显著性，检查对数似然值，检查残差的空间自相关是否消除。

# 检查 SAR 残差的空间自相关 sar_resid = sar.u moran_sar = (n / W.sum()) * (sar_resid.T @ W @ sar_resid) / (sar_resid.T @ sar_resid) print(f"SAR 残差 Moran's I: {moran_sar:.4f}") # 对比三个模型的 AIC print(f"SAR AIC: {sar.aic}, SEM AIC: {sem.aic}, SDM AIC: {sdm.aic}")

如果 SAR 残差的 Moran's I 仍然显著，说明 SAR 没抓干净空间效应，换 SEM 或 SDM。AIC 越小越好，但别只看 AIC，还要看系数是否符合经济含义。ρ 大于 1 或小于 -1 直接判定模型有问题。

4. 避坑与排查：空间计量截面数据最常见的 5 个翻车现场

4.1 现象：ρ 估计出来是负的，且绝对值很大

原因：W 矩阵行标准化之前没有把对角线设为 0，或者行顺序和 df 错位。对角线不为 0 会导致自己影响自己，ρ 的符号完全乱掉。

解决：构造 W 之后立刻np.fill_diagonal(W, 0)，然后检查df.index和 W 的行索引是否一一对应。如果 df 做过sort_values，W 必须重新构造。

4.2 现象：模型不收敛，报“Matrix is singular”

原因：W 矩阵有孤立点，某一行全为 0，导致 (I - ρW) 不可逆。或者自变量之间存在完全共线性。

解决：检查W.sum(axis=1)是否有 0，有就放宽距离阈值或手动指定最近邻。自变量共线性用 VIF 检查，VIF 大于 10 的变量删掉。

4.3 现象：SDM 的 θ 和 ρ 都不显著，但 SAR 的 ρ 显著

原因：SDM 引入了太多参数，截面数据样本量不够，自由度损失严重。

解决：如果 LR 检验接受 SDM 退化成 SAR，就直接用 SAR。别硬上 SDM，截面数据一般建议样本量至少 100 以上再考虑 SDM。

4.4 现象：直接效应和间接效应的分解结果和 β、θ 差很远

原因：ρ 较大时，(I - ρW)^{-1} 会放大 β 和 θ 的效应，直接效应不等于 β 是正常的。

解决：解释结果时以分解后的直接效应和间接效应为准，不要直接拿 β 说事。如果 ρ 接近 1，分解结果会很不稳定，这时候要考虑模型是否过度拟合。

4.5 现象：换了距离阈值，ρ 的显著性变了

原因：空间权重矩阵的构造对结果影响很大，阈值太小导致邻居太少，阈值太大导致空间效应被稀释。

解决：做敏感性分析，试 3 到 5 个阈值，看 ρ 的符号和显著性是否稳定。如果只在某个阈值下显著，结论不可靠。

5. 进阶技巧：用空间计量做政策评估的边界与验证

截面数据做空间计量，最容易被人质疑的就是内生性。ρWy 里的 y 和误差项相关，OLS 估计有偏，所以必须用 ML 或 IV。ML 是默认选择，但样本量小于 50 时 ML 的小样本性质不好，这时候可以考虑 GMM。不过截面数据做 GMM 需要找工具变量，一般用 WX 或 W 的高阶矩，实操中不太好找。

另一个进阶方向是空间断点回归。如果你有一个政策在某个边界上实施，边界两侧的样本可以看作准自然实验，这时候空间权重矩阵可以构造为“是否在边界同一侧”，ρ 的解释就变成了政策溢出效应。这个做法在区域经济评估里越来越常见，但对数据要求很高，边界两侧的样本要足够多。

验证模型是否靠谱，我一般会做两个检查。第一，把 W 矩阵随机置换 1000 次，看 ρ 的分布是否集中在 0 附近。如果随机置换后的 ρ 仍然显著，说明模型抓到的不是空间效应，而是某种全局趋势。第二，留出法交叉验证：随机删掉 10% 的样本，用剩下的估计模型，预测删掉的样本，看预测误差和 OLS 比是否更小。

# 随机置换检验 rho_perm = [] for _ in range(1000): perm = np.random.permutation(n) W_perm = W_std[perm][:, perm] sar_perm = ML_Lag(y, X, w=W_perm) rho_perm.append(sar_perm.betas[0]) # 看真实 rho 在置换分布中的位置 p_value = np.mean(np.abs(rho_perm) >= np.abs(sar.betas[0])) print(f"置换检验 p 值: {p_value:.4f}")

这个检验跑起来慢，1000 次大概要十几分钟，但能有效排除伪空间效应。我一般只在论文投稿前跑一次，日常分析用 Moran's I 就够了。

最后说个血泪经验：空间计量的结果对 W 矩阵极其敏感，同一份数据换一种 W 构造方式，ρ 可能从 0.3 跳到 0.7。所以别只报一个 W 的结果，至少报邻接矩阵和反距离矩阵两套，让读者自己判断稳健性。截面数据没有时间维度，你没法用固定效应去吸收空间异质性，W 的选择就是你的核心假设，写论文时要把构造理由讲清楚。希望帮到你。

本文还有配套的精品资源，点击获取