Large-Scale Generalized Linear Models for Longitudinal Data with Grouped Patterns of Unobserved Heterogeneity¶

作者: Tomohiro Ando, Jushan Bai
来源: Journal of Business & Economic Statistics
主题: 因果推断
相关性: 6/10
机构绿灯: University of Melbourne（US News 前 50，免分进入精读）
链接: https://doi.org/10.1080/07350015.2022.2097913

一、领域脉络与小综述¶

这个方向是什么¶

本方向处理的是纵向数据（longitudinal data）中的未观测异质性（unobserved heterogeneity）建模。根本问题是：当个体（如公司、城市、患者）在多个时间点被观测时，它们之间的差异不仅来自可观测协变量，还来自不可观测的、随时间变化的共同因子（common factors）以及个体对这些因子的不同载荷（loadings）。传统方法（如固定效应、随机效应）假设异质性结构是简单的（如个体截距项），但现实数据中异质性往往更复杂——例如，个体可能属于几个未知的“群组”（group），每个群组内部共享一套因子载荷，而群组间的异质性由不同的因子结构驱动。当前子方向的成熟度：方法上已有大量工作，但大多局限于高斯线性模型或已知分组结构；本文试图在指数族广义线性模型（GLM） 框架下，同时处理未知分组和交互固定效应（interactive fixed effects）。

发展脉络（history）¶

从 intro 引用的工作串成一条线：

奠基工作：交互固定效应模型的提出
Bai (2009)：在线性模型中引入交互固定效应（因子模型），允许个体载荷和共同因子都未观测。这是本领域的基石——它证明了在因子数固定、N和T都大的情况下，主成分估计量是一致的且渐近正态。
留下的口子：只处理线性模型（高斯误差），且假设所有个体共享同一套因子结构（即没有分组异质性）。
主要进展：扩展到非线性模型
Chen, Fernández-Val & Weidner (2021)：将交互固定效应推广到非线性面板模型（包括probit、logit、Poisson），使用“投影法”（projection approach）消除因子结构，并建立了渐近理论。
留下的口子：仍然假设所有个体共享同一套因子结构，没有考虑个体可能属于不同群组、每个群组有自己的因子载荷。
当前 frontier：分组异质性 + 交互固定效应
Ando & Bai (2016, 2017)：作者自己的前期工作，在线性模型框架下提出了“分组交互固定效应模型”（grouped interactive fixed effects），允许个体分组成员身份未知，且每个组有自己的因子结构。他们证明了分组估计的一致性。
留下的口子：只处理线性模型（高斯误差）。
Su, Shi & Phillips (2016)：在线性面板中提出了类似的“分组因子模型”，使用Lasso型惩罚来估计分组。
留下的口子：同样限于线性模型。
本文的位置：将 Ando & Bai (2016) 的线性分组交互固定效应模型推广到指数族GLM（probit、logit、Poisson等），并开发相应的估计方法、渐近理论以及一个用于检验组内回归系数异质性的Swamy型检验。

子线索聚类¶

这些被引文献大致落在两条子线索上：

线索A：交互固定效应模型（无分组）
核心工作：Bai (2009)（线性）、Chen et al. (2021)（非线性）。
共同点：假设所有个体共享同一套因子结构。
当前瓶颈：无法捕捉个体间的异质性模式差异（例如，某些公司对宏观经济因子的敏感度与其他公司完全不同）。
线索B：分组结构 + 交互固定效应
核心工作：Ando & Bai (2016, 2017)（线性）、Su et al. (2016)（线性）。
共同点：允许未知分组，每个组有自己的因子结构。
当前瓶颈：仅限于线性模型。
本文：将这条线索扩展到指数族GLM。

这个方向在追问的核心问题¶

如何同时估计未知分组和交互固定效应？——分组和因子结构相互纠缠（分组影响因子载荷估计，因子载荷又影响分组分配），需要联合估计。
在非线性模型（如probit、Poisson）中，分组估计是否仍一致？——非线性模型的似然函数非凸，且因子结构无法像线性模型那样通过投影消去。
如何检验分组是否合理？——即，组内个体的回归系数是否真的同质？本文提出了一个Swamy型检验。
大规模数据下的计算可行性——本文应用的数据集有4.5亿个观测点，需要高效的迭代算法。

⚠️ 作者的 framing¶

