Improved rates for prediction and identification for partially observed linear dynamical systems
Authors: Holden Lee
Venue: NeurIPS 2021 Submitted (2021)
Keywords: linear dynamical system
system identification
control theory
Hankel matrix
SVD
low-rank approximation

Abstract: 
Identification of a linear time-invariant dynamical system from partial observations is a fundamental problem in control theory. Particularly challenging are systems exhibiting long-term memory. A natural question is how learn such systems with non-asymptotic statistical rates depending on the inherent dimensionality (order) $d$ of the system, rather than on the possibly much larger memory length. We propose an algorithm that given a single trajectory of length $T$ with gaussian observation noise, learns the system with a near-optimal rate of $\widetilde O\left(\sqrt\frac{d}{T}\right)$ in $\mathcal{H}_2$ error, with only logarithmic, rather than polynomial dependence on memory length. We also give bounds under process noise and improved bounds for learning a realization of the system. Our algorithm is based on multi-scale low-rank approximation: SVD applied to Hankel matrices of geometrically increasing sizes. Our analysis relies on careful application of concentration bounds on the Fourier domain -- we give sharper concentration bounds for sample covariance of correlated inputs and for $\cal H_\infty$ norm estimation, which may be of independent interest.

23 As simple and tractable family of dynamical systems. LTI systems are a central object of study for control theory and time series analysis. The problem of prediction and filtering for known system dates back to [Ka660]. However, in many machine learning applications, the system is unknown and Improved rates for prediction and identification for partially observed linear dynamical systems
Anonymous Author(s)
Affiliation
Address
email
Abstract
Identification of a linear dine-invariant dynamical system from partial observa- tions is fundamental problem in control theory. Particularly challenging are sys- tems exhibiting long-term memory. A natural question is how learn such systems with con-asympiotic statistical rates depending on the inherent dimensionality (or- der) d of the system, rather than on the possibly much larger memory length. We propose an algorithm that given a single trajectory of length T with gaussian ob- servation noise, learns the system with near-optimal rate of Ō 4 in H2 er- ror, with only logarithmic, rather than polynomial dependence on memory length. We also give bounds under process noise and improved bounds for learning a realization of the system. Our algorithm is based on multi-scale low-rank approx- imation: SVD applied to Hankel matrices geometrically increasing sizes. Our analysis relies on careful application of concentration bounds on the Fourier do- main we give sharper concentration bounds for sample covariance of correlated inputs and for Hoc estimation, which may be of independent interest.
10
11
12
13
14
Introduction
15
We consider the problem of prediction and identification of unknown partially-observed linear Lime-invariant (LTI) dynamical system with stochastic noise,
axtt) = Ar(t - 1) Bu(! 1) + E(t)
(1)
yt) = Cx(()) + Du(t) + 17()).
(2)
18 with single trajectory length T, given access only input and output data. Here, u(t) Rd inputs, xr(t Rdd are the hidden states, y(t) Rd observations (or outputs). {(t) NO0 Eg) 20 and (() N(0 Eny) are iid gaussian noise, and A E Rdxd B E Rdxd., C E Rd, xd D E Rdd xd..
21 are matrices. Partial observability refers to the fact that we do not observe the state r(t), but rather
22 a noisy linear observation yt)).
26 must be learned from input and output data. Identification of an unknown system is often necessary first step for robust control [DMM* 19 BMR18 18). In a long line of recent work, the interplay between
27
28 machine learning and control theory has borne fruit an improved understanding of the statistical
29 and online learning guarantees for prediction, identification, and control for unknown systems. In machine learning, LTI systems also serve as a simple model problem for learning from correlated
30
Submitted to 35th Conference Neural Information Processing Systems (NeurIPS 2021). Do distribute.

---

