SDS-3
Some Properties of Estimation Methods
for Structural Relationships in Non-stationary Errors-inVariables Models
Naoto Kunitomo, Naoki Awaya
and Daisuke Kurisu

November 2017

Statistics & Data Science Series back numbers:
http://www.mims.meiji.ac.jp/publications/datascience.html

Some Properties of Estimation Methods for
Structural Relationships in Non-stationary
Errors-in-Variables Models ∗Naoto Kunitomo †Naoki Awaya ‡andDaisuke Kurisu §September 20, 2017

Abstract
Kunitomo and Sato (2017) developed a new statistical method called the Separating Information
Maximum Likelihood (SIML) estimation for the multivariate non-stationary errors variables models. In this paper we compare the SIML estimation and the maximum likelihood (ML) estimation
when there are non-stationary trends and noise components in the underlying multivariate timeseries. We have found that the Gaussian likelihood function can have some peculiar shape andthe computation of ML estimation can be unstable in some cases. The SIML estimation has theasymptotic robust properties under general conditions of the existence of fourth order moments
and it does not have any computational problem. We show the asymptotic properties of the ML
estimation and SIML estimation as well as their finite sample properties.
Key Words
Non-stationary economic time series, Errors-variables models, Co-integrated trends, Kn −transfomation,
Likelihood Function, Maximum Likelihood (ML), Separating Information Maximum Likelihood(SIML), Asymptotic Robustness∗

This is a preliminary memorandum and the research has been supported by JSPS Grant-in-Aid
for Scientific Research JP17H02513.
†
School of Potitical Science and Economics, Meiji University, Kanda-Surugadai 1-1, Chiyoda-ku
101-8301, Tokyo, JAPAN, naoto.kunitomo@gmail.com‡
Graduate School of Economics, University of Tokyo, Bunkyo-ku, Hongo 7-3-1, Tokyo 113-0033§
Graduate School of Economics, University of Tokyo, Bunkyo-ku, Hongo 7-3-1, Tokyo 113-0033,
JAPAN

1

1. Introduction
There is a vast amount of published research on the use of statistical time series
analysis of macroeconomic time series. One important distinction of macroeconomic
time series from the standard time series analysis in other areas has been the mixtureof non-stationarity and measurement errors including apparent seasonality. The
measurement errors in economic times series are often essential in the publishedoﬃcial macro-economic data, but enough attention has not been in the analysis of
non-stationary multivariate time series analysis.
In this regard, Kunitomo and Sato (2017) have developed a new statisticalmethod called the Separating Information Maximum Likelihood (SIML) estimation for the multivariate non-stationary errors-in-variables models. Earlier and related literature on the non-stationary economic time series are Engle and Granger(1987) and Johansen (1995), which have dealt with the econometric modeling ofmultivariate non-stationary and stationary time series and developed the notion of
co-integration, but they did not have paid much attention on measurement errors.
The problem discussed by Kunitomo and Sato (2017) is related to their work, but it
has diﬀerent aspects due to the fact that its main focus was on the non-stationarity,seasonality and measurement error in the non-stationary errors-in-variables models. Many economists tend to use the oﬃcial macro-economic data, but then theyshould have taken more attention on the measurement errors and oﬃcial seasonal
adjustments. There are important consequences by using economic data with themeasurement errors, which are common in many published macro-economic timeseries by central governments.
On the other hand, in the literature of statistical time series analysis, the state
space modeling for non-stationary time series has been developed by Akaike (1989)
and Kitagawa (2010) and there have been many applications in many fields including
control engineering and statistical seismology already reported (see Ohtsu, Pengand Kitagawa (2015) for instance). Their statistical method may look diﬀerent
from time series econometrics at the first glance. By looking their method carefully
and comparing with the econometric analysis of non-stationary multivariate timeseries, the underlying statistical problems are quite similar with respect to the nonstationarity and measurement errors. Then the study of Kunitomo and Sato (2017)
can be regarded as an investigation of statistical inference problem of state spacemodeling, which is related to an earlier work by Chiang, Jiang and Park (2008).
The main purpose of this paper is to compare the the SIML estimation and the
maximum likelihood (ML) estimation, which are two diﬀerent methods to estimate
the non-stationary errors-in-variables models when there are non-stationary trends
and noise components. In this paper we will report that the Gaussian likelihood
function show some peculiar shape and the computation of the ML estimation canbe unstable in some cases. Also when there are co-integrated relations among trends
2with the rank being smaller than the dimension of observations, there could be an
additional computational problem. On the other hand, the separating information
maximum likelihood (SIML) method gives an alternative way to overcome the underlying diﬃculty of computation while it has reasonable statistical properties. The
SIML estimation has the asymptotic robust properties under general conditions ofexistence of second order moments. We show the asymptotic properties as well as
its finite sample properties of the ML estimation and SIML estimation as well as
simulations.
In Section 2 we will present a general formulation of the problem and give simple
examples to illustrate the problem in this study. Then In Section 3, we will discuss
the simple one-dimensional case and the non-stationary model has the random walkplus noise and we will develop one common (non-stationary) factor case in multivariate errors-in-variables models. Section 4 we investigate the Gaussian likelihood
function and its shape. We give the consistency result of the ML estimation. It
seems to be new although it is not surprising result. In Section 5, we will discuss
some extensions and then present some concluding remarks in Section 6. Somemathematical derivations will be given in Appendix.

2. The Model
Let yij be the i−th observation of the j−th time series at i for i = 1, · · · , n; j =′′1, · · · , p. We set yi = (y1i , · · · , ypi ) be a p × 1 vector and Yn = (yi ) (= (yij )) be an
n×p matrix of observations and we denote y0 as the initial p×1 vector. We consider
the situation when the underlying non-stationary trends xi (= (xji )) (i = 1, · · · , n)′
are not necessarily the same as the observed time series and let vi = (v1i , · · · , vpi )be the vectors of the noise components, which are independent of xi . Then we use
the additive decomposition form
(2.1)

