IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

1

Continuous-time DC kernel —
a stable generalized first-order spline kernel

arXiv:1710.03932v1 [] 11 Oct 2017

Tianshi Chen, Member, IEEE,

Abstract
The stable spline (SS) kernel and the diagonal correlated (DC) kernel are two kernels that have been
applied and studied extensively for kernel-based regularized LTI system identification. In this note, we
show that similar to the derivation of the SS kernel, the continuous-time DC kernel can be derived by
applying the same “stable” coordinate change to a “generalized” first-order spline kernel, and thus can
be interpreted as a stable generalized first-order spline kernel. This interpretation provides new facets to
understand the properties of the DC kernel. In particular, we derive a new orthonormal basis expansion
of the DC kernel, and the explicit expression of the norm of the RKHS associated with the DC kernel.
Moreover, for the non-uniformly sampled DC kernel, we derive its maximum entropy property and show
that its kernel matrix has tridiagonal inverse.

Index Terms
System identification, regularization methods, kernel, orthonormal basis expansion, MaxEnt, tridiagonal inverse.

I. INTRODUCTION
Linear time invariant (LTI) system identification is a classical topic in system identification. The current
standard solution to this topic is the maximum likelihood/prediction error method(ML/PEM), see e.g.,
T. Chen is with the School of Science and Engineering, The Chinese University of Hong Kong, Shenzhen, 518172, Shenzhen,
China. e-mail: tschen@cuhk.edu.cn. He is also holding a visiting position at the Department of Electrical Engineering,
Linköping University, Linköping, Sweden.
This work was supported by the Thousand Youth Talents Plan funded by the central government of China, the general
project funded by NSFC under contract No. 61773329, the Shenzhen Research Projects Ji-20170189 and Ji-20160207 funded by
Shenzhen Science and Technology Innovation Council, the President’s grant under contract No. PF. 01.000249 and the start-up
grant under contract No. 2014.0003.23 funded by the Chinese University of Hong Kong, Shenzhen, as well as by a research
grant for junior researchers under contract No. 2014-5894, funded by Swedish Research Council.

October 12, 2017

DRAFT

2

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

[1]. A new solution is the kernel-based regularization method that is first proposed in [2] and further
studied in [3], [4], [5], see [6] for a survey of this method. Recent progress for this method include, e.g.,
the kernel design [7], [8], [9], analysis of the hyper-parameter estimators [10], [11], input design [12],
[13] and development of the dual theory in frequency domain [14].
This method uses the impulse response model and solves a regularized least squares problem with a
suitably designed and tuned kernel. The kernel plays a similar role as the parametric model structure
in ML/PEM: on the one hand, it decides in what space the estimated impulse response is searched and
on the other hand, it encodes the prior knowledge regarding the underlying system to be identified. The
first two kernels for this method are the stable spline (SS) kernel [2] and the diagonal correlated (DC)
kernel [4]. On the one hand, these two kernels are derived in different ways. The SS kernel is obtained
by applying a “stable” coordinate change to the second-order spline kernel [cf. (10)], and the DC kernel
is obtained by mimicing the behavior of the optimal kernel [4, eq. (56)] for this method. On the other
hand, these two kernels share some common features [15], [7], which are fundamental for developing
systematic methods to design kernels for this method. In particular, we have shown in [7] that both kernels
belong to the class of amplitude modulated locally stationary kernels, and simulation-induced kernels,
leading to a machine learning perspective and a system theory perspective to design more general kernels,
respectively.
The SS kernel still has some other features inherited from its mother kernel – the spline kernel [cf.
(8)], which were previously regarded as unique for the SS kernel. For example, the orthonormal basis
expansion of the SS kernel with respect to a suitably chosen measure can be simply derived by applying
the “stable” coordinate change to that of the second-order spline kernel [2], which is key for developing
efficient implementation algorithms for this method [16]. Moreover, the same technique applies when
deriving the norm of the reproducing kernel Hilbert space (RKHS) induced by the SS kernel and the
maximum entropy property of the SS kernel [17]. Interestingly, as will be shown shortly, these features
in fact also hold for the continuous-time DC kernel. The key lies in to show that the DC kernel can
actually be derived by applying the same “stable” coordinate change to a “generalized” first-order spline
kernel. The DC kernel can thus be interpreted as a stable generalized first-order spline kernel. This
interpretation provides new facets to understand the properties of the DC kernel. In particular, we derive
a new orthonormal basis expansion of the DC kernel and the explicit expression of the norm of the
RKHS associated with the DC kernel. Moreover, for the non-uniformly sampled DC kernel, we derive
its maximum entropy property and show that its kernel matrix has tridiagonal inverse, which extend the
corresponding results in [17], [18].

DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

II. S YSTEM

IDENTIFICATION WITH

3

K ERNEL - BASED R EGULARIZATION M ETHOD

A. Problem Statement
We consider continuous-time LTI stable and causal systems, which are described by
y(t) = (g ∗ u)(t) + v(t),

(1)

