序列线性模型中的估计与变量选择：基于SCAD惩罚的方法及其应用 袁艺文（Yiwen Yuan）， 尚俊峰（Junfeng Shang）， 顾超（Chao Gu）

《Mathematics》：Estimation and Variable Selection in Sequential Linear Models: SCAD-Penalized Method with Applications Yiwen Yuan, Junfeng Shang and Chao Gu

　　摘要 序列线性模型可以用来描述响应变量依赖于滞后结果和固定效应变量的数据。为了进行估计、变量选择和高精度的响应预测，我们提出了一种基于平滑剪切绝对偏差惩罚（SCAD）的序列线性模型惩罚方法。我们通过模拟比较了SCAD惩罚方法与普通最小二乘（OLS）、Lasso和自适应Lasso在序列线性模型中的表现。模拟结果表明，SCAD惩罚方法在估计方面具有更好的准确性和精确度，在变量选择方面也有更好的预测能力。我们将所提出的方法应用于两个真实数据集，以进一步说明SCAD惩罚方法在序列线性建模中的性能。

1. 引言
在临床领域、生物学研究和计量经济学中，许多数据是动态的，并且会随时间变化。为了对时间序列数据中的变量之间的动态关系进行建模和分析，自回归分布滞后（ADL）模型是一个强大的工具，因为它能够捕捉短期和长期的依赖性，并有助于探索复杂的时间动态[2]。ADL模型结合了自回归项（因变量的滞后）和分布滞后项（自变量的滞后），使研究人员能够研究解释变量的过去值如何影响因变量的当前值。这在过去信息对未来结果有关键作用的情境中特别有用。
在许多情况下，感兴趣的时变因变量不仅依赖于其之前的状态，还与时变预测变量相关。对于这种情况，可以采用[3]中的序列建模方法。这种序列线性建模方法可以利用序列时间点的最后一个状态和预测变量有效地对因变量进行建模。序列线性模型旨在处理序列信息，非常适合通过时变的滞后结果和预测变量来预测因变量。例如，为了应对序列预测建模，参考文献[3]提出了在中间时间点顺序预测结果，然后将这些结果作为后续时间点的模型预测变量，因为涉及恢复时间或慢性疾病的干预措施需要在不同的中间随访点进行结果测量。
普通最小二乘（OLS）回归也可以通过纳入滞后因变量的过去值来建模动态数据，以展示动态变量对预测过程的影响。然而，参考文献[2]指出，由于滞后因变量导致自变量的系数偏向较低的值，从而产生了相关效应。为了解决序列时间序列数据中的这个问题，参考文献[2]在OLS回归中引入了滞后因变量，其中滞后因变量系数表示相对于自变量的时序效应。参考文献[4]将估计残差与SCAD惩罚和Yule-Walker方程结合起来，以确定自回归误差的阶数并估计自回归参数。参考文献[5]提出了一类新的层次滞后结构（HLag），它将滞后选择嵌入到凸正则化器中，提高了预测性能。
受到ADL和上述相关工作的启发，我们采用序列线性建模方法，并利用序列线性模型通过纳入滞后因变量和时变自变量来描述因变量。与ADL模型一样，序列线性建模在估计时也面临过拟合和多重共线性的挑战，特别是当潜在解释变量的数量很大时。为了解决这些挑战，最近在惩罚回归技术方面的进展，如Lasso及其变体，在动态模型（如ADL）中的变量选择和参数估计中越来越受欢迎。这些方法对回归系数施加惩罚，从而得到既可解释又可预测的稀疏模型。
我们提出使用SCAD惩罚方法在序列线性建模中进行变量估计和选择。SCAD方法由[7]引入，以其处理大规模变量选择问题的能力而闻名，同时在系数估计中比Lasso提供更少的偏差。由于其在保持模型准确性的同时进行变量选择的卓越能力，SCAD惩罚方法已成为高维回归模型中的首选。我们的目标是评估使用SCAD惩罚估计的序列线性模型的性能，并将其与Lasso和自适应Lasso[8,9]进行比较，后者因其在变量选择和正则化方面的有效性而被广泛应用。我们的研究通过比较这些惩罚回归方法在估计和变量选择方面的能力来检验它们的优缺点。具体来说，我们评估SCAD惩罚方法在预测准确性和系数估计方面是否优于Lasso和自适应Lasso。通过这样做，我们旨在为惩罚回归技术及其在动态建模中的应用文献做出贡献，提供关于动态和时间序列数据正则化方法选择的见解。
尽管自[7]以来SCAD已被广泛使用，但其在时间序列和动态设置中的应用仍然比在纯横截面回归中要有限。参考文献[4]研究了具有自回归误差的变系数模型中的SCAD。参考文献[10]在部分线性模型中建立了SCAD的一致性。参考文献[11]提出了通过BIC型标准进行数据驱动的SCAD调整参数选择。参考文献[12]分析了自回归模型中的SCAD型惩罚，以区分平稳和非平稳状态。参考文献[13]为高维相关面板数据中的自适应惩罚方法推导出了预言不等式。参考文献[14]使用SCAD进行动态面板模型中的变点检测。向量自回归文献中也看到了相关的层次滞后正则化[5]。
在这些工作中，[13,14]的动态面板模型与序列线性模型最为相关。然而，在这两种模型中，滞后维度是固定的，并且模型在每个时间点不会重新拟合。因此，[3]中的序列线性模型不是这两种框架的特例。它的定义特征是：（i）滞后结果系数向量的维度随时间索引t增长；（ii）模型在每个时间点都会根据更新的因变量测量历史重新拟合。这两个特征将这种公式化置于现有动态模型SCAD理论的范围之外。本文旨在扩展SCAD理论到这种设置。特别是在每个时间点t，我们为高维情况下的SCAD惩罚估计量建立了预言性质，其中名义滞后维度可能超过横截面样本大小。
本文的增量贡献可以总结如下：
我们将SCAD惩罚框架扩展到序列线性模型，在这种模型中，滞后结果系数向量的维度随时间索引增长，这是现有动态面板或自回归模型SCAD理论未涵盖的设置。
我们为所提出的估计量建立了预言性质（定理1），包括活动集上的选择一致性和渐近正态性，并推导出了滞后阶数选择一致性的推论（推论1）。我们的论点将[14,15]中的高维SCAD理论适应到模型（2）的序列滞后设置。
我们通过低、中、高维设置下的模拟以及两个真实的纵向应用证明，所提出的估计量在平均绝对预测误差、相对风险和经验覆盖率方面优于OLS、Lasso和自适应Lasso，并且能够恢复正确的有效滞后阶数。
本文的其余部分组织如下。第2节介绍序列线性模型和正则条件。第3节发展了SCAD惩罚方法，陈述了带有证明的预言性质定理，给出了滞后阶数选择的推论，描述了坐标下降算法和调整参数选择，并讨论了算法和统计收敛性。第4节介绍了模拟研究。第5节包含两个应用。第6节进行总结和讨论。

