测量误差模型方差多变点的估计及收敛速度

The Estimator and Convergence Rate of Variance Change Point of Measurement Error Model

摘要: 当已知测量误差模型中误差的方差存在变点时，对方差变点用特征函数构造了一个含有调节参数的“CUSUM型估计量”，研究了方差变点估计量的弱(强)相合性以及收敛速度，并结合“二元分割法”推广至多个方差变点的估计。利用基于数据驱动的调节参数选取方法选取适合的调节参数，并应用含有调节参数的“CUSUM型估计量”对黄金价格的涨跌幅的方差进行实证分析，结果表明基于调节参数“CUSUM型估计量”得到的方差变点与实际相符，且估计量更稳健。

Abstract: When the variance of the measurement error model is known to have change points, a “CUSUM type estimator” with tuning parameters is constructed by using a characteristic function for each change point. The weak (strong) consistency and convergence rate of the variance change point estimation are studied, and the “binary segmentation method” is extended to estimate multiple variance change points. The data-driven tuning parameter selection method is used to select suitable tuning parameters, and using the “CUSUM type estimate” with tuning parameters to make an empirical analysis of the variance of the rise and fall of gold prices. The results show that the variance change points obtained based on the tuning parameter “CUSUM type estimator” are consistent with the reality, and the estimators are more robust.

1. 引言

随着大数据技术的迅猛发展，数据的收集越来越简捷，在数据收集过程往往会发生预测变量不能被直接测量，或者由于测量仪器的不精确、测量环境的不稳定以及测量人员的记录误差等各方面因素使得测量值与真实值之间有一定偏差的现象，我们将这个偏差称为测量误差。测量误差数据普遍存在各个领域中且更贴近现实研究情况，因此引起了很多学者的广泛关注，Fuller (1987)，Cheng和Ness (1999)以及Carroll等(2002)分别对线性测量误差模型 [1] 、测量误差的回归模型 [2] 以及非线性测量误差模型 [3] 做了深入研究。

考虑一般的测量误差模型：

$Y_j=X_j+{\eta }_j,j=1,\cdots ,n$ (1)

其中： $X_j$ 是独立同分布的不可直接观察的真实变量； $Y_j$ 是变量 $X_j$ 在测量误差 ${\eta }_j$ 下的观测值；测量误差 ${\eta }_j$ 是独立同分布的，服从均值0，方差为 ${\sigma }^2$ 的正态分布，且 ${\eta }_j$ 和 $X_j$ 相互独立， $j=1,\cdots ,n$ 。模型(1)是最简单的线性测量误差模型，也是应用最广的一类测量误差模型。

在实际应用中，模型(1)中变量的分布可能在某个未知时刻发生变化，此时称模型(1)为含有变点的测量误差模型，若模型(1)中变量的分布在多个未知时刻发生变化，则称模型(1)为含多个变点的测量误差模型。Chang和Huang (1997)应用了似然比方法对测量误差模型中变量 $X_j$ 的均值是否存在变点进行了检验，并给出了均值变点的极大似然估计 [4] 。王黎明(2002)应用了似然比方法对测量误差回归模型中回归系数是否存在变点进行了检验，并给出了回归系数变点的极大似然估计 [5] 。You (2006)对测量误差模型中回归系数的变点进行了研究 [6] 。而关于测量误差模型中方差变点的研究很少，Dong (2013，2016)分别用特征函数法 [7] 和非参数的U统计量法 [8] 研究了测量误差模型(1)中的方差变点的估计量及其收敛速度。