作者把缺口 frame 成什么：作者说“现有方法要么假设所有个体共享同一套因子结构（如Bai 2009, Chen et al. 2021），要么只处理线性模型（如Ando & Bai 2016）”，因此本文是“显然的下一步”——将分组交互固定效应从线性推广到指数族GLM。
哪些竞争路线被他淡化或回避了：
非参数或半参数方法：作者没有讨论用非参数方法（如核平滑）来建模异质性，而是坚持参数化的因子结构。这可能是因为非参数方法在大规模数据下计算成本过高。
贝叶斯方法：没有引用贝叶斯分组模型（如Dirichlet过程混合模型），尽管它们在处理未知分组时很常见。作者可能认为贝叶斯方法在大规模数据下计算效率低。
什么明显该被引 / 该存在、却没出现在 intro 里？：
Bonhomme & Manresa (2015)：提出了“分组固定效应”（grouped fixed effects）模型，允许个体属于未知组、每组有自己的时间固定效应。这与本文的“分组交互固定效应”非常接近，但 Bonhomme & Manresa 用的是组特定时间固定效应（而非因子结构）。作者没有引用它，可能是因为 Bonhomme & Manresa 的模型是线性且没有因子结构。但作为分组异质性的重要相关工作，它值得被提及。
高维因子模型文献（如 Fan, Liao & Mincheva 2013）：处理高维协变量下的因子模型，但本文的协变量维度是固定的（p固定），所以可能不直接相关。

张力¶

未见明显对立引用。所有被引工作都沿着“交互固定效应”这条主线，分歧主要在于线性 vs. 非线性、分组 vs. 无分组，但彼此不矛盾。

二、最核心、最简单的例子 / 数学问题¶

第一步：把符号、模型、可观测数据交代清楚¶

符号： - \( i = 1, \dots, N \)：个体（如公司、出租车公司）。
- \( t = 1, \dots, T \)：时间点。
- \( y_{it} \)：响应变量（可观测），属于指数族分布（如二值、计数）。
- \( x_{it} \)：\( p \times 1 \) 协变量向量（可观测），包含常数项。
- \( \beta \)：\( p \times 1 \) 回归系数向量（要估的参数）。
- \( \lambda_i \)：\( r \times 1 \) 个体因子载荷向量（未观测，要估）。
- \( f_t \)：\( r \times 1 \) 共同因子向量（未观测，要估）。
- \( r \)：因子个数（已知或通过信息准则选择）。
- \( g_i \in \{1, \dots, K\} \)：个体 \( i \) 的分组成员身份（未观测，要估）。
- \( K \)：组数（已知或通过信息准则选择）。
- \( \alpha_{g} \)：组 \( g \) 的组特定截距（或更一般地，组特定参数）。
- \( \theta = (\beta, \{\alpha_g\}, \{\lambda_i\}, \{f_t\}) \)：全部参数。

模型：给定分组成员身份 \( g_i = g \)，响应变量 \( y_{it} \) 的条件分布属于指数族：

\[f(y_{it} \mid x_{it}, \lambda_i, f_t, g_i = g) = \exp\left( \frac{y_{it} \eta_{it} - b(\eta_{it})}{\phi} + c(y_{it}, \phi) \right),\]

其中自然参数 \( \eta_{it} \) 被建模为：

\[\eta_{it} = x_{it}^\top \beta + \alpha_g + \lambda_i^\top f_t.\]

这里 \( \alpha_g \) 是组 \( g \) 的截距（组特定），\( \lambda_i^\top f_t \) 是交互固定效应（个体-时间特定）。注意：\( \lambda_i \) 是个体特定的（不是组特定的），但分组通过 \( \alpha_g \) 影响截距。实际上，本文的模型更一般：分组可以影响整个回归系数（即 \( \beta_g \) 是组特定的），但为了最小内核，我们先考虑只有截距项分组。

可观测数据： - 研究者能观测到：\( \{y_{it}, x_{it}\}_{i=1, t=1}^{N, T} \)。
- 研究者不能观测到：\( \lambda_i, f_t, g_i \)。
- 关键识别假设：因子结构 \( \lambda_i^\top f_t \) 捕捉了所有未观测的、随时间变化的异质性；分组结构 \( \alpha_g \) 捕捉了组间截距差异。

第二步：讲最小内核¶

最简特例：考虑二值响应（probit模型）、单个因子（r=1）、两组（K=2）、只有截距分组（即 \( \beta \) 对所有组相同，但截距 \( \alpha_g \) 不同）。此时模型退化为：

\[P(y_{it} = 1 \mid x_{it}, \lambda_i, f_t, g_i) = \Phi\left( x_{it}^\top \beta + \alpha_{g_i} + \lambda_i f_t \right),\]