where t ≥ 0 is the time index, y(t), u(t), v(t) ∈ R are the measured output, input and disturbance of
the system at time t, respectively, g(t) is the impulse response of the LTI system, and (g ∗ u)(t) is the
convolution of the impulse response g(t) with the input u(t) (evaluated at t). Since g(t) = 0 for t < 0
due to the causality assumption, the convolution (g ∗ u)(t) takes the form of
Z ∞
g(s)u(t − s)ds,
(g ∗ u)(t) =

(2)

s=0

where the unknown input u(t) with t < 0 is set to zero. The problem is to estimate g(t) as well as
possible based on the measured data {y(t)}N
t=1 and u(t) with t ≥ 0.
B. Kernel-based Regularization Method
To estimate the impulse response g(t) from {y(t)}N
t=1 and u(t) with t ≥ 0 is an ill-conditioned
problem. To overcome this difficulty, the kernel-based regularization method first introduces a positive
semidefinite kernel k(t, s; θ) with t, s ≥ 0 and constrain the search for a suitable impulse response within
the reproducing kernel Hilbert space (RKHS) Hk induced by k(t, s; θ)1 , where θ is a hyper-parameter
vector that contains the parameters used to parameterize the kernel. In particular, the regularized least
squares criterion is used to estimate g(t):
ĝ(t) = arg min
g∈Hk

N
X
(y(t) − (g ∗ u)(t))2 + γkgk2Hk ,

(3)

t=1

where k·kHk is the norm of Hk and γ > 0 is a regularization parameter and controls the tradeoff between
the data fit and the regularization term.
The performance of the kernel-based regularization method depends on several factors, such as the
choice of the kernel k(t, s; θ), the estimation of the hyper-parameter θ , and the input design. Assume
1

The kernel k(t, s; θ) sometimes will be written as k(t, s) for simplicity. Recall that a function k : X × X → R with X
P
being a metric space is called a positive semidefinite kernel, if it is symmetric and satisfies m
i,j=1 ai aj k(xi , xj ) ≥ 0 for any
m ∈ N, {x1 , · · · , xm } ⊂ X and {a1 , ..., am } ⊂ R. According to the Moore-Aronszajn Theorem, see e.g., [19], to every positive
semidefinite kernel k(x, x′ ) there exists a unique RKHS Hk with k(x, x′ ) as the reproducing kernel, i.e., with kx , k(x, ·),
k(x, x′ ) has the reproducing property hf, kx iHk = f (x), for f ∈ Hk , x ∈ X.

October 12, 2017

DRAFT

4

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

that a kernel k(t, s; θ) has been chosen, an estimate of θ has been found, and an input u(t) has been
designed. Then by setting γ = σ 2 and defining A ∈ RN ×N with its (t, s)th element At,s as follows
a(t, s) = (k(t, ·) ∗ u)(s), At,s = (a(·, s) ∗ u)(t),

(4)

we get the solution to (3) by the representer theorem [6, Theorem 3, page 671]:
ĝ(t) =

N
X

ĉs a(t, s),

(5)

s=1

where ĉs is the sth element of ĉ = (A + γIN )−1 Y .

C. Kernels for Regularized Impulse Response Estimation
Many kernels have been introduced, e.g., the stable spline (SS) kernel in [2] and the tuned correlated
(TC) kernel and the diagonal correlated (DC) kernel in [4]:
kSS (t, s; α) =

e−α(t+s) e−α max{t,s} e−3α max{t,s}
−
,
2
6

(6a)

kTC (t, s; β) = e−β(t+s) e−β|t−s| , β > 0,

(6b)

kDC (t, s; α, β) = e−α(t+s) e−β|t−s| , α > 0, β ≥ 0,

(6c)

where (6b) is a special case of (6c) with α = β [4] and is also called the first-order stable spline (SS-1)
kernel.
The SS kernel (6a) and the TC kernel (6b) can be derived [3] by applying a “stable” coordinate
change to the second-order and first-order spline kernel, respectively. To be specific, we first recall the
Sobolev space and its associated kernel [20], [21], [22]. The Sobolev space is the general term used for
0 defined on [0, 1] is defined
a functional space whose norm involves derivatives. The Sobolev space Wm

as the set of functions f : [0, 1] → R such that the next conditions hold:
1) for i = 0, · · · , m − 1, its ith order derivative f (i) is absolutely continuous2 and moreover, f (i) (0) =
0;

2) its mth order derivative f (m) ∈ L2 ([0, 1]).
0 over R is defined through the classical inner product in L2 ([0, 1]):
The inner product on Wm

hf, hiWm0 = hf (m) , h(m) iL2 ([0,1])
Z 1
f (m) (τ )h(m) (τ )dτ.
=

(7)

0

2

Let X be an interval of R. A function f : X → R is said to be absolutely continuous on X if for any ǫ > 0, there
P
is δ > 0 such that when any finite sequence of pairwise disjoint subinterval (ai , bi ) of X satisfies i (ai − bi ) < δ, then
P
(1)
i |f (ai ) − f (bi )| < ǫ. Moreover, f is absolutely continuous on a compact interval [a, b] if and only if f has a derivative f
R
x
almost everywhere and f (1) is Lebesgue integrable, and f (x) = f (a) + a f (1) (τ )dτ .

DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

5

0 is a RKHS with the reproducing kernel, often called the spline
It can be verified [20], [21], [22] that Wm

kernel,
S
wm
(τ, ν)

=

Z

0

1

Gm (τ, s)Gm (ν, s)ds, τ, ν ∈ [0, 1],

(8)

/(m − 1)! with (x)+ = x for x ≥ 0 and (x)− = 0
where the Green’s function Gm (τ, s) = (τ − s)m−1
+

for x < 0. In particular, the spline kernels (8) with m = 1, 2 are called the first-order and second-order
spline kernel, respectively:
w1S (τ, ν) = min{τ, ν}, τ, ν ∈ [0, 1],

(9)

1
1
w2S (τ, ν) = τ ν min{τ, ν} − (min{τ, ν})3 , τ, ν ∈ [0, 1].
2
6

(10)

Then it is easy to verify the following result.
Proposition 2.1 ([3]): Consider the SS kernel (6a) and the TC kernel (6b). Then the following result
holds
kSS (t, s; α) = w2S (e−αt , e−αs ),

(11)

kTC (t, s; β) = w1S (e−2βt , e−2βs ).

(12)

It follows from Proposition 2.1 that the SS and TC kernels have several properties inherited from
their mother kernel (8), such as the orthonormal basis expansion, the explicit expression of the norm,
and the maximum entropy property. Interestingly, the DC kernel can be derived by applying the same
“stable” coordinate change to a “generalized” first-order spline kernel so that the DC kernel have similar
properties as the SS and TC kernels inherited from the spline kernel (9).
III. S TABLE G ENERALIZED

FIRST- ORDER

S PLINE K ERNEL

To show this, we recall the generalized Sobolev space and its associated kernel [20], [21], [22]. The
generalized Sobolev space is derived by replacing the ith order derivatives i = 1, · · · , m, in the definition

0 , by more general derivatives. Specifically, the generalized mth order derivative of f are defined
of Wm

as follows:
Dm f =

d 1 d 1
d 1
···
f,
dτ a1 dτ a2
dτ am

(13)

where the functions ai (τ ) with τ ∈ [0, 1], i = 1, · · · , m are functions such that Dm f is well-defined,
and moreover,
d 1 d 1
d 1
f, M2 f =
f,
dτ am
dτ am−1 dτ am
d 1
d 1
· · · , Mm−1 f =
···
f.
dτ a2
dτ am

M0 f = f, M1 f =

October 12, 2017

(14)

DRAFT

6

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

0 defined on [0, 1] is defined as the set of functions f : [0, 1] → R
Then the generalized Sobolev space W̃m

such that the following conditions hold:
1) for i = 0, · · · , m − 1, Mi f is absolutely continuous and moreover, Mi f (0) = 0;

2) Dm f ∈ L2 ([0, 1]).

Note that the functions f : [0, 1] → R such that Mi f (0) = 0, i = 0, · · · , m − 1 can be represented in
the following form
f (τ ) = am (τ )
f (τ ) = am (τ )

Z

Z

τ

(Dm f )(s)ds, m = 1,

(15a)

0
τ

am−1 (tm−1 )dtm−1
0

am−2 (tm−2 )dtm−2 · · ·

Z

t1

Z

tm−1
0

(Dm f )(s)ds, m > 1.

(15b)

0

By interchanging the integration order in (15), we have
Z τ
G̃m (τ, s)(Dm f )(s)ds
f (τ ) =

(16)

0

where



0

τ ≤s

, m = 1,
 a (τ ) τ > s
m
Z τ
Z τ
a1 (t1 )dt1
G̃m (τ, s) = am (τ )
t1
s
Z τ
am−1 (tm−1 )dtm−1 , m > 1.
a2 (t2 )dt2 · · ·
G̃m (τ, s) =

(17a)

(17b)

tm−2

0 over R is defined through the classical inner product in L2 ([0, 1]):
Moreover, the inner product on W̃m

hf, hiW̃ 0 = hDm f, Dm hiL2 ([0,1])
m
Z 1
(Dm f )(τ )(Dm h)(τ )dτ.
=

(18)

0

0 is a RKHS with the reproducing
It can be verified [20], [21], [22] that the generalized Sobolev space W̃m

kernel
GS
wm
(τ, ν)

=

Z

1

G̃m (τ, s)G̃m (ν, s)ds,

(19)

0

which is called the generalized spline kernel here.
Then the following result can be proved for the DC kernel3 .
3

All proofs of the propositions are deferred to the Appendix.

DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

7

Proposition 3.1: Consider the DC kernel (6b) and the generalized spline kernel (19). If we take
m = 1, a1 (τ ) = τ ρ , ρ > −0.5,

(20)

then we have
α−β
),
2β

(21)

w1GS (τ, ν; ρ) = τ ρ ν ρ min{τ, ν}, τ, ν ∈ [0, 1].

(22)