自Page (1954)首次提出累积和(Cumulative Sum，简记CUSUM)方法对变点进行连续性检验后 [9] ，CUSUM方法被许多统计学家改进并应用于变点的检测与估计。Horváth和Kokoszka (1997)使用CUSUM方法对单变量序列中的变点进行估计时，在数值模拟过程中发现调节参数对变点估计的精确度有显著影响 [10] 。Fremdt (2013)对单变量序列是否存在变点提出了改进的CUSUM型检验统计量，数据模拟过程发现当调节参数 $\gamma \in \left(0,0.5\right)$ 时变点检验效果更稳健 [11] 。谭常春和江敏(2020)利用蒙特卡洛方法研究了CUSUM型统计量中调节参数对单变量序列均值变点估计效果的影响，研究发现当变点变化幅度较大时，无论变点真实位置如何，变点估计值基本不受调节参数取值的影响；但当变点变化幅度较小时，调节参数的不同取值对不同变点位置的估计结果有显著影响，并由此提出了基于数据驱动的调节参数的选取方法，使得CUSUM型变点估计量更稳健 [12] 。本文结合Dong (2013)的基础上，在已知测量误差模型(1)中测量误差的方差存在变点时，关于方差变点提出了基于调节参数“CUSUM型估计量”，并研究了变点估计量的强弱相合性以及收敛速度，进一步结合“二元分割法”推广至多个方差变点的估计。最后，对2004年1月~2023年9月黄金价格的涨跌幅周数据进行实证分析，应用谭常春和江敏(2020)数据驱动的调节参数选取方法选取了适合的调节参数，并应用含有调节参数的“CUSUM型估计量”得到黄金价格涨跌幅方差变点的估计量。

2. 模型与主要结果

假设模型(1)中测量误差 ${\eta }_j,j=1,\cdots ,n$ 服从均值为0的正态分布。且 $\left\{{\eta }_j\right\}$ 存在一个变点，即：

