Notebook

移动平均(MA)模型¶

原理讲解¶

模型介绍¶

另一类时间序列模型为“移动平均过程”（Moving Average Process，简记 MA）。记一阶移动平均过程为 MA（1）：

$y_{t}=\mu+\varepsilon_{t}+\theta \varepsilon_{t-1}\tag{1}$

其中，{ $\varepsilon_{t}$ } 为白噪声，而 $\varepsilon_{t}$ 的系数被标准化为1。由于 $y_t$ 可以被看成是白噪声的移动平均，故得名。考虑使用条件 MLE 估计。为此，假设 { $\varepsilon_{t}$ } 为独立同分布且服从正态分布 $N\left(0, \sigma_{\varepsilon}^{2}\right)$ 。如果已经知道 { $\varepsilon_{t-1}$ } ，则

$y_{t} | \varepsilon_{t-1} \sim N\left(\mu+\theta \varepsilon_{t-1}, \sigma_{\varepsilon}^{2}\right)\tag{2}$

假设 $\varepsilon_0=0$ ，则可以知道 $\varepsilon_{1}=y_{1}-\mu$ 。给定 $\varepsilon_{1}$ ，则可以知道 $\varepsilon_{2}=y_{2}-\mu-\theta \varepsilon_{1}$ 。以此类推，则可以使用递推公式“ $\varepsilon_{t}=y_{t}-\mu-\theta \varepsilon_{t-1}$ ”计算出全部 $\left\{\varepsilon_{1}, \varepsilon_{2}, \cdots, \varepsilon_{T}\right\}$ 。这样，在给定 $\varepsilon_0=0$ 的条件下，就可以写下样本数据 $\left\{y_{1}, y_{2}, \cdots, y_{T}\right\}$ 的似然函数，然后使用数值方法求解其最大化问题。

更一般地，对于 $q$ 阶移动平均过程，记为 MA $(q)$ ：

$y_{t}=\mu+\varepsilon_{t}+\theta_{1} \varepsilon_{t-1}+\theta_{2} \varepsilon_{t-2}+\cdots+\theta_{q} \varepsilon_{t-q}\tag{3}$

也可以进行条件 MLE 估计，即在给定“ $\varepsilon_{0}=\varepsilon_{-1}=\cdots=\varepsilon_{-q+1}=0$ ”的条件下，最大化样本数据的似然函数。

MA(1) 模型最大似然函数推导¶

$y_t=μ+ε_t+θε_{t-1}$

$y_1=μ+ε_1+θε_0=μ+ε_1 ⇒ε_1=y_1-μ ⇒y_1∼N(μ,σ_ε^2)$

$y_2=μ+ε_2+θε_1=(1-θ)μ+θy_1+ε_2 ⇒ε_2= y_2-θy_1-(1-θ)μ ⇒y_2∼N((1-θ)μ+θy_1,σ_ε^2)$

以此类推： $y_3∼N((1-θ+θ^2)μ+θy_2-θ^2 y_1,σ_ε^2)$

因此：

$f\left(y_{1}\right)=\frac{1}{\sqrt{2 \pi \sigma_{\varepsilon}^{2}}} e^{-\frac{\left(y_{1}-\mu\right)^{2}}{2 \sigma_{\varepsilon}^{2}}}, f_{y_{2} | y_{1}}\left(y_{2} | y_{1}\right)=\frac{1}{\sqrt{2 \pi \sigma_{\varepsilon}^{2}}} e^{-\frac{\left(y_{2}-(1-\theta) \mu-\theta_{1}\right)^{2}}{2 \sigma_{\varepsilon}^{2}}}\tag{4}$

$\ln L=-\frac{1}{2} \ln 2 \pi-\frac{1}{2} \ln \sigma_{\varepsilon}^{2}-\frac{\left(y_{1}-\mu\right)^{2}}{2 \sigma_{\varepsilon}^{2}}-\frac{T-1}{2} \ln 2 \pi-\frac{T-1}{2} \ln \sigma_{\varepsilon}^{2}-\sum_{t=2}^{T} \frac{\left(y_{t}-c_{1} \mu-c_{2}\right)}{2 \sigma_{\varepsilon}^{2}}\tag{5}$

