面向UAR模型及其高阶UAR模型的修正的SCAD惩罚估计方法

Modified SCAD Penalty Estimation for UAR Model and Higher-Order UAR Model

摘要: 平滑剪切绝对偏差(SCAD)惩罚估计虽已应用于不确定自回归(UAR)模型，但其依赖样本信息而未纳入先验信息；同时，能融合先验信息的修正的SCAD惩罚估计方法尚未拓展至UAR模型及其高阶UAR模型。为此，本文将修正的SCAD惩罚估计方法拓展至UAR模型及其高阶UAR模型，构建了面向UAR模型及其高阶UAR模型的修正的SCAD惩罚估计方法，并通过两个数值算例分别验证其可行性与有效性。

Abstract: While the Smoothly Clipped Absolute Deviation (SCAD) penalty estimation has been applied to the Uncertain Autoregressive (UAR) model, it relies solely on sample information and fails to incorporate prior knowledge. Meanwhile, the modified SCAD penalty estimation method, which can integrate such prior information, has not yet been extended to the UAR model or its higher-order counterparts. To address this gap, this paper extends the modified SCAD penalty estimation to both the UAR model and its higher-order forms. The corresponding estimation framework is constructed, and two numerical examples are presented to respectively validate its feasibility and effectiveness.

1. 引言

自回归(Autoregressive, AR)模型是传统时间序列模型的重要组成部分，在金融、环境等领域得到了广泛应用。Yule [1]首先给出了AR模型，之后，众多学者[2]-[5]对AR模型开展了进一步的完善和提升，并取得了丰硕的成果。然而，AR模型在处理具有不确定性的数据时存在局限，为此，Liu等[6] [7]将不精确数据视为不确定变量，构建不确定自回归(Uncertain Autoregressive, UAR)模型，并应用于二氧化碳浓度预测等方面。Lu和Peng [8]给出了高阶UAR模型。

在构建UAR模型和高阶UAR模型后，阶数确定和参数估计成为需要解决的两个重要问题[9]-[13]。在模型阶数确定方面，Liu和Yang [9]、Zhang等[14]、Chen和Yang [15]均采用假设方式确定UAR模型阶数。在模型参数估计方面，Yang等[10]提出了最小绝对偏差估计方法来估计UAR模型的未知参数，但该方法容易出现非唯一解的问题，Zhang等[14]、Chen与Yang [12]利用LASSO估计和岭估计来估计UAR模型的未知参数，通过收缩偏大的系数，避免了最小绝对偏差估计方法非唯一解问题，但是相对应的会丢失一些观测信息，而且往往难以收缩到0且不易达到变量选择的目的。Chen和Yang [15]提出了运用极大似然估计方法来确定UAR模型的未知参数，通过该方法既可以估计未知的自回归参数，也可以得到扰动项的不确定分布参数。Gao与Xin [16]进一步提出极大似然估计方法来确定UARMA模型的未知参数，并将其应用于金融市场数据分析。相关结果表明，该模型在包含异常波动的数据场景中具有更好的稳健性与拟合灵活性。

在上述研究中，研究者们通常预设UAR模型阶数已知，然后根据这个假设的阶数进行参数估计。实际上，模型阶数是未知的，如何客观地确定模型阶数来进行参数估计是十分重要的。为此，Fan和Li [17]提出了用于变量选择和参数估计的平滑剪切绝对偏差(Smoothly Clipped Absolute Deviation, SCAD)惩罚函数，SCAD惩罚函数主要是在高维数据中进行变量选择，进而达到参数估计的目的，而Li和Yang [18]通过引入SCAD惩罚项，对于给定的${\lambda }_n$ ，实现了对UAR模型阶数与未知参数的同步估计。但是该方法未考虑到先验信息的影响，Ng和Yu [19]提出了一种修正的SCAD惩罚函数，建立了约束变量选择问题的信息准则AIC和BIC，用来确定最优调谐参数，进而提升广义线性回归(Generalized Linear Regression, GLR)模型的拟合优度和稀疏性。不过，该修正的SCAD惩罚函数尚未应用于UAR模型及其高阶UAR模型的参数估计问题。

