# Parameter Estimation for p-Order Random Coefficient Autoregressive (RCA) Models Based on Kalman Filter

Authors: Mohammed Benmoumen; Jelloul Allal; Imane Salhi
Journal: Journal of Applied Mathematics
 (2019)

---

## Abstract

In this paper we elaborate an algorithm to estimate p-order Random Coefficient Autoregressive Model (RCA(p)) parameters. This algorithm combines quasi-maximum likelihood method, the Kalman filter, and the simulated annealing method. In the aim to generalize the results found for RCA(1), we have integrated a subalgorithm which calculate the theoretical autocorrelation. Simulation results demonstrate that the algorithm is viable and promising.


---

## Body

## 1. Introduction

Random Coefficient Autoregressive(
R
C
A
) processes have been widely studied in the literature for modeling time series exhibiting nonlinear behavior. The R
C
A process was introduced by Andel (1976) who also studied its properties. He had obtained conditions for the existence of a singly infinite process which is second-order stationary satisfying:(1)
X
t
=
∑
k
=
1
p
ϕ
k
+
β
k
t
X
t
-
k
+
ε
t where ϕ
k’s are fixed coefficients, (
β
k
(
t
)
,
t
=
1
,
…
,
p)is a sequence of i.i.d. random vectors with mean zero and constant covariance matrix C
=
(
σ
i
,
j
), ε
(
t
) is a white noise process with mean zero and variance σ
ϵ, and β
k
(
t
) and ε
(
t
) are assumed to be mutually independent.Multiple studies have emerged after Conlisk (1974, 1976) has derived conditions for the stability of RCA models, Robinson (1978) has considered statistical inference for the RCA model, and Nicholls and Quinn [1] have extended the results of Andel to the multivariate RCA model.Many authors have investigated the estimation of parameters in the random coefficient autoregressive. Among them, we cite Nicholls and Quinn who obtained the least squares and the maximum likelihood estimates. Under certain assumptions, they established the strong consistency as well as the asymptotic normality for both estimates. For a detailed early study, we refer to Nicholls and Quinn [1]. On their side, Thavaneswaran and Abraham [2] apply Godambe’s theorem (1985) to obtain optimal estimates for R
C
A models.Recently several authors are interested in theR
C
A model. Aue and Horvath (2011) propose a unified quasi-likelihood procedure for the estimation of the unknown parameters of RCA(1) models that works for both stationary and nonstationary processes. They also establish the weak consistency and the asymptotic normality for this procedure.Liang et al. (2012) have described moment properties forR
C
A processes and the corresponding squared processes; they have given also a joint prediction study of the mean and volatility.A new algorithm was proposed by Allal and Benmoumen [3] to estimate first-order RCA’s parameters. This algorithm combines quasi-maximum likelihood method, the Kalman filter, and the Powell’s method. Our contribution aims to extend previous method to higher-order of RCA process.As we shall see, the proposed algorithm exacts initial values. To provide the concerned values we are obliged to calculate the theoretical autocorrelation ofp-order RCA models. To ward off tedious computation we implement an algorithm, more details in Section 4, who calculate the autocorrelations in numerical way.This paper is organized as follows. In Section2, we present definition and some basic properties of RCA models. In Section 3 we recall Kalman’s algorithm and apply it to calculate likelihood function for RCA models. In Section 4, we describe our estimating algorithm. The performance of this algorithm is examined by Monte Carlo simulations and compared to the quasi-maximum likelihood method, in Section 5. Finally, we achieve with a conclusion.


## 2. Stationarity and Moment Properties