kDC (t, s; α, β) = w1GS (e−2βt , e−2βs ;

where

The kernel (22) is called the generalized first-order spline kernel here. Recall from Proposition 2.1 that
the SS kernel (6a) is called the stable spline kernel, because it can be obtained by applying a “stable”
coordinate change to the second-order spline kernel. Now Proposition 3.1 shows that the DC kernel (6c)
can be obtained by applying the same “stable” coordinate change to the generalized first-order spline
kernel (22). The DC kernel can thus be called a stable generalized first-order spline kernel. As will be
shown below, this finding provides new facets to understand the properties of the DC kernel.
IV. N ORM

AND

O RTHONORMAL BASIS E XPANSION

In this section, we derive the norm and the orthonormal basis expansion of the generalized first-order
spline kernel and the DC kernel from their mother kernel (9).
To this goal, we first recall that when X is compact and the positive semidefinite kernel k(x, x′ ) with
x, x′ ∈ X is continuous, Hk has the following property by Mercer’s Theorem, see e.g., [23], [24, Thm.

17, page 90], [25, Thm. 1, page 34].
Let µ be a Borel measure on X and L2 (X, µ)4 be the space of functions f for which
∞. Then the integral operator Lk φ(x) on L2 (X, µ):
Z
k(x, x′ )φ(x′ )dµ(x′ ), x ∈ X.
Lk φ(x) ,

2
X (f (x)) dµ(x)

R

<

(23)

X

∞ 5
has at most countably many positive eigenvalues {λi }∞
i=1 and orthonormal eigenfunctions {φi }i=1 , and

the positive semidefinite kernel k(x, x′ ) has a series expansion
k(x, x′ ) =

∞
X

λi φi (x)φi (x′ ),

(25)

i=1

4

When X is a subset of R and µ is the Lebesgue measure, L2 (X, µ) will be written as L2 (X) for simplicity.

5

If for some λ, the homogenous integral equation
Lk φ(x) = λφ(x), x ∈ X

(24)

has solutions other than φ(x) = 0, λ and the solutions of (24) are called the eigenvalues and eigenfunctions of the integral
operator Lk , respectively.

October 12, 2017

DRAFT

8

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

√
which converges uniformly and absolutely on X × X . Moreover, { λi φi }∞
i=0 forms an orthonormal basis

of Hk , which gives an alternative representation of Hk by [25, Thm. 4]:
Hk = {f : f =

∞
X

fi φi , with

∞
X
f2
i

i=1

i=1

λi

< ∞}.

(26)

Now we recall the eigenvalues and the orthonormal eigenfunctions of the first-order spline kernel (9)
on L2 ([0, 1]), which are well-known, see e.g., [26, equations (119)-(120), page 196]:
Z 1
min{τ, ν}ϕi (ν)dν = λi ϕ(τ ), i = 1, 2, · · ·
0

where λi =

1

(i −

1 2 2 , ϕi (τ )
2) π

1
1
= 2 2 sin((i − )πτ ).
2

(27)

Moreover, Mercer’s Theorem [23], [24, Thm. 17, page 90], [25, Thm. 1, page 34] guarantees that the
series
min{τ, ν} =

∞
X

λi ϕi (τ )ϕi (ν)

(28)

i=1

√
converges absolutely and uniformly on [0, 1] × [0, 1], and { λi ϕi }∞
i=1 forms an orthonormal basis of
W10 and W10 has an equivalent representation:
W10

= {h|h(τ ) =

∞
X
i=1

hi ϕi (τ ), τ ∈ [0, 1],

∞
X
h2
i

i=1

λi

< ∞},

(29)

and the norm of h can be computed according to
khk2W10


 Z 1
∞
X
h2i
(1)
2
(h (τ )) dτ
=
=
λi
0

(30)

i=1

Noticing the similarity between (9) and (22), we can prove the next result for the generalized first-order
spline kernel (22).
Proposition 4.1: Consider the generalized first-order spline kernel (22). Let the eigenvalues and orthonormal eigenfunctions of the first-order spline kernel (9) take the form of (27). Then the following
results hold:
∞
1) Let φi (τ ) = τ ρ ϕi (τ ), i = 1, 2, · · · . Then {λi }∞
i=1 and {φi }i=1 are the eigenvalues and orthonormal

eigenfunctions of the generalized first-order spline kernel (22) on L2 ([0, 1], µ(ν)) with the measure
µ(ν) such that dµ(ν) = ν −2ρ dν , respectively.

2) The series
τ ρ ν ρ min{τ, ν} =

∞
X

λi φi (τ )φi (ν)

(31)

i=1

converges absolutely and uniformly on [0, 1] × [0, 1].
DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

9

√
0
0
3) { λi φi }∞
i=1 forms an orthonormal basis of W̃1 , and W̃1 has an equivalent representation:
W̃10

= {f |f (τ ) =

∞
X

fi φi (τ ),

i=1

τ ∈ [0, 1],

∞
X
f2
i

i=1

λi

< ∞},

(32)

