N4SID-VAR Method for Multivariable Discrete Linear Time-variant
System Identification
Alexander E. Robles and Mateus Giesbrecht
School of Electrical and Computer Engineering, University of Campinas, Av. Albert Einstein 400, Campinas - SP, Brazil
Keywords:
System Identification, Subspace Methods, Time-variant System Identification.
Abstract:
In this paper, a method for multivariable discrete linear time-variant system identification is presented. This
work is focused on slowly multivariable time-variant systems, so that it is possible to define time intervals,
defined as windows, in which the system can be approximated by time-invariant models. In each window,
a variation of N4SID that uses Markov parameters is applied and a state space model is estimated. For that
reason the proposed method is defined as N4SID-VAR. After obtaining the models for all windows, the error
between system model outputs are calculated and compared to the system outputs. The N4SID-VAR was
tested with a time-variant multivariable benchmark and the results were accurate. The proposed method was
also compared to the MOESP-VAR method and, for the tested benchmark, the N4SID-VAR was faster and
more accurate than the MOESP-VAR algorithm.
1 INTRODUCTION
System identification consists in the search for a
mathematical model that can describe the behavior of
a dynamical system, from the observed input-output
signals (Ljung, 1999), (Katayama, 2005). A signifi-
cant part of activities and researches in system identi-
fication focuses on time-invariant dynamical systems.
However, there are innumerable systems in nature that
are multivariable with nonlinear and time-variant be-
havior. To deal with last problem, the time-variant
systems can be approximated by linear time-invariant
systems, as long as these systems vary slowly (Tama-
riz et al., 2005).
During the last two decades, subspace-based
methods have been extensively studied to address the
problem of identifying multivariable discrete linear
time-invariant systems (Katayama, 2005). From that
methods, the most popular are the MOESP (Verhae-
gen and Dewilde, 1992) and the N4SID (Overschee
and Moor., 1994). Both methods have a mathematical
support in linear matrix algebra.
The MOESP method is based on LQ decomposi-
tion of a matrix formed by input-output data, where:
L is a lower triangular matrix and Q is an orthogo-
nal matrix. From a block of the matrix L a singular
value decomposition (SVD) is performed, from that
it is possible to find out the system order and its ob-
servability matrix. With this last matrix it is possible
to obtain the matrices C and A corresponding to the
model in state space. The final step is to form a linear
equation and apply the least squares method and esti-
mate the matrices B and D of the model. The method
is detailed in the section 4 of this paper.
Another subspace method called N4SID (Numer-
ical Algorithms for Subspace State Space System
Identification), the same way as the MOESP, is based
on a LQ decomposition of data matrices. However,
in this case this decomposition is interpreted as the
oblique projection of future outputs in the subspace
of the past inputs and outputs, towards the future in-
puts. From these projections the system states are
estimated. With the states, inputs and outputs, the
matrices A, B,C and D can be determined using a
simple least squares method. This method and a
variant proposed in (Clavijo, 2008) are presented in
the section 5. A method inspired by MOESP and
called MOESP-VAR, was introduced and developed
in (Tamariz et al., 2005). The method is initialized
by splitting the total input and output data into data
groups that are associated with time intervals in which
the system exhibits a slow change and can be ap-
proximated by a time-invariant system. The MOESP
method is applied to the data of each interval, result-
ing in a linear time-invariant for the system in each of
the time intervals. Following the same concept, the
N4SID-VAR method is proposed in this article. The
first step of the proposed method is to split the data
502
Robles, A. and Giesbrecht, M.
N4SID-VAR Method for Multivariable Discrete Linear Time-variant System Identification.
DOI: 10.5220/0006907505020509
In Proceedings of the 15th International Conference on Informatics in Control, Automation and Robotics (ICINCO 2018) - Volume 1, pages 502-509
ISBN: 978-989-758-321-6
Copyright © 2018 by SCITEPRESS – Science and Technology Publications, Lda. All rights reserved
into intervals, defined as time windows and then, to
each window data a variation of the N4SID developed
in (Clavijo, 2008) is applied.
In addition to the development of the algorithm, in
this paper, a comparison between the MOESP-VAR
and the N4SID-VAR is made using the same time
window sizes and a slow time-variant benchmark.
The computational time to process the algorithms is
also evaluated.
The article is divided into eight sections. In the
next section the state space system identification is
presented. In the third section, the procedure to deter-
mine an extended state space model, which serves as
the basis for the MOESP and N4SID methods, is de-
tailed. In the following section the MOESP method is
presented. The next section is focused on the variant
of the N4SID method that uses the Markov parame-
ters. In the sixth section the method proposed in this
article is developed. In the seventh section, the re-
sults of the N4SID-VAR and MOESP-VAR methods
to identify a time-variant benchmark are presented.
Finally the last section, the conclusions of this article
are exposed.
2 STATE SPACE SYSTEM
IDENTIFICATION
The identification of discrete multivariable linear
time-invariant systems using subspace methods al-
lows to determine a causal and time invariant model,
estimated only from the system inputs and outputs.
The main advantage of this approach is that a multi-
variable system can be modeled, without any a priori
assumption about the system order (state vector di-
mension). A discrete time-invariant model can be de-
scribed by the following equation in the state space:
x(k +1) = Ax(k) + Bu(k)
y(k) = Cx(k) +Du(k)
(1)
where x(k) R
n
is defined as the system state at in-
stant k, A R
n×n
is the state transition or system
matrix, B R
n×m
is a matrix that relates the input
u(k) R
m
to the state, C R
l×n
is the matrix that re-
lates the output y(k) R
l
to the state and D R
l×m
is the matrix that relates the outputs to the inputs
(Katayama, 2005).
3 EXTENDED STATE SPACE
MODEL
In the same way that time-invariant systems can be
represented by the equation (1), there are also other
ways to represent the relation between the input, out-
put and state vectors. In this section an extended
model that is useful for the subspace methods is pre-
sented.
For a time instant t, it is defined that inputs before
that instant are null. With this, from the equation (1),
it is possible to substitute the relation between inputs,
outputs and states between the instants t and t +k 1,
where k is an integer, in the following way:
y
t|k1
=
C
CA
.
.
.
CA
k1
x(t)+ (2)
+
D 0 0 0
CB D 0 0
.
.
.
.
.
.
.
.
.
.
.
.
CA
K2
B ··· CB D
u
t|k1
where
y
t|k1
=
y(t)
y(t + 1)
.
.
.
y(t + k 1)
(3)
u
t|k1
=
u(t)
u(t + 1)
.
.
.
u(t + k 1)
(4)
the dimensions of the concatenated output vectors and
inputs presented in equations 3 and 4, are respec-
tively: y
t|k1
R
kp×1
and u
t|k1
R
km×1
.
Rewriting the equation 2, can be obtain the fol-
lowing relation between the data matrices:
y
t|k1
= O
k
x(t) + Ψ
k
u
t|k1
(5)
where O
k
is a observability matrix, Ψ
k
R
kp×km
is
the Toeplitz matrix, as detailed in equations 6 and 7.
O
k
=
C
CA
.
.
.
CA
k1
(6)
N4SID-VAR Method for Multivariable Discrete Linear Time-variant System Identification
503
Ψ
k
=
D 0 0 0
CB D 0 0
.
.
.
.
.
.
.
.
.
.
.
.
CA
K2
B ··· CB D
(7)
If the vectors u
t|k1
, y
t|k1
and x(t) for t =
0. .. N 1 are concatenated side by side, the following
matrices can be defined:
U
0|k1
= [ u
0|k1
u
1|k
.. . u
N1|k+N2
] R
km×N
(8)
Y
0|k1
= [ y
0|k1
y
1|k
.. . y
N1|k+N2
] R
kp×N
(9)
X
N1
= [ x(0) x(1) .. . x(N 1) ] R
n×N
(10)
where X
N1
is a state matrix. With the concatenated
matrices the following extended model is written:
Y
0|k1
= O
k
X
N1
+ Ψ
k
U
0|k1
(11)
from the extended model the MOESP method is de-
veloped in the following section.
4 MOESP METHOD
The MOESP (Verhaegen and Dewilde, 1992) method
is based on the LQ decomposition of a matrix formed
by input and output data into two matrices: a matrix
L, which is lower triangular, and a matrix Q, which is
formed by linearly independent columns. The input
and output data are concatenated in the Hankel matri-
ces, presented in the equations 8 and 9.
Then an extended state space model is formed
with the input and output Hankel matrices (11). With
the extended model, U
0|k1
and Y
0|k1
, is applied the
LQ decomposition.
U
0|k1
Y
0|k1
=
L
11
0
L
21
L
22
Q
T
1
Q
T
2
(12)
where L
11
R
km×km
e L
22
R
kp×kp
are lower trian-
gular matrices, Q
T
1
R
km×N
e Q
T
2
R
kp×N
are or-
thogonal and L
21
R
kp×km
.
From this decomposition, the equation (11) can be
rewritten as shown below:
L
21
Q
T
1
+ L
22
Q
T
2
= O
k
X
N1
+ Ψ
k
L
11
Q
T
1
(13)
Post-multiplying 13 by Q
2
yields
L
22
= O
k
X
N1
Q
2
(14)
where Q
T
1
Q
2
= 0, Q
T
2
Q
2
= I
kp
, due to the orthonor-
mality between vectors that form these matrices. The
observability matrix O
k
can be obtained and the di-
mension of the system n after applying a SVD to
L
22
R
kp×kp
.
Let SVD of L
22
be given by
L
22
=
U
1
U
2
Σ
1
0
0 0
V
T
1
V
T
2
(15)
where U
1
R
kp×n
and U
2
R
kp×(kpn)
. Then, from
(15) and (14) it is possible to write:
O
k
X
N1
Q
2
= U
1
Σ
1
V
T
1
(16)
so the extended observability matrix can be defined as
follow
O
k
= U
1
Σ
1/2
1
(17)
Define O
k
as a matrix O
k
with an offset of one
row block, it is possible to write the following rela-
tion:
O
k
=
CA
CA
2
.
.
.
=
C
CA
.
.
.
A = O
k
A (18)
After applying the pseudoinverse of O
k
on both
sides of the equation 18 , the matrix A is estimated.
A = O
k
O
k
(19)
the matrix C is readily given by
C = O
k
(1 : p, 1 : n) (20)
The matrices B and D can be estimated using the
following relationships: First taking advantage of the
orthogonality between U
1
and U
2
.
U
T
2
L
22
= U
T
2
U
1
Σ
1
V
T
1
= 0 (21)
U
T
2
O
k
= U
T
2
U
1
Σ
1/2
1
= 0
Next step, multiplying both sides of equation (13)
by U
T
2
the following relationship is found:
U
T
2
L
21
L
1
11
= U
T
2
Ψ
k
(22)
Splitting U
T
2
in blocks with l columns defined as
L
i
and splitting the matrix U
T
2
L
21
L
1
11
in blocks with m
columns defined as M
i
, so the following relationship
is valid.
M
1
M
2
.. . M
k
=
L
1
L
2
.. . L
k
Ψ
k
(23)
ICINCO 2018 - 15th International Conference on Informatics in Control, Automation and Robotics
504
Afterwards, replacing Ψ
k
, matrix defined in the
equation 7, a linear relation 24 can be found. In
fact 24 is a linear equation that has as variables the
matrices B and D. Define
˜
L
i
=
L
i
.. . L
k
R
(kpn)×(k+1i)p
,i = 2, .. .,k, and substituting in
equation 23 the following overdetermined linear sys-
tem is obtained.
L
1
˜
L
2
O
k1
L
2
˜
L
3
O
k2
.
.
.
.
.
.
L
k1
˜
L
k
O
1
L
k
0
D
B
=
M
1
M
2
.
.
.
M
k1
M
k
(24)
A solution for the overdetermined system 24 can
be found with the least squares method, resulting on
estimates for the matrices B and D.
5 N4SID-VAR USING MARKOV
PARAMETERS
In this section, another method to identify multivari-
able discrete time-invariant systems using subspace
methods is presented. This method is called N4SID
and was developed by Van Overschee and De Moor
(Overschee and Moor., 1994). The algorithm is ini-
tialized by calculating the oblique projection of future
outputs Y
f
, on the vector W
p
which is the concatena-
tion of the inputs U
p
and outputs Y
p
passed in the di-
rection of the future inputs U
f
(Katayama, 2005). The
classical N4SID estimates the system model using the
least squares method, but in this work is considered a
variation that was introduced and developed in (Clav-
ijo, 2008). By definition:
W
P
=
U
P
Y
P
(25)
and by the extended state space model presented in
section 3.
Y
f
= O
k
X
f
+ Ψ
k
U
f
(26)
Consider dthe following LQ decomposition:
U
f
U
p
Y
p
Y
f
=
L
11
0 0 0
L
21
L
22
0 0
L
31
L
32
L
33
0
L
41
L
42
L
43
L
44
Q
T
1
Q
T
2
Q
T
3
Q
T
4
(27)
where L
44
= 0.
Then it is possible to rewrite the equation 27 as
follows.
U
f
W
p
Y
f
=
R
11
0 0
R
21
R
22
0
R
31
R
32
0
Q
T
1
Q
T
2
Q
T
3
(28)
the relations between the matrices presented in equa-
tion 27 and 28, are as follows:
R
11
= L
11
R
21
=
L
21
L
31
R
22
=
L
22
0
L
32
L
33
R
31
= L
41
R
32
=
L
42
L
43
Q
T
1
= Q
T
1
Q
T
2
=
Q
T
2
Q
T
3
Q
T
3
= Q
T
4
W
p
=
U
p
Y
p
From 28
U
f
= R
11
Q
T
1
(29)
Q
T
1
= R
1
11
U
f
The matrix Q
T
2
can be written as:
W
p
= R
21
Q
T
1
+ R
22
Q
T
2
(30)
R
22
Q
T
2
= W
p
R
21
Q
T
1
Q
T
2
= R
22
W
p
R
21
Q
T
1
the matrix of future outputs Y
f
, can be rewritten from
28
Y
f
= R
31
Q
T
1
+ R
32
Q
T
2
(31)
Using the equations 29 and 30 in 31, the following
relation is given for Y
f
.
Y
f
= R
31
R
1
11
U
f
+ R
32
R
22
W
p
R
21
Q
T
1
(32)
= R
32
R
22
W
p
+
R
31
R
32
R
22
R
21
R
1
11
U
f
The extended state space model 11, considering
only future data, is given by:
Y
f
= O
k
X
f
+ Ψ
k
U
f
(33)
Comparing 32 and 33, can be obtained two impor-
tant relations to N4SID
R
32
R
22
W
p
= O
k
X
f
(34)
Ψ
k
=
R
31
R
32
R
22
R
21
R
1
11
(35)
N4SID-VAR Method for Multivariable Discrete Linear Time-variant System Identification
505
applying the SVD to the equation 34
R
32
R
22
W
p
=
U
1
U
2
Σ
1
0
0 0
V
T
1
V
T
2
(36)
where U
1
R
kp×n
and U
2
R
kp×(kpn)
. So from 34
and 36 we have:
O
k
= U
1
Σ
1/2
1
(37)
We also have the Toeplitz matrix Ψ
k
R
kp×km
of
the equation 35,
Ψ
k
=
R
31
R
32
R
22
R
21
R
1
11
In the original N4SID method the future states ma-
trix X
f
is also obtained from the equations 34 and 36
and with the inputs, outputs and states, the matrices
A, B, C and D from the state space model are esti-
mated using the least squares method. Alternatively,
the matrices can be estimated from the Toeplitz ma-
trix Ψ
k
, as proposed by (Clavijo, 2008) and detailed
in the sequence.
The first column block of the matrix defined above
represents the impulse responses, also called Markov
parameters, given by:
G
(k)
:=
D k = 0
CA
k1
B k 6= 0
D
CB
.
.
.
CA
k1
B
=
G
0
G
1
.
.
.
G
k
(38)
With these impulse responses, the following Han-
kel matrix can be formed
H
k
= O
k
C
k
=
CB CAB CA
2
B . ..
CAB CA
2
B CA
3
B . ..
.
.
.
.
.
.
.
.
.
.
.
.
(39)
The reachability matrix C
k
R
nXkm
can be esti-
mated using the observability matrix and the equation
39
C
k
= O
k
H
k
=
B AB A
2
B .. . A
k1
B
(40)
Finally, the matrix B is the first n × m block from
C
k
and the matrix D is G
0
B = C
k
(1 : n,1 : m) (41)
D = G
0
From the observability matrix 37, the matrices A
and C are estimated, the same way as is done in the
MOESP method, according to the presented in the
equations 19 and 20.
6 N4SID-VAR METHOD
In real life systems are not time-invariant, increasing
the complexity of the identification problem. The sys-
tem can be modeled in the state space, represented in
equation 42, where the state space system matrices
A
(k)
,B
(k)
,C
(k)
,D
(k)
change over the time k.
x(k +1) = A
(k)
x(k) +B
(k)
u(k)
y(k) = C
(k)
x(k) +D
(k)
u(k)
(42)
There is a version of MOESP for identify slowly
time-variant systems, called MOESP-VAR (Tamariz
et al., 2005). The principle of the method is as fol-
lows: given a time-variant system, time intervals, also
called windows, are defined with a quantity of data
(inputs and outputs). These windows are set so that
the system does not undergo significant changes dur-
ing each window. The MOESP method is applied to
the data of each window, and with this the state space
model matrices A,B,C and D that represent the sys-
tem in that window are estimated. The process is re-
peated until all the data windows are modeled.
The study of subspace methods led to the develop-
ment of a version of N4SID for the problem of identi-
fication of multivariable linear time-variant systems,
named in this work as N4SID-VAR. The method
works as follows: The first step is to define intervals
or time windows (where system variations are slow).
Each window contains a subset of input and output
data. The next step is to apply the N4SID to each
of the windows. With the algorithm the quadruple of
matrices A
w j
,B
w j
,C
w j
,D
w j
is estimated for each win-
dow.
The total number of windows j can be determined
by dividing the total number of available data (inputs
and outputs) N and the number of data per each time
window N
w
, as follows in the equation.
j =
N
N
w
(43)
Each quadruple represents the state space model
within a time interval, as defined in the equation 44.
A
w j
= A
(k)
0 k jN
w
1 (44)
The subsets of input data (U
w1
,. ..,U
w j
) and out-
puts (Y
w1
,. ..,Y
w j
) of the system, in each of the win-
dows are defined as:
U
w1
= [u(0) u(1) .. . u(N
w
1)] (45)
U
w2
= [u(N
w
) u(N
w
+ 1) .. . u(2N
w
1)]
.
.
. =
.
.
.
U
w j
= [u( jN
w
N
w
) u( jN
w
N
w
+ 1) .. .
ICINCO 2018 - 15th International Conference on Informatics in Control, Automation and Robotics
506
.. . u( jN
w
1)]
Y
w1
= [y(0) y(1) .. . y(N
w
1)] (46)
Y
w2
= [y(N
w
) y(N
w
+ 1) .. . y(2N
w
1)]
.
.
. =
.
.
.
Y
w j
= [y( jN
w
N
w
) y( jN
w
N
w
+ 1) .. .
.. . y( jN
w
1)]
The extended models for outputs and inputs (past
and future) in each window are given by 47 and 48:
Y
w j|p
= O
k
X
w j|p
+ Ψ
k
U
w j|p
(47)
Y
w j| f
= O
k
X
w j| f
+ Ψ
k
U
w j| f
(48)
From these extended models it is possible to start
the N4SID-VAR algorithm using the variation pre-
sented in the previous section. The flowchart of the
algorithm is presented in the figure 1.
Figure 1: Steps of N4SID-VAR algorithm.
7 RESULTS
In order to test the proposed N4SID-VAR algorithm,
a benchmark from the reference (Tamariz et al., 2005)
was selected. In this benchmark, the matrix A is
slowly time-variant.
A
k
=
a(k) b(k)
1 1
(49)
where
a(k) =
1
2
(k/2500)
2
+
1
2
(k/2500)
b(k) =
1
16
(k/2500)
4
1
8
(k/2500)
3
13
16
(k/2500)
2
3
4
(k/2500)
1
2
The matrices B
k
,C
k
and D
k
are constant over the
time
B
k
=
2 1
1 1
; C
k
=
1 3
1 2
D
k
=
1 3
1 1
(50)
The system was excited by a two-dimensional
white noise with N = 1000 samples. In addition,
a noise equivalent to the 30% of the output gener-
ated after excitation of the system with the input was
added. This was done to verify the robustness of the
algorithms when applied to noisy data. The variation
of the parameters a(k) and b(k) is plotted in the figure
2.
Figure 2: Variation of A
k
elements.
At first a size of each data window N
w
= 50 (in-
puts and outputs) was defined. In the time interval for
each window the system varies slowly, and can be ap-
proximated by a time-invariant model. The MOESP-
VAR and N4SID-VAR, were executed 2000 times.
N4SID-VAR Method for Multivariable Discrete Linear Time-variant System Identification
507
For each of the runs, a new set of white noise in-
puts was generated, always following the same mean
and covariance characteristics. This is done to en-
sure greater consistency in the comparison between
the two methods. In the figures 3 and 4 the system
and model outputs are presented for each of the meth-
ods for one of the executions.
Figure 3: In the solid red line the system outputs are shown
and in blue stars the outputs of the model obtained with the
MOESP-VAR are shown for one of the executions with a
window of N
w
= 50 data.
Figure 4: In the solid red line the system outputs are shown
and in blue stars the outputs of the model obtained with the
proposed method N4SID-VAR are shown, for one of the
executions with window of N
w
= 50 data.
To compare the quality of the methods, the mean
square error was calculated between the actual and es-
timated outputs, which is defined by the equation 51,
the variable m is the number of runs of each algo-
rithm. The mean of the results of the mean square
error calculation for the two methods after 2000 runs
are shown in the table 1.
e =
m
r=1
1
2N
N
i=1
(y
r(i)
ˆy
r(i)
)
2
m
(51)
In the figure 5 the mean squared er-
rors are shown after 2000 runs, for the
following data window sizes: N
w
=
Table 1: Mean square error.
N
w
= 50 e
N4SID-VAR 1.1521
MOESP-VAR 2.8896
Figure 5: Mean squared errors for the N4SID-VAR and
MOESP-VAR methods, for windows varying between 50
and 500.
Figure 6: Detail of mean square errors for windows ranging
from 50 to 100.
50,60, 70,80,90,100,150,200,250,300,350, 400,500.
Another graph showing more about the performance
of N4SID-VAR is figure 6, where windows in the
range of N
w
= 50,60,70,80, 90,100, are presented in
detail.
The mean square error obtained, presented in the
table 1, with the N4SID-VAR was lower than with
the MOESP-VAR. For all window sizes evaluated, the
proposed method presents the smallest mean square
error. One intuition about this is that the algorithm
is based on a version of N4SID that does not use
the method of least squares, one only has to find the
Toeplitz matrix formed with the Markov parameters.
In the case of MOESP, the least squares method is part
of the algorithm to estimate the matrix B and D of the
model.
ICINCO 2018 - 15th International Conference on Informatics in Control, Automation and Robotics
508
Thus, the small number of data in a window
causes the MOESP algorithm to have a less accu-
rate result, an important observation is that, for both
algorithms, the error tends to increase with window
size, because, in this way, the system or benchmark
presents larger variations within the same window and
its representation by a system invariant in time be-
comes less exact.
The average execution time T
m
for all the execu-
tions for the window size 50 was calculated. The
N4SID-VAR method was 5 times faster than MOESP-
VAR, in other words the time compilation was 5 times
smaller if compared to the MOESP-VAR execution
time. The results are displayed in the table 2. This re-
sult is also expected due to the lower need for matrix
inversions in the N4SID method.
Table 2: Average runtime.
N
w
= 50 T
m
(segundos)
N4SID-VAR 2.87
MOESP-VAR 14.35
8 CONCLUSIONS
From the results presented in this article it is con-
cluded that the deterministic method N4SID-VAR has
a good performance based on the results, a smaller
quadratic error in comparison to MOESP-VAR. The
algorithm also presents a shorter execution time be-
cause it uses the impulse responses (Markov parame-
ters) in the estimation of the matrices B and D and not
the least squares method.
For future work, it is possible to study the accu-
racy of the proposed algorithm for higher order sys-
tems. In addition, other evolutionary heuristic tech-
niques to identify time-variant systems developed in
(Giesbrecht and Bottura, 2015) and (Robles and Gies-
brecht, 2017) can be compared with the proposed
method.
REFERENCES
Clavijo, D. G. (2008). Metodos de Subespaos para Identifi-
cao de Sistemas:Propostas de Alteraes, Implementaes
e Avaliaes. Tese Mestrado- UNICAMP.
Giesbrecht, M. and Bottura, C. P. (2015). Recursive
immuno-inspired algorithm for time variant discrete
multivariable dynamic system state space identifica-
tion. International Journal of Natural Computing Re-
search,5(2), April-June 69-100.
Katayama, T. (2005). Subspace Methods for System Identi-
fication. Springer.
Ljung, L. (1999). System Identification-Theory for the User.
Prentice-Hall, Upper Saddle River, N.J., 2nd edition
edition.
Overschee, P. V. and Moor., B. D. (1994). N4SID:
Supspace Algorithms for the Identification of Com-
bined Deterministic-Stochastic Systems. Automatica,
30(1):75-93.
Robles, A. E. and Giesbrecht, M. (2017). Mtodo co-
evolutivo para identificao de sistemas variantes no
tempo. XIII Simpsio Brasileiro de Automao In-
teligente (SBAI), Porto Alegre.
Tamariz, A. D. R., Bottura, C. P., and Barreto, G. (2005).
Iterative MOESP Type Algorithm for Discrete Time
Variant System Identification. Proceedings of the 13th
Mediterranean Conference on Control and Automa-
tion.
Verhaegen, M. and Dewilde, P. (1992). Subspace model
identification - part 1 : The output-error state-space
model identification class of algorithms. International
Journal of Control 56, 5, 11871210.
N4SID-VAR Method for Multivariable Discrete Linear Time-variant System Identification
509