${\eta }_j~\{\begin{array}{l}N\left(0,{\sigma }_1^2\right),j=1,\cdots ,k^{*},\\ N\left(0,{\sigma }_2^2\right),j=k^{*}+1,\cdots ,n.\end{array}$ (2)

其中： ${\sigma }_1^2\ne {\sigma }_2^2$ 和 $k^{*}$ 均未知。不失一般性，不妨假设 $\exists \delta ,0<\delta <\frac{1}{2}$ ，使得 $n\delta \le k^{*}\le n\left(1-\delta \right)$ ，即变点不在端点处。为了方便起见， $\hat{k}$ 和 $\hat{\tau }=\frac{\hat{k}}{n}$ 均为变点的估计量。记 $z_j={\text{e}}^{it{Y}_j}$ ，其中 $i=\sqrt{-1}$ 为复数单位；复数z的

模记为 $\left|z\right|$ ；随机变量V的特征函数为 $f_V^{Ft}\left(t\right)=E{\text{e}}^{itv}$ ； $N\left(0,{\sigma }_j^2\right)$ 的特征函数为 $f_j^{Ft}\left(t\right)$ ， $j=1,2$ ； $a=f_X^{Ft}\left(t\right)f_1^{Ft}\left(t\right)$ ， $b=f_X^{Ft}\left(t\right)f_2^{Ft}\left(t\right)$ ， $\left|t\right|\le 1$ 。

记：

$U_k\left(t\right)={\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }\left[\frac{1}{k}\underset{j=1}{\overset{k}{{\displaystyle \sum }}}{\text{e}}^{it{Y}_j}-\frac{1}{n-k}\underset{j=k+1}{\overset{n}{{\displaystyle \sum }}}{\text{e}}^{it{Y}_j}\right]$ (3)

其中： $\gamma$ 为调节参数， $\gamma \in \left[0,1\right)$ ， $\gamma =0$ 表示无调节参数。

方差变点的估计量定义为：

$\hat{\tau }=\frac{\hat{k}}{n}=\frac{1}{n}\arg \underset{\begin{array}{c}\delta n\le k\le \left(1-\delta \right)n\\ 0\le t\le 1\end{array}}{\max }\left|U_k\left(t\right)\right|$ (4)

定理1 设 $Y_j,j=1,\cdots ,n$ 是模型(1)的独立的观测值，且该模型中的测量误差 $\left\{{\eta }_j\right\}$ 的方差在未知时刻 $k^{*}$ 发生改变，则：

(i) 由(4)定义的方差变点的估计量 $\hat{\tau }$ 是 ${\tau }^{*}$ 的弱相合估计，且 $\underset{n\to \infty }{\lim }\left|\hat{\tau }-{\tau }^{*}\right|g_1\left(n\right)=0$ 依概率成立，其中： $g_1\left(n\right)=n^{\frac{1}{2}}{\mathcal{l}}^{-1}\left(n\right)$ 。

(ii) $\hat{\tau }$ 是 ${\tau }^{*}$ 的强相合估计，且 $\underset{n\to \infty }{\lim }\left|\hat{\tau }-{\tau }^{*}\right|g_2\left(n\right)=0$ a.s. 其中： $g_2\left(n\right)=n^{\frac{1}{2}}{\left(\log n\right)}^{-\frac{1}{2}}{\mathcal{l}}^{-1}\left(n\right)$ 。

其中： $\mathcal{l}\left(n\right)$ 是满足 $\underset{n\to \infty }{\lim }\mathcal{l}\left(n\right)=\infty$ 的慢变函数， $g_1\left(n\right)$ 与 $g_2\left(n\right)$ 是满足 $\underset{n\to \infty }{\lim }{g}_j\left(n\right)=\infty$ 的非负函数。

要证明这一点需要一下引理：

引理1 [13] (Hoeffding不等式)假设 $0\le X_j\le 1$ ， ${\mu }_j=E\left(X_j\right)$ ， $\bar{X}=\frac{1}{n}{\displaystyle {\sum }_{j=1}^n{X}_j}$ ， $\mu =E\left(\bar{X}\right)$ ， $S_j={\displaystyle {\sum }_{j=1}^j{X}_i}$ ，则对于 $0<x<1-\mu$ ，

$P\left(\bar{X}-\mu \ge x\right)\le P\left(\underset{1\le j\le n}{\max }\left(S_j-E{S}_j\right)\ge nx\right)\le \exp \left\{-2n{x}^2\right\}$ ，

如果存在着 $a_j\le b_j$ ，使得 $a_j\le X_j\le b_j$ ， $j=1,2,\cdots ,n$ ，那么对任意的 $x>0$ ，

$P\left(\bar{X}-\mu \ge x\right)\le P\left(\underset{1\le j\le n}{\max }\left(S_j-E\left(S_j\right)\right)\ge nx\right)\le \exp \left\{\frac{-2n{x}^2}{{\displaystyle {\sum }_{j=1}^n{\left(b_j-a_j\right)}^2}}\right\}.$

下面给出定理1的证明：

(i) 由(3)计算容易得到：

$E{z}_j=E{\text{e}}^{it{Y}_j}=\{\begin{array}{l}{f}_X^{Ft}\left(t\right)f_1^{Ft}\left(t\right),j=1,\cdots ,k^{*};\\ f_X^{Ft}\left(t\right)f_2^{Ft}\left(t\right),j=k^{*}+1,\cdots ,n,\end{array}$

经简单计算，可知：

$E{U}_k\left(t\right)=\{\begin{array}{l}{\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }E\left(\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}\right),1\le k\le k^{*};\\ {\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }E\left(\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}\right),k^{*}<k\le n,\end{array}$

$=\{\begin{array}{l}{\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }E\left(\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-\frac{1}{n-k}\left({\displaystyle \underset{j=k+1}{\overset{k^{*}}{\sum }}{\text{e}}^{it{Y}_j}}+{\displaystyle \underset{j=k^{*}+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}\right)\right),1\le k\le k^{*};\\ {\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }E\left(\frac{1}{k}\left({\displaystyle \underset{j=1}{\overset{k^{*}}{\sum }}{\text{e}}^{it{Y}_j}}+{\displaystyle \underset{j=k^{*}+1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}\right)-\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}\right),k^{*}<k\le n,\end{array}$

$=\{\begin{array}{l}{\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }{f}_X^{Ft}\left(t\right)\left[f_1^{Ft}\left(t\right)-\frac{k^{*}-k}{n-k}{f}_1^{Ft}\left(t\right)-\frac{n-k^{*}}{n-k}{f}_2^{Ft}\left(t\right)\right],1\le k\le k^{*};\\ {\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }{f}_X^{Ft}\left(t\right)\left[\frac{k^{*}}{k}{f}_1^{Ft}\left(t\right)-\frac{k-k^{*}}{n-k}{f}_2^{Ft}\left(t\right)-f_2^{Ft}\left(t\right)\right],k^{*}<k\le n,\end{array}$

$=\{\begin{array}{l}{\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }{f}_X^{Ft}\left(t\right)\frac{n-k^{*}}{n-k}\left[f_1^{Ft}\left(t\right)-f_2^{Ft}\left(t\right)\right],1\le k\le k^{*};\\ {\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }{f}_X^{Ft}\left(t\right)\frac{k^{*}}{k}\left[f_1^{Ft}\left(t\right)-f_2^{Ft}\left(t\right)\right],k^{*}<k\le n,\end{array}$

$=\{\begin{array}{l}{\tau }^{1-\gamma }{\left(1-\tau \right)}^{-\gamma }\left(1-{\tau }^{*}\right)\left(a-b\right),1\le k\le k^{*};\\ {\tau }^{-\gamma }{\left(1-\tau \right)}^{1-\gamma }{\tau }^{*}\left(a-b\right),k^{*}<k\le n.\end{array}$ (5)

注意到，

$E{U}_{k^{*}}\left(t\right)={\left({\tau }^{*}\right)}^{1-\gamma }{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(a-b\right).$ (6)

所以，

$\begin{array}{l}\left|E{U}_{k^{*}}\right|-\left|E{U}_k\right|\\ =\{\begin{array}{l}{\left({\tau }^{*}\right)}^{1-\gamma }{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(a-b\right){\tau }^{1-\gamma }-{\left(1-\tau \right)}^{-\gamma }\left(1-{\tau }^{*}\right)\left(a-b\right),1\le k\le k^{*};\\ {\left({\tau }^{*}\right)}^{1-\gamma }{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(a-b\right){\tau }^{1-\gamma }-{\tau }^{-\gamma }{\left(1-\tau \right)}^{1-\gamma }{\tau }^{*}\left(a-b\right),k^{*}<k\le n,\end{array}\\ =\{\begin{array}{l}{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(a-b\right)\left[{\left({\tau }^{*}\right)}^{1-\gamma }-{\tau }^{1-\gamma }{\left(\frac{1-{\tau }^{*}}{1-\tau }\right)}^{\gamma }\right],1\le k\le k^{*};\\ {\left({\tau }^{*}\right)}^{1-\gamma }\left(a-b\right)\left[{\left(1-{\tau }^{*}\right)}^{1-\gamma }-{\left(1-\tau \right)}^{1-\gamma }{\left(\frac{{\tau }^{*}}{\tau }\right)}^{\gamma }\right],k^{*}<k\le n,\end{array}\\ \ge \{\begin{array}{l}{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(a-b\right)\left[{\left({\tau }^{*}\right)}^{1-\gamma }-{\tau }^{1-\gamma }\right],1\le k\le k^{*};\\ {\left({\tau }^{*}\right)}^{1-\gamma }\left(a-b\right)\left[{\left(1-{\tau }^{*}\right)}^{1-\gamma }-{\left(1-\tau \right)}^{1-\gamma }\right],k^{*}<k\le n.\end{array}\end{array}$

若令 $f_1\left(\tau \right)={\tau }^{1-\gamma },\left(0<\gamma <1\right)$ ， $f_2\left(\tau \right)={\left(1-\tau \right)}^{1-\gamma },\left(0<\gamma <1\right)$ 由拉格朗日中值定理知，

$\left|E{U}_{k^{*}}\right|-\left|E{U}_k\right|\ge \{\begin{array}{l}{\left(1-{\tau }^{*}\right)}^{1-\gamma }\left(1-\gamma \right){\left({\tau }^{*}\right)}^{-\gamma }\left({\tau }^{*}-\tau \right)\left(a-b\right),1\le k\le k^{*};\\ {\left({\tau }^{*}\right)}^{1-\gamma }\left(1-\gamma \right){\left(1-{\tau }^{*}\right)}^{-\gamma }\left({\tau }^{*}-\tau \right)\left(a-b\right),k^{*}<k\le n.\end{array}$

记： $\bar{\tau }=\left(1-\gamma \right){\left[{\tau }^{*}\left(1-{\tau }^{*}\right)\right]}^{-\gamma }$ ， $\beta =\min \left\{{\tau }^{*},1-{\tau }^{*}\right\}$ ， $\rho =\left|b-a\right|$ 。

则 $\left|E{U}_{k^{*}}\right|-\left|E{U}_k\right|\ge \bar{\tau }\left|\tau -{\tau }^{*}\right|\beta \rho$ ，由三角不等式，

$\left|E{U}_{k^{*}}\right|-\left|E{U}_k\right|\le \left|U_k-E{U}_k\right|+\left|U_{k^{*}}-E{U}_{k^{*}}\right|+\left|U_{k^{*}}\right|-\left|U_k\right|.$

由 $\hat{k}$ 的定义知 $\left|U_{\hat{k}}\right|\ge \left|U_{k^{*}}\right|$ ，由(5)和(6)可得：

$\begin{array}{l}\bar{\tau }\left|\hat{\tau }-{\tau }^{*}\right|\beta \rho \le \left|E{U}_{k^{*}}\right|-\left|E{U}_{\hat{k}}\right|\le 2\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k-E{U}_k\right|\\ \begin{array}{c}\le 2\underset{k}{\max }\left\{\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k\left(t\right)-E{U}_k\left(t\right)\right|,\underset{\begin{array}{c}{k}^{*}\le k\le \left(1-\delta \right)n\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k\left(t\right)-EU{\left(t\right)}_k\right|\right\}.\end{array}\end{array}$ (7)

不失一般性，不妨设 $\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k\left(t\right)-E{U}_k\left(t\right)\right|$ ，当最大值为 $\underset{\begin{array}{c}{k}^{*}\le k\le \left(1-\delta \right)n\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k\left(t\right)-EU{\left(t\right)}_k\right|$ 的情形类似可证。由(7)

$\begin{array}{l}\bar{\tau }\left|\hat{\tau }-{\tau }^{*}\right|\beta \rho \le 2\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k-E{U}_k\right|\\ =2\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }{\left(\frac{k\left(n-k\right)}{n^2}\right)}^{1-\gamma }\left|\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}-\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}\right|\\ \le 2\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}\right|+2\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}\right|,\end{array}$

则有：

$P\left(\left|\hat{\tau }-{\tau }^{*}\right|>\epsilon \right)\le P\left(\frac{2}{\bar{\tau }\beta \rho }\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|U_k-E{U}_k\right|>\epsilon \right).$

所以，

$\begin{array}{c}P\left(\left|\hat{\tau }-{\tau }^{*}\right|>\epsilon \right)\le P\left(\frac{2}{\bar{\tau }\beta \rho }\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{k}{\displaystyle \underset{j=1}{\overset{k}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}\right|>\frac{\epsilon }{2}\right)\\ +P\left(\frac{2}{\bar{\tau }\beta \rho }\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{n-k}{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}\left({\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right)}\right|>\frac{\epsilon }{2}\right)\\ \triangleq \text{I}+\text{II}\text{.}\end{array}$ (8)

注意到： $\left|{\text{e}}^{it{Y}_j}-E{\text{e}}^{it{Y}_j}\right|\le 2$ ，且此刻 $\delta n\le k\le k^{*}$ ， $\delta n\le n-k^{*}<n-k$ ，由Hoeffding不等式，(8)的第一项

$\begin{array}{c}\text{I}=P\left(\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{k}\left[{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-E{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}\right]\right|>\frac{k\bar{\tau }\beta \rho \epsilon }{4}\right)\\ \le 2\exp \left\{-2n\delta {\left(\frac{\bar{\tau }\beta \rho \epsilon }{4}\right)}^2\right\}\to 0n\to \infty .\end{array}$ (9)

第二项：

$\begin{array}{c}\text{II}=P\left(\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|\frac{1}{n-k}\left[{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-E{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}\right]\right|>\frac{k\bar{\tau }\beta \rho \epsilon }{4}\right)\\ \le 2\exp \left\{-2n\delta {\left(\frac{\bar{\tau }\beta \rho \epsilon }{4}\right)}^2\right\}\to 0n\to \infty .\end{array}$ (10)

同理，当 $n\to \infty$ 时，(8)的第二项也 $\to 0$ ，这意味着变点估计 $\hat{\tau }$ 是 ${\tau }^{*}$ 的弱相合估计。

由(8)、(9)和(10)进一步可得：

$\begin{array}{c}P\left(g_1\left(n\right)\left|\hat{\tau }-{\tau }^{*}\right|>\epsilon \right)\le P\left(\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}-E{\displaystyle \underset{j=1}{\overset{k}{\sum }}{\text{e}}^{it{Y}_j}}\right|>\frac{\bar{\tau }\beta \rho \epsilon }{4g_1\left(n\right)}\right)\\ +P\left(\underset{\begin{array}{c}\delta n\le k\le k^{*}\\ \left|t\right|\le 1\end{array}}{\max }\left|{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}-E{\displaystyle \underset{j=k+1}{\overset{n}{\sum }}{\text{e}}^{it{Y}_j}}\right|>\frac{k\bar{\tau }\beta \rho \epsilon }{4g_1\left(n\right)}\right)\\ \le 4\exp \left\{-2n\delta {\left(\frac{\bar{\tau }\beta \rho \epsilon }{4g_1\left(n\right)}\right)}^2\right\}.\end{array}$ (11)

若 $g_1\left(n\right)=n^{\frac{1}{2}}{\mathcal{l}}^{-1}\left(n\right)$ 时，则有： $\underset{x\to \infty }{\lim }\left|\hat{\tau }-{\tau }^{*}\right|g_1\left(n\right)=0$ 依概率成立，其中： $\mathcal{l}\left(n\right)$ 是一个趋于∞的慢变函数。

对于 $k^{*}\le k\le \left(1-\delta \right)n$ ，同理可证。

(i) ：由(11)知：

${\displaystyle {\sum }_{j=1}^{n}P}\left(g_2\left(n\right)\left|\hat{\tau }-{\tau }^{*}\right|>\epsilon \right)\le 4{\displaystyle {\sum }_{j=1}^n\exp }\left\{-2n\delta {\left(\frac{\bar{\tau }\beta \rho \epsilon }{4g_2\left(n\right)}\right)}^2\right\}$ ，

若取 $g_2\left(n\right)=n^{\frac{1}{2}}{\left(\log n\right)}^{-\frac{1}{2}}{\mathcal{l}}^{-1}\left(n\right)$ ，则有 ${\displaystyle {\sum }_{j=1}^n\exp \left\{-2n\delta {\left(\frac{\bar{\tau }\beta \rho \epsilon }{4g_2\left(n\right)}\right)}^2\right\}}<\infty$ ，由Borel-Cantelli引理知 $\hat{\tau }$ 几乎处处收敛到 ${\tau }^{*}$ 。且 $\underset{x\to \infty }{\lim }\left|\hat{\tau }-{\tau }^{*}\right|g_2\left(n\right)=0$ a.s.成立。对于 $k^{*}\le k\le \left(1-\delta \right)n$ ，同理可证。

3. 测量误差模型中方差多变点的估计

假设模型(1)中测量误差 ${\eta }_j,j=1,\cdots ,n$ 服从均值为0的正态分布。且 $\left\{{\eta }_j\right\}$ 存在m个变点(m已知)，即

${\eta }_j=\{\begin{array}{l}N\left(0,{\sigma }_1^2\right),j=1,2,\cdots ,k_1;\\ N\left(0,{\sigma }_2^2\right),j=k_1+1,\cdots ,k_2;\\ ⋮\\ N\left(0,{\sigma }_m^2\right),j=k_{m-1}+1,\cdots ,k_m;\\ N\left(0,{\sigma }_{m+1}^2\right),j=k_m+1,\cdots ,n.\end{array}$ (12)

“二元分割法”最早是由Sen和Srivastava (1975)在检测均值变点时提出的 [14] 。Bai (1997)也对二元分割的理论进行了介绍 [15] ，与其他检测方法相比，二元分割法由于方法简单，易操作且准确率高等优势被广泛地应用于多变点的检测与估计问题中。应用二元分割法对(12)中的多个方差变点估计具体步骤如下：

步骤一：对全体数据利用(4)式估计出第一个变点 ${\hat{k}}_1$ ；

步骤二：变点 ${\hat{k}}_1$ 将数据分成两个子序列，第一部分为 $Y_1,Y_2,\cdots ,Y_{{\hat{k}}_1}$ ，第二部分为 $Y_{{\hat{k}}_1+1},Y_{{\hat{k}}_1+2},\cdots ,Y_n$ ，在这两段样本内分别计算 ${\max }_{1\le k\le {\hat{k}}_1}\left|U_k\left(t\right)\right|$ 和 ${\max }_{\hat{k}\le k\le n}\left|U_k\left(t\right)\right|$ ；

步骤三：比较 ${\max }_{1\le k\le {\hat{k}}_1}\left|U_k\left(t\right)\right|$ 和 ${\max }_{\hat{k}\le k\le n}\left|U_k\left(t\right)\right|$ 的大小，选择两者中较大的那个，并选择较大的所对应的k值，记为 ${\hat{k}}_2$ ；

步骤四：重复步骤2和步骤3，直到m个变点被估计出来。

4. 基于数据驱动调节参数的选择

本文采用谭常春和江敏(2020)提出的基于数据的调节参数选取算法：

第一步：设 $Y_1,Y_2,\cdots ,Y_{k_0},Y_{k_0+1},\cdots ,Y_n$ 是带有变点 $k_0$ 的样本序列，记调节参数 $\gamma$ 可能取值构成的集合为 $\Gamma =\left\{{\gamma }_1,{\gamma }_2,\cdots ,{\gamma }_l\right\}$ 为集。对给定 ${\gamma }_i$ ，可得到一个变点估计量，记为 ${\hat{k}}_{{\gamma }_i}$ ，有

${\hat{k}}_{{\gamma }_i}=\arg \underset{1<k<n}{\max }\left|U_k\right|,1\le i\le l$

第二步：计算 ${\hat{k}}_{{\gamma }_1},{\hat{k}}_{{\gamma }_2},\cdots ,{\hat{k}}_{{\gamma }_l}$ 的均值k：

$\overline{\stackrel{^}{k}}=\frac{1}{l}{\displaystyle \underset{i=1}{\overset{l}{\sum }}{\hat{k}}_{{\gamma }_i}}$ (13)

第三步：定义

$r_s=\arg \underset{r_i}{\min }\left|{\hat{k}}_{{\gamma }_i}-\overline{\stackrel{^}{k}}\right|$ (14)

即为基于数据驱动方法选取的合适的调节参数。

5. 实证分析

在统计学中，方差常用于刻画风险，因此变点的检测及估计不仅应用于质量控制中也常应用于金融学、经济学以及风险管理等领域。在金融市场中，投资者不仅关注于金融产品收益率的均值变化，也关心所投资金融产品的风险变化(即方差的波动)。

本文选取2004年1月~2023年9月黄金价格的涨跌幅周数据(数据来源于https://cn.investing.com/)作为研究对象，结合“二元分割法”并利用基于数据驱动的调节参数选取法选取适合的调节参数，进一步应用含有调节参数的“CUSUM型估计量”对黄金价格涨跌幅的方差变点进行估计，样本量 $n=1042$ 。首先应用(4)式含有调节参数的“CUSUM型估计量”对这1042个黄金价格涨跌幅数据进行方差变点估计，取调节参数 $\gamma$ 可能取值的集合 $\Gamma =\left\{0.1,0.2,0.25,0.3,0.5,0.75,0.9\right\}$ ，得到变点的估计值见表1所示。此时变点估计值不同，应用基于数据驱动的调节参数的选取方法中的第二步，得到 $\overline{\stackrel{^}{k}}=613.16$ ，结合第三步，发现适应于这一过程的调节参数 ${\gamma }_s$ 取0.1，0.2，0.25，0.3，0.5中任何一个都是合适的，得到的变点估计值都是401，所对应的实际时间是2011年8月28日，造成这种现象的原因主要是2011年的欧债危机增加了欧元国家及世界经济的不确定性，股指大幅度下跌。2011年8月6日国际评级机构标准普尔宣布将美国主权信用评级由AAA调降到AA+，这次信用评级的下调引起全球性的恐慌，也是金价大幅暴涨及波动的直接推手。

Table 1. The estimated value of the first variance change point estimated by the CUSUM estimator under different adjustment parameters γ and the corresponding actual time

表1. 不同的调节参数 $\gamma$ 下CUSUM型估计量估计的第一个方差变点的估计值及对应的实际时间

应用“二元分割法”，变点401 (即2011年8月28日)将2004年1月~2023年9月黄金价格涨跌幅的1042个周数据一分为二，得到2004年1月~2011年8月和2011年9月~2023年9月两个子序列，应用(4)式含有调节参数的“CUSUM型估计量”分别对上述两个子序列进行方差变点估计，样本量分别为401和640。调节参数 $\gamma$ 可能的取值集为适合上一过程中的调节参数组成，即 $\Gamma =\left\{0.1,0.2,0.25,0.3,0.5\right\}$ ，得到上述两个子序列的方差变点估计值见表2和表3所示。

Table 2. The second variance change point estimate of the CUSUM-type estimator under different adjustment parameters γ and the corresponding actual time

表2. 不同的调节参数 $\gamma$ 下CUSUM型估计量估计的第二个方差变点估计值及对应的实际时间

对于表2，此时变点估计值不同，应用基于数据驱动的调节参数的选取方法中的第二步，得到 $\overline{\stackrel{^}{k}}=280$ ，结合第三步，发现适应于这一过程的调节参数 ${\gamma }_s$ 取0.1，0.2，0.25，0.3中任何一个都是合适的，变点估计值都是252，所对应的实际时间为2008年10月26日，出现变点的原因主要是2008年8月~10月金融危机全面爆发，雷曼兄弟申请破产保护，次贷危机迅速演化为金融危机，黄金遭遇抛售，金价出现波动。

同理，表3中变点估计值也不全相同，应用基于数据驱动的调节参数的选取方法中的第二步，得到 $\overline{\stackrel{^}{k}}=583.6$ ，结合第三步，发现适应于这一过程的调节参数 ${\gamma }_s$ 取0.1，0.2，0.25，0.3中任何一个都是合适的，变点估计值都是628，对应的实际时间为2015年11月29日，出现变点的原因主要是2015年9月美联储加息落空，货币贬值，投资性实物黄金销售回暖，助推金价上行及波动。

Table 3. The third variance change point estimate of the CUSUM-type estimator under different adjustment parameters γ and the corresponding actual time

表3. 不同的调节参数 $\gamma$ 下CUSUM型估计量估计的第三个方差变点估计值及对应的实际时间

综合分析，适用于上述全过程的调节参数 ${\gamma }_s$ 取0.1，0.2，0.25，0.3中任一个都是合适的，即为最优调节参数，并且由最优调节参数得到的2004年1月~2023年9月黄金价格的涨跌幅的方差变点(如图1中红线所示)更符合实际情况。

Figure 1. Weekly data on gold price increases and decreases from January 2004 to September 2023

图1. 2004年1月~2023年9月黄金价格涨跌幅周数据

6. 结论

本文对测量误差模型中的方差变点用特征函数构造了一个含有调节参数的“CUSUM型估计量”，并研究了方差变点估计量的相合性及收敛速度，结合“二元分割法”推广至多个方差变点的情况。利用基于数据驱动的调节参数选取法选取适合的调节参数，并进一步应用含有调节参数的“CUSUM型估计量”对黄金价格的涨跌幅的方差进行实证分析，结果表明基于调节参数“CUSUM型估计量”得到的方差变点与实际相符且更稳健。该理论丰富了方差变点估计的研究方法，并可以推广至其它线性模型的变点估计问题，使得变点估计量更稳健。

基金项目

新疆维吾尔自治区自然科学基金项目(2023D01A37，2022D01A219)。

NOTES

^*通讯作者。