2. 序列线性模型
我们有一组具有重复结果的数据和一系列时间点。在每个时间点t，我们有一个基于预测变量顺序预测结果的模型，其中有一个参数/系数向量，p是预测变量的数量，以及一个随t变化的参数/系数向量，表示为。我们也可以在任何时间点t获得后续结果，表示为。随着时间索引t每增加一个单位，结果可以基于协变量或预测变量以及之前的预测结果来预测。预测结果的分布被认为是正态的，表示为。我们重复之前的模型，直到终止时间点。
在继续进行序列建模之前，我们可以通过一个例子来可视化这个设置。可以在不同的时间点收集患者的测量数据，这些测量数据包括检测疾病的结果以及所有随时间变化的预测变量，例如影响这种疾病的许多因素。然后我们可以自然地假设当前时间点的结果是由该时间点的预测变量以及之前的结果状态共同预测的。
序列线性模型如下：
在第一个时间点，线性模型是
其中误差项呈正态分布。
在第二个时间点，后续线性模型是
其中，结果向量是。
在第三个时间点，后续线性模型是
其中。
对于任何时间点t，序列线性模型可以写为
（1）
其中。这种递归一直应用到终止时间点。
我们注意到，序列线性模型涉及跨时间点的滞后因变量，因此在拟合模型和进行预测之前需要验证关键假设。下面给出了概念性描述。在引入模型（2）之后，会有精确的数学陈述。
因变量和自变量之间的关系是线性的。
每个时间点的误差项是独立的，并且遵循具有恒定方差的正态分布。
每个时间点的自变量之间没有完全共线性。
序列线性模型允许对每个观测值加入固定效应。
在引入模型（2）之后，这些关键假设被以数学形式表述，用于证明。对于方程（1）中的序列建模，我们看到这个模型是线性预测变量和自回归时间序列回归的组合。在第3节中，我们提出了惩罚方法来进行估计。这些方法具有稀疏性质，意味着采用适当的惩罚项可以将一些参数估计值缩小到零。依靠这种稀疏性质，可以同时实现变量选择和估计。在模型（1）中估计参数向量 和 目标是在 和 中选择预测变量/特征以及自回归部分的滞后数量。关于滞后数量，例如，如果估计发现只有 和 是显著的（即不会缩小到零），我们得出结论，当前时间点的响应变量可以由最后两个点的预测变量预测，因此滞后阶数等于2。随着时间的推移，响应变量通常与远期的值变得越来越不相关，相应的滞后系数会缩小到零。从理论上讲，当前时间的响应变量只依赖于最近的时间点，因此有效滞后阶数不能是一个很大的整数。因此，我们的目的是提出一种有效的惩罚方法，选择 和 中的预测变量和 中的滞后阶数。
模型（1）作为一个结合了线性回归和自回归时间序列的规范，滞后阶数由所提出的方法确定。为了简洁起见，我们将所有预测变量排列成一个矩阵，其中n是样本大小，p是预测变量的数量。
直到时间t的响应/结果变量由 表示，这是一个 矩阵。
为了表示模型（1），我们使用时间t的结果作为预测变量，并将所有预测变量组合成一个矩阵。为了以简洁的形式表示方程（1）中的序列线性模型，可以将其重写为
（2）
其中 是时间t的预测变量矩阵， 是时间t的响应/结果向量， 是时间t的误差项向量。我们假设 是独立同分布（i.i.d.）的，具有方差 。同样，模型（2）重复进行，直到 。在初始时间点，模型不包含该项；而在后续时间点则包含。具体来说，在序列线性模型（2）中，我们将设计矩阵划分为两部分：一部分收集对应于非零回归系数的列，另一部分收集对应于零回归系数的列。我们还用符号表示第t时刻的第j个参数对应的列。我们还对模型施加了弱自相关性（或无自相关性）和异方差性不存在的假设。设s为一个正整数，用于标记滞后项。不失一般性，与[14,16]中的假设类似，这些假设如下：

A1. 随着时间的推移，非零回归系数和零回归系数逐渐变得不相关。
A2. 该过程是弱（协方差）平稳的。即，均值向量在t时刻是常数；方差在t时刻是有限的且恒定的；自协方差仅依赖于滞后项s，而不依赖于t。
A3. 对于所有足够大的n和每个固定的滞后项s，活跃集格拉姆矩阵（active-set Gram matrix）几乎肯定是正定的。
A4. 活跃集由n个独立同分布（i.i.d.）的随机变量组成，其分布不随t变化。

我们注意到A1假设误差项与预测变量和滞后项相互独立。A2表达了随时间推移的弱线性依赖性，而不是要求完全独立。除非施加额外的分布假设，否则A2并不意味着完全的统计独立性。A3是线性时间序列模型中一个温和且常用的假设。任何具有有限二阶矩的严格平稳过程都自动满足A3，因此下面备注1中探讨的严格平稳构造足以满足这一要求。A4是活跃集参数的可识别性条件，它保证了活跃集格拉姆矩阵在每个n时刻都是可逆的，从而使得在活跃集上的回归是适定的。在高维情况下，完整的格拉姆矩阵必然是秩不足的，因此不能对完整参数向量施加可识别性要求。在第3.2节中开发SCAD惩罚方法时，我们进一步假设最小特征值有一个统一的下界（A7）。A4和A5共同保证了模型的适定性和其活跃集参数的可估计性。

