Time Series Analysis and Examples |
Consider the univariate AR(
![p](images/timeseriesexpls_timeseriesexplseq13.gif)
) process
![y_t = \alpha_0 + \sum_{i=1}^p \alpha_i y_{t-i} + \epsilon_t](images/timeseriesexpls_timeseriesexplseq280.gif)
Define the design matrix
![{x}](images/timeseriesexpls_timeseriesexplseq281.gif)
.
![{x}= [1 & y_p & ... & y_1 \ \vdots & \vdots & \ddots & \vdots \ 1 & y_{t-1} & ... & y_{t-p} ]](images/timeseriesexpls_timeseriesexplseq282.gif)
Let
![{y}= (y_{p+1}, ... ,y_n)^'](images/timeseriesexpls_timeseriesexplseq283.gif)
.
The least squares estimate,
![\hat{{a}}=({x}^'{x})^{-1}{x}^'{y}](images/timeseriesexpls_timeseriesexplseq284.gif)
,
is the approximation to the maximum likelihood
estimate of
![{a}=(\alpha_0,\alpha_1, ... ,\alpha_p)](images/timeseriesexpls_timeseriesexplseq285.gif)
if
![\epsilon_t](images/timeseriesexpls_timeseriesexplseq104.gif)
is assumed to be Gaussian
error disturbances. Combining
![{x}](images/timeseriesexpls_timeseriesexplseq281.gif)
and
![{y}](images/timeseriesexpls_timeseriesexplseq8.gif)
as
![{z}= [{x}\,\vdots\,{y}]](images/timeseriesexpls_timeseriesexplseq286.gif)
the
![{z}](images/timeseriesexpls_timeseriesexplseq287.gif)
matrix can be decomposed as
![{z}= {q}{u}= {q}[{{r}}& {w}_1 \ 0 & {w}_2 ]](images/timeseriesexpls_timeseriesexplseq288.gif)
where
![{q}](images/timeseriesexpls_timeseriesexplseq289.gif)
is an orthogonal matrix and
![{{r}}](images/timeseriesexpls_timeseriesexplseq290.gif)
is an upper
triangular matrix,
![{w}_1 = (w_1, ... ,w_{p+1})^'](images/timeseriesexpls_timeseriesexplseq291.gif)
, and
![{w}_2 = (w_{p+2},0, ... ,0)^'](images/timeseriesexpls_timeseriesexplseq292.gif)
.
![{q}^'{y}= [w_1 \ w_2 \ \vdots \ w_{t-p} ]](images/timeseriesexpls_timeseriesexplseq293.gif)
The least squares estimate that uses Householder transformation
is computed by solving the linear system
![{{r}}{a}= {w}_1](images/timeseriesexpls_timeseriesexplseq294.gif)
The unbiased residual variance estimate is
![\hat{\sigma}^2 = \frac{1}{t-p} \sum_{i=p+2}^{t-p} w_i^2 = \frac{w_{p+2}^2}{t-p}](images/timeseriesexpls_timeseriesexplseq295.gif)
and
![{\rm aic}=(t-p)\log(\hat{\sigma}^2) + 2(p+1)](images/timeseriesexpls_timeseriesexplseq296.gif)
In practice, least squares estimation does not require the
orthogonal matrix
![{q}](images/timeseriesexpls_timeseriesexplseq289.gif)
. The TIMSAC subroutines compute
the upper triangular matrix without computing the matrix
![{q}](images/timeseriesexpls_timeseriesexplseq289.gif)
.
Consider the additive time series model
![y_t = t_t + s_t + \epsilon_t,\hspace*{.25in}\epsilon_t\sim n(0,\sigma^2)](images/timeseriesexpls_timeseriesexplseq297.gif)
Practically, it is not possible to estimate parameters
, since the number of
parameters exceeds the number of available observations.
Let
denote the seasonal difference operator with
seasons and degree of
; that is,
.
Suppose that
. Some constraints on the trend and seasonal
components need to be imposed such that the sum of squares of
,
, and
is
small. The constrained least squares estimates are obtained by
minimizing
![\sum_{t=1}^t \{(y_t-t_t-s_t)^2 + d^2[s^2(\nabla^k t_t)^2 + (\nabla_l^m s_t)^2 + z^2(s_t+ ... +s_{t-l+1})^2]\}](images/timeseriesexpls_timeseriesexplseq305.gif)
Using matrix notation,
![({y}-{m}{a})^'({y}-{m}{a}) + ({a}-{a}_0)^'{d}^'{d}({a}-{a}_0)](images/timeseriesexpls_timeseriesexplseq306.gif)
where
![{m}= [{i}_t\,\vdots\,{i}_t]](images/timeseriesexpls_timeseriesexplseq307.gif)
,
![{y}= (y_1, ... , y_t)^'](images/timeseriesexpls_timeseriesexplseq58.gif)
,
and
![{a}_0](images/timeseriesexpls_timeseriesexplseq308.gif)
is the initial guess of
![{a}](images/timeseriesexpls_timeseriesexplseq309.gif)
. The matrix
![{d}](images/timeseriesexpls_timeseriesexplseq197.gif)
is a
![3tx 2t](images/timeseriesexpls_timeseriesexplseq310.gif)
control matrix in which structure varies
according to the order of differencing in trend and season.
![{d}= d [{e}_m & 0 \ z{f}& 0 \ 0 & s{g}_k ]](images/timeseriesexpls_timeseriesexplseq311.gif)
where
![{e}_m & = & {c}_m\otimes{i}_l,\hspace*{.25in} m=1,2,3 \ {f}& = & [1 & 0 & ...... ...s & \ddots & \ddots & \ddots & 0 \ 0 & ... & 0 & -1 & 3 & -3 & 1 ]_{tx t}](images/timeseriesexpls_timeseriesexplseq312.gif)
The
![nx n](images/timeseriesexpls_timeseriesexplseq313.gif)
matrix
![{c}_m](images/timeseriesexpls_timeseriesexplseq314.gif)
has the same structure
as the matrix
![{g}_m](images/timeseriesexpls_timeseriesexplseq315.gif)
, and
![{i}_l](images/timeseriesexpls_timeseriesexplseq316.gif)
is the
![lx l](images/timeseriesexpls_timeseriesexplseq317.gif)
identity matrix.
The solution of the constrained least squares method is equivalent
to that of maximizing the function
![l({a}) = \exp\{-\frac{1}{2\sigma^2}({y}-{m}{a})^'({y}-{m}{a})\} \exp\{-\frac{1}{2\sigma^2}({a}-{a}_0)^'{d}^'{d}({a}-{a}_0)\}](images/timeseriesexpls_timeseriesexplseq318.gif)
Therefore, the PDF of the data
![{y}](images/timeseriesexpls_timeseriesexplseq8.gif)
is
![f({y}|\sigma^2,{a}) = (\frac{1}{2\pi})^{t/2} (\frac{1}{\sigma})^t \exp\{-\frac{1}{2\sigma^2}({y}-{m}{a})^'({y}-{m}{a})\}](images/timeseriesexpls_timeseriesexplseq319.gif)
The prior PDF of the parameter vector
![{a}](images/timeseriesexpls_timeseriesexplseq309.gif)
is
![\pi({a}|{d},\sigma^2,{a}_0) = (\frac{1}{2\pi})^t (\frac{1}{\sigma})^{2t}|{d}^'{d}| \exp\{-\frac{1}{2\sigma^2}({a}-{a}_0)^'{d}^'{d}({a}-{a}_0)\}](images/timeseriesexpls_timeseriesexplseq320.gif)
When the constant
![d](images/timeseriesexpls_timeseriesexplseq220.gif)
is known, the estimate
![\hat{{a}}](images/timeseriesexpls_timeseriesexplseq321.gif)
of
![{a}](images/timeseriesexpls_timeseriesexplseq309.gif)
is the mean of the posterior distribution, where
the posterior PDF of the parameter
![{a}](images/timeseriesexpls_timeseriesexplseq309.gif)
is proportional to
the function
![l({a})](images/timeseriesexpls_timeseriesexplseq322.gif)
.
It is obvious that
![\hat{{a}}](images/timeseriesexpls_timeseriesexplseq321.gif)
is the minimizer of
![\vert{g}({a}| d)\vert^2 = (\tilde{{y}}-\tilde{{d}}{a})^'(\tilde{{y}}-\tilde{{d}}{a})](images/timeseriesexpls_timeseriesexplseq323.gif)
,
where
![\tilde{{y}} = [{y}\ {d}{a}_0 ]](images/timeseriesexpls_timeseriesexplseq324.gif)
![\tilde{{d}} = [{m}\ {d} ]](images/timeseriesexpls_timeseriesexplseq325.gif)
The value of
![d](images/timeseriesexpls_timeseriesexplseq220.gif)
is determined by the minimum ABIC
procedure. The ABIC is defined as
![{\rm abic} = t\log[\frac{1}t\vert{g}({a}| d)\vert^2] + 2\{\log[\det({d}^'{d}+{m}^'{m})] - \log[\det({d}^'{d})]\}](images/timeseriesexpls_timeseriesexplseq326.gif)
In this section, the mathematical formulas for state space
modeling are introduced. The Kalman filter algorithms are
derived from the state space model. As an example, the state
space model of the TSDECOMP subroutine is formulated.
Define the following state space model:
![{x}_t &=& {f}{x}_{t-1} + {g}{w}_t \ y_t &=& {h}_t{x}_t + \epsilon_t](images/timeseriesexpls_timeseriesexplseq327.gif)
where
![\epsilon_t\sim n(0,\sigma^2)](images/timeseriesexpls_timeseriesexplseq328.gif)
and
![{w}_t\sim n(0,{q})](images/timeseriesexpls_timeseriesexplseq329.gif)
.
If the observations,
![(y_1, ... ,y_t)](images/timeseriesexpls_timeseriesexplseq330.gif)
, and the initial
conditions,
![{x}_{0|}](images/timeseriesexpls_timeseriesexplseq331.gif)
and
![{p}_{0|}](images/timeseriesexpls_timeseriesexplseq332.gif)
, are
available, the one-step predictor
![({x}_{t| t-1})](images/timeseriesexpls_timeseriesexplseq333.gif)
of the state vector
![{x}_t](images/timeseriesexpls_timeseriesexplseq334.gif)
and its mean square error (MSE)
matrix
![{p}_{t| t-1}](images/timeseriesexpls_timeseriesexplseq335.gif)
are
written as
![{x}_{t| t-1} = {f}{x}_{t-1| t-1}](images/timeseriesexpls_timeseriesexplseq336.gif)
![{p}_{t| t-1} = {f}{p}_{t-1| t-1}{f}^' + {g}{q}{g}^'](images/timeseriesexpls_timeseriesexplseq337.gif)
Using the current observation, the filtered value of
![{x}_t](images/timeseriesexpls_timeseriesexplseq334.gif)
and its variance
![{p}_{t| t}](images/timeseriesexpls_timeseriesexplseq338.gif)
are updated.
![{x}_{t| t} = {x}_{t| t-1} + {k}_t e_t](images/timeseriesexpls_timeseriesexplseq339.gif)
![{p}_{t| t} = ({i}- {k}_t {h}_t){p}_{t| t-1}](images/timeseriesexpls_timeseriesexplseq340.gif)
where
![e_t = y_t - {h}_t{x}_{t| t-1}](images/timeseriesexpls_timeseriesexplseq341.gif)
and
![{k}_t = {p}_{t| t-1}{h}^'_t[{h}_t {p}_{t| t-1}{h}^'_t + \sigma^2{i}]^{-1}](images/timeseriesexpls_timeseriesexplseq342.gif)
.
The log-likelihood function is computed as
![\ell = -\frac{1}2\sum_{t=1}^t \log(2\pi v_{t| t-1}) - \sum_{t=1}^t \frac{e_t^2}{2v_{t| t-1}}](images/timeseriesexpls_timeseriesexplseq343.gif)
where
![v_{t| t-1}](images/timeseriesexpls_timeseriesexplseq344.gif)
is the conditional variance of
the one-step prediction error
![e_t](images/timeseriesexpls_timeseriesexplseq345.gif)
.
Consider the additive time series decomposition
![y_t = t_t + s_t + t\!d_t + u_t + {x}^'_t\beta_t + \epsilon_t](images/timeseriesexpls_timeseriesexplseq346.gif)
where
![{x}_t](images/timeseriesexpls_timeseriesexplseq334.gif)
is a
![(kx 1)](images/timeseriesexpls_timeseriesexplseq347.gif)
regressor vector and
![\beta_t](images/timeseriesexpls_timeseriesexplseq348.gif)
is a
![(kx 1)](images/timeseriesexpls_timeseriesexplseq347.gif)
time-varying coefficient vector.
Each component has the following constraints:
![\nabla^k t_t & = & w_{1t},\hspace*{.25in}w_{1t}\sim n(0,\tau_1^2) \ \nabla_l^m... ...m_{i=1}^6 \gamma_{it}(t\!d_t(i)-t\!d_t(7)) \ \gamma_{it} & = & \gamma_{i,t-1}](images/timeseriesexpls_timeseriesexplseq349.gif)
where
![\nabla^k = (1-b)^k](images/timeseriesexpls_timeseriesexplseq350.gif)
and
![\nabla_l^m = (1-b^l)^m](images/timeseriesexpls_timeseriesexplseq300.gif)
.
The AR component
![u_t](images/timeseriesexpls_timeseriesexplseq21.gif)
is assumed to be stationary. The
trading-day component
![t\!d_t(i)](images/timeseriesexpls_timeseriesexplseq351.gif)
represents the number of the
![i](images/timeseriesexpls_timeseriesexplseq216.gif)
th day of the week in time
![t](images/timeseriesexpls_timeseriesexplseq205.gif)
.
If
![k=3, p=3, m=1](images/timeseriesexpls_timeseriesexplseq352.gif)
, and
![l=12](images/timeseriesexpls_timeseriesexplseq353.gif)
(monthly data),
![t_t & = & 3t_{t-1} - 3t_{t-2} + t_{t-3} + w_{1t} \ \sum_{i=0}^{11} s_{t-i} & = & w_{2t} \ u_t & = & \sum_{i=1}^3 \alpha_i u_{t-i} + w_{3t}](images/timeseriesexpls_timeseriesexplseq354.gif)
The state vector is defined as
![{x}_t = (t_t,t_{t-1},t_{t-2},s_t, ... ,s_{t-11},u_t,u_{t-1},u_{t-2}, \gamma_{1t}, ... ,\gamma_{6t})^'](images/timeseriesexpls_timeseriesexplseq355.gif)
The matrix
![{f}](images/timeseriesexpls_timeseriesexplseq356.gif)
is
![{f}= [{f}_1 & 0 & 0 & 0 \ 0 & {f}_2 & 0 & 0 \ 0 & 0 & {f}_3 & 0 \ 0 & 0 & 0 & {f}_4 ]](images/timeseriesexpls_timeseriesexplseq357.gif)
where
![{f}_1 = [3 & -3 & \phantom{-}1 \ 1 & 0 & 0 \ 0 & 1 & 0 ]](images/timeseriesexpls_timeseriesexplseq358.gif)
![{f}_2 = [-1^' & -1 \ {i}_{10} & 0 ]](images/timeseriesexpls_timeseriesexplseq359.gif)
![{f}_3 = [\alpha_1 & \alpha_2 & \alpha_3 \ 1 & 0 & 0 \ 0 & 1 & 0 ]](images/timeseriesexpls_timeseriesexplseq360.gif)
![{f}_4 = {i}_6](images/timeseriesexpls_timeseriesexplseq361.gif)
![1^' = (1,1, ... ,1)](images/timeseriesexpls_timeseriesexplseq362.gif)
The matrix
can be denoted as
![g = [{g}_1 & 0 & 0 \ 0 & {g}_2 & 0 \ 0 & 0 & {g}_3 \ 0 & 0 & 0 \ 0 & 0 & 0 \ 0 & 0 & 0 \ 0 & 0 & 0 \ 0 & 0 & 0 \ 0 & 0 & 0 ]](images/timeseriesexpls_timeseriesexplseq364.gif)
where
![{g}_1 = {g}_3 = [1 & 0 & 0 ]^'](images/timeseriesexpls_timeseriesexplseq365.gif)
![{g}_2 = [1 & 0 & 0 & 0 & 0 & 0 ]^'](images/timeseriesexpls_timeseriesexplseq366.gif)
Finally, the matrix
![{h}_t](images/timeseriesexpls_timeseriesexplseq367.gif)
is time-varying,
![{h}_t = [1 & 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 & {h}^'_t ]](images/timeseriesexpls_timeseriesexplseq368.gif)
where
![{h}_t & = & [d_t(1) & d_t(2) & d_t(3) & d_t(4) & d_t(5) & d_t(6) ]^' \ d_t(i) & = & t\!d_t(i)-t\!d_t(7),\hspace*{.5in}i=1, ... ,6](images/timeseriesexpls_timeseriesexplseq369.gif)