Let(
X
(
t
)
,
t
∈
Z
) be a R
C
A
(
p
) model. In order to derive stationarity conditions for the processes X
(
t
), Nicholls and Quinn proposed the following vectorial presentation.Define thep random vector Y
(
t
) by Y
(
t
)
=
X
t
+
1
-
p
,
…
,
X
t
′, equation (1) may be rewritten in terms of Y
(
t
) by(2)
Y
t
=
M
+
D
t
Y
t
-
1
+
η
twhere the p
×
p matrix M is given by M
=
0
1
⋯
0
0
0
⋱
0
⋮
0
⋯
1
ϕ
p
ϕ
p
-
1
⋯
ϕ
1, D
(
t
)
=
L
′
⊗
B
(
t
), where B
(
t
)
=
[
β
p
(
t
)
,
…
,
β
1
(
t
)
] and L
=
[
0
,
…
,
0,1
]
∈
M
p
,
1 and η
(
t
)
=
L
′
⊗
ε
(
t
).Theorem 1.
IfF
t is the σ
-field generated by (
ϵ
(
s
)
,
β
(
s
)
)
,
s
≤
t and C
>
0 then there exists a unique F
t
-measurable second-order stationary solution to (1) if and only if M has all its eigenvalues within the unit circle and v
e
c
C
′
A
<
1, where C
=
E
[
B
t
B
′
(
t
)
] and A is the last column of the matrix (
I
d
-
M
⨂
M
)
-
1 (see Andel (1976) or Nicholls and Quinn [1]).In the next theorem, Lian et al. [4] give us autocorrelation structure and marginal variance of RCA(p) models.Theorem 2.
Consider the stationaryR
C
A
(
p
) in (1). (i)
The processes have the same autocorrelation structure as theA
R
(
p
) processes:(3)
ρ
k
=
∑
j
=
1
p
ϕ
j
ρ
k
-
j
,
k
≥
1
(ii)
The marginal variance of a stationaryR
C
A
(
p
) process is given by(4)
γ
0
=
σ
ε
2
1
-
∑
j
=
1
p
∑
k
=
1
p
ϕ
j
ϕ
k
+
σ
j
k
ρ
j
-
kIn practice, this theorem will enable us to deduce the starting conditions for Kalman filter. Hence, we develop in Section3 a recursive algorithm called variance to calculate the autocovariance function for RCA models.


## 3. Quasi-Maximum Likelihood and Kalman Filter

Assuming the joint normality of(
ϵ
(
t
)
,
t
∈
Z
) and (
B
(
t
)
,
t
∈
Z
), the conditional log-likelihood is given by(5)
L
x
;
θ
=
log
⁡
∏
t
=
1
n
f
x
t
∣
χ
t
-
1
=
∑
t
=
1
n
log
⁡
f
x
t
χ
t
-
1whereθ is the vector of unknown parameters, χ
t
=
(
x
t
,
…
,
x
t
) is a sample of observations at t
=
1
,
…
,
n, and f
(
x
t
∣
χ
t
-
1
) indicate the normal density function of x
t given χ
t
-
1, with mean x
^
t
∣
t
-
1
=
E
[
X
t
∣
χ
t
-
1
] and variance M
^
t
∣
t
-
1
=
E
[
(
x
t
-
x
^
t
∣
t
-
1
)
2
]. The quasi-maximum likelihood may be put in the following form:(6)
L
x
;
θ
=
-
n
2
l
o
g
2
π
-
1
2
∑
t
=
1
n
l
o
g
M
^
t
∣
t
-
1
-
1
2
∑
t
=
1
n
x
t
-
x
^
t
∣
t
-
1
2
M
^
t
∣
t
-
1Given this formula the quasi-maximum likelihood can be calculated using the Kalman filter; see Hamilton [5].In order to apply Kalman filter, we consider the appropriate state-space representation forR
C
A model in (1).(7)
Z
t
+
1
=
F
Z
t
+
G
ν
t
+
e
t
;
s
t
a
t
e
e
q
u
a
t
i
o
n
,
X
t
=
H
Z
t
;
o
b
s
e
r
v
a
t
i
o
n
e
q
u
a
t
i
o
nwhereZ
t
=
X
t
,
…
,
X
t
-
p
+
1
′
∈
M
p
,
1
(
R
), H
=
(
1,0
,
…
,
0
)
∈
M
1
,
p
(
R
), e
t
=
ε
t
,
0
,
…
,
0
′
∈
M
p
,
1
(
R
), ν
t
-
1
=
(
X
t
-
1
β
1
(
t
)
,
…
,
X
t
-
p
β
p
(
t
)
)
∈
M
p
,
1
(
R
)
(8)
F
=
ϕ
1
ϕ
2
⋯
ϕ
p
-
1
ϕ
p
1
0
⋯
0
0
⋮
⋮
⋱
⋮
⋮
0
0
⋯
1
0
∈
M
p
,
p
R
,
G
=
1
⋯
1
0
⋯
0
⋮
⋱
⋮
0
⋯
0
∈
M
p
,
p
R
,
e
t
~
N
0
,
Q
w
h
e
r
e
Q
=
σ
ϵ
2
0
⋯
0
0
0
⋯
0
⋮
⋮
⋱
⋮
0
0
⋯
0We should point out that the previous representation was proposed by Benmoumen in his M.Sc. Thesis [6].Now, we describe the Kalman filter in the aim to build the log-likelihood function. The Kalman filter is a recursive algorithm derived by Kalman [7] to provide an optimal forecast Z
^
t
+
1
∣
t
=
E
[
Z
t
+
1
∣
X
1
,
…
,
X
t
] of Z
t
+
1 given X
1
,
…
,
X
t, with the mean square error P
t
+
1
∣
t
=
V
[
Z
t
+
1
-
Z
^
t
+
1
∣
t
]
,
t
=
1
,
·
,
n. Given starting values Z
^
1
∣
0 and P
1
∣
0 which are derived from Theorem 2, the recursive procedure is as follows:(i)
Calculate the forecastingX
^
t
∣
t
-
1 of the observation X
t, and the error M
^
t
∣
t
-
1 of this forecast.
(ii)
Update the state vectorZ
^
t
∣
t. Compute P
t
∣
t the MSE of this updated projection.
(iii)
Calculate the forecastingZ
^
t
+
1
∣
t and the MSE P
t
+
1
∣
t of this forecast.Thus, we could construct the log-likelihood function using Kalman filter, so as to obtain the maximum likelihood estimators and in order to avoid the fastidious computation of partial derivatives ofL
(
x
;
θ
) we used the simulated annealing method (see Corana et al. [8]) which is a global optimization algorithm.