对于方程（1）和（2）中的序列模型，我们的目标是进行变量选择和估计，以确定描述时间变化数据的最佳模型，预测结果，并进行进一步的推断。在分别用和估计出和之后，下一个时间点的因变量预测值为（3），其中是预测公式。每个时间点的预测输出是当前时刻预测变量和之前时间点滞后响应的函数；当滞后响应可用时，使用真实的值，否则使用预测值。在稳定的时间序列环境中，滞后系数向量的元素通常随着滞后指数的增加而减小[17]，这进一步激发了第3节中惩罚方法进行的滞后阶数选择。

3. SCAD惩罚方法用于估计和变量选择
在前一节介绍了序列线性模型（1）和（2）之后，我们现在开发这些模型的参数估计和变量选择方法。为了使用最合适的模型进行预测，我们采用惩罚方法来改进估计和变量选择，并减少过拟合。一种典型的惩罚方法是Lasso[6]，它引入了一个惩罚项。Lasso估计量具有一致性和稀疏性，这使得模型在序列建模中更具解释性和效率。从Lasso开始，本节接着介绍了自适应Lasso，并提出了用于估计和变量选择的SCAD惩罚方法。

3.1 SCAD惩罚方法用于估计和变量选择
在接下来的内容中，我们将介绍本研究中用于估计和变量选择的方法。基本思想是对方程（2）进行惩罚，以获得估计量。回想一下，设计矩阵包含了所有预测变量和之前的结果。用表示的第i行向量。对于模型（2）中的回归，用表示结合了和的参数向量。的维度是，的维度是，因此的长度是。我们用表示的每个参数元素。

为了估计参数向量，我们从Lasso方法开始。根据模型（2），Lasso估计量是通过最小化惩罚平方误差损失来获得的：
（4）
其中是时间点t的第j个参数元素，是调整参数。

Lasso方法由[6]首次提出。正如文献中所建立的，Lasso是一种同时进行参数估计和变量选择的正则化方法，其估计量具有一致性和稀疏性。具体来说，惩罚项将不显著变量的系数缩减为零，赋予Lasso稀疏性，这是变量选择的基础。

在序列线性建模中，每个时间点都会对观测数据或生成的数据进行建模，从而在连续迭代之间产生依赖性。在这个框架中的一个常见挑战是多重共线性，它源于预测变量之间的显著相关性。Lasso回归方法通过自动选择相关预测变量并可能排除不相关的变量来处理多重共线性。它还通过考虑观测值的顺序排列并选择有助于预测未来值的滞后特征来解决时间依赖性问题。

当对中的不同参数赋予不同的权重时，该过程称为自适应Lasso[8]，本质上是一种加权惩罚方法。在某些规则条件下，的估计量为：
（5）
其中是参数元素的权重（例如，初步估计的倒数），是调整参数。

现在我们应用带有SCAD惩罚[7]的惩罚方法来估计序列线性模型（1）和（2）中的参数。估计量定义为平方误差和SCAD惩罚之和的最小化者：
（6）
其中表示绝对值，惩罚函数（7）按照[7,18]明确写为：
（7）
其中是控制惩罚平滑度的额外调整参数，控制惩罚的强度。方程（7）是一个在和处有节点的二次样条函数。它是连续的，其导数在区间外为零。

3.2 神谕性质
现在我们探讨序列线性模型中SCAD惩罚估计量的神谕性质。为了清晰起见，我们用表示。设为时间t的真实活跃集，其基数为。设为显著SCAD估计量的索引集，其基数为。我们将真实参数向量组织为，其中收集非零系数，收集零系数。类似地，我们将组织为，并设置。

根据[14]，为了在高维和稀疏框架内建立神谕性质，我们在[19]中施加了以下补充规则条件，与A1–A5结合使用：
A6. 存在一个正常数，使得对于所有和所有n，在时间t时成立。
A7. 存在一个正常数，使得对于所有和成立。
A8. 对于某些成立。
A9. 存在正常数和，使得并且当预测变量和滞后结果标准化时，限制了的列范数的平方值。A7是限制特征值条件，确保活跃集格拉姆矩阵在n中始终是严格正定的，从而即使在预测变量和滞后的数量远大于n时，也能正确选择变量。A8控制有效维数的增长速率，确保显著预测变量的数量K相对于n以受控的多项式速率增长，防止过拟合。A9确保最小非零系数与零之间有一个有序的间隔，确保信号不会被随机误差掩盖。更多背景信息见[14,15]。

在转向神谕定理之前，我们验证A1–A9是否可以共同满足。

备注1（假设A1–A9的兼容性）。假设A1–A9是相互一致的。一个满足所有这些假设的具体构造如下：设是一个双索引的独立同分布（i.i.d.）向量序列，其在i和t上都是独立的。为了简单起见，固定一个最大滞后和时不变的滞后系数，满足。定义的非零元素为对于，以及对于。然后过程承认一个严格平稳的因果解[17]，因此也是弱平稳的。对于任何稀疏指数，A1–A9都成立。注意，类似的构造也适用于任何非零滞后支持（例如，具有非连续滞后的稀疏支持），前提是相应的多项式满足相同的根条件。

尽管参数向量的名义维度随t增长，但A8要求其大部分元素为零。由此产生的有限有效滞后阶数正是使A3中的弱平稳性与增长的名义维度兼容的原因。第4节中的模拟研究使用了具有有效滞后阶数的，这是这种构造的一个明确实例。

现在我们将神谕性质作为一个正式定理陈述。
定理1（SCAD惩罚估计量的神谕性质）。假设模型（2）成立，假设A1–A9得到满足，并且调整参数λ满足、和。固定。进一步假设存在一个正定矩阵，使得随着的减少。设表示在（6）中SCAD准则的局部最小化器。那么，随着的减少，以下结论成立：
（i）选择一致性。
（ii）活跃集上的渐近正态性。设是一个固定的整数，设是一系列满足的矩阵，其中是一个正定矩阵。那么