and the norm of f can be computed according to
2 !
Z 1
∞
2
X
f
f
(τ
)
d
i
kf k2W̃ 0 =
dτ .
=
1
λi
dτ τ ρ
0

(33)

i=1

By Proposition 4.1 and noting (21), we can prove the next result for the DC kernel (6c).
Proposition 4.2: Consider the DC kernel (6c). Let the eigenvalues and orthonormal eigenfunctions of
the generalized first-order spline kernel (22) take the form of (31). Then the following results hold:
∞
1) Let ψi (t) = φi (e−2βt ), i = 1, 2, · · · . Then {λi }∞
i=1 and {ψi }i=1 are the eigenvalues and the

orthonormal eigenfunctions of the DC kernel (6c) on L2 ([0, ∞), ι(t)) with the measure ι(t) such

that dι(t) = 2βe2β(2ρ−1)t dt, respectively.
2) The series

e−α(t+s) e−β|t−s| =

∞
X

λi ψi (t)ψi (s)

(34)

i=1

converges absolutely and uniformly on Z × Z ′ with Z, Z ′ being any compact subsets of [0, ∞).
√
DC
3) { λi ψi }∞
induced by the DC kernel (6c), and
i=1 forms an orthonormal basis of the RKHS H
HDC has an equivalent representation:
HDC = {g|g(t) =

∞
X
i=1

gi ψi (t), t ≥ 0,

and the norm of g can be computed according to
kgk2HDC =

∞
X
g2
i

i=1

λi

∞
X
g2
i

i=1

λi

< ∞},

.

4) HDC and W̃10 are isometrically isomorphic and
2

Z ∞
1 (1)
2βe(4ρ+2)βt
g (t) + ρg(t) dt.
kgk2HDC =
2β
0

(35)

(36)

(37)

Remark 4.1: It should be noted that the orthonormal basis expansion of the DC kernel with respect to
the Lebesgue measure has been derived before in [27] and is different from (34). Then a natural question
is ”what is their difference?”. It turns out that these orthonormal basis expansions are optimal in some
sense with respect to certain kinds of inputs. Due to the limitation of the space, the details cannot be put
here but will be in an independent paper.

October 12, 2017

DRAFT

10

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

Remark 4.2: When ρ = 0, the DC kernel (6c) reduces to the TC kernel (6b). In this case, we have
Z ∞
1 2βt  (1) 2
kgk2HTC =
g (t) dt
e
(38)
2β
0
Comparing (38) with (37) shows from another perspective that the DC kernel is more flexible than
the TC kernel for regularize impulse response estimation. For illustration, if g(t) is constrained to be
g(t) = e−γt , then a necessary condition for g ∈ HTC is γ > β , but a necessary condition for g ∈ HDC

is γ > (2ρ + 1)β , where ρ > −0.5.
V. M AXIMUM E NTROPY I NTERPRETATION

OF

N ON - UNIFORMLY S AMPLED DC K ERNEL

In this section, for the non-uniformly sampled DC kernel (6c), we first derive its maximum entropy
(MaxEnt) property from its mother kernel (9) using arguments similar to [17] and then show that its
kernel matrix has tridiagonal inverse.
First, recall that a real-valued stochastic process w(i) with i = 0, 1, 2, · · · , is called a white Gaussian
noise if the w(i)’s are independent identically Gaussian distributed with zero mean and constant variance6.
Then we construct a Gaussian process f (τ ) defined on an ordered index set Γ = {τi |0 ≤ τi < τi+1 ≤
1, i = 0, 1, 2, · · · , } as follows:
f (τ0 ) = 0 with τ0 = 0,
f (τk ) =

τkρ

k
X
i=1

√
w(i − 1) τi − τi−1 , k = 1, 2, · · · .

(39)

It is easy to verify that f (τ ) with τ ∈ Γ has the generalized first-order spline kernel (22) with τ, ν ∈ Γ
as its covariance function. Then we can prove the next result for the kernel (22).
Proposition 5.1: Let h(τ ) with τ ∈ Γ be any stochastic process with h(τ0 ) = 0 for τ0 = 0. For any
n ∈ N, the stochastic process f (τ ) in (39) is the optimal solution to the MaxEnt7 problem
maximize
h(·)

H(h(τ1 ), h(τ2 ), · · · , h(τn ))

subject to E(h(τi )) = 0, i = 1, · · · , n


h(τ1 )
V
= τ1 ,
τ1ρ


h(τi ) h(τi−1 )
−
= τi − τi−1 , i = 2, · · · , n
V
ρ
τiρ
τi−1
6

(40)

Without loss of generality, we assume the variance is 1.

7

Recall that for a real-valued random variable X, the differential entropy H(X) of X is defined as H(X) =
R
− S p(x) log p(x)dx, where p(x) is the probability density function of X and S is the support set of X.

DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

11

where E(·) and V (·) represent the expectation and variance, respectively, and for simplicity, H(h(t1 ), h(t2 ), · · · , h(tn ))
denotes the differential entropy of [h(t1 ) h(t2 ) · · · h(tn )]T .