## 4. Estimating Algorithm for p-Order Random Coefficient Autoregressive Model parameters

The algorithm proposed here is a generalization of a procedure developed by Benmoumen and al. (2013) to estimate first-order RCA’s parameters.Recently, the same idea has been developed for parameters estimation in GARCH(1,1), ARCH(1), and ARCH(p) models by Benmoumen et al. (2011, 2014, and 2015).Before describing our algorithm MLKF (quasi-maximum likelihood and Kalman filter estimation), it is worthwhile to provide a subalgorithm which tests if parameters fulfill the conditions of stationarity; we will denote it by Test. The second subalgorithm, which we must provide, concerns the computation ofL
(
x
;
θ
) by Kalman filter; we will denote it by KF. These two subalgorithms will be implemented in our global estimating algorithm.Herein, we are interested in minimizing-
L
(
x
;
θ
)
=
l
(
x
;
θ
).Subalgorithm Test(θ)
if The eigenvalues of the matrix M have modul less than unity and C
A
<
1
then
Then go to next
else
Take the last point as starting point
end if
End SubalgorithmSubalgorithm variance(θ)
Solve the equations:A
ρ
1
,
⋯
,
ρ
p
′
=
0
,
⋯
,
0
′
Where(9)
A
i
j
=
ϕ
i
-
1
s
i
j
=
1
;
i
=
2
,
⋯
,
p
+
1
,
ϕ
i
-
j
+
ϕ
i
+
j
-
2
s
i
j
=
2
,
⋯
,
p
+
1
;
i
=
2
,
⋯
,
p
+
1
.
And(10)
ϕ
k
=
-
1
k
=
0
,
0
k
>
p
Calculate:(11)
γ
0
=
σ
ε
2
1
-
∑
j
=
1
p
∑
k
=
1
p
ϕ
j
ϕ
k
+
σ
j
k
ρ
j
-
k
End SubalgorithmSubalgorithm KF(θ)
Given the starting conditionsZ
^
0
∣
1 and P
0
∣
1
Calculate:X
^
0
∣
1
=
H
Z
^
0
∣
1 and M
^
0
∣
1
=
H
P
0
∣
1
H
′
for
t
=
1 to n
do
K
t, Z
^
t
∣
t, P
t
∣
t, Z
^
t
+
1
∣
t and P
t
+
1
∣
t
X
^
t
+
1
∣
t and M
^
t
+
1
∣
t
end for
(12)
s
o
m
=
1
2
l
o
g
2
π
for
t
=
1 to n
do
(13)
s
o
m
=
s
o
m
+
1
2
l
o
g
M
^
t
+
1
∣
t
+
1
2
X
t
-
X
^
t
+
1
∣
t
2
M
^
t
+
1
∣
t
end for
(14)
l
x
;
θ
=
s
o
m
End SubalgorithmIndeed, as we shall see later our algorithm is an iterative process requiring initial estimates of the parameters to commence the iterations. The consistent least squares estimates are suitable for this purpose.MLKF Algorithm
Step 1: Initialize: the vector parametersθ the step vector ν and the temperature T.
Step 2: Starting from the pointθ
i, generate a random point θ along the direction h: θ
=
θ
i
+
r
ν
m
h
e
h, where r is a random number generated in the range [
-
1,1
] by a pseudorandom generator; e
h is the vector of the hth coordinate direction; and ν
m
h is the component of the step vector ν along the same direction.
Step 3: Call sub algorithmT
e
s
t
(
θ
)
Step 4: Call sub algorithmK
F
(
θ
)
ComputeK
F
(
θ
i
) and K
F
(
θ
)
If
K
F
(
θ
)
≤
K
F
(
θ
i
) accept the new point
Else accept or reject the new point with acceptance probability p:(15)
p
=
e
x
p
K
F
θ
i
-
K
F
θ
T
generate a uniformly distributed random numberp
′ in the range [
0,1
]
If
p
′
<
p, the point is accepted otherwise it is rejected.
Step 5: Steps 2 to 4 are repeated for each coordinate directioni, i
=
1, …, m (m is the dimension of the vector parameter).
Step 6: Steps 2 to 5 are repeatedN
s times (N
s is the number of step variation) and the step vector ν is adjusted.
Step 7: Steps 2 to 6 are repeatedN
T times (N
T is the number of temperature reduction) the temperature is reduced following the rule: T
′
=
r
T
T with r
T
∈
[
0,1
].
Step 8: Steps 2 to 7 are repeated until a termination criterion is satisfied.
End Algorithm