yi = xi + vi (i = 1, · · · , n),where the trend components xi (i = 1, · · · , n) satisfies
(2.2)

(x)

∆xi = (1 − L)xi = wi
(x)

(x)

(x)′

with Lxi = xi−1 , ∆ = 1 − L, E(wi ) = 0, E(wi wi ) = Σx , and a sequence of′stationary components satisfies vi (i = 1, · · · , n) with E(vi vi ) = Σv and
(2.3)

vi =

∞
∑

Cj ei−j ,j=−∞with absolutely summable coeﬃcients Cj and a sequence of i.i.d. random vectors′
with E(ei ) = 0, E(ei ei ) = Σe .
(x)
We assume that wi and ei are the sequence of i.i.d. random vectors with Σe being3

(x)positive-semi-definite, and the random vectors wi and ei are mutually independent.
When vi = ei , we can interpret that it is a sequence of independent measurement
errors. The present additive decomposition is similar to the one given by Kitagawa
and Gersch (1984) and Kitagawa (2010).
In order to develop the SIML estimation, we first consider the situation when ∆xi
and vi (i = 1, · · · , n) are mutually independent and each of the component vectors
are independently, identically, and normally distributed as Np (0, Σx ) and Np (0, Σv ),′respectively. We use an n × p matrix Yn = (yi ) and consider the distribution of′′ ′np × 1 random vector (y1 , · · · , yn ) . Given the initial condition y0 , we have
(
)
′
′
(2.4)
vec(Yn ) ∼ Nn×p 1n · y0 , In ⊗ Σv + Cn Cn ⊗ Σx ,′where 1n = (1, · · · , 1) and


1 0
 1 1

Cn = 
 1 1
 1 ···
1 ···

(2.5)

···
0
1
1
1

0
···
···
1
1


0
0 

0 
.


0
1 n×n

Then, given the initial condition y0 the conditional maximum likelihood (ML) estimator can be defined as the solution of maximizing the conditional log-likelihood
function except a constant as′

L∗n = log |In ⊗ Σv + Cn Cn ⊗ Σx |−1/2
1
′ ′
′
′
− [vec(Yn − Ȳ0 ) ] [In ⊗ Σv + Cn Cn ⊗ Σx ]−1 [vec(Yn − Ȳ0 ) ] ,
2
where
′

Ȳ0 = 1n · y0 .

(2.6)

′

We use the Kn −transformation that from Yn to Zn (= (zk )) by
(
)
(2.7)
Zn = Kn Yn − Ȳ0 , Kn = Pn C−1
n ,
where

(2.8)




1
0 ··· 00
 −1 10 ··· 0 −1

0 ··· Cn =  0 −1 1
,

 0

0 −1 10000 −1 1
n×n
4

and
(2.9)