为了有效解决上述存在的问题，本文将修正的SCAD惩罚估计方法拓展至不确定时间序列分析中的UAR模型及其高阶UAR模型，以实现其阶数与参数的客观估计。

2. 面向UAR模型和高阶UAR模型的修正SCAD惩罚估计方法

Li和Yang [18]提出了一种面向UAR模型的传统的SCAD惩罚估计方法。然而，传统的SCAD惩罚估计方法依赖于固定的调谐参数，缺乏选择惩罚权重的灵活性，难以在实践中充分利用信息准则(AIC和BIC)进行优化。为了避免传统的SCAD惩罚函数的缺点，Ng和Yu [19]提出了一种修正的SCAD惩罚函数，该函数不仅保留了线性约束下的oracle属性，而且允许通过先验信息灵活选择调谐参数，他们将其应用于广义线性回归模型。本节将修正的SCAD惩罚函数应用于UAR模型和高阶UAR模型的参数估计。

2.1. 面向UAR模型和高阶UAR模型的修正SCAD惩罚估计方法

用于预测未来值的UAR模型定义如下：

$X_t={\displaystyle \sum _{i=1}^p{c}_i}{X}_{t-i}+{\epsilon }_t$ (1)

其中${\epsilon }_t$ 为扰动项，$c_i\left(i=1,2,\cdots ,p\right)$ 为UAR模型的不确定自回归系数，设$c={\left(c_1,c_2,\cdots ,c_p\right)}^{\text{T}}$ 为不确定自回归系数向量。设$p=n_0$ ，当$n_0$ 充分大时，(1)可以表示为高阶UAR模型[8]。

下面以高阶UAR模型为例，给出面向高阶UAR模型的修正SCAD惩罚估计方法。为确保模型参数估计满足实际情况，并提升估计精度，我们引入线性约束$Rc=r$ ，其中，$R$ 是$q\times n_0$ 约束矩阵，$r$ 是$q$ 维向量。该线性约束要求满足。此约束减少了参数空间的自由度，从而增强了模型的稳定性和稀疏性，尤其在样本量较小的情况下效果显著。

基于修正的SCAD惩罚估计函数，我们提出了面向高阶UAR模型的修正的SCAD惩罚估计方法，该估计通过在约束$Rc=r$ 下最小化以下目标函数，将先验信息和样本信息结合起来：

$Q\left(c|X_n\right)+n{\displaystyle \sum _{i=1}^{n_0}{P}_{{\lambda }_n}}\left(\left|c_i\right|\right)$ (2)

其中$Q\left(c|X_n\right)$ 是最小二乘目标函数，$P_{{\lambda }_n}\left(\left|c_i\right|\right)$ 是约束条件$Rc=r$ 下的修正SCAD惩罚函数。

面向高阶UAR模型的修正的SCAD惩罚估计具体如下：

$\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)=\underset{c_1,c_2,\cdots c_{n_0}}{\arg \min }\left\{{\displaystyle \sum _{t=n_0+1}^{n}E\left[{\left(X_t-{\displaystyle \sum _{i=1}^{n_0}{c}_i{X}_{t-i}}\right)}^2\right]}+\left(n-n_0\right){\displaystyle \sum _{i=1}^{n_0}{P}_{{\lambda }_n}\left(c_i\right)}\right\}$

s.t. $Rc=r$ ，其中修正的SCAD惩罚函数为：