## 5. Simulations

To examine the performance of our algorithm, we have carried out series of simulation experiments. In this study, we consider two examples of models RCA.(1)
(16)
X
t
=
β
1
t
X
t
-
1
+
0.36
+
β
2
t
X
t
-
2
+
ε
t
,
t
∈
Rwith (
ε
(
t
)
,
β
(
t
)
)
~
N
(
(
0,0
,
0
)
,
Ω
), where Ω
=
1
0
0
0
0.2176
0
0
0
0.2176.
(2)
(17)
X
t
=
0.8
+
β
1
t
X
t
-
1
+
-
0.15
+
β
2
t
X
t
-
2
+
ε
t
,
t
∈
Rwith (
ε
(
t
)
,
β
(
t
)
)
~
N
(
(
0,0
,
0
)
,
Ω
), where Ω
=
1
0
0
0
0.0919
0.0919
0
0.0919
0.1838For the models mentioned earlier, we generated 1000 replications of sample sizesn
=
100 and 200.The results of this experiment are displayed in Tables1–4 where for each estimator we give the mean and MSE, where we used notation QMLE for the quasi-maximum likelihood estimators and MLKF for the estimation by our algorithm.Table 1
Mean and MSE of estimated parameters for Example 1;n
=
100.
Parameters
MLKF Mean
MLKF MSE
QMLE Mean
QMLE MSE
ϕ
1
=
0
-0.00032
0.01491
-0.00051
0.03021
ϕ
2
=
0.36
0.32874
0.01265
0.29532
0.03127
σ
11
=
0.2176
0.20847
0.01506
0.17879
0.03414
σ
12
=
0
0.00416
0.00635
0.00631
0.00827
σ
22
=
0.2176
0.23316
0.01197
0.18413
0.03698
σ
ε
=
1
1.05658
0.05597
1.10734
0.07859Table 2
Mean and MSE of estimated parameters for Example 1;n
=
200.
Parameters
MLKF Mean
MLKF MSE
QMLE Mean
QMLE MSE
ϕ
1
=
0
0.00021
0.00523
0.00038
0.02003
ϕ
2
=
0.36
0.34857
0.00556
0.31857
0.02114
σ
11
=
0.2176
0.19859
0.00902
0.17016
0.02487
σ
12
=
0
0.00363
0.00412
0.00578
0.00723
σ
22
=
0.2176
0.21072
0.01041
0.19189
0.02865
σ
ε
=
1
1.02176
0.03397
1.09754
0.06682Table 3
Mean and MSE of estimated parameters for Example 2;n
=
100.
Parameters
MLKF Mean
MLKF MSE
QMLE Mean
QMLE MSE
ϕ
1
=
0.8
0.82356
0.01564
0.75590
0.03258
ϕ
2
=
-
0.15
-0.15886
0.00982
-0.18164
0.02042
σ
11
=
0.0919
0.08326
0.01013
0.13584
0.04782
σ
12
=
0.0919
0.03735
0.01126
0.02104
0.03152
σ
22
=
0.1838
0.21214
0.01342
0.24065
0.03356
σ
ε
=
1
0.93012
0.05263
0.89702
0.08237Table 4
Mean and MSE of estimated parameters for Example 2;n
=
200.
Parameters
MLKF Mean
MLKF MSE
QMLE Mean
QMLE MSE
ϕ
1
=
0.8
0.79036
0.00523
0.76137
0.02278
ϕ
2
=
-
0.15
-0.17777
0.00358
-0.17776
0.00935
σ
11
=
0.0919
0.10248
0.00647
0.08250
0.00838
σ
12
=
0.0919
0.06556
0.00512
0.05217
0.00912
σ
22
=
0.1838
0.21036
0.00788
0.23481
0.01014
σ
ε
=
1
0.97054
0.01463
0.94957
0.03327Remark 3.
(i)
The least squares estimate ofθ
=
(
ϕ
,
γ
,
σ
ε
2
), where ϕ
=
ϕ
1
,
…
,
ϕ
p
′ and γ
=
v
e
c
h
(
C
) is given by (18)
ϕ
^
=
∑
t
=
1
N
Y
t
-
1
Y
′
t
-
1
-
1
∑
t
=
1
N
Y
t
-
1
X
t
σ
^
ε
2
=
1
N
∑
t
=
1
N
u
^
2
t
-
γ
^
z
¯
γ
^
=
∑
t
=
1
N
z
t
-
z
¯
z
t
-
z
¯
′
-
1
∑
t
=
1
N
u
^
2
t
z
t
-
z
¯where u
^
(
t
)
=
X
(
t
)
-
ϕ
^
′
Y
(
t
-
1
)
,
t
=
1
,
…
,
N, z
(
t
)
=
K
n
v
e
c
(
Y
t
-
1
Y
′
(
t
-
1
)
), and z
¯
(
t
)
=
1
/
n
∑
t
=
1
n
‍
z
(
t
).
(ii)
The maximum likelihood estimate ofθ is defined by minimizing the following function: (19)
l
^
n
θ
=
1
n
∑
t
=
1
n
l
n
σ
ϵ
2
+
γ
′
z
t
+
1
n
∑
t
=
1
n
X
t
-
β
′
Y
t
-
1
2
σ
ε
2
+
γ
′
z
tUnder certain assumptions, both of the estimates, the least squares and the maximum likelihood, are strongly consistent and obey a central limit theorem. In fact, condition E
(
X
4
(
t
)
)
<
∞ is required for the strong consistency of the least squares estimates, and E
(
X
8
(
t
)
)
<
∞ for the asymptotic normality of these estimators. Respectively, in order that a central limit theorem exist for maximum likelihood estimators, conditions E
(
ε
4
(
t
)
)
<
∞ and E
(
β
k
4
(
t
)
)
<
∞, k
=
1
,
…
,
p will be required. For more details see Nicholls and Quinn [1].