Based on the stochastic process f (τ ) in (39), we define another Gaussian process g(t) defined on an
ordered index set T = {ti |0 ≤ ti < ti+1 ≤ ∞, i = 0, 1, 2, · · · , } as follows:
g(tk ) = e−2βρtk

n−1
X
i=k

p
w(n − 1 − i) e−2βti − e−2βti+1 ,

k = 0, · · · , n − 1, g(tn ) = 0 with tn = ∞.

(41)

It is easy to verify that g(t) with t ∈ T has the DC kernel (6c) with t, s ∈ T as its covariance function.
Then we can prove the next result for the DC kernel (6c).
Proposition 5.2: Let h(t) be any stochastic process with h(tn ) = 0 for tn = ∞. For any n ∈ N, the
stochastic process g(t) in (41) is the optimal solution to the MaxEnt problem
maximize
h(·)

H(h(t0 ), h(t1 ), · · · , h(tn−1 ))

subject to E(h(ti )) = 0, i = 0, · · · , n − 1


h(ti+1 )
h(ti )
V −2βρti+1 − −2βρti = e−2βti − e−2βti+1 ,
e
e

V



(42)

i = 0, 1, · · · , n − 2

h(tn−1 )
= e−2βtn−1 .
e−2βρtn−1

Proposition 5.2 leads to an interesting result that, the kernel matrix of the DC kernel (6c) with t, s ∈ T
has tridiagonal inverse, which is an extension of the result of [17] from the TC kernel (6b) to the DC
kernel (6c) and an extension of the result of [18] from the uniformly sampled case to the non-uniformly
sampled case.
Proposition 5.3: Consider the DC kernel (6c) with t, s ∈ T . Then its kernel matrix K DC ∈ Rn̄×n̄
with n̄ ≥ 3 has tridiagonal inverse, where the (i, j)th element of K DC is equal to kDC (ti , tj ; α, β) with
ti , tj ∈ T .

For illustration, we consider an example.
Example 5.1: We consider the inverse of K DC ∈ Rn̄×n̄ with n̄ = 10, α = −0.2, β = −0.3 and
ti , tj with i, j = 1, . . . , 10 take values from the set generated in Matlab with sort(rand(10,1)). By

Proposition 5.3, the inverse of K DC should be tridiagonal, which is confirmed by Fig. 1.
VI. C ONCLUSION
In this note, we have shown that the continuous-time diagonal correlated (DC) kernel can be interpreted
as a stable generalized first-order spline kernel. This interpretation provides new facets to understand the

October 12, 2017

DRAFT

12

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

1
300

2
3

200

4
100

5
0

6
7

-100

8
-200

9
10

-300

2

4

6

8

10

Fig. 1. Scaled image of (K DC )−1 in Example 5.1. The image is drawn by using imagesc in MATLAB, where the colder the
color the smaller the element of (K DC )−1 .

properties of the DC kernel. In particular, we derive a new orthonormal basis expansion of the DC kernel
and the explicit expression of the norm of the RKHS associated with the DC kernel. Moreover, for the
non-uniformly sampled DC kernel, we derive its maximum entropy property and show that its kernel
matrix has tridiagonal inverse. There are several interesting works that can be done in the future. For
example, it is interesting and possible to derive the explicit expression of the norm of the sampled DC
kernel for discrete time linear time invariant system identification, and compare the derived orthonormal
basis expansion with the one in [27].
A PPENDIX
A. Proof of Proposition 3.1
On the one hand, it is easy to check that under the assumption (20) the generalized spline kernel (19)
takes the form of (22). On the other hand, the DC kernel (6c) can be rewritten as follows:
kDC (t, s; α, β) = e(β−α)(t+s) kTC (t, s; β)
− β−α 

− β−α
2β
2β
−2βt
−2βs
= e
e
min{e−2βt , e−2βs }.

Then noting the above equation and (22) gives the result.
B. Proof of Proposition 4.1

We first consider the proof for 1). For any i ∈ N, we have
Z 1
τ ρ ν ρ min{τ, ν}φi (ν)dµ(ν)
0

=

Z

1

τ ρ ν ρ min{τ, ν}ν ρ ϕi (ν)ν −2ρ dν

0

= λi τ ρ ϕi (τ ) = λi φi (τ )

DRAFT

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

13

∞
and moreover, note that the orthonormality of {φi }∞
i=1 in L2 ([0, 1], µ(ν)) follows from that of {ϕi }i=1

in L2 ([0, 1]).
For 2) and 3), we first show that the kernel (22) is continuous. By Schwartz inequality, we have
|w1GS (τ + h, ν + k) − w1GS (τ, ν)|
≤ {w1GS (τ + h, τ + h)×
1

[w1GS (ν + k, ν + k) − 2w1GS (ν + k, ν) + w1GS (ν, ν)]} 2
+ {w1GS (ν, ν)×
1

[w1GS (τ + h, τ + h) − 2w1GS (τ + h, τ ) + w1GS (τ, τ )]} 2 .