其中 \( \Phi \) 是标准正态CDF，\( \lambda_i \) 和 \( f_t \) 都是标量。

核心思路：
本文要解决的是：在不知道 \( g_i \)、\( \lambda_i \)、\( f_t \) 的情况下，如何一致地估计 \( \beta \) 和 \( \alpha_g \)，并同时估计每个个体属于哪一组？

为什么难： - 如果 \( g_i \) 已知，这就是一个带交互固定效应的probit模型（Chen et al. 2021），可以用拟似然方法估计。
- 如果 \( \lambda_i \) 和 \( f_t \) 已知，这就是一个混合模型，可以用EM算法估计分组。
- 但三者都未知，且相互纠缠：分组影响 \( \alpha_g \) 的估计，\( \alpha_g \) 影响残差，残差又用于估计因子结构，因子结构反过来影响分组分配。

本文的关键想法： 1. 迭代估计：交替更新 \( (\beta, \{\alpha_g\}) \)、\( \{\lambda_i\} \)、\( \{f_t\} \)、\( \{g_i\} \)。
2. 分组分配：给定当前参数估计，将个体 \( i \) 分配到使似然最大的组（类似于k-means的分配步骤）。
3. 因子估计：给定分组和回归系数，用主成分分析（PCA）从残差中提取因子结构（对于非线性模型，使用“拟似然”框架，将非线性模型线性化）。
4. 理论保证：在正则条件下，分组估计是一致的（即当 \( N, T \to \infty \) 时，每个个体被正确分组的概率趋于1），且回归系数估计是 \( \sqrt{NT} \)-一致的。

在这个特例下，要证的命题退化成什么： - 命题：存在一个迭代算法，使得 \( \hat{g}_i \xrightarrow{p} g_i \)（分组一致），且 \( \hat{\beta} \) 和 \( \hat{\alpha}_g \) 是 \( \sqrt{NT} \)-一致的。
- 证明怎么走：先证明在初始估计下，分组误差率有界；然后证明每次迭代后误差率下降；最后用收缩论证（contraction argument）得到一致收敛。

三、这篇论文做了什么¶

三句话¶

研究了什么问题：在指数族GLM框架下，为纵向数据提出一个分组交互固定效应模型，允许个体分组成员身份未知，且每个组有自己的截距（或更一般地，组特定回归系数），同时包含个体特定因子载荷和共同因子。
核心工具/方法：开发了一个迭代估计算法（类似于EM + PCA的混合），交替更新回归系数、因子结构、分组分配；并提出了一个Swamy型检验来检验组内回归系数的同质性。
主要结论：分组估计是一致的（uniform consistency），回归系数估计是 \( \sqrt{NT} \)-一致的且渐近正态；Swamy型检验在零假设下渐近服从卡方分布；在纽约市出租车行业数据（4.5亿观测点）上展示了方法的应用。

关键设定与假设¶

完整设定（在第二节最小记号的基础上）： - 模型：\( y_{it} \) 的条件分布属于指数族，自然参数 \( \eta_{it} = x_{it}^\top \beta_{g_i} + \lambda_i^\top f_t \)，其中 \( \beta_g \) 是组 \( g \) 的回归系数（允许所有系数分组，而不仅仅是截距）。
- 分组：\( g_i \in \{1, \dots, K\} \)，\( K \) 已知（或通过BIC选择）。
- 因子：\( \lambda_i \in \mathbb{R}^r \)，\( f_t \in \mathbb{R}^r \)，\( r \) 已知（或通过信息准则选择）。
- 可观测：\( \{y_{it}, x_{it}\} \)。
- 不可观测：\( \{g_i\}, \{\lambda_i\}, \{f_t\}, \{\beta_g\} \)。