$P_{{\lambda }_n}\left(c_i\right)=\{\begin{array}{l}\begin{array}{cc}{\lambda }_n{k}^{-1}\left[{\left(1+\left|c_i\right|\right)}^k-1\right],& \begin{array}{cc}\begin{array}{cc}& \end{array}& \left|c_i\right|\le {\lambda }_n,\end{array}\end{array}\hfill \\ \begin{array}{cc}\frac{-{\left(1+{\lambda }_n\right)}^{k-1}}{2\left(a-1\right)}\left[c_i^2-2a{\lambda }_n\left|c_i\right|+\left(2a-1\right){\lambda }_n^2\right]+\frac{{\lambda }_n}{k}\left[{\left(1+{\lambda }_n\right)}^k-1\right],& \end{array}{\lambda }_n<\left|c_i\right|\le a{\lambda }_n,\hfill \\ \begin{array}{cc}\frac{a-1}{2}{\left(1+{\lambda }_n\right)}^{k-1}{\lambda }_n^2+\frac{{\lambda }_n}{k}\left[{\left(1+{\lambda }_n\right)}^k-1\right],& \begin{array}{cc}& \left|c_i\right|>a{\lambda }_n,\end{array}\end{array}\hfill \end{array}$ (3)

其中，调节参数$a=3.7,{\lambda }_n>0$ 。

注：当$k=1$ 时，修正的SCAD惩罚函数退化为Fan和Li [18]提出的SCAD惩罚函数。

在选择最优调谐参数时，需要综合考虑其与样本量的关系。除了使用来控制惩罚强度外，还可以采用AIC、BIC等信息准则来确定最优调谐参数。这些信息准则基于修正的SCAD惩罚函数和有效参数的自由度。面向高阶UAR模型的修正的SCAD惩罚估计方法通过引入非凸惩罚项，提出了新的AIC和BIC信息准则，有效解决了传统准则在高惩罚设置下的局限性(此时AIC/BIC曲线可能出现单调性)，从而提升了模型的稀疏性和适应性。

下面通过计算有效自由度$\text{eff}\left({\lambda }_{\text{n}}\right)$ 确定最佳调谐参数${\lambda }_n$ ，其中$\text{eff}\left({\lambda }_{\text{n}}\right)$ 定义如下

$\text{eff}\left({\lambda }_n\right)=\text{trace}\left(H^R\left({\lambda }_n\right)\right)$ ，

其中$H^R\left({\lambda }_n\right)$ 是考虑惩罚和约束的帽子矩阵，$H^R\left({\lambda }_n\right)$ 定义如下：

$H^R\left({\lambda }_n\right)=X_n\left(A_n^{-1}-A_n^{-1}{R}^T{\left(R{A}_n^{-1}{R}^T\right)}^{-1}R{A}_n^{-1}\right)X_n^T$ ，

其中$X_n$ 为设计矩阵，$A_n=X_n^T{X}_n+n{\lambda }_n{I}_{n_0}$ 为带惩罚项的设计矩阵，$R$ 为约束矩阵，$I_{n_0}$ 是$n_0\times n_0$ 的单位矩阵。帽子矩阵$H^R\left({\lambda }_n\right)$ 将观测数据投影到拟合值上同时整合了惩罚项和约束条件。有效自由度$\text{eff}\left({\lambda }_n\right)$ 反映了模型在这些约束条件下的复杂程度。

为解决上述AIC和BIC在高惩罚或约束条件下的局限性，我们提出了一种约束条件下的信息准则，其中AIC和BIC的表达式如下：

$\text{AIC}=-2\log L\left(\widehat{c}|X\right)+2\text{eff}\left({\lambda }_n\right)$，

$\text{BIC}=-2\log L\left(\widehat{c}|X\right)+\text{eff}\left({\lambda }_n\right)\log n$，

在给定约束条件下，合理选择AIC和BIC能在模型拟合度与稀疏性之间取得平衡。该准则能有效惩罚小系数，同时尽量减少对大系数的影响，从而避免AIC和BIC在高惩罚下出现的单调性问题。

2.2. 求解未知参数向量c

在本节中，我们讨论如何利用修正的SCAD惩罚估计确定UAR模型及其高阶UAR模型的阶数并估计其未知参数。首先，提出一种基于AIC和BIC信息准则的算法，用于从候选参数集中筛选最优调参组合。最后，通过计算均方绝对误差(MAD)，我们设定最优调谐参数，从最优候选参数集中选取最佳调谐参数，从而确定UAR模型及其高阶UAR模型的参数估计。

定理[18]：设不确定观测值$\left\{X_1,X_2,\cdots ,X_n\right\}$ 是独立同分布的不确定变量，其标准不确定分布为${\Phi }_1,{\Phi }_2,\cdots ,{\Phi }_n$ ，满足高阶UAR模型。为计算参数向量$\left(c_1,c_2,\cdots ,c_{n_0}\right)$ ，利用修正的SCAD惩罚估计方法的最优解$\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)$ 也是下列目标函数的最优解：