√
(n)
(n)
Pn = (pjk ) , pjk =

]
[
2
2π
1
1
cos
(k − )(j − ) .
2n + 1
2
2
n + 12′′−1
By using the spectral decomposition C−1= Pn Dn Pn and Dn is a diagonaln Cnmatrix with the k-th element

dk = 2[1 − cos(π(2k − 1))] (k = 1, · · · , n) .
2n + 1

Then the conditional likelihood function given the initial condition is proportionalto(2.10)

Ln =

n
∑

1
log |a∗kn Σv + Σx |−1/2 −

2 k=1

k=1

where
(2.11)a∗kn (= dk )n∑′zk [a∗kn Σv + Σx ]−1 zk ,)][ (π 2k − 1(k = 1, · · · , n) .
= 4 sin
2 2n + 1
2

We have used two transformations on the non-stationary time series into the sequence of independent random variables zk (k = 1, · · · , n) which follows Np (0, Σx +a∗kn Σv ), and the coeﬃcients a∗kn is a dense sample of 4 sin2 (x) in (0, π/2) 1 .
′It may be natural to use zk zk to estimate a∗kn Σv + Σx since it is the variancecovariance matrix of zk . We notice that a∗kn → 0 as n → ∞ for a fixed k. Whenk is small, a∗kn is small and we can∑expect that k = kn depending n is still smalln∗when n is large. However, (1/mn ) m
k=1 akn is not small if mn is near to n, whichsuggests the condition mn /n → 0 as n → ∞. The separating information maximum
likelihood (SIML) estimator of Σ̂x can be defined byn
1 ∑′zk zk .
Σ̂x,SIM L =mn k=1m(2.12)

This estimator of the variance-covariance of non-stationary trends is trying to use
the information on trends in the frequency domain, which corresponds to only thetrend parts without measurement errors from the time series observations. ForΣ̂x , the number of terms mn should be dependent on n. Then we need the order
requirement that mn = O(nα ) and 0 < α < 1, which is the first property of themacro-SIML estimation.
From our construction of the SIML estimation the essential features of estimation
do not depend on the noise terms as long as the noise terms are stationary. It has
been pointed out by Kunitomo and Sato (2017).
We have used a∗kn , which is slightly diﬀerent ak in Kunitomo and Sato (2013) and the lattercorresponds to akn = na∗kn .
1

5

3. Simple Cases
3.1

An Illustration

Let yi be the i−th observation of time series for i = 1, · · · , n and yn = (yi ) be an
n×1 vector of observations. (y0 is the initial observation.) We consider the situation
when the underlying non-stationary trends xi at i = 1, · · · , n are not necessarily thesame as the observed time series and followed by
(x)

(3.1)

xi = xi−1 + vi
(x)where vi are the independent random variables followed by N (0, σx2 ) and x0 is theinitial variable. Let vi are the independent random variables followed by N (0, σv2 ),which are independent of xi (i = 1, · · · , n). Then we consider the additive modelyi = xi + vi (i = 1, · · · , n) .

(3.2)

Then the log-likelihood function is proportional to
(3.3)

Ln =

n
∑
k=1

where
(3.4)

zk2
1∑
,
2 k=1 a∗kn σv2 + σx2
nlog |a∗kn σv2 + σx2 |−1/2 −a∗kn = 4 sin2)][ (π 2k − 1(k = 1, · · · , n) .
2 2n + 1

By defining τ = σx2 /σv2 (≥ 0), we rewrite −(1/2)Ln as
(3.5)

n
n
∑
1 ∑ zk22∗[log σv + log(akn + τ )] + 2
L1n =.
σv k=1 a∗kn + τ
k=1

Since zk ∼ N (0, a∗kn σv2 + σx2 ) (k = 1, · · · , n), then the maximum likelihood estimatorof σv2 is given by
(3.6)

n
∑

zk2
n k=1 a∗kn + τ

1
2
σ̂v.M
L =

and the concentrated (normalized) likelihood function) in this simple case is proportional to
(3.7)

n
n
1 ∑ zk2
1∑
L1n (τ ) = log[]+
log[a∗kn + τ ] + 1 .
n k=1 a∗kn + τ
n k=1

6

It may not be straightforward to obtain the maximum likelihood estimator of τ
because the likelihood function may not be a simple function. Also Hirotsugu Akaike(1985) argued that we should impose the restriction 0 < τ ≤ 1 in our setting becauseof the usefulness of models. He advocated to estimate the models by using the ABIC
maximization procedure.
On the other hand, we approximate (1.9) as
LSI
n =

(3.8)

mn
∑

n
1∑
zk2
,
2 k=1 a∗kn σv2 + σx2

mlog |a∗kn σv2 + σx2 |−1/2 −k=1

and we have some requirement on mn /n = o(1). We notice that a∗kn = o(1) whenk = 1, · · · , mn and n → ∞. Then the SIML estimator σx2 can be given by
n
1 ∑
z2 .
mn k=1 km

2σ̂x.SIM
L =(3.9)

However, the macro-SIML estimation of σv2 is not the same as the original (finance)
SIML estimation because a∗kn = O(1) for k = n−mn +1, · · · , n. One way to estimate
σv2 is to use the fact that
(3.10)E[nn
1∑ ∗ 2
1∑ 2
zk ] = σx2 + (a )σ −→ σx2 + 2σv2 (n → ∞).
n i=1
n i=1 kn v

It is because a∗kn −→ 0 when kn /n → 0 (n → ∞) and a∗kn → 4 (kn /n → 1, n → ∞).
Then a possible SIML estimator σv2 can be given by(3.11)1 1
2
σ̂v.SIM
L (1) = [n∑2z 2 − σ̂x.SIM
L] ,2 n k=1 kwith the restriction of the resulting positivity.
For the estimation problem of high-frequency financial data, Kunitomo and Sato(2013) have suggested to use(3.12)

2σ̂v.SIML =n
1 ∑ −1 2a z ,ln k=n−l kn knwhere akn = na∗kn and ln = o(n) for the hight frequency asymptotics. However,in ∑the present case of Macro-SIML with a∗kn it is straightforward to show thatn∗−1
1k=n−ln +1 akn → 1/4 as n → ∞, the macro-SIML estimator is given byln(3.13)

2σ̂v.SIM
L (2) =

n
1 ∑
1 2
2
a∗−1
kn zk − σ̂x.SIM L
ln k=n−l +14nwith the restriction of the resulting positivity.
73.2

A Non-stationary Common Trend CaseLet yi be the i−th observation of p-dimensional time series for i = 1, · · · , n and′yi = xi + vi . Also let Yn = (yi ) be an n × p (p ≥ 1) matrix of observations. We
assume that the vectors xi satisfy
xi = xi−1 + πµ∗i ,(3.14)

where π is a p × 1 vector, µ∗i is a sequence of i.i.d. (one-dimensional) randomvariables 2 following N (0, σµ2 ) and vi are i.i.d. (p-dimensional) random variables
following Np (0, Σv ) with the variance-covariance matrix Σv , which is a non-singularmatrix. We take b = σµ π, A = a∗kn Σv and apply the matrix formulas such that fora positive definite A we have′′

|A + bb | = |A|[1 + b A−1 b](3.15)
and′′′[A + bb ]−1 = A−1 − A−1 b[1 + b A−1 b]−1 b A−1(3.16)′for Σx = bb .
Then the likelihood function Ln is proportional to (−1/2) times
L1n =
=

n [
∑

]
′
−1
2
a∗−1
−1
∗−1 ′ −1∗−1 ′∗kn (zk Σv b)
log |akn Σv | + log(1 + akn b Σv b) + akn zk Σv zk − ∗akn + b′ Σ−1
v b
k=1
[
]
′
n
n
n
−1
2
∑
∑
∑
a∗−1
−1
∗−1 ′
∗−1
∗
kn (zk Σv b)
log |akn Σv | +akn zk Σv zk +log(1 + akn c) −
,
∗
a
+
c
kn
k=1
k=1
k=1

where we set′

c = σµ2 π Σ−1v π .

(3.17)

We need a normalization for the vector π. If we take a simple normalization, themaximum likelihood estimator of π could be a quite complicated solution of the′likelihood equation even when p = 2. One possible normalization is to take β =′(1, −β 2 ) and then the maximization is not a trivial task.
As an alternative way to solve the present problem is to use the conditions that′E[zk zk ] = Σx + o(1) fork = 1, · · · , mn
and

2

1
′
E[a∗−1kn zk zk ] = Σv + Σx + o(1) fork = n + 1 − mn , · · · , n .
4

The notation µ∗i is diﬀerent from µi and µ∗i = ∆µi (i = 2, · · · , n) in Kunitomo and Sato (2017).

8

The rank of matrix Σx is one while the matrix Σv has a full rank. Then it may benatural to consider the characteristic equation[](3.18)
Σ̂x.SIM L − λΣ̂v.SIM L β̂ SIM L = 0 ,
and
n
1 ∑
′
zk zk ,
mn k=1m(3.19)(3.20)

Σ̂x.SIM L =

n
1 1∑
′
Σ̂v.SIM L (1) = [zk zk − Σ̂x.SIM L ] ,2 n k=1or(3.21)n1
1 ∑′Σ̂v.SIM L (2) =a∗−1
kn zk zk − Σ̂x.SIM L ,ln k=n+1−l4nwhere we have the restrictions of the resulting positive definiteness of the estimatedvariance-covariance matrix.
We use the notation that
)
(
′
′
,
Y
−
1
ȳ
(3.22)
Zn = (zk ) = Pn C−1nn0nΣ̂v.SIM L is a SIML estimator of Σv , and λ is the (scalar) eigen value. Because therank of Σx is degenerated and it is one, we need to take the smallest eigenvalue
λ1 . We have the β̂ SIM L , which is called the SIML estimator of β. A simplified(consistent) estimation may be given by(3.23)

Σ̂x.SIM L × β̂ SIL = 0 ,that is,(3.24)

Σ̂x.SIM L × [1−β̂ 2.SIL]=0.

We can solve as(3.25)−1β̂ 2.SIL = Σ̂22x.SIM L Σ̂21x.SIM L ,which can be the least squares estimator for the transformed variables and it iscalled the SILS estimator.

9

3.3

On Autocorrelations of Noise terms

The main statistical problem of the time series decomposition models in Kunitomo
and Sato (2017) was the estimation of non-stationary trend components and seasonal components. The results do not much depend on the specification of noise
components. When there are autocorrelations in the noise terms as (2.3), the analysis of frequency domain on the underlying errors-in-variables models may give someinsight on the issue.
When the vector sequence of noise terms vi (i = 1, · · · , n) follows the stationary
stochastic process as (2.3), we can represent the spectral density p × p matrix of vi
as(3.26)

∞
∞
∑
1 ∑
2iλj
fv (λ) = (Cj e )Σe (Cj e−2iλj )π j=−∞
j=−∞

(−

π
π
≤ λ ≤ ),
22

where i2 = −1 (see Chapter 7 of Anderson (1971) for instance.)
Then the p × p spectral density matrix of the transformed vector process ∆yi (=yi − yi−1 ) can be also represented as(3.27)f∆y (λ) =]1[
Σx + (1 − e2iλ )fv (λ)(1 − e−2iλ ) .
πSince the transformed random vectors z k (k = 1, · · · , n) correspond to the Fouriertransformed vectors of ∆yi (i = 1, · · · , n) except the initial condition y0 , it is possible to estimate the spectral density matrix fv (λ) from the sequence of observationsof zk (k = 1, · · · , n) in principle.
The eﬃcient estimation on fv (λ) could be certainly possible, but we need toinvestigate further related problems such as the kernel estimation and the eﬀects
of small sample sizes usually encounterd in analysing macro-economic time series
by the multivariate economic time series, which are quite diﬀerent from the initial
purrpose of the present investigation,it is beyond scope of this paper and we willtake have another occasion.

4. Gaussian Likelihood Function and Maximum Likelihood
Estimation
It may be natural to think that we could apply the general principle of the maximum
likelihood (ML) method. For the resulting simplicity, we use the case of i.i.d. noises
for the ML estimation. One of interesting aspects of the present problem is the fact
the ML method does not necessary give a satisfactory solution. We illustrate thisproblem by using Example 1 in Section 2 of Kunitomo and Sato (2017). We set thetrue parameter values as σµ2 = 0.4, β2 = 1.0 and
(
)
(
)
′
0.45 0.23
12
Σv =, Σx = σµ ππ , π =.
0.23 0.4−β210−100000−200000−150000

loglikelihood−50000

0

beta2

−4

−2

0

2

4

value of beta2

Figure 1 : Likelihood Function of β2 (n = 1, 000)
Then we generate a set of simulated observations as a typical example and we havedrawn two Gaussian likelihood functions of β2 in Figures 1-3 when 1, 000, given thetrue values for other parameters. We have found that the likelihood function could
have some peculiar forms in some cases as illustrated by Figures 1 and 2. Also
we have found that the Gaussian likelihood function is rather flat with respect tothe correlation coeﬃcient parameter. It seems that these are some of importantconsequences in the non-stationary errors-in-variables models.
We have drawn one wrong likelihood function in Figure 4 as an illustration on theassumption of Gaussian distributions. We generated the random variables followedby the uniform distribution on [−2, +2]. We have found that the variance-covarianceestimation crucially depends on the assumption of Gaussianity as we had expected.

Now we investigate the analytical properties of the likelihood function in a detail.

11−6000−12000−10000−8000

loglikelihood

−4000

−2000

0

beta2

−40

−20

0

20

40

value of beta2

Figure 2 : Likelihood Function of β2 (n = 1, 000)

−20000−40000−60000loglikelihood

0

rho

−1.0

−0.5

0.0

0.5

1.0

value of rho

Figure 3 : Likelihood Function of ρ (n = 1, 000)12

loglikelihood

−3500

−3000

−2500

−2000

rho

0.0

0.2

0.4

0.6

0.8

1.0

value of rho

Figure 4 : Wrong Likelihood Function of ρ (n = 1, 000)
We normalize the Gaussian likelihood function by n and rewrite
L∗1n

n
n
1∑ ∗
1∑
=
|a Σv | +log(1 + a∗−1
kn c)n k=1 kn
n k=1
]
[
n
1
1 ∑ ∗−1
′
′
′
−1
−1
∗
a tr Σv (Σv − ∗bb )Σv (zk zk − (akn Σv (θ0 ) + b(θ0 )b(θ0 ) ))+n k=1 kn
akn + c[]n
1 ∑ ∗−1
1
′
′
−1
−1 ∗
+
a tr Σv (Σv −bb )Σv (akn Σv (θ0 ) + b(θ0 )b(θ0 ) ))n k=1 kn
akn + c1∑ ∗∗|a Σv | + L12n
+ L∗13n + L∗14n (say),=n k=1 knnwhere Σv (θ0 ) andf b(θ0 ) are evaluated at the true parameter values.
We prepare the next lemma.
Lemma 4.1 : Let a p × p random vector zk follows Np (0, Q). Then for any p × p
matric Ak ,
(4.1)