Now we let h, k → 0, then w1GS (ν + k, ν + k) → w1GS (ν, ν) because w1GS (τ, ν) is continuous along the
diagonal τ = ν and w1GS (ν + k, ν) → w1GS (ν, ν) because w1GS (τ, ν) is continuous in τ . This means that
the right-hand side of the above inequality converges to zero as h, k → 0 and thus the kernel (22) is
continuous. Then noting that the measure µ(ν) is a Borel measure, the results follow from [25, Thm. 1,
page 34] and [25, Thm. 4, page 37], respectively. This completes the proof.
C. Proof of Proposition 4.2
The proof 1) to 3) can be done by applying the coordinate change τ = e−2βt to 1) to 3) of Proposition
4.1. So we only prove 4) below.
Given any f ∈ W̃10 , it follows from (32) that there exist fi , i = 1, 2, · · · , such that f (τ ) =
P∞
P∞
−2βt , we get a unique g(t) , f (e−2βt ) =
−2βt ) =
i=1 fi φi (τ ) with τ ∈ [0, 1]. Taking τ = e
i=1 fi φi (e
P∞
P
P∞ 2
DC
2
2
2
since ∞
i=1 fi ψi (t) ∈ H
i=1 fi /λi < ∞, and moreover, kf kW̃ 0 = kgkHDC =
i=1 fi /λi . Sim1

ilarly, given any g ∈ HDC , it follows from (35) that there exist gi , i = 1, 2, · · · such that g(t) =
P∞
i=1 gi ψi (t) with t ≥ 0. Taking t = log τ /(−2β), we get a unique f (τ ) , g(log τ /(−2β)) =
P∞ 2
P∞
P∞
2
0
i=1 gi /λi < ∞, and moreover, kgkHDC =
i=1 gi φi (τ ) ∈ W̃1 , since
i=1 gi ψi (log τ /(−2β)) =
P
2
DC
and W̃10 are isometrically isomorphic. Given any g ∈ HDC , we
kf k2W̃ 0 = ∞
i=1 gi /λi . Therefore, H
1

have

kgk2HDC

October 12, 2017


d f (τ ) 2
=
dτ
dτ τ ρ
0
!2
log τ
d g( −2β )
dτ
dτ τ ρ


1
1 (1) log τ
log τ 2
g (
) + ρg(
) dτ
τ 2ρ+2 2β
−2β
−2β

=

kf k2W̃ 0
1

=

Z

=

Z

1
0
1
0

Z

1

DRAFT

14

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

Taking the coordinate change τ = e−2βt in the above equation and simple calculation yields (37). This
completes the proof.

D. Proof of Propositions 5.1 and 5.2
The proofs of Propositions 5.1 and 5.2 are similar. Due to the limitation of space, we only give the
proof for the latter.
Define
l(ti ) =

h(ti )
,i
−2βρt
i
e

= 0, · · · , n − 1

(43)

Then let L = [l(t0 ) l(t1 ) · · · l(tn−1 )]T , V = [h(t0 ) h(t1 )

· · · h(tn−1 )]T , and B = diag([e−2βρt0 e−2βρt1 · · · e−2βρtn−1 ]). Then we have V = BL. Since B is

nonsingular, it holds that H(V ) = H(L) + log det(B) according to [28, Corollary to Thm. 8.6.4]).
Moreover, since B is independent of l(t) and h(t), the MaxEnt problem (42) is equivalent to
maximize

H(l(t0 ), l(t1 ), · · · , l(tn−1 ))

subject to

E(l(ti )) = 0, i = 0, · · · , n − 1,

l(·)

V (l(ti+1 ) − l(ti )) = e−2βti − e−2βti+1 , i = 0, · · · , n − 2
V (l(tn−1 )) = e−2βtn−1 ,
By [17, Thm. 1], the optimal solution to the above MaxEnt problem is the Gaussian process
l(tk ) =

n−1
X
i=k

p
w(n − 1 − i) e−2βti − e−2βti+1 ,

(44)

k = 0, · · · , n − 1, l(tn ) = 0 with tn = ∞.

Clearly, (44) and (43) implies that the Gaussian process g(t) in (41) is the optimal solution to the MaxEnt
problem (42). This completes the proof.

E. Proof of Proposition 5.3
It follows from (42) that the Gaussian process g(t) with t ∈ T defined in (41) can be rewritten as
follows:
g(ti+1 ) = e−2βρ(ti+1 −ti ) g(ti )
p
+ e−2βρti+1 e−2βti − e−2βti+1 w̄(i), i = 0, . . . , n − 2,
g(tn−1 ) = e−β(2ρ+1)tn−1 w̄(n − 1)

DRAFT

(45)

October 12, 2017

SHELL et al.: BARE DEMO OF IEEETRAN.CLS FOR JOURNALS

15