$\underset{c_1,c_2,\cdots ,c_{n_0}}{\min }\left\{{\displaystyle \sum _{t=n_0+1}^n{\displaystyle {\int }_0^1{\left({\Phi }_t^{-1}\left(\alpha \right)-{\displaystyle \sum _{i=1}^{n_0}{c}_i{\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)}\right)}^{2}d\alpha +\left(n-n_0\right){\displaystyle \sum _{i=1}^{n_0}{P}_{{\lambda }_n}\left(\left|c_i\right|\right)}}}\right\}$ ， (4)

其中${\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)=\{\begin{array}{l}{\Phi }_{t-i}^{-1}\left(\alpha \right),c_i\ge 0\hfill \\ {\Phi }_{t-i}^{-1}\left(1-\alpha \right),c_i<0\hfill \end{array}$ ，$i=1,2,\cdots ,n_0$ 。

为了初始化修正的SCAD惩罚估计算法，我们定义以下目标函数：

$f\left(c_1,c_2,\cdots ,c_{n_0}\right)=\left\{{\displaystyle \sum _{t=n_0+1}^n{\displaystyle {\int }_0^1{\left({\Phi }_t^{-1}\left(\alpha \right)-{\displaystyle \sum _{i=1}^{n_0}{c}_i{\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)}\right)}^{2}d\alpha +\left(n-n_0\right){\displaystyle \sum _{i=1}^{n_0}{P}_{{\lambda }_n}\left(\left|c_i\right|\right)}}}\right\}$ ， (5)

$\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)=\underset{c_1,c_2,\cdots ,c_{n_0}}{\arg \min }f\left(c_1,c_2,\cdots ,c_{n_0}\right)$ (6)

其中${\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)=\{\begin{array}{l}{\Phi }_{t-i}^{-1}\left(\alpha \right),c_i\ge 0,\hfill \\ {\Phi }_{t-i}^{-1}\left(1-\alpha \right),c_i<0.\hfill \end{array}$

为求解(4)所对应的非凸惩罚最小化问题，采用局部二次近似思想对惩罚项进行二次化，从而在每次迭代中转化为带权二次规划问题，便于使用牛顿–拉弗森(Newton-Raphson)算法求解。具体地，设第m次迭代的当前估计为

${\widehat{c}}^{\left(m\right)}=\left({\widehat{c}}_1{}^{\left(m\right)},{\widehat{c}}_2{}^{\left(m\right)},\dots ,{\widehat{c}}_{n_0}{}^{\left(m\right)}\right)$

当$\left|{\widehat{c}}_i{}^{\left(m\right)}\right|$ 非常接近0时，为增强稀疏性可直接令$\left|{\widehat{c}}_i{}^{\left(m\right)}\right|=0$ ；否则在$\left|{\widehat{c}}_i{}^{\left(m\right)}\right|$ 的邻域内，对$p_{{\lambda }_n}\left(\left|c_i\right|\right)$ 作二阶泰勒展开并利用$\left|c_i\right|$ 的局部线性化，可得

$\begin{array}{l}{p}_{{\lambda }_n}\left(\left|c_i\right|\right)\approx p_{{\lambda }_n}\left(\left|{\widehat{c}}_i^{\left(m\right)}\right|\right)+{p^{\prime }}_{{\lambda }_n}\left(\left|{\widehat{c}}_i^{\left(m\right)}\right|\right)\text{sgn}\left(c_i^{\left(m\right)}\right)\left(c_i^{\left(m\right)}-{\widehat{c}}_i^{\left(m\right)}\right)+\frac{\text{1}}{\text{2}}p"\left(\left|{\widehat{c}}_i^{\left(m\right)}\right|\right){\left(c_i^{\left(m\right)}-{\widehat{c}}_i^{\left(m\right)}\right)}^{\text{2}}\\ =p_{{\lambda }_n}\left(\left|{\widehat{c}}_i^{\left(m\right)}\right|\right)+\frac{{p^{\prime }}_{{\lambda }_n}(\left|{\widehat{c}}_i^{\left(m\right)}\right|)}{2\left|{\widehat{c}}_i^{\left(m\right)}\right|}\left(c_i^{\left(m\right)}-{\left({\widehat{c}}_i^{\left(m\right)}\right)}^{\text{2}}\right),c_i\ne 0\end{array}$

将上述二次化近似代入(4)，并略去与$c$ 无关的常数项，可得到第$m$ 次迭代的近似目标函数

$f^{\left(m\right)}\left(c\right)={\displaystyle \sum _{t=n_0+1}^n{\displaystyle {\int }_0^1{\gamma }_{\text{t}}{\left(\alpha ;c\right)}^2}}d\alpha +\left(n-n_0\right){\displaystyle \sum _{i=1}^{n_0}\left|c_i^2\right|}$

因此(4)每次迭代中等价于一个带权ridge型的约束/无约束最小二乘子问题。进一步，在将积分项离散化后，$f^{\left(m\right)}\left(c\right)$ 的梯度与海森矩阵可显式计算，据此使用牛顿–拉弗森迭代更新。若同时存在线性约束$Rc=r$ ，则可在每次更新后采用KKT系统求解或将无约束更新结果投影到可行域，以保证迭代点满足约束，重复迭代直至收敛。由此(4)的积分型目标函数可与修正SCAD惩罚项共同转化为可计算的数值优化流程。

设${\tilde{X}}_1,{\tilde{X}}_2,\cdots ,{\tilde{X}}_n$ 是一个不确定的时间序列，$c_1^{*},c_2^{*},\cdots ,c_{n_0}^{*}$ 分别是$c_1,c_2,\cdots ,c_{n_0}$ 的估计值。定义信息准则为

$\text{MAD}=\frac{1}{n-n_0+1}{\displaystyle \sum _{t=n_0+1}^n{\displaystyle {\int }_0^1\left|{\Phi }_t^{-1}\left(\alpha \right)-{\displaystyle \sum _{i=1}^{n_0}{c}_i^{*}{\gamma }_{t-i}^{-1}\left(\alpha ,c_i^{*}\right)}\right|d\alpha }}$ ， (7)

其中MAD是评估预测模型性能的指标，用于衡量实际观测值与模型预测值之间的平均绝对差异，此时，具有最小信息准则的调谐参数${\lambda }_n$ 对应的UAR模型阶数即为最优模型阶数。

因此，未知参数求解算法如下所示：

步骤1：通过(6)初始化参数$\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)$ 。