### 5.1. Comparison with Quasi-Maximum Likelihood Estimators

In this series of simulation we compare our algorithm (MLKF) versus the quasi-maximum likelihood (QMLE). In each method we use simulated annealing algorithm for optimization and the least squares estimators for initiation.As is to be expected, the MLKF estimation procedure has performed better, as is seen from the fact that the sample mean square errors (MSE) are generally smaller than for the quasi-maximum likelihood estimators (QMLE). Hence, we can conclude that the performance of our estimation procedure is promising.


## 5.1. Comparison with Quasi-Maximum Likelihood Estimators

In this series of simulation we compare our algorithm (MLKF) versus the quasi-maximum likelihood (QMLE). In each method we use simulated annealing algorithm for optimization and the least squares estimators for initiation.As is to be expected, the MLKF estimation procedure has performed better, as is seen from the fact that the sample mean square errors (MSE) are generally smaller than for the quasi-maximum likelihood estimators (QMLE). Hence, we can conclude that the performance of our estimation procedure is promising.


## 6. Conclusion

In this paper, we constructed an algorithm for calculating the covariance matrix of RCA(p) models and we have generalized a procedure developed by Benmoumen et al. (2013) to estimate first-order RCA’s parameters.The log-likelihood function is constructed using the Kalman filter and is numerically maximized applying simulated annealing method. The results of our simulation study show that our estimation approach succeeds and it performs better than the competitor.



---

*Source: 8479086-2019-05-13.xml*