where w̄(i) with i = 0, 1, . . . is a white Gaussian noise with zero mean and unit variance. It follows
from (45) that the Gaussian process g(t) with t ∈ T is also a Markov process with order 1. Then by [7,
Lem. A.3], we have the kernel matrix K DC has tridiagonal inverse. This completes the proof.
R EFERENCES
[1] L. Ljung, System Identification - Theory for the User. Upper Saddle River, N.J.: Prentice-Hall, 2nd ed., 1999.
[2] G. Pillonetto and G. D. Nicolao, “A new kernel-based approach for linear system identification,” Automatica, vol. 46, no. 1,
pp. 81–93, 2010.
[3] G. Pillonetto, A. Chiuso, and G. D. Nicolao, “Prediction error identification of linear systems: a nonparametric Gaussian
regression approach,” Automatica, vol. 47, no. 2, pp. 291–305, 2011.
[4] T. Chen, H. Ohlsson, and L. Ljung, “On the estimation of transfer functions, regularizations and Gaussian processes Revisited,” Automatica, vol. 48, pp. 1525–1535, 2012.
[5] T. Chen, M. S. Andersen, L. Ljung, A. Chiuso, and G. Pillonetto, “System identification via sparse multiple kernelbased regularization using sequential convex optimization techniques,” IEEE Transactions on Automatic Control, no. 11,
pp. 2933–2945, 2014.
[6] G. Pillonetto, F. Dinuzzo, T. Chen, G. De Nicolao, and L. Ljung, “Kernel methods in system identification, machine
learning and function estimation: A survey,” Automatica, vol. 50, no. 3, pp. 657–682, 2014.
[7] T. Chen, “On kernel design for regularized LTI system identification,” Automatica, vol. X, no. XX, p. to appear, 2017.
[8] M. Zorzi and A. Chiuso, “The harmonic analysis of kernel functions,” arXiv preprint arXiv:1703.05216, 2017.
[9] G. Prando, A. Chiuso, and G. Pillonetto, “Maximum entropy vector kernels for mimo system identification,” Automatica,
vol. 79, pp. 326–339, 2017.
[10] G. Pillonetto and A. Chiuso, “Tuning complexity in regularized kernel-based regression and linear system identification:
The robustness of the marginal likelihood estimator,” Automatica, vol. 58, pp. 106–117, 2015.
[11] B. Mu, T. Chen, and L. Ljung, “On asymptotic properties of hyperparameter estimators for kernel-based regularization
methods,” CoRR, vol. abs/1707.00407, 2017.
[12] Y. Fujimoto and T. Sugie, “Informative input design for kernel-based system identification,” Automatica, p. to appear, 2017.
[13] B. Mu and T. Chen, “On input design for regularized LTI system identification: Power-constrained input,” CoRR,
vol. abs/1708.05539, 2017.
[14] J. Latarie and T. Chen, “Transfer function and transient estimation by Gaussian process regression in frequency domain,”
Automatica, vol. 52, pp. 217–229, 2016.
[15] T. Chen and L. Ljung, “Constructive state-space model induced kernels for regularized system identification,” in 19th IFAC
World Congress, (Cape town, South Africa), pp. 1047–1052, 2014.
[16] F. Carli, A. Chiuso, and G. Pillonetto, “Efficient algorithms for large scale linear system identification using stable spline
estimators,” in IFAC symposium on system identification, 2012.
[17] T. Chen, T. Ardeshiri, F. P. Carli, A. Chiuso, L. Ljung, and G. Pillonetto, “Maximum entropy properties of discrete-time
first-order stable spline kernel,” Automatica, vol. 66, pp. 34 – 38, 2016.
[18] F. P. Carli, T. Chen, and L. Ljung, “Maximum entropy kernels for system identification,” IEEE Transactions on Automatic
Control, vol. 62, no. 3, pp. 1471–1477, 2017.
[19] N. Aronszajn, “Theory of reproducing kernels,” Transactions of the American mathematical society, pp. 337–404, 1950.
[20] R. A. Adams and J. J. Fournier, Sobolev Spaces. Academic press, 2003.

October 12, 2017

DRAFT

16

IEEE TRANSACTIONS ON AUTOMATIC CONTROL, VOL. XX, NO. X, X 2017

[21] V. G. Maz’ja, Sobolev Spaces. Springer, 2013.
[22] G. Wahba, Spline models for observational data. Philadelphia: SIAM, 1990.
[23] J. Mercer, “Functions of positive and negative type, and their connection with the theory of integral equations,” Philosophical
transactions of the royal society of London. Series A, containing papers of a mathematical or physical character, pp. 415–
446, 1909.
[24] H. Hochstadt, Integral equations. John Wiley & Sons, 1973.
[25] F. Cucker and S. Smale, “On the mathematical foundations of learning,” American Mathematical Society, vol. 39, no. 1,
pp. 1–49, 2002.
[26] H. L. V. Trees, Detection, Estimation, and Modulation Theory: Part I. John Wiley & Sons, 2nd ed., 2001.
[27] T. Chen, G. Pillonetto, A. Chiuso, and L. Ljung, “Spectral analysis of the dc kernel for regularized system identification,”
in Proc. 54th IEEE Conference on Decision and Control(CDC), (Osaka, Japan), pp. 4017–4022, 2015.
[28] T. M. Cover and J. A. Thomas, Elements of information theory. John Wiley & Sons, 2012.

DRAFT

October 12, 2017