31
data in stateful environments, and can give insight into understanding the successes of reinforcement learning [Rec19] TR19 and recurrent neural networks (HMRT8).
32
33
Partial observability poses a significant challenge to system identification: In the fully observed setting, given access to æt!), there is no obstacle to learning the matrices directly through linear regression. However, in the partially observed setting, the most natural form of the optimization problem is non-convex.
34
36
37
Systems exhibiting long-term memory are particularly challenging to learn. Restricting to marginally stable systems, this occurs when the spectral radius of A, p(A), close to I, and it implies that the output at a particular time cannot be accurately estimated without taking into account inputs over many previous time-steps-or oo the order of O (I-p(A] times steps. Such systems often arise in practice. A particular class of such systems are those exhibiting mulltiscale behavior, with different state variables that change on vastly different timescales [CRIO]. . For example, the body's 's pH level is affected both by long-term changes on a timescale of days or weeks, as well as breathing rate which changes over timescale of seconds. For such systems, it makes sense discretize the scale of the fastest changing variable, which leads long memory for the slowest-changing variable. With few exceptions, existing guarantees for learning partially observed LTI systems degrade the memory length increases. However, counting the number of parameters in the model right measure of statistical complexity is the intrinsic dimensionality of the system, not the memory length. This leads to the following natural question.
38
39
40
41
42
43
44
45
47
[ 2 suggests that the
48
49
50
Question: How can we learn partially observed LTI systems with (non-asymptotic) statistical rates that depend on the intrinsic dimensionality of the system, rather than the memory length?
51
52
Despite the simplicity of the question, little in the way of theoretical results are known. We focus on the particular problem of learning the impulse response functtin of the system-whic whhch fully determines its input-output behavior-in H2 norm. This is a natural norm for prediction problems 55 as it measures the expected prediction error under random input. Known guarantees for learning the IR depend on the memory length. One particularly undesirable consequence is that for a continuous system with time discretization A going to O, the memory scales as 1/A (while the system order stays constant), leading to suboptimal estimation by arbitrarily large factor.
53
54
56
58
59
Our key contribution is an algorithm and analysis that gives statistical rates that are optimal up to logarithmic factors. Unlike previous works, our rates depend on the system order d-the natural dimensionality of the problem- and only logarithmically on the memory length of the system. Our algorithm is based taking low-rank approximation (SVD) of the Hankel matrix, which is widely used technique in system identification. We consider a multiscale version of this algorithm, where we repeat this process for geometric sequence of sizes of the Hankel matrix. This is essential for obtaining stronger theoretical guarantee. In the setting of zero process noise, we prove that our algorithm achieves near-optimal ① d(d. d,) rates in H2 error for the learned system.
60
61
64
85
66
Our analysis relies careful application of concentration bounds the Fourier domain to give 68 sharper concentration bounds for sample covariance and Hox norm estimation, which may be of independent interest. While we consider our algorithm in a simple setting, we hope that this is first step to understanding and improving more complex subspace identification algorithms. Indeed, SVD and related spectral methods are a standard step used in subspace identification algorithms
89
70
71
72 such as N4SID; our analysis suggests that SVD has an important "de-noising effect".
73 We also give improved bounds for system identification, that is, learning the matrices A, B, C.D
74 using the Ho- Kalman algorithm [HK66), with Ō Lddd, +dz) rates.
1.1 Notation
Norms. We use ||:|| to denote the 2-norm of vector. For matrix A, let |A|| |||1/ denote its operator norm, p(A) denote its spectral radius (maximum absolute value of cigen-

---

value), and | denote its Frobenius norm. For a matrix-valuee function M(t) E Cd1 78
|| IF Σ. [M Let (A) denote the rth singular value of A.
80 Fourier transform. Given a matrix-valuee function F 7 Cmxn define the (discrete-time) Fourier transform as the function F: R Cm given by F(w) = Ei= F(t)e 2niwt -
Σ
81
82 Matrices. Given a sequence (F(t))==1 where each F(t) E Cm define Hankelaxt F)) as the 85 am x bn block matrix such that the (i, block is [Hankelaxo(F)lss = F(i + 1). Given 84 sequence (F(t)) where each F()) E Cm define the Toeplitz matrix as the block matrix such that the (ih block is (Toppax((F)s:
F(i 1)) For matrix A, let A ,A" A denote its transpose, Hermitian (conjugate transpose), and pseudoinverse, respectively. For vector-valued function {a,...,b) b) + R", let Varb R [[(a-b)+1)n denote the the vertical concatenation of (()).
87
88
Control theory. For matrix A [dxd define its resolvent as (A(2) (zI A)- For linear dynamical system D given by
-1
89
let TD T uu-yy denote the transfer function from to y (response to input). Then DD 1 Pu-ty C + D C( A) D. Let T { C: |=|=
91
= 1} be the unit circle in the complex plane. For matrix-valued function F: 93 T + Cdidd define the H2 and Hos norms by
[= dz | No = V 2T /
|
Ho. = sup [ =ET
For function F: No - Ch xdo define its Z-transform to be = = (m))=-. Con- sidered as a function T C, we can take its H2 and Hoc norms. Overloading notation, we will let ,, || 12F' for p 2, oc. The H2 and H. norms can be interpreted as the Frobenius and operator norms of the linear operator from input output, i..., they measure the average power of
94
35
96
98
the output signal under random or worst-case input, respectively. For background on control theory,
99
see e.g.. IZG++96
100
Constants.
proofs, C' may represent different constants from line to line.
2
Main results
101
102
We consider the problem of prediction and identification for an unknown linear dynamical sys- Our main goal is to obtain error guarantees in Hg norm, which determines prediction enror under random input [0019 Lemma 3.3].
103
tem I (
104
Problem 2.1. Consider the partially-observed LTI system D I 2 with gaussian inputs u(t) ~ N (0, Ia.) for S T. Suppose that the system stable, that is, p(A) < 1, and that we observe
105
106
107
single trajectory of length T started with x: (0)) = 0, that observe u(l) NO, Ia..) and y(1)
for 0 < T.
The goal is to learn LTI system D with the aim of minimizing | &Dll Equivalently, letting
(t) 100-(CA--P)
D. CCAL-IB,
t=0
121
denote the impulse response function (also called the Markov parameters) of the system, the goal is 111 to learn an impulse response F minimizing
F
H2
112 Note that learning F* is sufficient to fully understand the input-output behavior the system, but
113 we may also ask to recover the system parameters A, B, C, D up to similarity transformation (see
Theorem 2.3
115 Previous results 10019 SRD19
roughly depend polynomially on the "memory" T-p(aj which blows up the spectral norm of A approaches 1. In the setting of zero process noise, our goal is
3

---

to obtain rates that are Ō ( with only poly-logarithmic dependence on i-plaj See VT
Figure I for comparison.
118
We assume that ((A) because if D is not stable, it is in general impossible to learn DD with finite M2 error, system with infinite response can have arbitrarily small response any finite time interval. However, it may still be possible to learn the response up to time L < T in this
122 case [BRI9]] The marginally stable case (p(A) = 1) is an important case we leave to future work.
Method
Rollout type
Min samples
IR error
Least
squares
(IR)
Multi
T
[TBPR17]
Least squares (IR) [0019
Single
L
V
Nuclear norm minimization
Multi
min
0 V
(SF20))
Multi
dt.
rank-d SVD (Theorem 2.2)
Single
Here, L is the memory length for the system, which is T-p(A)) for well-conditioned systems. Rollout type refers to whether we have access to single trajectory multiple trajectories. Min samples refers to the minimum number samples (up to log factors) before bounds are operational. IR error refers to the error in the impulse response in Frobenius/Hg norm. Logarithmic factors are omitted.
Figure 1:
123 In our Algorithm | we first use linear regression to obtain a noisy estimate of the impulse re- sponse. Next, following standard system identification procedures, we form the Hankel matrix Hankel (F) with the entries of F on its diagonals. Because the true Hankel matrix
124
125
CB
CAB
CA
CAB
CA3B
Hankelzx1(FF)
CAL-IB B
CA
128
has rank d, we take low-rank SVD Rt of the Hankel matrix to "de-noise" the impulse response. We can then read off the estimated impulse response by averaging over the corresponding diagonal of Rr. For technical reasons, we need to repeat this process
127
128
a geometric sequence of sizes of the Hankel matrix: L L, L/2 L/2, L/4 x L/4, and so forth. This is because the low- rank approximation objective for a x é Hankel matrix encourages the diagonals that are ( to
129
130
be close-as those are the diagonals with the most entries-ant hence estimates F* well when
132
t = (()). In other words, low-rank estimation for Hankel (F) is only sensitive to the portion of the signal that is at timescale L. Repeating this process ensures that we cover all scales.
133
134
Our main theorem is the following.
135
Theorem 2.2. There a constant C₁ such that following holds. In the setting of Problem 2.1
suppose that F* the impulse response function, T is such that T 2 C₁ Ldd, log (Ldu)) Etrune
196
F 122,00))l Vdn + ||G 1121,00)11 7.
2//2
1/2
and M (0, C.C.A.A..... CAL-1)T
137
R [L Let < 8 S and o = v |
|F |Σ Llog 3
Ldu
M| Then with probabil-
138
ity at least 8, Algorithml ] learns an impulse response function F such that
139
F F Fr 0 ov (dy + du + log (5)) log L
Strune Vd+ d | (F"L(1,,0)) (LL.ool)le)
T
140 In the absence process noise (when Σg O), when L and T are chosen large enough, the first
141 term dominates, and ignoring log factors, the dependence O (v (Vldd We expect this to be

---

Algorithm Learning impulse response through multi-scale low-rank Hankel SVD
1: Input: Length L (power of 2), time T.
2: Part I: Linear regression to recover noisy impulse response
3: Let u(t) (0, la..) for 0 v T, and observe the outputs y(t) E &Rt,
4: Solve the least squares problem
(ly(t) - F u(())]?
mil FiSupp (.22
10
obtain the noisy impulse response F: [0, njnz 1 Xdu
IR
5: Part 2: Low-rank Hankel SVD to de-noise impulse response
6: Let F(0) = F(O).
7:
for = to logz L do
8:
Let l = 2*
9: Let Ry be the rank-d SVD Hankelexe(F)) (i.e., argmin || Hankelexe (FF)))
10:
For 1 < <<< C, let F(1) be the dy x d, matrix given by F(I) = EE+jje(Re)y where
(1)s denotes the (i, block of the matrix.
11:
end for
12: Output: Estimate impulse response F.
142
the optimal sample complexity up to logarithmic factors. However, in the presence of process noise, there is an undesirable factor of vI (M which (for well-conditioned matrices) is expected to
143
be O (1 (I))())) or OLL). We leave it an open problem to improve the guarantees in this setting.
144
145
Remark 1. The T.-factor dependence the process noise is unavoidable with the current algorithm: when the process noise has covariance Mn = / and decays after L steps, it can cause perturbations of size O (VT.) compared to the noiseless system. Even in the case d = 1, when the impulse response function is ae kt/L for a known k, the noise will cause the estimate of to be off by O(VL). and hence the H2 norm of the impulse response to be off by O(L). Our algorithm only regresses on previous inputs, but in the presence of process noise, a better approach is to regress on both the previous inputs u(t) and outputs yt) and then take a (weighted) SVD, as in N4SID (Qin06).
148
148
149
150
151
Remark 2. A burn-in time of (2(L) is information-theorrtieall required to get poly(d) rates. At- tempting to extrapolate an impulse response function from time o(L) to time L can magnify errors by (d), because the finite impulse response of a system of order d can approximate polynomial of degree d on j0, L].
154
155
156
We also show the following improved rates for learning the system matrices, by combining Hoc bounds for the learned impulse response with stability results for the Ho-Kalman algorithm (019). Because the input-output behavior is unchanged under a similarity transformation (A, B, C) 4
159
(W AW, W-1 B, CW), we can only learn the parameters up to similarity transformation.
Theorem 2.3. Keep the assumptions and notation of Theorem 2.2 suppose D is observable and controllable, and let
é' = L (d, + dh log (4))
+ Etrune-
T
Let H Hankel (x(t-1) (F*). Suppose that e' O min (H-)). Then with probability at least
163 1 1- 8, the Ho- Kalman algorithm (Algorithm in (0019) with Ti = L,, T2 = I T--11 I.-1) returns A, B,
184 such that there exists a unitary matrix W satisfying
max fl|c Tw|]
1---w-B -|-|-B-- w- [| o(va.e')
A- 1- -- =0 Omiu (H- V'. ( 1 + +1))
Omiu(H-)
5

---

165 As L can be chosen to make ttunne negligible, this gives Ō ((Laldd rates, however, with factors depending the minimum eigenvalue of H. This is improvement over the (Va rates in [0019].
188
We prove Theorem 2.2 in Section 4 and Theorem 2.3 in Appendix
3
Related work
189
170
We survey two classes of methods for learning partially observable LDS's, subspace identification and improper learning. With the exception of [RJR20], all guarantees have sample complexity depending on the memory length L, which we wish to avoid.
171
172
173
3.1
Subspace identification
The basic idea of subspace identification [Lju98 Qin06 VODMI2 is to learn certain structured matrix (such as Hankel matrix), take a best rank approximation (using SVD or another linear di- mensionality reduction method), and learn the system matrices A, B, C, D up to similarity transfor-
175
178
177
mation. Usage of spectral methods circumvents the fact that the most natural optimization problem
178
for A, B,C, C D is non-convex. However, classical guarantees for these methods are asymptotic.
Recently, various authors have given non-asymptotic guarantees for system identification algo- rithms. (0019] analyzed the Ho-Kalman algorithm (HK66) in this setting. [SRD19] consider the setting where system order is unknown and give an end-to-end result for prediction, while [TMP20) consider the problem of online filtering, that is, recovering x(t)'s to some linear transformation.
181
182
183
An alternate, empirically successful approach is that of nuclear norm minimization or regulariza- tion (PST133. (SOF20) (building on (CQXY16)) give explicit rates of convergence, and show that the algorithm has lower minimum sample complexity and is casier tune.
186
Our algorithm is based on the classical approach of taking a low-rank approximation of the Hankel matrix, but we repeat this process with Hankel matrices of sizes L LxL, L, x/2xL/2, L/4x and so forth; this is key modification that allows us to obtain better statistical rates. Our analysis builds on the analyses given in [0019 SOF20]. As essential part of the analysis is analyzing linear regression for correlated inputs, where extend the work of (DMR19 to MIMO systems, as explained below.
187
188
189
190
191
3.1.1
Linear regression with correlated inputs
192
An important step in obtaining non-asymptotic rates for system identification is analyzing linear regression for correlated inputs. The most challenging step is to lower-bound the sample covariance matrix of inputs to the linear regression. A lower bound, rather than matrix concentration result, is
199
194
195
sufficient [Men14] SMT 18 MT19;;
however, concentration result is obtainable in our setting.
1BPR17 give non-asymptotic bounds for learning the finite impulse response for SISO system in 8 Fourier norm; however, they require L rollouts of size O(L) and hence (2(L²) timesteps. Addressing the more challenging single-rollout setting, [0019 obtain bounds for single rollout of (L) timesteps, by using concentration bounds for random circulant matrices [4M14] to de-
197
198
199
200
rive concentration inequalities for the covariance matrix. These concentration inequalities for the covariance matrix were improved (by logarithmic factors) by [DMR19]. Although [DMR19) give analysis in the SISO case, we show in Theorem A.2 the results can be extended the MIMO case with an e-net argument.
201
204
3.2
Improper learning using autoregressive methods
205
Instead of solving the statistical problem of identifying parameters, another line of work develops algorithms for regret minimization in online learning. The goal is simply do well in predicting
208
future observations, with small loss (regret) compared to the best predictor in hindsight; the learned predictor is allowed to be improper, that is, take different functional form. In the stochastic case,

---

209
this allows prediction almost as well as if the actual system parameters were known; however, the framework also allows for adversarial noise.
210
211
One popular strategy for improperly learning the system is to learn a linear autoregressive filter over previous inputs and observations, ARMA model. Naturally, because larger hypothesis class, the statistical rates depend L rather than the system order d.
212
optimizing over
214
Theorem 4.7] consider the problem of online prediction for a fully or partially observed LDS, and give a regret bound that depends polynomially on the memory length L. Their approach works even for marginally stable systems, that is, systems with ( S 1. See also [AHMSI3
2LS+
215
216
217
HSZ17
HLS+18
KMTM19
TP20 RJR20] for previous work using autoregressive methods.
Of particular interest to is [JR200, which gives rates independent of spectral radius. Building on HSZI7), they observe that it suffices to regress on previous inputs and outputs projected to 219
220
lower-dimensee dimensional space. Their algorithm works in the setting of process noise and competes with the Kalman filter, but only when A KC has real eigenvalues, where K' the Kalman gain.
221
4
Proof of main theorem
222
In this section, we prove Theorem 2.2 The proof hinges the following lemma, which shows that if we observe a low-rank matrix plus noise, then taking a low-rank SVD can have a de-noising effect, producing a matrix that is closer to the true matrix.
224
225
228
Lemma 4.1 (De-noising effect of SVD). There exists a constant such that the following holds. Suppose that A E Cm a rank-r matrix, À A + E, and A, the rank-r SVD of A. Then
227
A A S CVT |E]|
Compare this with the original error | A 228
|| which can only be bounded by
V mmnnmm, n} |[]]] When applied to the d-SVD of the Hankel matrix, this gives a factor of Vd rather than VL for the error.
230
Proof. We have
|., A
2r
A, A
(4)
()
All + || 1l.))
2r
(5)
< 2r |E
(6)
where 4 follows from A, A having rank at most 2r, 5) follows from the triangle incquality,
232
233 and 6 follows from Weyl's Theorem: [ A AlB ar+1 S a (A) || ||
To prove Theorem 2.2 we will need to obtain bounds for F: 10. , 1} R.d. xd. learned from linear regression in Hox norm. The following is our main technical result.
234
235
Lemma 4.2. There are C1, C2 such that the following hold. Suppose y F * + G 1) where u(l) ~ N(0, Ia.). E(L) N(O, E,), ni()) NOO, for 0 < < T. and (
C [0) Let F argmin (0....... Rxd IT- ly(() (F (.) Mo- = (G*(0),
238
E (B++))dxd,, and Etrunc = F [1+ (1+1,00) Ha Vdg +
240 ||a ||
22/2 For 0 < 8 1, T Ld log (Idz), - 5 a < L 6-L, L/,
[|(F F²) (attottii
S SC² [V (VIBull L/ du + dy, + log ))) + V || L'Ld, log
Lddy ||
+ Strunc
with probability at least 1 1-6 6.
7

---

242 In the case Σg = O, this roughly says that the error in the learned impulse response, F F", over any interval of length L', has all Fourier coefficients bounded in spectral norm by Ldd,, what we expect if the error from linear regression uniformly distributed over all frequencies.
243
T
244
A complete proof is in Appendix A we give a brief sketch. First, because the errors are Gaussian, the error from linear regression, F F*, follows a Gaussian distribution. To bound its covariance, we lower-bound the smallest singular value of the sample covariance of the inputs (Lemma A.I Appendix A.Ip. Here, the difficulty is that the inputs are correlated - the input at time is u---L- Fortunately, the translation structure means it is close to a submatrix of an infinite block Toeplitz matrix, which becomes block diagonal in the Fourier domain. This "decoupling" allows us show concentration. Compared to the SISO setting in DMR 191, require extra s-net argument. Once we have a bound on the covariance, we can bound any (F- F* [w)) by matrix concentration 252
247
248
249
250
251
253
Appenddix A.2 to bound the Hoc norm it suffices bound this over grid of w's (Lemma
Bounding the error in Hoo norm of the impulse response allows us to bound the error in operator norm of the Hankel matrix, as the following lemma shows.
255
Lemma 4.3. For any F Z + Cuxnn we have ||Hankelaxo (F)II S || x
257
Proof. Note that when Z + C", Supp(u) S 10, 1), we have (F * v)bbbta-1 Hankelaxol F)m--11- Hence for any o set to outside of [0, 1]. using Parseval's Theorem
258
259
and the fact that the Fourier transform of a convolution is the product of the Fourier transforms, we have
260
[Haakelaxo(F)os--aul S
2ll2
Fr
S
sup we[0,1]
F(w) |0l = 1l ||l||2
This shows that (omnkolaxx((F)) S |l
Theorem 2.2 will follow from the following bound after application of the triangle inequality.
263
Lemma 4.4. There are C₁, C2 such that the following holds for the setting of Problem 2.7
Suppose L is power of 2, and T Z CjLdy log Ldx). Let ||F 1(L+1,,0)) Ho dd +
264
285
G 114,,ll
|
and Ma-10 = (O,C,CA....
E [R(++)dxdd, C.LL,IA
CAI-IIT
Then with
286
probability at least 8. the output given by Algorithm
satisfies
( (T-F")Auu Fu)ul sca ca(v/
Exll ( (dy + du + log ( (t)) log
|E Lddy log Idd
| + Serume vd
T
T
287 Proof. We are in the situation of Lemma
4.2 with (t) = CA Let He = Hankelexe
288 and H2 = Hankclexe(F")) Suppose < L even. Note that
He = (ankellx((F))) Hankelaxe(F F*)
289 where Hj = Hankelexe(F")) is rank-d matrix, with error term is bounded by
|Haanelexee F*)|| I [I( Firune) 11 (1 1122-11lxx
by Lemma 4.3
/ dy log ( 7. )
V
T
Ldku
+ Exll LLd log
| + Etrune
by Lemmal 4.2
270 with probability at least - 5. Let Rs be the rank-d SVD of He. Then by Lemma 4.1
RR H o [Vj]lHaaeeleeef -F=)))))
(8)
8

---

271 Now letting F(t)
(Rels when 2 < 272 minimizes the sum of squared errors)
E
t S l we have (using the fact that the mean
Ree Hill
Σ
Σ
(Re)) FF
t= 2 +HI
Σ t=
2 ( Σ (1=()F())) (F 1-t+1
273 Note that we only get a lower bound with factor of l if we restrict to that (()), i.e., restrict to diagonals that have many entries. This is the reason we will have to repeat this process for multiple sizes. Hence
275
1 8/2
R Hill..
Together with (8) and 7 this gives with probability 2 1 - δ that
|| F 14 ll. sc( | d (d, + du - log (+)) |E Ldd, log
Lddt..
Etrune Vd v
| || +
T
T
277
Replacing 8 by log, I. using a union bound over powers of 2, and summing gives
|| | Fu) -0 -o(v |2n |l (d, - d. log (1) log
E2 l| Ldd,, log ( Ldu log L
T M || + Etrune Etrune Vd
+
278
Proof of Theorem 2.2
We have the bound in Lemma 4.4 and also the same bound
280
for
(F F* (0)
after applying Lemma 4.2 to (F F Finally, note that
(F - -FF)
F ((,,,)) | and use the triangle inequality.
281
5
Experiments
We compared three algorithms for learning the impulse response function: least-squares, and low- rank Hankel SVD with and without the multi-scale repetition. We include details of the experimental setup in Appendix D Note that to reduce the number of scales, we consider use a slight modification of our Algorithm which triples the size at each iteration instead.
285
286
The plots show the error || 1 Flly: where F is the estimated impulse response on [1, T.). averaged over 10 randomly generated LDS's, as function of the time T elapsed. We consider systems of order d = 1,3, 5, 10, and memory lengths L L= 27, 81.
290
Using SVD significantly reduces the error, supporting our theory which shows that SVD has a "de- noising" effect. Additionally, multiscale SVD has better performance than naive SVD when is moderate, L is large, and data is limited, but the performance is similar in data-rich setting.
291
292
d=1, L=27
d=5, L=81
5417
ettor
14
01
300
1260
1500
2000
800
500
1200
7550
1800
time
time

---

References
293
294
[AHMS13] Oren Anava, Elad Hazan, Shie Mannor, and Ohad Shamir. Online learning for time series prediction. In COLT 2013 The 26th Annual Conference on Learning Theory, June 12-14, 2013, Princeton University, NJ, USA, pages 172-184, 2013.
295
296
[BMR18] Ross Boczar, Nikolai Matni, and Benjamin Recht. Finite-data performance guarantees for the output-feedback control of an unknown system. In 2018 IEEE Conference Decision and Control (CDC), pages 2994-2999 IEEE, 2018.
299
300
[BTR13] Badri Narayan Bhaskar, Gongguo Tang, and Benjamin Recht. Atomic norm denoising with applications to line spectral estimation. IEFE Transactions on Signal Processing. 61(23):587755999 2013.
301
302
[6QXY16] Jian-Feng Cai, Xiaobo Qu, Weiyu Xu, and Gui-Bo Ye. Robust recovery of complex exponential signals from random gaussian projections via low rank hankel matrix re- construction. Applied and computational harmonic analysis, 41(2):470-490, 2016.
304
305
306
[CR10] Shaunak Chatterjee and Stuart Russell. Why are dbns sparse? In Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics, pages 81 31-88. JMLR Workshop and Conference Proceedings, 2010.
307
309
[DMM 19] Sarah Dean, Horia Mania, Nikolai Matni, Benjamin Recht, and Stephen Tu. On the sample complexity of the linear quadratic regulator. Foundations of Computational Mathematics, pages 1-47, 2019.
310
311
312
[DMR19] Boualem Djehiche, Othmane Mazhar, and Cristian R Rojas. Finite impulse response models: A non-asymptotic analysis of the least squares estimator. arXiv preprint arXiv:1911112994 2019.
313
315
[FPST13] Maryam Fazel, Ting Kei Pong, Defeng Sun, and Paul Tseng. Hankel matrix rank minimization with applications to system identification and realization. SIAM Journal on Matrix Analysis and Applications, 34(3)944-977,, 2013.
316
317
318
[GLS+20] Udaya Ghai, Holden Lee, Karan Singh, Cyril Zhang, and Yi Zhang. No-regret predic- tion in marginally stable systems. In COLT 2020 The 33rd Annual Conference Learning Theory, July 9-12, 2020, pages 1-44, 2020.
321
[HK66] BL Ho and Rudolph E Kalman. Effective construction of linear state-variable models from input/output functions. dr-Auttomaisternngstehhiiii 14(1-12):55-548, 1966.
322
323
[HLS+18] 18] Blad Hazan, Holden Lee, Karan Singh, Cyril Zhang, and Yi Zhang. Spectral filtering for general linear dynamical systems. In Advances in Neural Information Processing Systems, pages 4634 4643, 2018.
325
326
[HMR18]
Moritz Hardt, Tengyu Ma, and Benjamin Recht. Gradient descent learns linear dynam- ical systems. The Journal of Machine Learning Research, 19(1):005--108, 2018.
327
328
[HSZ17] Blad Hazan, Karan Singh, and Cyril Zhang. Learning linear dynamical systems via spectral filtering. In Advances in Neural Information Processing Systems, pages 1-2, 2017.
331
[Ka660] Rudolph Emil Kalman. A new approach to linear filtering and prediction problems. Journal of Basic Engineering, 12.1:35-45, 1960.
332
333
[KMR14]
Felix Krahmer, Shahar Mendelson, and Holger Rauhut. Suprema of chaos processes and the restricted isometry property. Communications on Pure and Applied Mathemat- ics, 67(11))1877-1190 2014.
10

---

336
[KMTM19]
Mark Kozdoba, Jakub Marecek, Tigran T. Tchrakian, and Shie Mannor. On-line learn- ing of linear dynamical systems: Exponential forgetting in kalman filters. In Proceed- ings of the AAAI Conference on Artificial Intelligence, volume 33, pages 4098-44 4105, 2019.
337
338
339
340
[Lju98] Lennart Ljung. System identification: Theory for the User. Prentice Hall, Upper Saddle River, NJ, 2 edition, 1998.
341
[Men14] Shahar Mendelson. Learning without concentration. In Conference Learning The- ory, pages 25-39, 2014.
[MT19] Nikolai Matni and Stephen Tu. tutorial on concentration bounds for system identifi- cation. In 2019 IEFE 58th Conference Decision and Control (CDC), pages 3741- 3749. IEEE, 2019.
347
[0019] Samet Oymak and Necmiye Ozay. Non-asymptotic identification of Iti systems from single trajectory. In 2019 American Control Conference (ACC), pages 5655-5661. IEEE, 2019.
348
350
[Qin06] Joe Qin. An overview subspace identification. Computers chemical engineer- ing, 30(10-12):500--1513 2006.
351
352
19e Benjamin Recht. A tour of reinforcement learning: The view from continuous control. Annual Review of Control, Robotics, and Autonomous Systems, 2:253-2733 2019.
353
[RJR20] Paria Rashidinejad, Jiantao Jiao, and Stuart Russell. Slip: Learning to predict in un-
355
known dynamical systems with long-term memory. arXiv preprint arXiv:2010.05899, 2020.
356
RRV 13] Mark Rudelson, Roman Vershynin, al. Hanson-wright inequality and sub- sub-gaussian concentration. Electronic Communications in Probability, 18, 2013.
358
359
[SBR19] Max Simchowitz, Ross Boczar, and Benjamin Recht. Learning linear dynamical sys- tems with semi-parametric least squares. arXiv preprint arXiv:1902.00768, 2019.
361
[SGA20] Arnab Sarker, Joseph E Gaudio, and Anuradha M Annaswamy. Parameter estimation bounds based on the theory of spectral lines. arXiv preprint wrXiv:2006.12668 2020.
362
363
[SMT+18) 18] Max Simchowitz, Horia Mania, Stephen Tu, Michael I Jordan, and Benjamin Recht. Learning without mixing: Towards a sharp analysis of linear system identification. In Conference On Learning Theory, pages 439 473, 2018
364
365
386
[SOF20] Yue Sun, Samet Oymak, and Maryam Fazel. Finite sample system identification: Im- proved rates and the role of regularization. In 2nd L4DC 2020 - Learning for Dynamics & Control, June 1 12- 2020, 2020.
367
368
389
[SRD19] Tuhin Sarkar, Alexander Rakhlin, and Munther A Dahleh. Finite-time system identification for partially observed Iti systems of unknown order. arXiv preprint 1rX1v:190..01848, 2019.
370
371
[TBPR17) Stephen Tu, Ross Boczar, Andrew Packard, and Benjamin Recht. Non-asymptotic analysis of robust control from coarse-grained identification. arXiv preprint arXiv: 2017.
373
374
[TMP20] Anastasios Tsiamis, Nikolai Matni, and George Pappas. Sample complexity of kalman filtering for unknown systems. In Learning for Dynamics and Control, pages 435-444, PMLR, 2020.
376
377
[TP20] Anastasios Tsiamis and George Pappas. Online learning of the kalman filter with logarithmic regret, arXiv preprint arX1v:2002.05141, 2020.

---

[TR19] Stephen Tu and Benjamin Recht. The gap between model-based and model-free meth- ods on the linear quadratic regulator: An asymptotic viewpoint. In Conference on Learning Theory, pages 3036-3083. PMLR, 2019.
381
382
[Ver18] Roman Vershynin. High-dimensional probability: An introduction with applications in data science, volume 47. Cambridge University Press, 2018.
[VODMI2] Peter Van Overschee and BL De Moor. Subspace Identification for Linear Systems. Springer Science & Business Media, 2012.
387
[W]20] Andrew Wagenmaker and Kevin Jamieson. Active learning for identification of linear dynamical systems. In COLT 2020 The 33rd Annual Conference on Learning Theory, July 9-12, 2020, 2020.
388
389
IZDG* [D] Kemin Zhou, John Comstock Doyle, Keith Glover, et al. Robust and optimal control, volume 40. Prentice hall New Jersey, 1996.
391
12

---

Checklist
392
393
1. For all authors...
384
(a) Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope? [Yes]
395
(b) Did you describe the limitations of your work? [Yes] The main limitations are the
sensitivity to process noise (Remark I after Theorem 2.2 elaborated in Appendix F
398
the restriction to Gaussian noise, and the assumption p(A) < 1.
399
(c) Did you discuss any potential negative societal impacts of your work? [N/A] The work is primarily theory result.
400
401
(d) Have you read the ethics review guidelines and ensured that your paper conforms to them? [Yes]
402
403
2. If you are including theoretical results...
404
(a) Did you state the full set of assumptions of all theoretical results? [Yes]
405
(b) Did you include complete proofs of all theoretical results? [Yes]
406
3. If you ran experiments...
407
(a) Did you include the code, data, and instructions needed to reproduce the main experi- mental results (either in the supplemental material or as a URL)? [Yes] Code is in the supplement. Executing either the Julia file or notebook directly produces the graphs
408
409
410
(b) Did you specify all the training details (e.g., data splits, hyperparameters, how they were chosen)? [Yes] Appendix D
411
412
(c) Did you report error bars (e.g., with respect to the random seed after running exper- iments multiple times)? [No] Graphs are averages over 10 runs; visual inspection of raw data shows that the variability is small.
414
415
(d) Did you include the total amount of compute and the type of resources used (e.g., type of GPUs, internal cluster, or cloud provider)? [No] Code runs in 5 minutes on laptop.
416
417
4. If you are using existing assets (c.g., code, data, models) or curating/relcasing new assets...
418
(a) If your work uses existing assets, did you cite the creators? [N/A]
419
(b) Did you mention the license of the assets? [N/A]
420
(c) Did you include any new assets either in the supplemental material or as a URL? [N/A]
421
(d) Did you discuss whether and how consent was obtained from people whose data you' 'e using/curating? [N/A]
424
(e) Did you discuss whether the data you are using/curating contains personally identifi- able information or offensive content? [N/A]
425
426
5. If you used crowdsoureing or conducted research with human subjects...
427
(a) Did you include the full text of instructions given to participants and scrcenshots, applicable? [N/A]
428
429
(b) Did you describe any potential participant risks, with links to Institutional Review Board (IRB) approvals. if applicable? [N/A]
430
431
(c) Did you include the estimated hourly wage paid participants and the total amount spent participant compensation? [N/A]
432
13


Reviews:
[{"review_id": "naiUOIW1GDC", "reviewer": "Reviewer_BcJQ", "summary": "The paper proposes an algorithm for learning LTI systems with zero process noise. Given an input trajectory generated from independent Gaussian input, the algorithm estimates the impulse response function, which in turn can be used to estimate system parameters up-to similarity transformation.\n\nThe key idea is the multi-scale smoothing of the impulse response: A first estimate is obtained using linear regression. Denoising of the response function estimate is then performed by low-rank approximation of Hankel matrices of different sizes, with the insight being that the size of the matrix determines which parts of the response functions are preserved under the low-rank approximation.  \n\nThe proposed algorithm achieves an error rate that scales only logarithmically with the system memory. The efficacy of the algorithm is demonstrated on a synthetic problem.", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Overall, this is a well-written paper. The contribution is clear and the flow of ideas is logical. \n\nMy main concern is that the experiments were conducted exclusively on synthetic data. While It is beneficial to have control over the parameters of the problem. It is equally important to demonstrate the efficacy of the algorithm in real-world data.\n\nOther comments:\n* L96-98 gives the impression that Frobenius norm lower bounds the operator norm whereas the opposite should be true.\n* L130: It is not clear what \"diagonals that are $\\Theta(l)$\" means. Also, do you mean skew-diagonals?\n* L137: What is G? Also, it would be good to define the indicator function.\n* Algorithm 1 - Line 4: $(F * u)(t)$  (use parentheses)\n* L146: $\\Sigma_x$ instead of $\\Sigma_y$ for process noise.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "I do not see a potential negative societal impact in the paper.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "3"}, "extra_text": {}}, {"review_id": "fYQ4_XDhuwS", "reviewer": "Reviewer_iHne", "summary": "This work considers the problem of learning linear dynamical systems with partial observability both in terms of learning the impulse response and the system parameters upto similarity transformation (using Ho-Kalman algorithm). The main technical contribution is Algorithm 1 and Theorem 2.2 which presents the error bounds for Algorithm 1. This shows that the impulse response can be learned upto frobenius error $\\sqrt{d/T}$, without any dependence on the mixing time/memory length of the dynamics. The main technical tool used in Algorithm 1 is the least squares estimation of the Hankel matrix at various scales to cover dependences from inputs at various time instants and denoising it with SVD since the true Hankel matrix is constrained to have rank d. \n\nThe experiments illustrate the advantage of the algorithm compared to least squares regression. ", "questions": "", "limitations": "", "rating": 7, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "I think the algorithm and its analysis is novel on both technical and conceptual level. It solves an important problem with near optimal bounds which does not depend on the memory length by using the elegant idea of SVD to denoise the Hankel matrix estimate obtained through least squares regression. Rest of the algorithm seems to be fairly standard. \n\nI think this deserves to be published in NeurIPS 2021.\n\nSome comments:\n\n1. The improved rates only hold when there is no additional noise in the hidden dynamics of the system. (This is mentioned in the paper)\n2. The parameter recovery in Theorem 2.3 has an additional factor L, which is not clearly explained. \n3. It would be great if the magnitude of $\\epsilon_{\\mathsf{trunc}}$ is mentioned. \n4. It would be useful to mention the computational complexity of the algorithm. ", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "Yes.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "4 hours"}, "extra_text": {}}, {"review_id": "E-QCaQ0FtC", "reviewer": "Reviewer_ZPKh", "summary": "This work studies the estimation of linear dynamical systems with unobserved state. The authors propose a novel algorithm which relies on obtaining low-rank approximations to the Hankel matrix. They provide several results on the performance of their algorithm, including finite time bounds on the error it obtains in estimating the system’s impulse response. As compared to previous work, the main novelty of this paper is eliminating a factor of the horizon in the case where there is no process noise.", "questions": "", "limitations": "", "rating": 5, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "Obtaining finite-time bounds for the estimation of linear dynamical systems—in both the setting where the state is observed and the state is unobserved—has been a well-studied problem in the machine learning community over the last several years, and this works fits in that literature. A primary goal of existing works has been to obtain error bounds which do not scale with the “memory” of the system, 1/(1-\\rho(A)). While achievable in the case where the state is observed, this has been a more difficult goal to reach in the case when the state is unobserved. Several existing works have addressed this—for example, [1] introduces the notion of \"phase rank\" to avoid the memory dependence—but it is not clear that the issue has been completely solved (phase rank could be exponential in the dimension in some cases, which is also highly undesirable). This work claims to remove the memory dependence without incurring poor dependence on other quantities, but is only able to do so in the case with 0 process noise.\n\nPros:\n- As noted above, removing the 1/(1-\\rho(A)) dependence is an important goal, and as yet has not been completely achieved in the case when the state is unobserved. As such, the problem considered here is relevant and, if solved, would be a contribution. \n- The authors do solve the problem in the case when the process noise is 0, eliminating the dependence on the memory in this setting.\n\nCons:\n- The primary issue with this paper is that the memory dependence can only be removed in the case when there is no process noise, i.e. the results improve on the state of the art only when there is no process noise. This is a very restrictive assumption, is likely unrealistic in practice, and it is not clear how interesting this setting is to the learning theory community (existing works have all assumed non-zero process noise). \n- Even in the noiseless setting, the memory dependence is not removed for the recovery of the system parameters (i.e. A,B,C,D) but only for the recovery of the impulse response. \n- The authors restrict to the case when \\rho(A) < 1 while existing works also encompass the case when \\rho(A) = 1. \n- The algorithm takes as input the memory length of the system, but it is not clear how one should actually set this without knowing the system parameters. \n- More care should be taken in comparing to [1]. [1] introduces the idea of phase rank, and provides bounds that scale with the phase rank instead of 1/(1-\\rho(A)) in some cases (even with non-zero process noise). No discussion of this is given—it should be clarified what precisely the bounds in [1] work out to be, and what the actual memory dependence suffered is. At minimum the authors should not claim that no existing works obtain bounds that are memory-independent and should show explicit examples where their results improve over the phase rank dependence in [1]. \n\nWhile this work does make a contribution, the restrictive assumptions noted above and the fact that [1] partially solves the problem already are, in my opinion, significant enough that it does not warrant an accept. \n\n[1] Max Simchowitz, Ross Boczar, and Benjamin Recht. Learning linear dynamical systems with semi-parametric least squares. arXiv preprint arXiv:1902.00768, 2019.", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "The authors do not discuss potential societal impact. However, the work is primarily theoretical so it is difficult to determine what negative societal impacts it may have. Limitations are discussed.", "review_text": "", "extra_scores": {"needs_ethics_review": "No", "ethics_review_area": "", "time_spent_reviewing": "2"}, "extra_text": {}}]

Author rebuttal:

Meta-review: The paper studies the problem of learning a partially observed linear dynamic system and provides non-asymptotic system identification guarantees that scale with the dimensionality of the system, rather than the memory length. This is definitely an interesting direction, but the reviewers identified several limitations with the results in the current paper. Most notably is the fact that the paper only concerns systems with no process noise.

I am sympathetic to the fact that making progress on a hard problem can take some time and may require making assumptions to obtain preliminary results and develop new ideas. However, it would be worthwhile to explore the consequences of the assumptions more to provide a more complete story. In the context of the present paper, as suggested by the reviewer, it would help to understand if other algorithms have a memory dependence in the no-process-noise regime and also if there are lower bounds in the more general setting. Such auxiliary results would complete the picture and make the paper much stronger. 

Thus, we recommend rejecting the paper at this time.
Decision: Reject