步骤2：根据参数$\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)$ ，遍历${\lambda }_n$ ，步长为0.01利用2.1节的方法计算一组AIC和BIC，选取最小AIC值对应的${\lambda }_n$ ，记为${\lambda }_n^{\text{AIC}}$ ，选取最小BIC值对应的${\lambda }_n$ 记为${\lambda }_n^{\text{BIC}}$ 。

步骤3：将${\lambda }_n={\lambda }_n^{\text{AIC}}$ 和${\lambda }_n={\lambda }_n^{\text{BIC}}$ 分别带入(4)，使用牛顿–拉弗森算法通过迭代公式

${\widehat{c}}_{i\left(m+1\right)}={\widehat{c}}_{i\left(m\right)}-{\left(\frac{{\partial }^{2}f}{\partial c_i^2}\right)}^{-1}\left(\frac{\partial f}{\partial c_i}\right),i=1,2,\cdots ,n_0,m=1,2,\cdots$ ，

更新参数${\widehat{c}}_i$ 分别获得两组高阶UAR模型参数${\widehat{c}}^{\text{AIC}}$和${\widehat{c}}^{\text{BIC}}$。

步骤4：将两组系数向量${\widehat{c}}^{\text{AIC}}$和${\widehat{c}}^{\text{BIC}}$分别带入(7)，计算MAD，选取使MAD最小化的参数$\left(c_1^{*},c_2^{*},\cdots c_{n_0}^{*}\right)$ 作为当前最优参数估计，并记录对应的最优调谐参数${\lambda }_n^{*}$ ，以及高阶UAR模型参数向量${\widehat{c}}^{*}$ 。

步骤5：使用步骤4得出${\widehat{c}}^{*}={\left({\widehat{c}}_1^{*},{\widehat{c}}_2^{*},\cdots ,{\widehat{c}}_{n_0}^{*}\right)}^{\text{T}}$替换上一轮得到的$\widehat{c}={\left({\widehat{c}}_1,{\widehat{c}}_2,\cdots ,{\widehat{c}}_{n_0}\right)}^{\text{T}}$代入