证明。证明通过归约为[20]中的定理2来进行，该定理建立了在参数数量增加时的SCAD惩罚估计量的神谕性质。我们验证序列线性模型（2）在A1–A9和定理1的调整参数条件下满足该定理的假设。

渐近情况是固定的t。在每个固定的t时刻，模型（2）是对独立同分布（i.i.d.）受试者的n样本回归（由A5保证），时间依赖性仅通过设计矩阵中的回归量进入。这与[20]中假设的跨受试者独立同分布设置相匹配。

首先，模型（2）是具有参数向量的高斯线性模型。这是[20]第1.2节中一般似然设置的一个特例。在这种特殊情况下，Fisher信息简化为，而[20]中的似然规则条件（E）、（F）和（G）转化为模型（2）的设计和误差分布条件。

本文的A1–A9与[20]中的条件（A）–（H）之间的对应关系如下：惩罚条件（A）–（D）由定理1中的SCAD惩罚在调整参数速率下得到满足。特别是，在[20]的第3.1.1节中最终得到了相关结果。条件（E）（独立同分布的观测值具有共同的密度，可识别的模型，得分在真实值处的均值为零）是根据A5（具有高斯误差的独立同分布的受试者，提供共同的支持）、A4（活跃集的Gram矩阵是正定的，提供活跃集上的可识别性）和A1（条件均值为零的误差，在正确的规范下得分均值为零）得出的。条件（F）（Fisher信息是正定的，具有有界的特征值和得分分量的有界四阶矩）是根据A7（活跃集的Gram矩阵有下界）、A6（列范数有上界）和A5（的四阶矩是有限的）得出的。条件（G）（对数密度的三阶导数规则性）对于高斯线性模型来说是显而易见的，因为对数密度是二次的。条件（H）（最小信号条件）是根据A9结合调整参数的速率得出的，这两者共同得出了该条件。稀疏率A8提供了[20]第2定理中所需的参数数量的增长条件，前提是相对于那里施加的速率约束，该参数足够小。其余的条件A2（渐近不相关性）和A3（弱平稳性）描述了序列线性模型框架，并支持对假设的解释。在这些条件下，[20]第2定理得出了SCAD准则（6）的局部最小化器的存在性，该最小化器满足条件，并且在活跃集上具有一致性，其速率与A9相结合，从而得出了定理1的结论（i）。同一定理还得出了定理1（ii）中的结论。对于有限维SCAD预言器属性的类似论证出现在[7]的第1和第2定理中，并且已经被[15]扩展到高维设置中。□定理1（i）的选择一致性结论自动意味着一个滞后阶数选择保证，我们将其记录为以下推论。推论1（滞后阶数选择一致性）。定义时间t时的真实有效滞后阶数为，并将其基于SCAD的估计器定义为，约定为空集的最大值为零。在定理1的条件下，随着。证明。对于的滞后系数对应于的坐标，因为与。映射是在滞后指数和的滞后坐标之间的双射。在这种双射下，当且仅当，对于SCAD选择的对应物也是如此。因此，事件意味着，并且根据定理1（i）。□我们注意到推论1表明，SCAD惩罚方法可以在单次优化中同时且一致地执行预测器选择和滞后阶数选择，而无需对候选滞后阶数进行单独的信息准则搜索。在名义维度随t增长的序列设置中，这是一个有意义的计算和统计优势。3.3. 算法实现由于定理1中建立的稀疏性质，SCAD惩罚方法选择了真实模型，产生了一组稀疏的解和对于显著预测变量及其显著滞后的近似无偏系数。在正交的特殊情况下，SCAD目标（6）在坐标中是可分离的，具有封闭形式的最小化器（9），其中是符号函数，表示第j个坐标的未惩罚最小二乘系数。对于模型（2）中出现的非正交设计，方程（9）不再是问题（6）的精确封闭形式解。相反，我们通过坐标下降和SCAD阈值迭代计算：在每个坐标上，应用阈值规则（9），同时保持所有其他坐标不变（因此（9）中的是当前坐标j的偏残差单变量最小二乘更新），算法迭代直到目标的相对变化低于收敛容忍度（在我们的实验中我们使用了）。这是在ncvreg R包3.16.0中实现的求解器，我们在整个模拟和应用中都使用了它。该算法的收敛性质以及所得估计器的统计收敛速率在第3.5节中进行了讨论。最优对可以通过使用交叉验证的二维网格搜索获得，这在计算上代价较高。Fan和Li [7]建议将设置为跨不同场景的稳健选择，并表明的选择不会显著影响相对于的整体性能。我们在整个过程中采用，并通过交叉验证选择。在序列线性模型的背景下，我们使用方程（7）-（9）中的SCAD惩罚来增强变量选择和细化系数估计。稀疏性质结合SCAD惩罚的平滑惩罚，对于处理时间序列数据的序列和相关性质非常有益。我们注意到，SCAD惩罚[7]对于参数估计和变量选择更为有效，特别是在高维数据中，因为它促进了稀疏性，实现了渐近无偏性，并满足了预言器属性。与Lasso类似，SCAD惩罚方法移除了系数较小的变量。由于SCAD惩罚是非凹的，它对小系数施加了重罚，同时几乎保持了大系数的无偏性。随着样本量的增加，SCAD惩罚方法一致地识别出真正的非零系数集。根据推论1，这特别识别了序列线性建模中的真实滞后阶数，从而提高了模型（1）和（2）中的预测准确性。如下一节中的模拟结果所示，SCAD惩罚提供了比Lasso和自适应Lasso更好的惩罚估计和变量选择解决方案。与Lasso不同，SCAD不是对系数施加统一的惩罚，而是施加了一个非凹惩罚，有效地将小系数缩小到接近零，同时保留了大系数以进行更准确的估计。在假设A7（受限特征值）、A8（稀疏率）和A9（最小信号）下，即使预测变量的数量超过了样本量，SCAD惩罚方法也通过定理1实现了预言器属性。另见[21]，以获取高维中SCAD类型方法的更广泛概述。正则化参数控制了预测变量系数的收缩，系数向量结合了独立变量的参数和的滞后系数。请注意，序列线性建模通常涉及不同时间点观测值之间的依赖性，SCAD惩罚方法通过鼓励对潜在模式的简洁表示来帮助缓解多重共线性问题。Xie和Huang [10]使用SCAD惩罚在部分线性模型中促进稀疏性，同时使用多项式样条估计非参数函数。他们的工作表明，SCAD惩罚不仅实现了一致的变量选择，还准确识别了正确的模型结构，从而在真实模型下实现了有效的参数估计。3.4. 调整参数的选择在惩罚方法中选择正则化调整参数的过程至关重要。在存在时间依赖性的序列线性模型中，需要仔细考虑以平衡对系数的惩罚和重要特征的保留。因此，必须适当地选择调整参数。如前一小节所讨论的，SCAD惩罚方法通过将不显著的参数缩小到零并减少显著参数估计的偏差来估计序列线性模型中的参数。因此，SCAD惩罚方法一致地选择了显著参数。Fan和Li [7]表明，SCAD方法适用于参数化和高维非参数化设置。惩罚似然估计器的收敛速率取决于正则化参数，因此需要适当选择的调整参数以实现有效的变量选择和准确的参数估计。为了便于调整参数的选择，通常使用交叉验证。序列线性模型涉及额外的复杂性，例如构建具有时间效应系数的滞后变量，这需要在惩罚系数和保留重要特征之间进行仔细的平衡。在这种情况下，SCAD惩罚被证明是有效的，这也是为什么它的性能超过了这里采用的其他方法，如下一节的模拟所示。对于具有不同惩罚项的惩罚方法，我们使用交叉验证来选择调整参数，它控制了参数估计的收缩权重。在带有Lasso方法的序列线性模型中，惩罚函数是，并且Lasso估计器允许逐坐标阈值规则，其中表示未惩罚的最小二乘系数，是符号函数。对于自适应Lasso方法，逐坐标的惩罚是，其中是分配给系数的权重。自适应Lasso对每个系数的惩罚是根据相应的Lasso估计进行调整的。使用自适应权重的目的是为了对某些变量进行降权，对其他变量进行加权，从而实现更有效的变量选择。具有较大Lasso估计的变量获得较小的权重，而具有较小Lasso估计的变量获得较大的权重。在通过各种方法拟合模型之前，我们将数据按4:1的比例分为训练集和测试集。训练集用于估计参数、识别显著变量和调整正则化参数。测试集用于评估模型的预测性能。在拟合之前的一个步骤是调整正则化参数以最小化均方误差。交叉验证是一种强大的技术，以确保泛化性能的稳健性并避免序列线性建模中的过拟合。通过将训练数据分为K折交叉验证中的训练集和验证子集，该过程有助于SCAD惩罚方法选择显著变量并防止过度收缩。它还与Lasso和自适应Lasso结合使用，以评估性能并选择最佳的正则化调整参数。预测误差作为交叉验证中的标准。在训练集固定的情况下，交叉验证估计预期预测误差，验证集提供对预测模型的评估。在K折交叉验证中，使用可用的折来拟合模型，而保留的折用于评估性能。选择最佳调整参数的规则是最小化估计的预测误差。根据方程（7），SCAD惩罚方法如下选择具有系数的变量。当较小时，SCAD惩罚的行为类似于Lasso，将参数缩小到零。当较大时，SCAD惩罚是恒定的，避免了大参数的偏差。当适中时，SCAD惩罚对参数赋予不同的权重。3.5. 算法收敛性和收敛速率在定义了SCAD惩罚估计器、建立了其预言器属性、描述了用于计算它的坐标下降算法并指定了调整参数选择程序之后，我们通过解决两个进一步的问题来结束该方法的理论发展。首先，第3.3节中的坐标下降算法是否可靠地产生了定理1适用的估计器？其次，随着样本量的增加，这个估计器以什么速率接近真实参数？由于SCAD目标（6）是非凸的，全局优化通常是难以处理的，因此这两个问题都不是显而易见的。尽管如此，两者都有我们在这里总结的理论答案。算法收敛性。使用[22]中的坐标下降算法最小化（6）中的SCAD惩罚目标，在R包ncvreg中实现。在每一步中，算法通过软阈值化目标的单变量二次主要项来更新单个坐标。Breheny和Huang [22]在他们的命题1中表明，由此过程产生的迭代序列收敛到（6）的逐坐标最小值，这也是一个局部最小值，因为SCAD惩罚的方向导数在任何地方都存在。在我们的实现中，当目标的相对变化低于容忍度时，就宣布收敛。此外，他们[22]建立了一个凸性诊断，用于识别参数空间中SCAD目标局部凸的区域；在这样的区域内，局部最小值是唯一的全局最小值，坐标下降的输出与那里的全局优化器一致。选择预言器局部最小值。SCAD目标可能具有多个局部最小值，定理1的预言器属性适用于其中之一。因此，一个相关的问题是，我们的坐标下降求解器返回的特定局部最小值是否是实现预言器属性的那个。Fan、Xue和Zou [23]提供了最清晰的答案。他们表明，在对预言器估计器和SCAD惩罚的局部化条件下，由-一致估计器（如Lasso）初始化的局部线性近似（LLA）算法在有限次迭代中以趋于1的概率收敛到预言器估计器。我们的实现使用坐标下降而不是显式的LLA，但ncvreg通过热启动构建正则化路径，因此每个SCAD拟合都是从路径上前一个点获得的稀疏、接近Lasso的解决方案开始的。这与[23]中要求的-一致初始化非常相似，他们的理论保证对于求解器的实际行为提供了有价值的信息。统计收敛速率。定理1中声明的预言机属性直接意味着活跃集上估计器的收敛速率。具体来说，定理1(ii)中的渐近正态性声明产生了，其中是真实活跃集的基数。这是经典的参数速率，与事先知道活跃集的预言机最小二乘估计器所达到的速率相同。结合定理1(i)中的选择一致性结论，这表明SCAD惩罚的序列线性模型估计器在支持恢复和收敛速率上都与预言机估计器相匹配。综合这三个观察结果，为第4节中报告的数值性能提供了理论基础。该算法通过[22]收敛到一个局部最小值，当算法通过[23]适当初始化时，这个局部最小值与预言机估计器一致，并且根据定理1，预言机估计器达到了参数速率。4. 模拟 4.1. 模拟设置为了检验SCAD惩罚方法在序列线性模型(2)中的性能，我们进行了一个固定样本量和不同数量预测因子的模拟研究，以及基于五十个时间点的连续结果集。基于每个模拟设置，我们将SCAD惩罚方法与OLS、Lasso和自适应Lasso在传统线性和序列线性模型中进行比较。最优调整参数分别通过10折交叉验证在Lasso、SCAD惩罚和自适应Lasso方法中选择。考虑一个具有重复测量结果的数据集，这些结果在时间点上变化。训练一个模型来依次预测每个时间点的结果，基线数据矩阵为，参数为。在模拟研究中，预测因子的数量设置为低()、中()和高()维。在第一个时间点，模型为。对于任何后续时间点t，具有时间序列的序列线性模型为(10)，其中和。模型(10)迭代生成数据直到。在初始时间点，模型不涉及项。观测数量为。三个预测因子数量级别分别为和1000。我们独立地从均值为和协方差为的多变量正态分布中生成每个时间点的模拟数据矩阵。真实的回归系数为，其余的等于0。我们将有效滞后阶数设置为，真实的时间效应系数为。数据按比例分为训练集和测试集。这种设置是备注1中兼容构造的一个明确实例。生成的数据分别通过线性模型(LM)和序列线性模型(SLM)使用OLS、Lasso、自适应Lasso和SCAD惩罚方法进行拟合。4.2. 估计和变量选择的模拟结果估计系数的模拟结果和平均预测误差是在1000次重复实验上平均得到的。对于估计，使用平均预测误差作为衡量标准。表1展示了每种方法的平均预测误差。表1. 平均绝对预测误差(MAPE)。根据[24]，相对风险(RR)是比较方法的准确性指标。我们希望RR尽可能小，这表示更准确的系数估计。相对风险的公式为，其中表示的协方差矩阵。我们在终止时间点计算RR。完美得分为0，这意味着等于，而RR表示为零。表2展示了每种方法的RR值。表2. 相对风险(RR)。表1显示，对于所有模型和维度，SCAD惩罚方法的平均绝对预测误差在所有方法中是最小的。对于序列线性模型，表1表明SCAD惩罚方法的平均绝对预测误差相对较小，表明预测精度较高。例如，在具有时间序列的中等维度中，估计的为。Lasso、SCAD和自适应Lasso的平均预测误差分别为、和。因此，在具有时间序列的中等维度中，通过序列线性模型拟合的SCAD惩罚方法是最佳选择。同样，在高维度中，SCAD惩罚方法在SLM中获得了最小的MAPE。表2还表明，当通过序列线性模型拟合时，SCAD惩罚方法通常实现较低的平均RR值，证实了它在估计参数方面的最佳性能。此外，对于具有时间序列的生成数据，表2中的RR结果和表1中的MAPE结果表明序列线性建模比简单线性模型表现得更好。这证实了序列线性建模在模型(2)和(3)中的优势。当涉及固定效应和时间序列数据时，序列线性建模是一个好的选择，这反映了引言中讨论的ADL在动态数据中的作用。在生成的数据中，我们还进行了变量选择的模拟。我们设置了五个非零的真实值作为显著变量——第一个、第三个、第五个、第七个和第九个——而其余的都不显著。每个参数的个别经验区间计算为1000次重复实验中估计值的从到分位数的范围。为了检验使用OLS、Lasso、SCAD惩罚和自适应Lasso方法对显著变量估计的准确性，所有显著变量和显著变量数量的经验95%区间（以下简称“ACI”，遵循原始符号）在表3、表4和表5中为三个维度进行了总结。这些方法分别通过传统线性模型(LM)和序列线性模型(SLM)进行拟合。表3. 低维度中的个别95%经验区间(ACI)。表4. 中等维度中的个别95%经验区间(ACI)。表5. 高维度中的个别95%经验区间。表3给出了低维度设置下的这些区间。基于具有的模拟数据，我们预期显著变量数量的区间接近真实值（LM中为5；SLM中为6，包括滞后系数）。在LM中，自适应Lasso获得了最接近真实值的区间。然而，在SLM中，SCAD惩罚方法最好，显著变量数量（6到9.53）和个别参数的区间都最紧，滞后系数的置信区间为(0.59, 0.66)，这是所有方法中最紧的。在中等维度设置中，我们比较了三种惩罚方法。在SLM中，SCAD产生了显著变量数量的95%经验区间为(6, 11)和滞后系数的区间为(0.60, 0.68)，在三种方法中都是最紧的。详细信息见表4。表5显示，在高维度设置中，SCAD惩罚方法再次获得了最佳的95%经验区间。在具有预测因子和1000次重复实验的情况下，SCAD在LM中产生的区间为(5, 12.53)，在SLM中为(6, 17)，这是三种方法中最紧的——其个别系数的区间也接近真实值。其滞后效应系数的区间为(0.59, 0.68)，也比Lasso的更紧。5. 应用在前一节中，模拟结果表明SCAD惩罚方法比这里考虑的其他惩罚方法具有优势。为了进一步比较这些方法的有效性，我们将它们应用于两个纵向研究数据集——joineR包[18]中的Heart.valve和[25]中的Austin House Price——进行估计和变量选择。5.1. 心脏瓣膜数据 5.1.1. 数据描述这里应用的数据涉及主动脉瓣置换手术。它们来自R中的joineR包[18]，并在一项观察性研究中纵向测量手术后的心脏功能结果。对于每个患者，有效子集是患者观察时间点大于或等于5的部分。在这个子数据集中，观察时间点等于5的被用作训练集，时间点大于5的被用作测试集。在各种术前预测因子中，log.lvmi是响应变量。纵向列lvmi是随访访问时的心室质量指数。数据框架包括每个预测因子的988个观测值，包括性别（患者的性别，男性和女性）、年龄（手术当天的年龄，年）、状态（审查指标，死亡和在随访时丢失）、prenyha（术前纽约心脏协会分类）、大小（瓣膜的大小，毫米）、con.cabg（并发冠状动脉旁路移植）、creat（术前血清肌酐，mol/mL）、dm（术前糖尿病）、acei（术前使用ACE抑制剂）和emergence（手术紧急性，选择性、紧急、紧急）。5.1.2. 估计和变量选择的结果为了评估基于患者在五个随访访问时的左心室质量指数对数，图1显示了log.lvmi值的箱线图，显示了每个时间点的四分位数和平均值。手术后每次随访访问时标准化左心室质量指数的分布大致是对称的。只有在时间点5时，响应变量显示出右偏斜。中位数在时间上大致相同，表明连续时间点的log.lvmi值高度相关，当前响应强烈依赖于过去的时间点。图1. 5个时间点的响应变量(log.lvmi)的箱线图。从医学角度来看，lvmi值表示左心室质量指数(LVMI)，[26]讨论了LVMI，这是衡量高血压成人发病率和死亡率的重要指标。[18]显示LVMI与术后心脏恢复有很强的相关性。这里应用了两种类型的模型：简单线性模型和方程(3)的序列线性模型。每种方法比较的是OLS、Lasso、SCAD惩罚和自适应Lasso。表6报告了平均绝对预测误差。SCAD惩罚方法（连同自适应Lasso）在SLM中获得了0.001的MAPE，这是所有方法中最低的。表6. Heart.valve数据中的平均绝对预测误差(MAPE)。对于变量选择，LM中SCAD选择的变量的估计系数在表7中总结，SLM中的系数在表8中总结。表7显示了简单线性模型中log.lvmi的显著预测因子。表8显示，在序列线性模型中，作为响应变量的log.lvmi仅通过（第四个滞后变量）依赖于之前的时间点。因为连续时间点的滞后因变量高度相关，SCAD惩罚方法移除了其他滞后响应变量以提高预测准确性，因此只有一个滞后响应变量进入SLM。表7. LM（Heart.valve数据）中通过SCAD惩罚方法估计的显著变量系数。表8. SLM（Heart.valve数据）中通过SCAD惩罚方法估计的显著变量系数。为了直观展示这个数据集中的估计和变量选择，图2显示了线性模型中显著变量的系数估计。对于序列线性模型，只选择一个变量（，其估计系数已在表8中列出），因此单独的系数图没有信息量，因此省略了。LM中的SCAD方法检测到十三个显著变量，其估计系数显示在表7中，而SLM中的SCAD方法只检测到一个显著变量，因为其他滞后依赖变量高度相关。通过移除不必要的滞后依赖变量，SCAD方法提高了预测准确性。见表6，它显示SLM的预测行为显著改善。图2. LM中显著变量的系数估计。最后但同样重要的是，我们检查了SLM的回归假设，输出图显示在图3中。残差图显示误差项大致遵循正态分布，Q-Q残差图的两个尾部都有一些偏斜。误差项在尺度-位置图中的点范围内大致均匀分布，表明方差大致恒定且没有明显模式。显著的预测变量是相互独立的。根据这些图表，序列线性模型（SLM）的假设得到满足。图3. Heart.valve数据中的SLM诊断。