′

E[(tr(Ak zk zk ))2 ] = [tr(Ak Q)]2 + 2tr(Ak QAk Q) .

By using this lemma, as n −→ ∞ the second term converges to
(4.2)

p

L∗13n −→ 0 .

Then we can establish the next result and the proof is given in Appendix. Although
we expect that the ML estimation under the Gaussian assumption has good asymptotic properties, we could not find any proof in the present setting and then decided13

to give the proof. We can also expect the asymptotic normality of the ML estimator,
but we have omitted its discussion because we need more spaces.
Theorem 4.2 : Assume that vi (i = 1, · · · , n) are a squence of i.i.d. random
vectors and |Σv | ̸= 0. Then under the assumption of Gaussian distributions the
maximum likelihood estimator of β is consistent as n −→ ∞.
Remark 4.1 : It should be noted that in time series econometrics it has been
known that coeﬃcient parameter vector β can be estimated by using the standard
regression method if the observed variables are co-integrated. However, the ML
method and the SIML method estimate not only the coeﬃcient vector consistently,
but also the variance-covariance matrices of trend terms and the noise terms at thesame time. Johansen (1995) had developed the ML method without any noise terms
as well as seasonality. His main interested was the hypothesis testing problem of therank condition.
We should mention to the first part of Theorem 4.1 of Kunitomo and Sato (2017)
on the asymptotic property of the SIML estimator under general conditions. They
have given the asymptotic normality of the SIML estimator under general moment
conditions.
Theorem 4.3 (Kunitomo and Sato 2017) : Assume the non-stationary errors-invariables model given by (2.1)-(2.3) and |Σv | ̸= 0. Then under the assumption ofexistence of fourth order moments the SIML estimator of β is consistent as n −→ ∞.
Although the SIML estimator has asymptotic normality, we have omitted this aspectin this paper. (See Kunitomo and Sato (2017) for the details.)