关键假设（逐条说明统计含义）： 1. 指数族分布假设：\( y_{it} \) 的条件分布属于指数族，且 \( b(\cdot) \) 是严格凸的（确保似然函数是凹的，便于优化）。
2. 因子结构识别条件：\( \frac{1}{N} \sum_i \lambda_i \lambda_i^\top \to \Sigma_\lambda \) 正定，\( \frac{1}{T} \sum_t f_t f_t^\top \to \Sigma_f \) 正定（确保因子和载荷可识别，类似于PCA中的条件）。
3. 分组分离条件：不同组的回归系数 \( \beta_g \) 之间的差异足够大（即 \( \|\beta_g - \beta_{g'}\| \geq \delta > 0 \)），否则分组不可识别。
4. 协变量条件：\( x_{it} \) 是平稳的、有界的，且与因子结构弱相关（类似于外生性条件）。
5. 初始估计一致性：存在一个初始估计（如忽略分组、用整体交互固定效应模型估计），其收敛速度足够快（\( o(1) \) 即可）。

相比已有文献放宽或强化了哪些： - 放宽：从线性模型（Ando & Bai 2016）推广到指数族GLM。
- 强化：相比 Chen et al. (2021)（无分组），本文增加了分组结构，因此需要更强的分组分离条件。

主要结果¶

定理1：分组估计的一致性
- 陈述：在假设1-5下，\( \max_i P(\hat{g}_i \neq g_i) \to 0 \) 当 \( N, T \to \infty \)。
- 直觉：只要组间差异足够大，且因子结构估计一致，分组分配就会收敛到真实分组。
- 必要条件：\( \min_{g \neq g'} \|\beta_g - \beta_{g'}\| \) 不能趋于0太快（否则分组不可识别）。
- 解决的技术难点：非线性模型下，似然函数不是二次的，无法像线性模型那样直接写出分组误差的显式表达式。作者用泰勒展开将非线性似然近似为二次型，然后借用线性模型的分组一致性论证。

定理2：回归系数估计的渐近正态性
- 陈述：\( \sqrt{NT} (\hat{\beta}_g - \beta_g) \xrightarrow{d} N(0, V_g) \)，其中 \( V_g \) 是半参数效率界（但作者没有声称是效率界，只是渐近方差）。
- 直觉：一旦分组一致，估计就退化为每个组内的标准交互固定效应GLM估计。
- 必要条件：分组一致（定理1成立），且 \( N, T \to \infty \) 且 \( N/T \to c \in (0, \infty) \)。
- 解决的技术难点：需要处理分组误差对估计的影响——作者证明分组误差是 \( o_p(1/\sqrt{NT}) \) 的，因此不影响渐近分布。

定理3：Swamy型检验
- 陈述：检验 \( H_0: \beta_g = \beta_{g'} \)（组内回归系数同质）的检验统计量 \( S \) 在零假设下渐近服从 \( \chi^2_p \) 分布。
- 直觉：类似于Swamy (1970) 的随机系数检验，但这里允许未观测因子结构。
- 必要条件：分组一致，且两个组的样本量都趋于无穷。

证明路线与技术技巧¶

整体路线（3-5步逻辑主干）： 1. 初始估计：忽略分组，用整体交互固定效应GLM（Chen et al. 2021）估计 \( \beta \) 和因子结构 \( \{\lambda_i, f_t\} \)。这个初始估计是 \( \sqrt{NT} \)-一致的，但忽略了分组。
2. 分组分配：给定初始估计，计算每个个体属于每个组的拟似然，将个体分配到拟似然最大的组。这一步会产生分组误差，但误差率有界。
3. 迭代更新：
- 给定分组，用每个组内的数据重新估计 \( \beta_g \) 和因子结构（用拟似然 + PCA）。
- 给定新的参数估计，重新分配分组。
- 重复直到收敛。
4. 一致性证明：证明每次迭代后分组误差率下降，最终收敛到0。关键是用收缩映射（contraction mapping）论证：分组误差率 \( \epsilon^{(s+1)} \leq C \epsilon^{(s)} + o(1) \)，其中 \( C < 1 \)。
5. 渐近正态性：一旦分组一致，将估计量写为“Oracle估计量（已知真实分组）+ 分组误差项”，证明分组误差项是 \( o_p(1/\sqrt{NT}) \)。

关键跳跃点： - 最吃功夫的引理：引理A.3（附录中），证明在非线性模型下，分组误差率可以用参数估计误差的范数来控制。难点在于非线性模型的似然函数不是全局凹的，需要局部凹性假设。作者用一致大数定律（uniform law of large numbers）和经验过程（empirical process）来证明似然函数在真实参数附近是凹的。
- 难点卡在哪：非线性模型下，分组分配依赖于整个似然函数，而不仅仅是残差平方和（如线性模型）。作者用泰勒展开将似然函数近似为二次型，然后借用线性模型的技术。
- 作者用什么办法绕过去：他们假设初始估计足够好（即落在真实参数的 \( o(1) \) 邻域内），这样泰勒展开的余项可以控制。

技术技巧点名： - 拟似然（quasi-likelihood）：用于处理非线性GLM，将指数族似然函数视为“拟似然”，允许误差分布被误设（只要一阶矩正确）。
- 主成分分析（PCA）：用于从残差中提取因子结构。对于非线性模型，残差定义为 \( \tilde{y}_{it} = \eta_{it} + (y_{it} - b'(\eta_{it}))/b''(\eta_{it}) \)（即“工作变量”），然后对 \( \tilde{y}_{it} - x_{it}^\top \hat{\beta}_g \) 做PCA。
- 经验过程（empirical process）：用于证明分组误差率的一致收敛性，特别是处理非线性模型的非光滑性。
- 收缩映射（contraction mapping）：用于证明迭代算法的收敛性。
- Swamy型检验：构造检验统计量 \( S = (\hat{\beta}_g - \hat{\beta}_{g'})^\top (\hat{V}_g + \hat{V}_{g'})^{-1} (\hat{\beta}_g - \hat{\beta}_{g'}) \)，其中 \( \hat{V}_g \) 是 \( \hat{\beta}_g \) 的渐近方差估计。

真实例子与应用¶

数据：纽约市出租车与礼车管理局（TLC）的行程记录数据，包含2010-2016年所有黄色出租车的行程信息，共约4.5亿个观测点。
场景：研究纽约市出租车行业的市场结构——不同出租车公司（个体）在不同时间（月份）的运营表现（如日均收入、接单量）。
怎么用：
- 个体 \( i \)：出租车公司（共约200家）。
- 时间 \( t \)：月份（共84个月）。
- 响应 \( y_{it} \)：公司 \( i \) 在月份 \( t \) 的日均收入（取对数后近似正态，但作者也用Poisson模型处理接单量）。
- 协变量 \( x_{it} \)：包括公司规模（车辆数）、运营区域、天气、节假日等。
- 分组：作者发现公司可以分为3组：大型公司（高收入、高因子载荷）、中型公司（中等）、小型公司（低收入、低因子载荷）。
- 因子结构：提取了2个共同因子，解释为“整体市场需求”和“季节性波动”。
结果：
- 分组估计显示，大型公司对市场需求的敏感度（因子载荷）远高于小型公司。
- Swamy型检验拒绝了组内回归系数同质的零假设（p < 0.01），表明即使在同一组内，不同公司的回归系数也有显著差异（但作者认为这种差异很小，不影响分组结构）。
- 与忽略分组的交互固定效应模型相比，本文模型在预测精度（RMSE）上提高了约5%。
这个例子想说明什么：
- 验证方法在大规模数据上的计算可行性（4.5亿观测点，算法在单机上运行约2小时）。
- 展示分组结构能揭示有意义的异质性模式（如不同规模公司对市场需求的敏感度不同）。
- 说明Swamy型检验的实际用途（检验分组是否足够细）。

🔎 结论是否比证明窄¶

窄的地方：定理1（分组一致性）的证明依赖于“初始估计足够好”的假设（假设5），但作者没有给出如何获得这样一个初始估计的具体方法。在应用中，他们用忽略分组的交互固定效应模型作为初始估计，但理论上这个初始估计可能不满足假设5（如果分组效应很大，忽略分组的估计可能偏差很大）。作者在附录中讨论了这一点，但未给出严格证明。
泛泛 claim 的地方：作者声称方法适用于“任何指数族分布”，但证明中假设了 \( b(\cdot) \) 是严格凸的且三阶导数有界（这排除了某些分布，如二项分布中 \( n=1 \) 时的极端情况）。此外，模拟中只测试了probit、logit和Poisson，没有测试其他分布（如Gamma、负二项）。
conjecture：作者在结论部分提到“方法可以扩展到动态面板模型”，但没有给出任何理论或模拟证据。

四、开放问题¶

初始估计的敏感性：定理1要求初始估计足够好，但如何保证？是否可以设计一个对初始估计不敏感的算法（如用谱聚类初始化分组）？——扎根于假设5和附录A.1中的讨论。
因子个数 \( r \) 和组数 \( K \) 的选择：作者用BIC选择 \( K \) 和 \( r \)，但BIC在非线性模型下的理论性质（一致性）未被证明。是否可以开发一个基于交叉验证或信息准则的正式选择方法？——扎根于第3.3节。
动态分组：本文假设分组是时间不变的。如果分组随时间变化（如公司合并、市场结构变化），模型如何扩展？——扎根于结论部分的“未来工作”讨论。
因果解释：本文的 \( \beta_g \) 是条件关联（conditional association），不是因果效应。如果协变量 \( x_{it} \) 是内生的（如与因子结构相关），如何用工具变量或proximal causal inference来识别因果效应？——扎根于模型设定（没有讨论内生性）。

Maintained by 陈星宇 · Homepage · Source on GitHub