5.2. 奥斯汀房价数据
5.2.1. 数据描述
该数据集来自Kaggle [25]，最初来源于Eric Pierce的一个项目。2021年，奥斯汀的房地产市场是最活跃的市场之一，这些房源列表提供了近年来市场变化的见解。数据集包括影响房价的各种因素。在清除不相关的预测变量后，该应用中保留了38个变量。数据框包含每个预测变量的15,171个观测值，包括hasAssociation（房源是否与业主协会相关联）、hasSpa（房屋是否有水疗中心）、numOfPhotos（Zillow房源列表中的照片数量）、avgSchoolRating（Zillow房源列表中所有学校类型的平均评分）和numOfBathrooms（房产中的浴室数量）。

5.2.2. 估计和变量选择的结果
目标是预测每个年份建造的房产在数据收集时的最新可用价格。对于数据清洗，我们首先将字符变量转换为因子。然后，我们移除包含缺失值的观测值。接下来，我们根据建造年份识别重复的房产，选择建造年份在2000年或之后的房产。我们创建了一个新的数据框，其中包含两个变量“yearbuilt”和“repeated numbers”，这使我们能够准确地匹配原始数据集并选择每个数据段。我们将数据分为训练集和测试集：训练集包含重复的yearbuilt值，其计数小于或等于180，而剩余数据构成测试集。经过这些处理后，我们得到了37个独特的时间点。序列线性模型将第37个时间点作为结果，前36年的房产价格作为预测变量。我们在序列线性模型中应用OLS、Lasso、SCAD惩罚和自适应Lasso来比较这些方法。图4绘制了37个时间点上最新房产价格的平均值。平均值波动较大，并显示出随时间增长的整体趋势，因此房价在序列线性建模中应具有时间序列依赖变量和其他固定效应独立变量。