$\left|\frac{{\widehat{c}}_i^{*}-{\widehat{c}}_i}{{\widehat{c}}_i}\right|,i=1,2,\cdots ,n_0$ ，

如果全部${\widehat{c}}_i^{*}$ 满足条件

$\left|\frac{{\widehat{c}}_i^{*}-{\widehat{c}}_i}{{\widehat{c}}_i}\right|\le {10}^{-3},i=1,2,\cdots ,n_0$ ，

则高阶UAR模型的最优参数向量${\widehat{c}}^{*}$ 被获取，并记录对应的最终${\lambda }_n$ 。否则返回步骤2。

2.3. 残差分析与假设检验

高阶UAR模型残差如下所示：

${\widehat{\epsilon }}_t={\widehat{X}}_t-{\displaystyle \sum _{i=1}^{n_0}{c}_i^{*}{X}_{t-i}}$ .

因此，不确定残差的期望值计算方法如下：

$\widehat{e}=\frac{1}{n-n_0}{\displaystyle \sum _{t=n_0+1}^{n}E\left[{\widehat{\epsilon }}_t\right]}$ .

下面公式用于计算干扰的方差：

${\widehat{\sigma }}^2=\frac{1}{n-n_0}{\displaystyle \sum _{t=n_0+1}^{n}E\left[{\left({\widehat{\epsilon }}_t-\widehat{e}\right)}^2\right]}$ .

接下来，我们介绍线性不确定变量的逆不确定性分布

${\psi }^{-1}\left(\alpha \right)=\left(1-\alpha \right)x+\alpha y,x>y,x,y\in R$ .

假设${\epsilon }_t~N\left(e^{*},\sigma *\right)$ ，参数$e^{*}$ 和${\sigma }^{*}$ 由Yang和Liu [7]给出的方法估计出来，其中不确定假设检验采用

$G_0:e=e^{*}和\sigma ={\sigma }^{*}$ 与$G_1:e\ne e^{*}和\sigma \ne {\sigma }^{*}$ ，

在显著性水平$\theta$ 下，检验定义为

$H=\left\{\begin{array}{l}\left({\widehat{\epsilon }}_{n_0},{\widehat{\epsilon }}_{n_0+1},\cdots ,{\widehat{\epsilon }}_n\right):存在至少\theta 个满足n_0\le t\le n的t值,使得\\ 当{\psi }_t^{-1}\left(1-\theta \right)<{\Phi }_{\beta }^{-1}\left(\theta \right)或{\psi }_t^{-1}\left(\theta \right)<{\Phi }_{\beta }^{-1}\left(1-\theta \right)时{\widehat{\epsilon }}_t~{\psi }_t\end{array}\right\}$

其中${\Phi }_{\beta }^{-1}\left(\theta \right)=e^{*}+\frac{{\sigma }^{*}}{\pi }\sqrt{3}\ln \frac{\theta }{1-\theta }$ 。

若残差向量$\left({\widehat{\epsilon }}_{n_0},{\widehat{\epsilon }}_{n_0+1},\cdots ,{\widehat{\epsilon }}_n\right)$ 属于$H$ ，即$\left({\widehat{\epsilon }}_{n_0},{\widehat{\epsilon }}_{n_0+1},\cdots ,{\widehat{\epsilon }}_n\right)\in H$ ，若拒绝$G_0$ ，则认为模型不通过检验，否则不能拒绝$G_0$ 。

2.4. 数值实验

例1 在这个示例中，一组不精确的观测值$X_1,X_2,\cdots ,X_{30}$ 由以下UAR模型生成

$X_t=1.0129X_{t-1}+0.0326X_{t-2}+{\epsilon }_t$ ，

其中$X_1,X_2,\cdots ,X_{30}$ 是相互独立的线性不确定变量，下面将说明如何利用修正的SCAD惩罚估计方法确定模型阶数并估计未知参数。若不考虑惩罚函数，可采用最小二乘法进行估算，得到