5. Simulations
There has not been any simulation result on alternative estimation methods for the
non-stationary time series models with errors-in-variables as far as we know. It
may be because we need a careful treatment of the non-stationarity with errors-invariables. However, since there are many situations with macro-economics variables
that we observe the non-stationarity and measurement errors, it is worthwhile toinvestigate the related issues by using simulations.
In our simulations of the SIML method, we have set σµ2 = 1, σv2 = 0.5, 2 or 4,β2 = 1.5, n = 80 or 400 and mn = nα with α = 0.6 or 0.7 as the parametrization,and the number of Monte Carlo repetition is 1,500 in each case. We can summarize
our setting of simulations as(Σx =Σx,11 Σx,12
Σx,12 Σx,22)(=

2.25 1.5
1.5 1

)

(
, Σv =

14

Σv,11 Σv,12
Σv,12 Σv,22

)

(
=

σv2 0
0 σv2

)
.

Among many simulations we have summarize them as Table 1. In Table 1 thenumber inside the parentheses are the standard deviation of estimators calculated
by our simulations.
σv2

α

n
80

0.6
400
0.5
80
0.7
400
80
0.6
400
2
80
0.7
400
80
0.6
400
4
80
0.7
400

Σx,11
2.263
(0.904)2.239(0.542)2.294(0.706)2.296(0.393)2.378(0.948)2.318(0.534)2.629(0.833)2.410(0.534)2.639(1.040)2.377(0.558)3.118(1.005)2.601(0.451)

Σx,12
1.475
(0.602)1.484(0.362)1.454(0.469)1.499(0.263)1.438(0.617)1.500(0.352)1.452(0.543)1.479(0.440)1.469(0.669)1.503(0.368)1.457(0.636)1.483(0.298)

Σx,22
1.038
(0.415)1.006(0.244)1.086(0.337)1.044(0.182)1.163(0.455)1.060(0.245)1.438(0.471)
1.166(0.287)1.384(0.547)1.127(0.267)1.885(0.630)1.357(0.249)Σv,11 (2)0.495
(0.445)0.512(0.283)0.521(0.386)0.498(0.214)1.992(1.020)1.981(0.615)1.943(0.816)1.975(0.465)
3.927(1.837)
3.933(1.099)3.777(1.427)3.937(0.806)

Σv,12 (2)0.010
(0.305)0.009(0.185)0.030(0.250)0.007(0.137)0.006(0.699)0.006(0.411)
0.032(0.570)0.007(0.309)-0.010(1.214)-0.008(0.727)0.065(0.963)0.010(0.550)Σv,22 (2)0.499
(0.318)0.501(0.186)0.496(0.251)0.494(0.136)1.922(0.852)1.991(0.535)1.944(0.709)1.946(0.407)3.990(1.670)3.965(1.025)3.874(1.311)3.908(0.727)β2.SIL
1.420(0.114)1.475(0.038)1.339(0.131)1.436(0.047)1.233(0.231)
1.415(0.077)1.017(0.230)1.267(0.093)1.072(0.296)1.334(0.105)0.787(0.274)1.095(0.119)β2.SIM L (2)1.528(0.178)1.502(0.038)1.554(0.269)1.502(0.054)1.630(0.810)1.504(0.082)1.701(0.901)1.514(0.123)1.702(1.042)1.514(0.132)1.846(1.452)1.519(0.197)

Table 1: Sample mean of estimators