图4. 随时间点变化的最新房价平均值。对于估计和变量选择，SLM中通过SCAD选择的系数估计在表9中总结。当前房价取决于之前的时间点（），以及几个独立变量：hasAssociation、hasSpa、numOfPhotos、avgSchoolRating和numOfBathrooms。表9. SLM中通过SCAD惩罚方法选择的显著变量系数估计（奥斯汀房价数据）。图5中的地图显示了房产的分布情况。每个点代表一个房产，红色（最近价格最低）、黄色（第二低）、绿色（中等）和蓝色（最高）。地图揭示了哪些区域的房产价值相对较高。如表9所示，选定的变量——浴室数量、是否有景观、是否有水疗中心以及是否属于业主协会——以及一个滞后依赖变量（）出现了。由于高相关性，其他滞后依赖变量从序列线性模型中移除。这些选定的变量在近年来显著影响了德克萨斯州奥斯汀的房价。浴室数量和水疗中心的存在对房价有正面影响，而属于业主协会则对房价有负面影响。

表10报告了平均绝对预测误差。同样，SCAD惩罚方法在惩罚方法中实现了最小的值（89.089）。表10. SLM中平均（绝对）预测误差的比较（奥斯汀，德克萨斯州房价数据）。对于显著变量选择，序列线性模型中的Lasso系数在图6中显示。除了上述变量（浴室数量、景观、水疗中心、协会），Lasso选择的模型还包括几个滞后依赖变量。图6清楚地反映了稀疏性属性：不显著变量的估计系数被缩小到零。图6. 通过Lasso方法选择的显著变量系数估计。从图7中的自适应Lasso系数图中可以看出，选定的变量是Lasso选择的子集，并且没有保留任何滞后依赖变量：基于Lasso系数，自适应Lasso对每个系数应用不同的惩罚权重，并由于高相关性移除了所有滞后依赖变量。自适应Lasso选择的显著变量比Lasso少，这与Lasso相对于自适应Lasso的过拟合一致。由于显著变量较少，自适应Lasso在预测方面表现优于Lasso（见表10）。图7. 通过自适应Lasso方法选择的显著变量系数估计。在图8中，SCAD选择的变量包括自适应Lasso选择的变量，以及第25个时间点的滞后房产价格。浴室数量和水疗中心的存在对房价有正面影响，而属于协会则对房价有负面影响。