$\left(c_1^{*},c_2^{*},\cdots ,c_{10}^{*}\right)=\left(1.025112,0.41629,0.00051,0.0001,0,0,0,0,0,0\right)$

随后，将最小二乘法估计值作为初始值代入第2.2节面向UAR模型的修正的SCAD惩罚估计算法中。

本例中将不同惩罚强度k设为0.5和1，在算法的每次迭代中，生成两组${\lambda }_n$ ，具体结果如表1，选取对应最小MAD的${\lambda }_n=0.5722$ 作为最终的调参参数，因此高阶UAR模型的参数估计如下所示

$\left(c_1^{*},c_2^{*},\cdots ,c_{10}^{*}\right)=\left(1.0242,0.0325,0,0,0,0,0,0,0,0\right)$

因此，拟合的UAR模型为

$X_t=1.0242X_{t-1}+0.0325X_{t-2}+{\epsilon }_t$ .

该模型采用UAR(2)模型进行参数估计，实验证实修正的SCAD惩罚估计方法能有效应用于UAR模型的参数估计与阶数选择，实现参数稀疏性。与最小二乘法相比，该方法对真实参数值的估计更为接近。通过精准估计模型参数，显著提升了模型预测精度。基于修正的SCAD惩罚估计方法获得的参数估计值，我们进一步计算了下一阶段的预测值及预测区间。

Table 1. The model selection results for different values of k in Example 1

表1. 例1中不同k值的模型选择结果

	${\lambda }_n$	$\text{eff}\left({\lambda }_n\right)$	AIC/BIC	MAD
k = 1	0.5722	4.2471	−2.2129/−2.0146	5.2372
	0.4641	4.4251	−2.5260/−1.9360	6.0075
k = 0.5	0.6767	4.4111	−2.0615/−2.0012	7.0511
	0.5512	4.1767	−1.8813/−1.4513	7.1211

假设扰动项${\epsilon }_t$ 与观测值$X_{29},X_{30}$ 相互独立，采用修正的SCAD惩罚法估计的UAR模型中$\delta =0.6272$ ，由此得出$X_{31}$ 的95%置信区间为

$19.8190\pm 0.6272$ .

在上述模拟实验中，结果表明，修正的SCAD惩罚估计方法比最小二乘估计法更接近真实参数值，能正确估计模型参数$\left(c_1,c_2\right)=\left(1.0129,0.0326\right)$ ，提高预测精度，也验证了改进SCAD估计法的可行性。

例2 美国夏威夷莫纳罗亚天文台1997~2017年月度CO₂月均值排放量数据集表明，月度碳排放量水平是相互独立的线性不确定变量，如图1所示。

类似地，本例题采用$n_0=10$ 来预测后续未来值。将最小二乘估计值作为初始值代入改进的SCAD方法算法中。

Figure 1. The 95% confidence interval for ${\tilde{X}}_{22}$ in Example 2

图1. 例2中${\tilde{X}}_{22}$ 的95%置信区间

根据表2数据，当${\lambda }_n=3.34$ 时均方绝对误差达到最小值，此时得到高阶UAR的参数估计结果如下：

$\left(c_1^{*},c_2^{*},\cdots ,c_{10}^{*}\right)=\left(0.9922,0.0065,0,0,0,0.0016,0,0,0,0\right)$ .

Table 2. The model selection results for different values of k in Example 2

表2. 例2中不同k值的模型选择结果

	${\lambda }_n$	$\text{eff}\left({\lambda }_n\right)$	AIC/BIC	MAD
k = 1	3.2700	2.8981	−1.3954/−1.2513	9.0489
	3.3400	2.8978	−2.6560/−2.1360	8.7875
k = 0.5	3.2200	2.5978	−1.4676/−1.2314	10.9583
	2.9945	2.7150	−2.7954/−1.9420	8.8414

拟合的高阶UAR模型如下

$X_t=0.9922X_{t-1}+0.0065X_{t-2}+0.0016X_{t-6}+{\epsilon }_t$

第22年的碳排放量预测值计算如下：

$\mu ={\displaystyle \sum _{i=1}^{10}{c}_i^{*}E\left[X_{22-i}\right]}+\widehat{e}=407.9525$ .