We have done a number of simulations, but the results are similar to Table
1. There are several interesting findings. First, on the eﬀects of sample sizes theperformance of the estimators of the SIML estimation becomes as the sample size
increases as we had expected. Second, when the variances of noises are small,both the SILS estimator and the SIML estimator give reasonable estimates on thecoeﬃcient parameter, the former is slightly biased toward zero while the latter has
some correction of this bias. Third, when when the variances of noises are not small,the SILS estimator has a significant bias.
In addition to these simulations we have done a number of simulations when theunderlying distributions are not true and there are some stationary noises insteadof i.i.d. noises. In order to investigate the eﬀects of autocorrelations in the noise

15

terms, we consider the V AR2 (1) model given by
vi = Φvi−1 + v̄i ,
i.i.d.

v1 ∼ N (0, Σv ), v̄i ∼ N (0, Σ̄v ),vec(Σ̄v ) = (I4 − Φ ⊗ Φ)vec(Σv ),
)
(
√
0.4
0.3
0.4
×
0.45
√
Φ = diag(0.5, 0.5), Σv =0.3 0.4 × 0.45
0.45Then we summarize the simulation results as Table 2. Inside the parentheses inTable 2 are standard deviation of estimators.
n100
β2 = −1.5 200
400100
β2 = −0.5 200
400100
β2 = 0.5

200
400
100

β2 = 1.5

200
400

2
σx.SIML
0.234(0.083)0.225(0.065)0.216(0.045)0.251(0.086)0.221(0.052)0.208(0.048)0.264(0.105)0.224(0.061)
0.214(0.042)0.232(0.144)0.209(0.076)0.204(0.054)β2.SIL β2.SIM L (2)-1.154-1.512(0.275)(0.469)-1.208-1.493(0.176)(0.330)-1.281-1.497(0.112)(0.175)-0.258-0.456(0.161)(0.467)-0.337-0.488(0.102)(0.226)-0.403-0.494(0.074)(0.121)0.436
0.601(0.115)(0.345)0.465
0.572(0.102)(0.242)0.498
0.552(0.041)(0.113)
1.183
1.714(0.196)(0.773)1.274
1.585(0.162)(0.507)1.341
1.509(0.079)(0.147)

Table 2: Sample mean of estimators

As we had expected, the SIML estimation of the non-stationary part does notdepend on the autocorrelation structure of noise terms while the ML estimation
depends on the true structure of the underlying process. Although both the ML
16estimation and the SIML estimation give similar stable results on the coeﬃcient β2
in some cases, the ML estimation has sometimes computational diﬃculties when theabsolute size of coeﬃcients are large. (It is the case when β2 is relatively large as|β2 | = 1.5 in our simulations.) It may be due to the shape of the likelihood function
becomes some peculiar shapes in such cases 3 . But we have omitted the details of alarge number of simulations.
To summarize our simulations, the finite sample performance of the SIML estimation gives reasonable performances as the asymptotic theory has suggested as in
the previous sections.

6. Extensions
There can be several extensions of the problem we have been investigating and theresults obtained by Kunitomno and Sato (2017).
First, for the multivariate non-stationary (economic) time series, there are possibilities of co-integration in trends. In our framework, it may be interesting toconsider the general case of reduced rank cases whenrank[Σx ] = q , 1 ≤ q ≤ p ,

(6.1)
′

where Σx = BB and B is a p × q matrix.
Then the first example of Section 2 corresponds to the case when p = q = 1 and thesecond example of Section 3 corresponds to the case when q = 1 < p, p ≥ 2.
In the general case, however, there is a parametrization problem for the p × p matrix′
Σx , whose rank is q. We take a normalization as B Σ−1v B = (diag cii ). Then thealgebra of Section 3.2 can be extended by using the matrix formulae such that for apositive definite A we have′′

|A + BB | = |A||Iq + B A−1 B]|

(6.2)
and
′

′

′

[A + BB ]−1 = A−1 − A−1 B[Iq + B A−1 B]−1 B A−1

(6.3)
′

for Σx = BB .
Theorem 6.1 : Assume 1 ≤ q < p, p ≥ 2 and |Σv | ̸= 0.
(i) Assume that vi (i = 1, · · · , n) are a squence of i.i.d. random vectors and |Σv | ̸=
0. Under the assumptions of Gaussian distributions, the ML estimator of B isconsistent as n −→ ∞.
3

We are investigating this computational problem of the ML solutions in non-stationary multidimension cases.

17

(ii) Assume the non-stationary errors-in-variables model given by (2.1)-(2.3) and|Σv | ̸= 0. Under the assumption of existence of fourth order moments, the SIML
estimator of B is consistent as n −→ ∞.
Second, in some cases the second order (or higher order) diﬀerencing may beoften appropriate for economic time series. Since the likelihood function can becomplicated in general, we consider the ML estimation and the SIML estimation
when p ≥ 1 and d = 2, where(x)∆d xi = vi

(6.4)

,

(x) (x)′(x)

E[vi ] = 0, and E[vi vi ] = Σx .
(2)
We use the Kn −transformation that from the observation matrix Yn to Zn (=′(2)(zk )) by
(
)
(2)′
−2
Z(2)
n = (zk ) = Kn Yn − Ȳ0 , Kn = Pn Cn .

(6.5)

Then the separating information maximum likelihood (SIML) estimator of Σ̂x inthis case can be defined byn
1 ∑(2) (2)′
Σ̂x,SIM L =zk zk .
mn k=1m(6.6)

We prepare the next Lemma.
Lemma 6.2 : Let
Kn = (bij ) = Pn C−2
n .
(2)

(6.7)
′

Then for i, i = 1, · · · , mn , we have
(6.8)