我们得出结论，SCAD方法不仅在SLM中选择了预测变量，还选择了滞后依赖变量的顺序，这与推论1一致。这种联合选择属性使得SCAD方法特别适合序列线性模型。图8. 通过SCAD方法选择的显著变量系数估计。在惩罚方法中，SCAD方法在表10中实现了最小的平均绝对预测误差89.089。

最后但同样重要的是，我们检查了SLM的回归假设。输出图表显示在图9中。第一个残差图显示了非线性模式，表明在建模响应变量时可能需要一些非线性预测变量。Q-Q图显示了近似正态分布，但两个尾部略有偏斜。在尺度-位置图中，异方差性似乎存在，应在方差估计中加以考虑。杠杆图中可以看到一些异常值。总体而言，与SLM假设的偏差并不严重，这里使用的序列线性模型很好地捕捉了房价的趋势。估计和选定的变量对模型解释很有信息量。图9. 奥斯汀，德克萨斯州房价数据中的SLM诊断。

6. 结论性评论和讨论
在序列线性建模中，每个时间点的测量值使得可以预测中间时间点的数据，从而提高预测准确性，因为正确捕捉了相关结构。由于关注序列线性回归模型，滞后依赖变量被反复作为输入变量包含在下一个预测模型中。因此，我们探讨了序列线性建模及其在估计和变量选择中的性能。