对第22年碳排放量$X_{22}$ 进行计算置信区间。根据该方法，第22年碳排放量的95%置信区间为

407.9525 ± 3.6172.

通过与Chen和Yang [15]提出的最大似然估计(MLE)方法对比，修正的SCAD惩罚估计方法展现出更窄的置信区间估计，其区间宽度比MLE方法约窄0.65%。

在例2中，取$n_0=10$ 并将最小二乘估计值作为初始值代入改进的SCAD算法，进一步地，为展示本文方法在存在先验信息时的优势，现考虑高阶系数向量$c=\left(c_1,c_2,\dots ,c_{10}\right)$ 满足线性等式约束$Rc=r$ 结合例2的稀疏识别结果可见，模型主要由$$ $c_1,c_2,c_6$ 三个滞后项决定。基于短期演化可由关键滞后项的加权平均刻画，给出如下约束$c_1+c_2+c_6=1$ 。

在该先验约束下，本文2.2节的修正SCAD估计(见式(5)~(6))可自然推广为如下约束优化问题

$\underset{c_1,c_2,\cdots ,c_{n_0}}{\min }\left\{{\displaystyle \sum _{t=n_0+1}^n{\displaystyle {\int }_0^1{\left({\Phi }_t^{-1}\left(\alpha \right)-{\displaystyle \sum _{i=1}^{10}{c}_i{\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)}\right)}^{2}d\alpha +\left(n-10\right){\displaystyle \sum _{i=1}^{n_0}{P}_{{\lambda }_n}\left(\left|c_i\right|\right)}}}\right\},\text{s}\text{.t Rc}=\text{r}$ ，

其中${\gamma }_{t-i}^{-1}\left(\alpha ,c_i\right)=\{\begin{array}{l}{\Phi }_{t-i}^{-1}\left(\alpha \right),c_i\ge 0\hfill \\ {\Phi }_{t-i}^{-1}\left(1-\alpha \right),c_i<0\hfill \end{array}$ ，$i=1,2,\cdots ,10$ 。

计算上仍采用与例2一致的参数选取策略：在候选集合中比较不同k与${\lambda }_n$ 的模型选择结果，并以MAD最小原则确定最优惩罚参数，从而得到满足$\text{R}\widehat{\text{c}}=\text{r}$系数估计如下：

$\left(c_1^{*},c_2^{*},\cdots ,c_{10}^{*}\right)=\left(\text{0}\text{.9919},\text{0}\text{.0064},0,0,0,0.0016,0,0,0,0\right)$

在引入线性先验约束$Rc=r$ 后，所得到的约束修正SCAD估计仍保持例2中的稀疏识别结构，且关键滞后项系数之和满足先验约束从而使模型解释更符合实际背景。相较于标准SCAD方法的无约束估计，其结果一般不保证严格满足$Rc=r$ 因此可引入约束违背度$CV={\Vert \text{R}\widehat{\text{c}}-\text{r}\Vert }_2$进一步量化先验满足程度；本方法对应$CV\approx 0$ 在后续预测与置信区间计算中，仅需将例2中的参数估计值替换为满足约束的得到第22年的预测结果。

3. 结论

本文提出了一种面向UAR模型及其高阶UAR模型修正的SCAD惩罚估计方法，用于实现模型的客观定阶与参数估计。相较于仅依赖样本信息的最小二乘法与极大似然估计法，该方法通过融合先验信息并自适应选取调谐参数，不仅克服了传统方法忽视先验知识的局限，还实现了更高的估计精度。

两个数值算例表明，所提方法能有效获取UAR及高阶UAR模型的预测值与置信区间，且其估计误差率显著低于对比方法，从而验证了该方法的可行性与有效性。本研究为金融、环境等领域的不确定性数据预测建立了一个稳健的建模框架，也为合同节水管理等实际问题的预测提供了新思路。

需要指出，基于AIC与BIC准则的调谐参数选取方式，在处理UAR与高阶UAR模型时，难以协同优化模型的稀疏性与拟合优度。未来研究将把该方法拓展至UARMA模型及其参数估计问题中，以进一步提升其适用性。

基金项目

国家重点研发计划子课题(项目编号：2023YFC3209403-04-05)。