其中： $c_1=1-θ+θ^2-θ^3+......; c_2=θy_(t-1)-θ^2 y_(t-2)+θ^3 y_(t-3)+......$

statsmodels 库实现¶

详情请参考：ARMA 模型

In [1]:

from statsmodels.tsa.arima_model import ARMA
import pandas as pd

data = pd.read_excel('../数据/上证指数与沪深300.xlsx')
res = ARMA(data['sz'], order=(0,1)).fit()
res.summary()

Out[1]:

ARMA Model Results
Dep. Variable:	sz	No. Observations:	460
Model:	ARMA(0, 1)	Log Likelihood	-2897.310
Method:	css-mle	S.D. of innovations	131.267
Date:	Fri, 29 May 2020	AIC	5800.620
Time:	19:15:02	BIC	5813.013
Sample:	0	HQIC	5805.500

	coef	std err	z	P>\|z\|	[0.025	0.975]
const	2930.6519	11.858	247.141	0.000	2907.410	2953.894
ma.L1.sz	0.9396	0.013	73.789	0.000	0.915	0.965

Roots
	Real	Imaginary	Modulus	Frequency
MA.1	-1.0643	+0.0000j	1.0643	0.5000

matlab实现¶

MAlnL.m

function lnL = MAlnL(Beta,Y)
%Beta(1) = u ;Beta(2) = sita;Beta(3) = sigma2;
T = length(Y);
lnL1 = -0.5*log(2*pi)-0.5*log(Beta(3)^2)-(Y(1)-Beta(1))^2/(2*Beta(3)^2)-0.5*log(2*pi)*(T-1)-0.5*(T-1)*log(Beta(3)^2);
for i = 2:T
    for j = 1:i
        c11(j) = (-Beta(2))^(j-1);
    end
    c1 = sum(c11)*Beta(1);
    for j = 1:i-1
        c22(j) = Y(i-j)*(Beta(2)^j)*(-1)^(j+1);
    end
    c2 = sum(c22);
    lnL2(i) = (Y(i)- c1 - c2)^2;
end
lnL = lnL1-0.5*sum(lnL2)/Beta(3)^2;
lnL = -lnL;

% %导入数据
% data = xlsread('../数据/上证指数与沪深300.xlsx');
% f1 = data(:,1); f2 = data(:,2);
% 
% % 2.初始参数设定
% maxsize         = 2000;         % 生成均匀随机数的个数(用于赋初值)
% REP			    = 100;          % 若发散则继续进行迭代的次数
% nInitialVectors    = [maxsize, 3];    % 生成随机数向量
% MaxFunEvals    = 5000;         % 函数评价的最大次数
% MaxIter         = 5000;         % 允许迭代的最大次数
% options = optimset('LargeScale', 'off','HessUpdate', 'dfp','MaxFunEvals', ...
% MaxFunEvals, 'display', 'on', 'MaxIter', MaxIter, 'TolFun', 1e-6, 'TolX', 1e-6,'TolCon',10^-12);
% 
% % 3.寻找最优初值
% initialTargetVectors = unifrnd(0,10, nInitialVectors);
% RQfval = zeros(nInitialVectors(1), 1);
% for i = 1:nInitialVectors(1)
%     RQfval(i) = MAlnL (initialTargetVectors(i,:), f1);
% end
% Results          = [RQfval, initialTargetVectors];
% SortedResults    = sortrows(Results,1);
% BestInitialCond  = SortedResults(1,2: size(Results,2));    
% 
% % 4.迭代求出最优估计值Beta
% [Beta, fval exitflag] = fminsearch(' MAlnL ', BestInitialCond,options,f1);
% for it = 1:REP
% if exitflag == 1, break, end
% [Beta, fval exitflag] = fminsearch(' MAlnL ', BestInitialCond,options,f1);
% end
% if exitflag~=1, warning('警告：迭代并没有完成'), end

Eviews 估计¶

估计方程 y c ma(p)