对于本研究中提出的序列线性模型，我们引入了SCAD惩罚方法来改进估计和变量选择。我们建立了估计量的神谕属性（定理1），并展示了滞后阶数选择是一致的（推论1）。预计SCAD惩罚在估计和变量选择中的表现将优于Lasso和自适应Lasso。Lasso通过将系数缩小到零来进行变量选择，这在中等或高维度中很有用。自适应Lasso对参数估计应用自适应惩罚权重，改进了Lasso的表现。SCAD惩罚方法在不同系数大小下有不同的收缩方式：对于小系数，它的行为类似于Lasso，将其缩小到零；对于大系数，惩罚是恒定的，减少了偏差并保留了显著信号；对于中等系数，惩罚的增长速度比Lasso慢，减少了收缩。

我们进行模拟，以比较所提出的SCAD惩罚方法与Lasso和自适应Lasso的性能。在低维度设置中使用OLS作为参考。在低、中、高维度中，我们生成了具有1000个重复次的基数据集，误差项遵循标准正态分布。所提出估计器的算法和统计收敛性在第3.5节中通过引用[22,23]的结果进行了理论讨论。通过将滞后依赖变量作为每个时间点的预测变量，我们比较了OLS、Lasso、SCAD和自适应Lasso。OLS仅适用于低维度设置，当p超过n时表现不佳。

模拟结果显示，SCAD惩罚方法通过减少参数估计的偏差和降低预测误差而优于Lasso。SCAD惩罚提高了估计准确性和变异性，在中等和高维度设置中产生了更好的预测结果。结果通过MAPE、RR和经验95%区间（ACI）输出在各种方法和场景中进行了总结。95%的经验区间显示，SCAD惩罚方法得到的结果最接近真实设置，无论是显著变量的数量还是真实参数值。值得注意的是，在小数据集中OLS表现良好，而在中到高维度数据集中，SCAD惩罚方法在参数估计和变量选择方面更为优越。

为了进一步探索所提出的SCAD惩罚方法的性能，我们将其应用于两个真实数据集，并与OLS、Lasso和自适应Lasso进行了比较。第一个数据集是joineR包[18]中的心脏瓣膜置换手术数据。SCAD在线性和序列模型中都实现了最低的平均绝对预测误差。第二个数据集是德克萨斯州奥斯汀的最新房价[25]。尽管建模假设没有完全满足，但移除极端异常值或对响应变量取对数加强了预测变量与对数房价之间的线性关系。SCAD惩罚方法在变量选择中表现出更稳定的性能，包含了自适应Lasso也识别出的大多数显著变量。

在未来的工作中，我们将考虑分析具有强多重共线性的数据集，并在该设置中比较OLS、Lasso、SCAD和自适应Lasso。许多真实数据集涉及变量之间的高相关性，这激发了在从多变量时间序列中选择显著滞后依赖变量时应用各种惩罚方法到序列线性建模中的动机。我们还将探讨假设检验，以确定序列线性建模中可靠模型预测所需的最小样本量，并估计使模型对显著变量具有鲁棒性并产生良好预测性能所需的样本量。