[
]4
n
∑
π 2i − 1
1′(2) (2)) + O( ) .
bij bi′ ,j = δ(i, i ) 2 sin(
2 2n + 1nj=1

By using the above lemma, we have the next result, which is an extension of Kunitomo and Sato (2017) for the case of d ≥ 1. The result holds for an arbitrary integerd.
Theorem 6.3 : Assume p ≥ 1, d = 2 and mn /n −→ 0 as n −→ ∞.
(i) Assume that vi (i = 1, · · · , n) are a squence of i.i.d. random vectors and |Σv | ̸= 0.
Under the assumption of Gaussian distributions, the ML estimator of Σx is consistent as n −→ ∞.
(ii) Assume the non-stationary errors-in-variables model given by (2.1)-(2.3) and18

|Σv | ̸= 0. Under the assumption of existence of fourth order moments, the SIML
estimator of Σx is consistent as n −→ ∞.
It should be important to note that the diagonal elements akn (k = 1, · · · , n) should
be modified to
π
1
(2)
(6.9)
akn = [2 sin(k − )]4
2n + 1
2in the present case and we need the corresponding bias correction for estimating thevariance-covariance matrix Σv .
Remark 6.1 : The SIML should be a useful tool for the state space modeling of
non-stationary multivariate time series because it does not have any computational
problem and it has an asymptotic robustness. Also it is possible to show that theasymptotic normality of the SIML estimator under the general condition of (2.3)
with additional arguments ouitlined in Kunitomo and Sato (2017). For instance, theSIML method may give reasonable estimates not only for the coeﬃcients parameters,
but also the variance-covariance matrices when d = 1 and d = 2 by taking mn
appropriately.

7. Concluding Remarks
In this study, we have reported the asymptotic properties and finite sample properties of two estimation methods for the non-stationary errors-in-variables models.
We have compared the SIML estimation and the maximum likelihood (ML) estimation when there are non-stationary trends and noise components. We have foundthat the Gaussian likelihood function shows some peculiar shape and the computation of the ML estimation can be unstable in some cases. On the other hand,
the SIML estimation has the asymptotic robust properties under general conditionsof existence of second order moments. We have investigated the conditions for thegood asymptotic properties and the finite sample properties of the ML estimationand SIML estimation by both the theoretical analysis and simulations. We havefound that the SIML method gives reasonable estimates not only for the coeﬃcientsparameters, but also the variance-covariance matrices when d = 1 and d = 2 bytaking mn appropriately.
There are several possible extensions and discussions. First, it may be interesting to see to what extent the results reported in this paper are still valid when there
are seasonality and stationary components although the issue has been investigatedby Kunitomo and Sato (2017) in a simple setting. Second, there are several important issues on modeling the non-stationary trends and noises. Although we havediscussed some of them briefly in Section 6, but we need more investigations. The
determination of the number of non-stationary factors is one of examples, which
is discussed by Kunitomo (2017). Third, there can be many applications of themethods we have discussed as indicated in Kunitomo and Sato (2017).
19

The results of current investigations on these issues will be reported in otheroccasions.

References[1] Akaike, H. (1985), ”Prediction and entropy. A Celebration of Statistics,” A. C.
Atkinson and E. Fienberg, eds., Springer-Verlag, New York, 1-24.
[2] Akaike, H. (1989), ”Choosing priors and its applications,” in Bayesian Statistics
and Applications, edited by Suzuki and Kunitomo, Tokyo University Press, inJapanese.
[3] Amemiya, T. (1986), Advanced Econometrics, Blackwell.
[4] Anderson, T.W. (1971), Statistical Time Series Analysis, John-Wiley.
[5] Anderson, T.W. (1984), ”Estimating Linear Statistical Relationships,” Annalsof Statistics, 12, 1-45.
[6] Anderson, T.W. (2003), An Introduction to Statistical Multivariate Analysis,
3rd Edition, John-Wiley.
[7] Chang, Y., B. Jiang and J. Park (2008), ”Using Kalman Filter to Extract andTest for Common Stochastic Trends,” Unpublished Manuscript.
[8] Fuller, W. (1987), Measurement Error Models, John-Wiley.
[9] Johansen, S. (1995), Likelihood Based Inference in Cointegrated Vector Autoregressive Models, Oxford UP.
[10] Kitagawa, G. (2010), Introduction to Time Series Analysis, CRC Press.
[11] Kunitomo, N. (2017), “On Estimating the Number of Factors in Non-stationary
Errors-in-Variables Models,” In preparation.
[12] Kunitomo, N. and S. Sato (2013), “Separating Information Maximum Likelihood Estimation of Realized Volatility and Covariance with Micro-Market
Noise,” North American Journal of Economics and Finance (2013), 26, 282309.
[13] Kunitomo, N. (2017), “On Detecting the Number of Factors in Non-stationary
Errors-in-Variables Models,” Unpublished Manuscript.
[14] Kunitomo, N. and S. Sato (2017), ”Trend, Seasonality and Economic Time Series : the Nonstationary Errors-in-variables Models,” Discussion Paper CIRJEF-977, Graduate School of Economics, University of Tokyo,(http://www.e.utokyo.ac.jp/cirje/research/dp/2015), revised.
20[15] Ohtsu, K., H. Pend and G. Kitagawa, (2015), Time Series Modeling for Analysis and Control Advanced Autopilot and Monitoring Systems, Springer-Brief in
Statistics, Springer.

APPENDIX : Mathematical Derivations
In this Appendix, we give some details of the derivations omitted in the previous
sections.
Proof of Lemma 4.1 : When we have zk = (zik ) ∼ Np (0, Q), we can use therelationE[zik zjk zi′ k zj ′ k ] = qij qi′ j ′ + qii′ qjj ′ + qij ′ qi′ j ,(A.1)

where Q = (qij ). Then it is straightforward to obtain the results.
(Q.E.D.)
Proof of Theorem 4.2 : By using Lemma 4.1, it is possible to obtain the varianceof L∗13n , which converges to 0 as n −→ ∞. We set
1′−1bb ]Σ−1
Ak = a∗−1v
kn Σv [Σv − ∗akn + cand′Q0 = a∗kn Σv (θ0 ) + b(θ0 )b(θ0 ) .
Then−1
tr[Ak Q0 ] = tr(a∗−1kn Σv [Σv −1′a∗kn + c′

bb ]Σ−1v [akn Σv (θ0 ) + b(θ0 )b(θ0 ) ]) .

If Q = Q0 , then we havetr(Ak Q0 ) = tr(Ip ) = p .
Also we find
(A.2)

= E[
= E[
= [

n
∑

′

trAk (zk zk − Q0 )]2

k=1
n
∑

k,k′ =1n∑′′

E[(zk Ak zk )2 − (trAk Q0 )2 ]

k=1

=

n
∑

′

trAk (zk zk − Q0 )Ak′ (zk′ zk′ − Q0 )]2tr(Ak Q0 Ak Q0 ) .

k=1

21

If Q = Q0 , then it is 2pn.
Next, we investigate the last two terms of L∗1n in details. After simple abgebra, they
can be reexpressed as
1 ∑ ∗−1′
akn b(θ0 )Σ−1v b(θ0 )n k=1nL∗14 = trΣ−1v Σ(θ0 ) +

1∑
1
′
′
−1
−1
∗−1 −1−b [Σ−1v Σ(θ)Σv + akn Σv b(θ0 )b(θ0 ) Σv ]b∗n k=1 akn + cn= trΣ−1v Σ(θ0 )
n
[
1∑
1
′
−1
∗−1 ′ −1+(a∗kn + c)[b(θ0 ) Σ−1v b(θ0 ) − akn b Σv Σv (θ0 )Σv b∗∗
n k=1 akn (akn + c)
]
′
2
b(θ))
−(b Σ−1v= trΣ−1v Σ(θ0 )
n
[ ′
1∑
1
−1
+
b Σ−1
v (Σv − Σv (θ0 ))Σv b(θ0 )n k=1 a∗kn (a∗kn + c)
]
′
′
′
−1
−1
−1
2
+b(θ0 ) Σv b(θ0 )b Σv b − (b Σv b(θ0 )) .
The last term is non-negative because of the Cauchy-Schwarz’s inequality and its
minimum occurs at b = b(θ0 ) because Σv is positive definite. Then we need toevaluate the sum of non-degenerate terms as
(A.3)

n
n
1∑
1∑
′
−1
∗
log |akn Σv | +log[1 + a∗−1kn b Σv b]
n k=1
n k=1

and
(A.4)

1
trΣ−1
v Σv (θ0 ) +′n∑b Σ−1 (Σv − Σv (θ0 ))Σ−1 b

v

n k=1

a∗kn + b′ Σ−1
v b

v

.,

By using the inversion formula and its determinant in (3.1) and (3.2), the sum ofabove two terms can be written as
1∑′
log |a∗kn Σv + b b|
n k=1
n

1∑
′
′
tr[a∗kn Σv + b b]−1 [a∗kn Σv (θ0 ) + b b] .
n k=1
n

+

By using Lemma 3.2.2 of Anderson (2003), the minimum of each terms occurs at′′a∗kn Σv + b(θ0 )b(θ0 ) = a∗kn Σv (θ0 ) + b(θ0 )b(θ0 ) ,22that is, Σv = Σv (θ0 ). Hence the global minimum of the likelihood function occurs
iﬀ b = b(θ0 ) and Σv = Σv (θ0 ).
The rest of our arguments for the consistency follows from the general arguments(see Theorem 4.1.1 of Amemiya (1985), for instance) and we have the result.
(Q.E.D.)
Proof of Lemma 6.2 : We set
2π
11(j − )(k − ) ,2n + 1222π
1=(k − ) ,2n + 1
2

(A.5)

θkj =

(A.6)

θk

and
bkj = √
(2)

(A.7)

1
[(1 − eiθk )2 eiθkj + (1 − e−iθk )2 e−iθkj ]2n + 1for j, k = 1, · · · , n. Thenn−iθk∑
1 + eiθk(2)−iθk 4 1 + e
+
(1
−
e
)
+ 2n(1 − eiθk )2 (1 − e−iθk )2(2n + 1)[bkj ]2 = (1 − eiθk )4
i2θk
−i2θk
1
−
e
1
−
e
j=1
θk

θk

= (1 − eiθk )3 + (1 − e−iθk )3 + 2n(ei 2 1 − e−i 2 )4θk
= 2n[4 sin2 ]2 + o(1) .
2′

Also for k ̸= k we find thatn∑(2) (2)(2n + 1)[bkj bk′ j ]j=1= (1 − e ) (1 − e
iθk 2

n
n
∑
∑
i(θkj +θ ′ )i(θ +θ )
−iθ
iθ
2
2
′
k
k j + (1 − e)e) (1 − e k )e kj k′ j−iθk 2

j=1

j=1

n
n
∑
∑
i(θkj −θ ′ )−i(θ +θ )iθ ′ 2−iθ
−iθk 2−iθ22′kk j + (1 − e+(1 − e) (1 − e k )e) (1 − e k )e kj k′ j
j=1

j=1= (I) + (II) + (III) + (IV ) (say) .
Then after some algebra, we have
θk

θk

θ ′
k

θ ′
k

θk

θ ′
−i k2

θ ′
i k2

θk

θk

(I) + (IV ) = (−1)[(e−i 2 − ei 2 )2 (e−i 2 − ei 2 )2 ]eiθk(II) + (III) = (−1)[(e−i 2 − ei 2 )2 (e−eθk +θ ′k2,)2 ][1 − (e−i(θk −θ ′ )
k
2

+ (ei

(θk −θ ′ )
k
2

)] .

Hence we finnally find that
θ ′
k

θ ′
k

(I)+(II)+(III)+(IV ) = (−1)(e−i 2 −ei 2 )2 (e−i 2 −ei 2 )2 [1−ei23θk −θ ′k2−e−i

θk −θ ′k2+e−i

θk −θ ′
k2]

and then it is proportional to(I) + (II) + (III) + (IV ) = 4[sin2when k/n −→ 0 as n → ∞.
(Q.E.D.)24θkθ′][sin2 k ] + o(1)22