Distributionally Robust Linear Quadratic Control
Authors: Bahar Taskesen
Dan Andrei Iancu
Çağıl Koçyiğit
Daniel Kuhn
Venue: NeurIPS 2023 spotlight (2023)
Keywords: linear quadratic control
distributionally robust optimization
optimal transport
Wasserstein distance

Abstract: Linear-Quadratic-Gaussian (LQG) control is a fundamental control paradigm that is studied in various fields such as engineering, computer science, economics, and neuroscience. It involves controlling a system with linear dynamics and imperfect observations, subject to additive noise, with the goal of minimizing a quadratic cost function for the state and control variables. In this work, we consider a generalization of the discrete-time, finite-horizon LQG problem, where the noise distributions are unknown and belong to Wasserstein ambiguity sets centered at nominal (Gaussian) distributions. The objective is to minimize a worst-case cost across all distributions in the ambiguity set, including non-Gaussian distributions. Despite the added complexity, we prove that a control policy that is linear in the observations is optimal for this problem, as in the classic LQG problem. We propose a numerical solution method that efficiently characterizes this optimal control policy. Our method uses the Frank-Wolfe algorithm to identify the least-favorable distributions within the Wasserstein ambiguity sets and computes the controller's optimal policy using Kalman filter estimation under these distributions.

Distributionally Robust Linear Quadratic Control
Bahar Taşkesen
Dan A. Iancu Stanford University daniancu@stanford.
Çağal Koçyiğit
EPFL
University of Luxembourg
bahar taskesen@epfl..cl ch
cagil herpll.kccyiiittee lu
Daniel Kuhn
EPFL
daniel kwwn@epfl.nccc
Abstract
Linear-Quadenn Quadratic-Gausian (LQG) control is fundamental control paradigm that has been studied and applied in various fields such as engineering, computer science, economics, and neuroscience. It involves controlling system with linear dynamics and imperfect observations, subject to additive noise, with the goal of minimizing quadratic cost function depending the state and control variables. In this work, we consider a generalization of the discrete-ting time, finite-horizon LQG problem, where the noise distributions ambiguity sets centered at nominal (Gaussian) distributions. The objective is to minimize worst-case
unknown and belong to Wasserstein
across all distributions in the ambiguity set, including Non-Gaussian distributions. Despite the added complexity, we prove that a control policy that linear in the observations is optimal, in the classic LQG problem. We propose numerical solution method that efficiently characterizes this optimal control policy. Our method uses the Frolk algorithm identify the least- favorable distributions within the Wasserstein ambiguity sets and computes the controller's 's optimal policy using Kalman filter estimation under these distributions.
1. Introduction
The Linear Quadratic Regulator (LQR) is a classic control problem that has served as building block for numerous applications in engineering and computer science [3, 12], economics [29], or neuroscience [47]. It involves controlling system with linear dynamics and imperfect observations affected by additive noise, with the goal of minimizing a quadratic state and control cost. Under the assumption that noise terms are independent and normally distributed (a case referred to as or LQG), is well known that the optimal control policy depends linearly on the observations and can be obtained efficiently by using the Kalman filtering procedure and dynamic programming [8].
Linear- Quadratic-Gausiie
Motivated by practical settings where noise distributions may not be readily available or may not be Gaussian, this paper considers a discrete-time lnite-horizon generalization the LQG setting where noise distributions are unknown and are chosen adversarially from ambiguity sets characterized by Wasserstein distance and centered around nominal (Gaussian) distributions.
We show that, even under distributional ambiguity, the optimal control policy remains linear in the system's ''s observations. Our proof novel and does not rely on traditional recursive dynamic programming arguments. Instead, we re-parametrize the control policy in terms of the purified state observations and we derive an upper bound for the resulting minimax formulation by relaxing the ambiguity set (from Wasserstein ball into Gelbrich ball) while simultaneously restricting the controller linear dependencies. We then use convex duality prove that this upper bound matches
37th Conference on Neural Information Processing Systems (NcurIPS 2023).

---

a lower bound obtained by restricting the ambiguity set in the dual of the minimax formulation. This implies the optimality of linear output feedback controllers, thus generalizing the classic results to a distributionally robust setting.
We also find that the worst-case distribution is actually Gaussian, which leads to very efficient algorithm for finding optimal controllers. Specifically, we propose an algorithm based on the Frank- Wolfe first-order method that at every step solves sub-problems corresponding to classic LQG control problems, using Kalman filtering and dynamic programming. We show that this algorithm enjoys a sublinear convergence rate and is susceptible to parallelization. Lastly, we implement the algorithm leveraging Py PyTorch's 's automatic differentiation module and we find that it yields uniformly lower runtimes than direct method (based on solving semidefinite programs) across all problem horizons.
1.1. Literature Review
This paper is related the ample literature in control theory and engineering aimed designing controllers that are robust to noise. The classic LQR/LQG theory, developed in the 1960s, examined linear dynamical systems in either time or frequency domain, seeking to minimize combination of quadratic state and control costs (in time-domain) or the H2 norm of the system's 's transfer function (in frequency domain). Motivated by findings that LQG controllers do not provide the guaranteed robust stability properties of LQR controllers [15], much effort has been devoted subsequently to designing controllers that are robust to worst-case perturbations, typically evaluated in terms of the Hos norm of the system's 's transfer function (see. c.g., [16, 53] for a comprehensive review of Hos and H2 controllers). Because Hos controllers tend be overly conservative [32], various approaches have been proposed to balance the performance of nominal and robust controllers, e.g., by combining H2 and Ho approaches [7, 17]. A parallel stream of literature has considered risk-sensitive control [51], which minimizes an entropic risk measure instead of the expected quadratic cost. Although risk- sensitive control has distributinnally robust flavor (as the entropic risk measure equivalent to a distributionally robust quadratic objective penalized via Kullback-Leible divergence), risk-sensitive control models do not admit a distributionally robust formulation because the entropic risk measure is convex, but not coherent [22]. In contrast, our distributionally robust model provides a direct interpretation of the exact set of noise distributions against which the controller provides safeguards, and leads to computationally tractable framework for finding the optimal controller.
In this sense, our work is more directly related the literature on distributionally robust control, which seeks controllers that minimize expected costs under worst-case noise distributions [11, 33, 34, 41, 50, 52]. Closest to our work are [28, 33]. [33] proves the optimality of linear feedbackk control policies for related minimax LQR model with a Wasserstein distance but with perfect state observations. With perfect observations, the optimal policies the classic LQR formulation are independent of the noise distribution and are thus inherently already robust, so considering imperfect observations is what makes the problem significantly more challenging in our case. [28] studies minimax formulation based the Wasserstein distance with both state and observation noise but without any control policy, and focuses solely on the problem of estimating the states. Several papers have considered robust formulations with imperfect observations but for constrained systems [5, 34], which are more challenging, the common approach is to restrict attention to linear feedback policies for computational tractability, and without proving their optimality.
Also related is the recent literature stream on distributionally robust optimization using the Wasserstein distance [36]. Within this stream, the closest work is [38, 44], which consider the problem of minimax mean-squared-eroo estimation when ambiguity is modeled with Wasserstein distance from a nominal Gaussian distribution. Our proof builds on some ideas from these papers (e.g., relying on the Gelbrich distance in the construction of the upper bound), which it combines with ideas from control theory purified intput-feedback to obtain the overall construction. Also related is [2], which studies multistage distributionally robust problems with ambiguity sets given by a nested Wasserstein distance for stochastic processes and identifies computationally tractable cases. For broader overview of developments related to optimal transport and Wasserstein distance with an emphasis on computational tractability and applications in machine learning, we refer to [42].
Finally, our paper is also related to literature that documents the optimality of linear/affine policies in (distributionally) robust dynamic optimization models. [10, 30] prove optimality for one-dimensional linear systems affected by additive noise and with perfect state observations, but with general (convex) state and/or control costs, [27, 49] provide computationally tractable approaches to quantifying the

---

suboptimality of affine controllers in finite or infinite-horizon settings, and [9, 21, 25] characterize the performance of affine policies in two-stage (distributionally) robust dynamic models.
Notation. All random objects are defined on a probability space (L,F, P). Thus, the distribution of any random vector & S2 + R is given by the pushforward distribution P& = P { of P with respect to &. The expectation under IP is denoted by [x][]]. For any t E Σ4, we set [t] . .J.
2. Problem Definition
We consider a discrete-time linear dynamical system
1t A + Biu, wh VI E [1
(I)
with states Te E R", control imputs u, E R"", process noise w, E IR" and system matrices A₁ E Raxx and B, R=xm, The controller only has access to imperfect state measurements
Y/ C. + 02 At E 11]
corrupted by observation noise 2, E IR,, where C, E R and usually p S n (so that observ- ing Yt would not allow reconstructing Iy even there were no observation noise). The control inputs u, are causal, i.e., depend on the past observations yo..... .e but not on the future observa- tions ye+1...... 9r-1- More precisely, the set of feasible control inputs My is the set of random vectors (ug: U,,..., "T-1) where for every t there exists a measurable control policy 9:: Re((+1) -> I such that U, = (Pr (yo....... : ve). Controlling the system generates costs that depend quadratically on the states and the controls:
T-1 J=
Σ 1=
(IT Qiz, + uu, Rewe) + 2TQrTT,
where Qe E St and RG ST+ represent the state and input cost matrices, respectively, The exoge- nous random vectors zo, [w) and {vz 1=0 are mutually independent and follow probability distributions given by PP Pzo: De and {Pe,, ) respectively. As the control inputs are causal, the system equations (2) imply that 22, Ue and yr can be expressed as measurable functions of the exogenous uncertainties zo as well as w. and v.. 8 E [0], for every t. From now on we may thus assume without loss of generality that l = 1= x RRxX x RPx1 is the space of realizations of the exogenous uncertainties, F is the Borel a algeebra on l2 and TP = Pr, 6 (== w. 6 (o) where PP; s IP2 denotes the independent coupling of the distributions iP, and 1P2.
In this context, the classic LQG model assumes that P is known and Gaussian, and seeks E My that minimizes Wep [J]. Appendix SA reviews the standard approach for computing optimal control inputs by estimating states through Kalman filtering techniques and using dynamic programming.
In contrast, we assume that P is only known to belong to an ambiguity set W, and formulate distributionally robust LQG problem that seeks u E Uy to minimize the worst-case expected cost:
T-11 Σ (xTQQe Qiz, + R.u) - QQQTT
max EF PEW
We construct the ambiguity set W ball based on the Wasserstein distance. Specifically, assume that nominal Gaussian distribution P P... @ (oT (P. w, @ (fT v.) is available so
that PP 20 (0, Xo), Pw,, (O, Wt), and P.., (o, V.) for all t [T 1). and W is given by:
W Wa (S @ (o Wh.),
where
Wzo (P2.0 E P(R"): W(P... Pao) S Pro: Er [vo] 0}
Wm. = =PP,, E P(K") W(P.., Pw.,) < Pw.: Ep,, [w:] = 0}
Wr. (Po, P(R") W(f., P..) Po., Ez., [vo] =):
and W is the 2-Wasserstein distance Thus, by construction, all exogenous random vari- ables 2o, wo...., .., WT -1, Vo ...., UT-1 8T-1 independent under every distribution in W.
3

---

Definition (2- C-Wasserstein distance). The 2- Wasserstein distance between two distributions P1 and F2 on R with finite second moments is given by
inf W WPP,,P2) [EIIIFF
111 - -2-2211: (
where 11 P2) denotes the set of all couplings, that is, all joint distributions of the random variables E and &2 with marginal distributions P1 and P2, respectively.
Our model strictly generalizes the classic LQG setting,' which can be recovered by choosing Peo Pur, Pr.. O. The parameters p thus allow quantifying the uncertainty about the nominal model and building robustness to mis-specification. We emphasize that the Wasserstein ambiguity W contains many non-Gaussian distributions and it is not readily obvious that the worst-case distribution in (4) is in fact Gaussian. However, the set W is also non-convex, as it contains only distributions under which the exogenous uncertainties are independent, which makes the distributionally robust LQG problem potentially difficult to solve.
3. Nash Equilibrium and Optimality of Linear Output Feedback Controllers
We henceforth view the distributionally robust LQG problem as a zero-sum game between the controller, who chooses causal control inputs, and nature, who chooses a distribution I' e W. In this section we show that this game admits Nash equilibrium, where nature's 's Nash strategy is Gaussian distribution P W and the controller's Nash strategy is a linear output feedback policy based on the Kalman filter evaluated under IP*.
Purified Observations. Before outlining our proof strategy, we first simplify the problem formula- tion by re-parametrizing the control inputs in a more convenient form (following [5, 6, 27]). Note that the control inputs the IQQG formulation are subject to eyclic dependencies, as u4 depends on Vt. while Y/ depends on z; through (2), and depends again on ", through (1), etc. Because these dependencies make the problem hard to analyze, it preferable instead consider the controls as functions of a new set of so-called purified observations instead of the actual observations yt.
Specifically, we first introduce a fictitious noise-free system
2++1 = Acâ, Bru Vt TT and Dr Cjtt
AA [1]
with states à E R" and outputs út E IRP, which is initialized by u. = 0 and controlled by the same inputs u, as the original system (2). We then define the purified observation at time as 1) = ye- Úr and we use (To...,7rr- 1 to denote the trajectory of all purified observations.
-1)
As the inputs u, are causal, the controller can compute the fictitious state म; and output from the observations yo...., Vt- Thus, Te is representable as a function of yo. yt- Conversely, one show by induction that Y/ can also be represented as function of fo, Th- Moreover, any measurable function of Yoo
1/2 can be expressed as a measurable function of .,. ?h and vice- versa [27, Proposition II. So if we define U, as the set of all control inputs (uo, uI...., uT-1) so that de (/ ... Tr) for some measurable function t/ RP((+1) - R for every E [T 1). the above reasoning implies that U, = Uy.
In view of this, rewrite the distributionally robust LQG problem equivalently as
min 2,0.
max PEW
Ee [u] Ru + x TT@]]
p²
s.t.
u E My, x Hu + Gw,
Cr +
min n.u
max
Ep [u Ru + Qx]
PEW
(5)
s.t.
u E Uny
x = Hu + Gw,
where (2o, ... TT),
(vo, ... UT 11). y (yo. ........ -1)-1),
(No, wo..... WT 1).
V = (Uo..... -1), 77 = (mo..... 11), and R. Q. H, G and C' are suitable block matrices
lOur assumption that noise terms are aero-mean is consistent with the standard LQG model [8]. Requiring Erzg [zo] = 0 is assumed for clarity and without loss of generality.

---

(see Appendix $B for their precise definitions). The latter reformulation involving the purified observations 7) is useful because these are independent the inputs. Indeed, by recursively combining the equations of the original and the noise-free systems, one can show that // Dr + v for some block triangular matrix D) (see Appendix $B for its construction). This shows that the purified observations depend (linearly) on the exogenous uncertainties but not on the control inputs. Hence, the cyclic dependencies complicating the original system are eliminated in (5).
Subsequently, we also study the dual of (5), defined as
min Ep [u Ru - x Qx]
max PEW Uny
d'
S.L.
Hu + Gw.
The classic minimax inequality implies that p* 2 d' If can prove that p² d', that (5) has solution u* and that (6) has a solution IP*, then (u", ,P=) must be a Nash equilibrium of the zero-sum game at hand [43, Theorem 2]. However, because Ue is an infinitedimmesionall function space and W is an oninite-dimensionall non-convex set of non-parametric distributions, the existence of Nash equilibrium (in pure strategies) not all evident. Instead, our proof strategy will rely on constructing an upper bound for p* and lower bound for d', and showing that these match.
Upper Bound for p". We obtain an upper bound for p" by suitably enlarging the ambiguity set W and restricting the controllers u, to linear dependencies. We enlarge W by ignoring all information about the distributions W except for their covariance matrices, and by replacing the Wasserstein distance with the Gelbrich distance. To that end, we first define the Gelbrich distance on the space of covariance matrices.
Definition 2 (Gelbrich distance). The Gelbrich distance between the two covariance matrices Σ,, D2 E $ is given by
- TT Ei Ee -2 (2)x)=))*)
G(E) = /Tx
We are interested in the Gelbrich distance because of its close connection to the 2- Wasserstein distance. Indeed, it is known that the 2- Wasserstein distance between two distributions with zero means is bounded below by the Gelbrich distance between the respective covariance matrices.
Proposition 3.1 (Gelbrich bound [24, Theorem 2.1]). For any two distributions P1 and P2 on Rd with zero means and covariance matrices E1, Σg S4, respectively, we have W (FI,, IP2) 2 G(E,, E2).
Recalling that Xo, We and V. respectively denote the covariance matrices for Ip, We and 2 under the nominal distribution PP, we can then define the following Gelbrich ambiguity set for the exogenous uncertainties:
g = Gzo @ (8 gw.) @ (8T-)e))
where
Gmp (Pao ug P(RR) Eer, [2o] 0, Ep Ex Xo. G(Xo, Xo) S Pzo
Gw. = {Pmm E P(R"): Er.,, [wm]] = 0, Er Br(wew:] Wr, G(W,, W) S pmo,)
Gr, = (Pr, E P(R"): EE,, [ve] O, EF(ur ] = V., G(V,, V.) p..).
By construction, the random vectors To. {wi {w:} ITo and {uz) T-1 are thus mutually independent under any IP E G. In addition and as direct consequence of Proposition 3. G constitutes an outer approximation for the Wasserstein ambiguity set W, as summarized in the next result.
Corollary (Gelbrich hull). We have W C G.
Because g covers W. we henceforth refer to it the Gelbrich hull of the Wasserstein ambiguity set W. To finalize our construction of the upper bound on p², we focus on linear policies² of the form
Technically, the policies are affine because they include constant term, but we retain the more common terminology that focuses the dependencies.
5

---

u == v), where q = .1- and U is a block lower triangular matrix
Un.0
U1
11,1
U
UT--0
UT-1TT-1
The block lower triangularity of U ensures that the corresponding controller is causal, which in turn ensures that u E Up. In the following, we denote by U the set of all block lower triangular matrices of the form (7). An upper bound on problem (5) can now be obtained by restricting the controller's feasible set to causal controllers that are linear in the purified observations 1) and by relaxing nature's feasible set to the Gelbrich hull 9 of W. The resulting bounding problem given by
min U.q.zau S.t.
max Ep [u T Ru + Qu]
p' {
PEG
(8)
UEu, a- U (Dw u), I Hu + Gw.
As we obtained (8) by restricting the feasible set of the outer minimization problem and relaxing the feasible set of the inner maximization problem in (5), is clear that 5" ž p=. Recall also that problem (5) constitutes an infinitediimensional zero-sum game, where the agents optimize over measurable policies and non parammeetriic distributions, respectively. In contrast, the next proposition shows that problem (8) is equivalent to a mnite-dimensionll zero-sum game.
Proposition 3.2. Problem (8) is equivalent to the optimization problem
min max Tr ((DTU (R+ HHTHHTHHUU++ 2GT QHUD +GTDD) gERPT Veu
D"
W VERv
(9)
where
W = diag( Xo, Wo, Wr-1), Xo E W. We [T [ }
sM(T-1);
Gw W E s
G(Xo, xo)2 S Pr. G(W, W.) Pe, Vt [T Gv = EH = diag(Vo. Vr-1), E E G(V,, vi) po We T 1]] 1] }
We emphasize that Proposition 3.2 remains valid even if the nominal distribution P fails to be normal.
Note also that, while nature's ' feasible set in (8) is non-convex due to the independence conditions, the sets Gw and Gv convex and even semidefinite representable thanks to the properties of the squared Gelbrich distance." 3 By dualizing the inner maximization problem, one can therefore reformulate the minimax problem (9) as convex semidefinite program (SDP). Even though this SDP is computationally tractable in theory, it involves (T (mp + + p² )²) decision variables. For practically interesting problem dimensions, it thus quickly exceeds the capabilities of existing solvers.
Lower Bound for d'. To derive a tractable lower bound on d", we restrict nature's 's feasible set to the family Wy of all normal distributions in the Wasserstein ambiguity set W. The resulting bounding problem is thus given by
max PEWN
min
Ep [u Ru + x TTQz] s.t. u E Un,
d'
(10)
Hu + Gw.
As obtained (10) by restricting the feasible set of the outer maximization problem in (6), it is clear that d' < d' Next, we show that (10) can be recast as a finite-dimensional zero-sum game. This result critically relies the following known fact regarding the 2- Wasserstein distance between two normal distributions, which coincides with the Gelbrich distance between their covariance matrices.
Proposition 3.3 (Tightness for normal distributions [26, Proposition 7]). For any two normal distri- butions P, = =0 En) and P2- (0, Σz) with zero means we have W(P,, P2) G(E,, E2).
With this, we can provide a finite-dimensional reformulation, as summarized in the next result.
Note that the ambiguity sets Gw and Gv appearing in (9) involve the squared Gelbrich distance, G(\,, ). The reason is that G(2,, )²)² is known to be jointly convex in Mr, Va and semidefinite representable [38, Proposition 2.3], unlike the Gelbrich distance G(21, 22) itself, which is generally non-convex.

---

Proposition 3.4. Problem (10) is equivalent to the optimization problem
min Tr UT((TTUT (R+H' QH) QHUD - GEReT
o' {
max WEOW VEGV
(1I)
DEu
+Tr +Tr((U (R+ H H QHYU)))+ (R+ H HHQTHH
where Gw and Gv are defined exactly as in Proposition 3.2.
Proposition 3.4 relies on Proposition 3.3 and thus fails to hold unless PP is normal. Also, one can again reformulate (11) as tractable SDP by dualizing the inner minimization problem.
Conclusions. Propositions 3.2 and 3.4 reveal that problems (9) and (1I) are dual to each other, that is, they can be transformed into one another by interchanging minimization and maximization. The following main theorem shows that strong duality holds irrespective of the problem data.
Theorem 3.5 (Strong duality of (9) and (11)). We have D" = =d"
Theorem 3.5 follows immediately from Sion's classic minimax theorem [45], which applies be- cause Gw and Gy are convex well as compact thanks [38, Lemma A.6].
By weak duality and the construction of the bounding problems (9) and (11), we trivially have d' < d' S p¹ < p". Theorem 3.5 reveals that all of these inequalities are in fact equalities, each of which gives rise to non-trivial insight. The first key insight that (5) and (6) are strong duals.
Corollary (Strong duality of (5) and (6)). We have p" = d
We stress that, unlike Theorem 3.5, Corollary 2 establishes strong duality between two infinite- dimensional zero-sum games. The second key implication of Theorem 3.5 is that the distributionally robust LQG problem (5) is solved by linear output-feedback controller.
Corollary 3 (The controller's Nash strategy is linear in the observations). There exist U U and q* E Th such that the distributionally robust LQG problem (5) solved by u² = q² + U'y.
The identity = readily implies that (5) is solved by causal controller that is linear in the purified observations. However, any causal controller that is linear in the purified observations be reformulated exactly causal controller that is linear in the original observations y and vice-versa [6, Proposition 3]. Thus, Corollary 3 follows. The third key implication of Theorem 3.5 is that the dual distributionally robust LQG problem solved by normal distribution.
Corollary 4 (Nature's 's Nash strategy is normal distribution). The dual distributionally robust LQG problem (6) solved by a distribution P E WWy.
Corollary 4 is a direct consequence of the identity d' = d' Note that the optimal normal distribu- tion IP* is uniquely determined by the covariance matrices W and V* of the exogenous uncertain parameters, which can be computed by solving problem (11). That the worst-case distribution is actually Gaussian is not a-priori expected and is surprising given that the Wasserstein ball contains many non-Gaussian distributions.
4. Efficient Numerical Solution of Distributionally Robust LQG Problems
Having proven these structural results, we next turn attention to the problem of finding the optimal strategies. Our next result shows that, under a mild regularity condition, the optimal controller u* of the distributionally robust LQG problem (5) can be computed efficiently from P*.
Proposition 4.1 (Optimality of Kalman fitterbased feedback controllers). 18 > for all L E [T- then problem (6) is solved by a Gaussian distribation F. under which v; has a covariance matrix V > 0 for every t E [T 11, and (5) is solved by the optimal LQG controller corresponding to Ip* Additionally, the optimal value of problem (9) and its strong dual (11) does not change if we restrict Gw and GV to Giv and Gt, respectively, where
gi =[We Iw Xo Z Amin (8)) W, Z Amin
g4 = VV E Gv V/2 Amin Amin (V)) E (T- 11.
7

---

This implies that the optimal controller can be computed by solving a classic LQG problem conre sponding to nature's 's optimal strategy P", which can be done very efficiently through Kalman filtering and dynamic programming (see Appendix $A for details). It thus suffices to design an efficient algorithm for computing IP' which is uniquely determined by the covariance matrices (W+V*) that solve problem (11). To this end, we first reformulate (11) as f(W, V),
max W WEGW.
(12)
where restrict Gw and Gv to Sity and G respectively, due to Proposition 4.1, and where f(W, V) denotes the optimal value function of the inner minimization problem in (11). As (11) is a reformula- tion of (10) and the family of all causal purified output-feedback controllers matches the family of
causal output-feedback controllers, f(W V) can also be viewed as the optimal value of the classic LQG problem corresponding to the normal distribution P determined by the covariance matrices W and V. These insights lead to the following structural result.
Proposition 4.2. f(W, V) is concave and 3-smooth in (W,V) @W G for some B 0.
By Proposition 4.2, it is possible address problem (12) with a Frank- Wolfe algorithm [13, 18, 19, 20, 23, 35]. Each iteration of this algorithm solves a direction-findee subproblem, that is, a variant of problem (12) that maximizes the first-order Taylor expansion of f(W, V) around the current iterates. (W, V), Lw W) (Vv/(W, V), Ta - -V)
max Lw EBw Lo EGy
(13)
The next iterates are then obtained by moving towards a maximizer (Lw, LV) of (13), i.e. update (W.V) - (W, V) + (Li W, -.
where a is an appropriate step size. The proposed Frank- algorithm enjoys very low per- iteration complexity because problem (13) is separable. To see this, we reformulate (13) as
T-1
(Vx,f(W.V), -xx Σ (W, -W)) + (W, V), Lv, - VV)
max Tw
s.t. G(Ixxo,x1)) 2 Pa:. G(Iw., W.)² S Pu.. G(Iv, 1i)2 S pe, Lx, = Amin (Xo)I. Lw, 2 Amin (W.). Lv, F Amin (V)) Vt 11.
WE [ -
Hence, (13) decomposes into 27' + 1 separate subproblems that can be solved in parallel. That is, for any matrix Z (Xo, Wo...., Wr-1, Vo...., 1r we solve a separate subproblem of the form
{ (V2 L2 - 2) G(Lz, 2)² < p2
max
(14)
Amin (2)
These subproblems can be reformulated as tractable SDPs and are thus amenable to efficient off- the-shelf solvers. By [38, Theorem 6.2]. however, one can exploit the structure the Gelbrich distance in order to reduce (14) to a univariate algebraic equation that can be solved to any desired accuracy 8 >> 0 by a highly efficient bisection algorithm. We say that L22 is 5-approximate solution of problem (14) for some E (0. 1) if L is feasible in (14) and if
2 V), 2 Z) 2 8(V2f(W, (W, LZ - Z),
where Lt is an exact maximizer of (14). Note that, by the concavity of f(W, V), the inner product on the right-hand side is nonnegative and vanishes if and only if Z maximizes f(W, V) over the feasible set of (14). For further details we refer to Appendix ŞE in the supplementary material.
Remark (Automatic differentiation). Recall that f(W, V) coincides with the optimal value of the LQG problem corresponding to the normal distribution 12 determined by the covariance matrices W and V. By using the underlying dynamic programming equations, f(W, V) can thus be expressed in closed form as serial composition of O(T) rational functions functiol. (see Appendix $A for details). Hence, 2 V) can be calculated symbolically for any Z (Xp. Wo...., Wr--, Vo..... .V by repeatedly applying the chain and product rules. However, the resulting formulas are lengthy and cumbersome. We thus compute the gradients numerically using backpropagation. The cost of evaluating Vz](W.V) V) is then of the same order of magnitude as the cost of evaluating I JW, V).
A detailed description of the proposed Frank- method is given in Algorithm 1 below.
By [31, Theorem and Lemma 7], which applies thanks to the structural properties of S(W, V) established in Proposition 4.2, Algorithm attains a suboptimality gap of within ((1/e) iterations.
8

---

Algorithm Frank-Wolf algorithm for solving (12)
Input: initial iterates W, V, nominal covariance matrices W, V. oracle precision d (0,
1: set initial iteration counter k
2:
while stopping criterion is not met do for Z E (Xo, Wo; Wr-1, Vo.... Vr-1) do in parallel
3:
4:
compute Vzf(W.)) find a d-approximate solution L3 (14)
5:
6:
end
7:
f (Vw V), L. W)+ (Vv V). Li - V) (Li W. -
8:
(W.V) (
9:
end while
10:
Output: W and V
10a
10a
Execution time (6)
10
Ostimalty
10-1
10
10
10
MOSEK
10
10-2
Frank-Wolfe
10-5
10a
10a
102
10a
10a
102
Time horizon
# Iterations
(a)
(b)
Figure 1: (a) Execution time for MOSEK and Frank-Wolfe algorithm over 10 simulation runs as function of the horizon T (solid lines show the mean and the shaded areas correspond to 1 standard deviation). (b) Convergence of optimality gap for Frank-Wold algorithm with horizon = 10.
5. Numerical Experiments
All experiments are run on an Intel i7-8700 CPU (3.2 GHz) machine with 16GB RAM. All linear SDP problems modeled in Python 3.8.6 using CVXPY [1, 14] and solved with MOSEK [37]. The gradients of (W, V)) are computed via Pymanopt [48] with PyTorch''s automated differentiation module [39, 40).
Consider a class of distributionally robust LQG problems with n = m = = 10. We set A = 0.1 A to have ones the main diagonal and the superdiagonal and zeroes everywhere else (Ai,) 1 if = j or i = j and Ar,; = otherwise), and the other matrices to B, = C, = Q = R₁ = la- The Wasserstein radii are set to P..a = Pv. = Pe, 10m! ! The nominal covariance matrices of the exogenous uncertainties are constructed randomly and with eigenvalues in the interval [1, 2] (so as to ensure they are positive definite). The code publicly available the Github repository https: Control.
//gi ligithub... cemtub..oo//RR -coo//AO--oo
The optimal value of the distributionally robust LQG problem (5) can be computed by directly solving the SDP reformulation of (1 with MOSEK or by solving the nonlinear SDP (12) with our Frank- Wolfe method detailed in Algorithm 1. We next compare these two approaches in 10 independent simulation runs, where we set a stopping criterion corresponding to an optimality gap below 10-3 and we run the Frank Wolffe method with 0.95. Figure la illustrates the execution time for both approaches as function of the planning horizon T: runs where MOSEK exceeds 100s are not reported. Figure visualizes the empirical convergence behavior of the Frank- Wolfe algorithm. The results highlight that the Trank-Wolff algorithm achieves running times that are uniformly lower than MOSEK across all problem horizons and is able to find highly accurate solutions already after small number of iterations (50 iterations for problem instances of time horizon T = 10).

---

6. Concluding Remarks and Limitations
In view of the popularity of LQG models, the results in this work carry important theoretical and practical implications. Despite considering a generalization of the classic LQG setting where the noise affecting the system dynamics and the observations follows unknown (and potentially non- Gaussian) distributions, and that highly efficient methods can be adapted to tackle this more realistic (and more challenging) problem. Specifically, that control policies depending linearly on observations continue to be optimal and that the worst-case distribution turns out to be Gaussian is surprising from a theoretical angle and also has direct practical implications, because it allows leveraging the highly efficient Kalman filter in conjunction with dynamic programming and a Woolfe method to design an efficient computational procedure for solving the problem.
findings suggest that certain classic structural results continue to hold
The results also raise several important questions that warrant future exploration. First, it would be highly relevant to consider extensions where the system matrices are also affected by uncertainty, as this captures many applications of practical interest in, e.g., reinforcement learning or revenue management. Second, it would be worth exploring infinite horizon setting relaxing the assumption that the nominal distribution is Gaussian, as both assumptions may be limiting the practical appeal of the framework. Third, could also attempt to prove structural optimality results or design novel algorithms for generating high-quality suboptimal solutions for the more general setting involving constraints on states and/or control inputs. Lastly, one could improve the present algorithmic proposal by exploiting topological properties of the objective so as to guarantee linear convergence rates in the Frank-Wolf procedure.
This research was supported by the Swiss National Science Foundation under the NCCR Automation, grant agreement 51NF40_I 180505. Dan A. Iancu would like acknowledge INSEAD for financial support during the duration of the project.
Acknowledgements.
References
[1] Akshay Agrawal, Robin Verschueren, Steven Diamond, and Stephen Boyd. A rewriting system for convex optimization problems. Journal of Control and Decision, 5(1):42-60, 2018.
[2] Rohit Arora and Rui Gao. Data-driven multistage distributionally robust optimization with nested distance: Time consistency and tractable dynamic reformulations. Available at Optimiza- tion Online, 2022.
[3] François Auger, Mickael Hilairet, Josep M. Guerrero, Eric Monmasson, Teresa Orlowska- Kowalska, and Seichhio Katsura. Industrial applications of the Kalman filter: A review. IEEE Transactions on Industrial Electronics, 60(12):5458-547 2013.
[4] Richard Bellman, Some inequalities for the square root of a positive definite matrix. Linear Algebra and its Applications, 1(3):321-32,, 1968.
[5] Aharon Ben-Tal, Stephen Boyd, and Arkadi Nemirovski. Control of uncertainty-affecter discrete time linear systems via convex programming. Available at Optimization Online, 2005.
[6] Aharon Ben-Tal, Stephen Boyd, and Arkadi Nemirovski. Extending scope of robust optimiza- tion: Comprehensive robust counterparts of uncertain problems. Mathematical Programming, 107(1):63-89, 89, 2006.
[7] Denis S. Bernstein and Wassim M. Haddad. LQG control with an Hos performance bound: A Riccati equation approach. In American Control Conference, pages 796-802, 1988.
[8] Dimitri Bertsekas. Dynamic Programming and Optimal Control, volume I. Athena Scientific, 2017.
[9] Dimitris Bertsimas and Vineet Goyal. On the power and limitations of affine policies in two-stage adaptive optimization. Mathematical Programming, 134(2):411531,, 2012.
[10] Dimitris Bertsimas, Dan A. Iancu, and Pablo Parrilo. Optimality of affine policies in multistage robust optimization. Mathematics of Operations Research, 35(2):363-394, 2010.
10

---

[11] Dimitris Bertsimas, Dan A. lancu, and Pablo A. Parrilo. A hierarchy of near-optimal policies for multistage adaptive optimization. IEEE Transactions on Automatic Control, 56(12):2809-288) 2011.
[12] Shen- Yong Chen. Kalman filter for robot vision: A survey, IEEE Transactions on Industrial Electronics, 59(11):449-4420, 2012.
[13] Vladimir F. Demyanov and Aleksandr M. Rubinov. Approximate Methods in Optimization Problems. Elsevier, 1970.
[14] Steven Diamond and Stephen Boyd. CVXPY: A Python-embeddded modeling language for convex optimization. Journal of Machine Learning Research, 17(83)11-5, 2016.
[15] John Doyle. Guaranteed margins for LQG regulators. IEEE Transactions Automatic Control, (544):756-757 1978.
[16] John Doyle, Keith Glover, Pramod Khargonekar, and Bruce Francis. State-space solutions to standard H2 and Hoo control problems. American Control Conference, pages 1691-1696, 1988.
[17] John Doyle, Kemin Zhou, and Bobby Bodenheimer. Optimal control with mixed H2 and Hox performance objectives. In American Control Conference, pages 2065-2070, 1989.
[18] Joseph C. Dunn. Rates of convergence for conditional gradient algorithms near singular and nonsingular extremals. SIAM Journal on Control and Optimization, 17(2):1872211, 1979.
[19] Joseph C. Dunn. Convergence rates for conditional gradient sequences generated by implicit step length rules. SIAM Journal on Control and Optimization, 18(5):473-487, 1980.
[20] Joseph Dunn and S. Harshbarger. Conditional gradient algorithms with open loop step size rules, Journal of Mathematical Analysis and Applications, 62(2):432-444, 1978.
[21] Omar EI Housni and Vineet Goyal. On the optimality of affine policies for budgeted uncertainty sets. Mathematics of Operations Research, 46(2):674-711, 2021.
[22] Hans Föllmer and Alexander Schied. Stochastic Finance: An Introduction in Discrete Time. de Gruyter, 2011.
[23] Marguerite Frank and Philip Wolfe. An algorithm for quadratic programming. Naval Research Logistics, 111-2):95-110, 1956.
[24] Matthias Gelbrich. On formula for the T,2 Wasserstein metric between measures on Euclidean and Hilbert spaces. Mathematische Nachrichten, 147(1):185-20) 1990.
[25] Angelos Georghiou, Angelos Tsoukalas, and Wolfram Wiesemann. On the optimality of affine decision rules in robust and distributionally robust optimization. Available at Optimization Online, 2021.
[26] Clark R. Givens and Rac M. Shortt. A class of Wasserstein metries for probability distributions. Michigan Mathematical Journal, 31(2):331-240, 1984.
[27] Michael J. Hadjiyiannis, Paul J. Goulart, and Daniel Kuhn. An efficient method to estimate the suboptimality of affine controllers. IEEE Transactions on Automatic Control, 56(12):2841 22852283,, 2011.
[28) Bingyan Han. Distributionally robust Kalman filtering with volatility uncertainty. arXiv preprint 202.::332..55993 2023.
[29] Lars Peter Hansen and Thomas J. Sargent, Robust estimation and control under commitment. Journal of Economic Theory, (24(2):258-300 2005.
[30] Dan A. Iancu, Mayank Sharma, and Maxim Sviridenko. Supermodularity and affine policies in dynamic robust optimization. Operations Research, 61(4)9941-956, 2013.
11

---

[31] Martin Jaggi. Revisiting Frank-Wolfe: Projection-free free sparse convex optimization. In Interna- tional Conference on Machine Learning pages 427-435, 2013.
[32] Aren Karapetyan, Andrea lannelli, and John Lygeros. On the regret of Hox control. In IEEE Conference on Decision and Control, pages 6181-6186, 2022.
[33] Kihyun Kim and Insoon Yang. Distributional robustness in minimax linear quadratic control with Wasserstein distance. SIAM Journal on Control and Optimization, 61(2):458-445 2023.
[34] Georgios Kotsalís, Guanghui Lan, and Arkadi S. Nemirovski. Convex optimization for finite- horizon robust covariance control of linear stochastic systems. SIAM Journal on Control and Optimization, 59(1):296-315 2021.
[35] Evgenii Levitin and Boris Polyak. Constrained minimization methods. USSR Computa- tional Mathematics and Mathematical Physics, 6(5):1-50, 1966.
[36] Peyman Mohajerin Esfahani and Daniel Kuhn. Data-driven distributionally robust optimization using the Wasserstein metric: Performance guarantees and tractable reformulations, Mathemati- cal Programming, (71(12))115-166 2018.
[37] MOSEK ApS. The MOSEK Optimization Toolbox. Version 2019.
[38] Viet Anh Nguyen, Soroosh Shafieezadeh-Abadeh, Daniel Kuhn, and Peyman Mohajerin Esfa- hani. Bridging Bayesian and minimax mean square error estimation via Wasserstein distribu- tionally robust optimization. Mathematics of Operations Research, 48(1):1-37, 2023.
[39] Adam Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary De Vito, Zeming Lin, Alban Desmaison, Luca Antiga, and Adam Lerer. Automatic differentiation in PyTorch. In NIPS 2017 Autodiff Workshop, 2017.
[40]
Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Kopf, Edward Yang. Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala. Pytorch: An imperative style, pigh-performance deep learning library. In Advances in Neural Information Processing Systems, pages 8037, 2019.
[41] Ian R. Petersen, Matthiew R. James, and Paul Dupuis. Minimax optimal control of stochastic uncertain systems with relative entropy constraints. IEEE Transactions on Automatic Control, 45(3):398 412, 2000.
[42] Gabriel Peyré and Marco Cuturi. Computational optimal transport: With applications to data science. Foundations and Trends in Machine Learning, 1165-6):355-677 2019.
[43] R. Tyrrell Rockafellar. Conjugate Duality and Optimization. SIAM, 1974.
[44] Soroosh Setileezzaee-Abbadee Viet Anh Nguyen, Daniel Kuhn, and Peyman Mohajerin Esfa- hani. Wasserstein distributionally robust Kalman filtering. In Advances in Neural Information Processing Systems, pages 8483-8492, 2018.
[45] Maurice Sion. On general minimax theorems. Pacific Journal of Mathematics, 8(1):171-176. 1958.
[46] Joelle Skaf and Stephen P. Boyd. Design of affine controllers via convex optimization. 1EEE Transactions on Automatic Control, (5(1)),5(111)))))))1) 1):24)) 2010.
[47] Emanuel Todorov and Michael I. Jordan. Optimal feedback control as theory of motor coordination. Nature Neuroscience, 5(11):1226-123,, 2002.
[48] James Townsend, Niklas Koep, and Sebastian Weichwald. Pymanopt: A Python toolbox for optimization on manifolds using automatic differentiation. Journal of Machine Learning Research, 17(137):1-5, 2016.
12

---

[49] Bart P. G. Van Parys, Paul J. Goulart, and Manfred Morari. Infinite horizon performance bounds for uncertain constrained systems. IEEE Transactions on Automatic Control, 58111)2803-2817 2013.
[50] Bart P. G. Van Parys, Daniel Kuhn, Paul J. Goulart, and Manfred Morari. Distributionally robust control of constrained stochastic systems. IEEE Transactions on Automatic Control, 61(2):430-44 2016.
[51] Peter Whittle. Risk-sensitive Innaa//uaaraaii/Gaussiin control. Advances in Applied Probability, 13(4):764-777, 1981.
[52] Insoon Yang. Wasserstein distributionally robust stochastic control: A data-driven approach. IEEE Transactions on Automatic Control, 66(5):3863-3870, 2021.
[53] Kemin Zhou and John C. Doyle. Essentials of Robust Control. Prentice Hall, 1998.
13

---

Appendix
The supplementary material is structured as follows. Appendix ŞA presents the well-known solution to the classic LQG problem using dynamic programming and Kalman Filter estimation. Appendix $B provides the definitions of the stacked system matrices utilized in the compact formulation (5) of the distributionally robust LQG problem. Appendix &C contains the proofs of the formal statements in the main text and provides additional technical results. Appendix $D derives the SDP reformulation of the dual problem (11). Appendix $E, finally, elaborates on the bisection algorithm used for solving the lincarization oracle of the Frank-Won Wolfe algorithm.
A. Solution of the LQG Problem
The classic LQG problem can be solved efficiently via dynamic programming; sec, c.g., [8]. That is, the unique optimal control inputs satisfy uz K, , for every E 1]. where K, = is the optimal feedback gain matrix, and ã = 1p [ (v0. ................................ yr] is the minimum mean-squared-eroo estimator of 2z given the observation history up time Thanks to the celebrated separation principle, K2 can be computed by pretending that the system is deterministic and allows for perfect state observations, and &, can be computed while ignoring the control problem.
To compute K₁, one first solves the deterministic I. LQR problem corresponding to the LQG problem hand. Its value function x. Prsx time is quadratic in 2t, and obeys the backward recursion
P. A P1 + Q A, P++ (R, + B, P.+ BB, Pe114 E
(A.15a)
initialized by Pr Qr- The optimal feedback gain matrix K₁ then be computed from P.+1
Kr (R, + B, P.1 B, P441A Vt -T]
1].
(A.15b)
Since 2y and (yo..... . Mt) are jointly normally distributed, the minimum mean-squaree muaa-sqquaed-errrr estima- tor &n can be calculated directly using the formula for the mean of a conditional normal distribution. Alternatively, however, one can use the Kalman filter to compute 2; recursively, which is significantly more insightful and efficient. The Kalman filter also recursively computes the covariance matrix Σ,
of conditional on ....... Y/ and the covariance matrix Σ,+1je of rt+1 conditional on yu..... y/ evaluated under P. Specifically, these covariance matrices obey the forward recursion
Σ, Emp=11 - Eeje-IC (C + Vo) Vt ET Σ A,E,A + W.
(A.16)
initialized by Moj-1 Xo- Using D then define the Kalman filter gain
L. = N, V Vt W
which allows us to compute the minimum mean-squared-erroo estimator via the forward recursion
$c+1 A, Be + T4+1 (Vt+1 C41(AA, + B,,)) At [1-
initialized by 2o Lowo. One can also show that the optimal value of the LQG problem amounts to
Tr((Q. P.E)) + Σ T(PP(AAEE Wi-1)) Tr(PoXoo)
(A.17)
B. Definitions of Stacked System Matrices
The stacked system matrices appearing in the distributionally robust LQG problem (5) are defined as follows. First, the stacked state and input cost matrices Q (T(++1) and R E SmT are set to
Qu
Ro
Q1
R1
and
R
QT
RT
14

---

respectively, Similarly, the stacked matrices appearing in the linear dynamics and the measurement
equations C E I xx(++), G E (=(++))x((+1)) and H E p=([+1)xmT are defined as
Co
0 G
A
A
A
C
G
Cr-1 0
A
AT
A7.
and
Bo
0
ABBo
A
H =
Af Bo
A B
AZ.
respectively, where A I1 Aş for every t and A Ia for t.
Using the stacked system matrices, we can now express the purified observation process 1/ as a linear function of the exogenous uncertainties and that is not impacted by also [5, 46]
Lemma B.1. We have 7 Dw + where D = CG.
Proof of Lemma B.I. The purified observation process is defined 7) = ÿ. Recall now that the observations of the original system satisfy y = Cr + v. Similarly, one readily verifies that the observations of the fictitious noise-free system satisfy Next, recall that the state of the original system satisfies x = Hu u + Gw, and note that the state the fictitious noise-free system satisfies & Hu. Combining all of these linear equations finally shows that t cancels out and that 1) CGw + tu= Dw t
C.. Thus, have 7/ C(a 4)
C. Proofs
C.1. Additional Technical Results
It is well known that every causal controller that is linear in the original observations can be reformulated as a causal controller that linear purified observations r) and vice versa [5, 46). Perhaps surprisingly, however, the one-to-one transformation between the respective coefficients of and rj not linear. To keep this paper self-contained, we review these insights in the next lemma.
Lemma C.1. If u = Un + q for some (T + UCH)-UU 'U and q' (T +T+ (-UCH) Conversely, if = U'y + q for some U E U and
u and q ReT then u. U'y + for U"
Lf E RT then u Un + q for U =I -1U and q = (T U'CH)-14
Proof of Lemma C.I. If u = Un + for some U E 21 and E R then we have
u = Un + = U 0) + - Uy UC + q= -y UCHUCH
where the second equality follows from the definition of r7. the third equality holds because y = Ca+t,
and the last equality exploits our earlier insight that ÿ C... The last expression depends only y and u. Solving for u yields u = U'y y + q', where U = (I UCH) and = (I+ Note that (I UCH) is indeed invertible because + UCH is lower triangular matrix with all diagonal entries equal to one, ensuring determinant of one.
Similarly, if U'y + o' for some U" E l and of E R then we have
U' + = + "0
Solving for v yields Un + q. where U (T- -C-UUC- U and (I-U(CC Note again that (T C) H)) is indeed invertible because (I U'CH) a lower triangular matrix with all diagonal entries equal to one.
D
15

---

C.2. Proofs of Section 3
Proof of Proposition 3.2. In problem (8), both u and z are linear in w and v, i.e. U+ Uv and x Hu + Gw Hq + HU Dw + HUv + Gw. By substituting the linear representations of
and x into the objective function of problem (8), we obtain the following equivalent reformulation.
min 9E3PT VEu
max FEG
Ep [w= U (R + H QHIUU + 2 UJ H' QG H'QG- G QG) w]
+ E [ (
(R + H QU)
(R+
For any fixed EG we can express the expectation in the objective function of the above problem in terms of the covariance matrices W Ep (ww' " and V = =Ellv'' 'J. Thus, the problem becomes
min 9ER"T
Tr ((DTU (R+ QU)UD+ 2G QHUU ())
W.V.P
UEM
+ To ((U ( Q)) +4
(A.18)
s.t. PEG, W Ep[www)] J, V -Elev'
Recall now the definition G, and note that the requirements G(Xo, Xo) S Pzo G(W,, Wi) Pwu and (V,, V.) S Pr; are equivalent the convex constraints G(Xo, Xo)² S P2. G(W,, W.)² pe. and G(V, V/)2 S p respectively, for all & E [T' 1j. The definition of G also implies that
W 1pp[ww 1 = diag(Xo, Wo...., Wr-1) and V = Ep(ev = diag (Vo.................. ..., Vr-1).
Problem (A.18) thus constitutes a relaxation of problem (9). Indeed, the feasible set of the inner maximization problem in (A.18) is subset of the feasible set of the inner maximization problem in (9). Moreover, for any W and V feasible the inner maximization problem in (9), the distribution T-
P = P2, @ (8) TD Pw.) @ (8 defined through Pz., = (N, Pw, = -V W,) and P.. (0, Vi). [T' 1]. is feasible in the inner maximization problem in (A.18) with the same objective value. The relaxation is thus exact, and the optimal values of (8), (9) and (A.18) coincide.
Proof of Proposition 3.4. Recall that the space U4 of all causal output-feedback controllers coincides with the space U1, of all causal purified oueput-feedback controllers. We can thus replace the feasible set Ut of the inner minimization problem in (10) with Ug. Hence, for any fixed E WAs, the inner minimization problem in (10) constitutes classic LQG problem. By standard LQG theory [8], it is solved by a linear output-feedback controller of the form u U'y-q' for some U' E U and o' E RIF"; see also Appendix $A. Lemma C.I shows, however, that any linear output-feedback controller can be equivalently expressed as linear purified-output feedback controller of the form u Un + for some U U and q E R² In summary, the above reasoning shows that the feasible set of the inner minimization problem in (10) can be reduced to the family of all linear purified-output feedback controllers without sacrificing optimality. Thus, problem (10) equivalent
max PEWN
min q,U,a,u s.t.
Ee [u' Ru +x Qz]]
Ueu, u u= -= Un, x= Hu + Gw.
Using a similar reasoning as in the proof of Proposition 3.2, we can now substitute the linear representations of and x into the objective function and reformulate the above problem as
max W.V.P
min GERPT Vell
Tr ((DTU UT (R+HT QH)UD++ T QHUD +G QG) )))))
+ (R+ +q (R+HTO Q11)g
PEWN; W = Ep|ww]] J] V Ep[uoT].
As Wx contains only normal distributions, Proposition 3.3 implies that W(P.g mg Pao) = Xo). W(P.., we, ((,,,,,,.,)>> = G(W,, Wi) and W(F., v., PP..) = G(V,, V3) for all E [T' - 1]. We may thus replace the requirement W(Pro (0.. P.n) S Pang in the definition of WA by G(Xo, Xo) S Pzo, which is equivalent the convex constraint G G(Xo, Xo)² 2= Pa.. The conditions the marginal distributions of w, and u, E [T' 1], admit similar reformulations. The definition of Ww also implies that
W = [ww ] diag(Xo, Wo..., ..., Wr-1) and V Ep(ve' J diag (V.................... VT-1).
16

---

Thus, the feasible set of the outer maximization problem in (11) constitutes a relaxation of that in (10). One readily verifies that the relaxation is exact by using similar arguments as in the proof of Proposition 3.2. Thus, the claim follows.
Proof of Theorem 3.5. By Proposition 3.2, pT coincides with the minimum of (9). Similarly, by Proposition 3.4 d" coincides with the maximum of (11). Note that problems (9) and (11) only differ by the order of minimization and maximization. Note also that 24 is convex and closed, Gw and Gy are convex and compact by virtue of [38, Lemma A.6], and the (identical) trace terms in (9) and (11) are bilinear in (W, V) and (U, The claim thus follows from Sion's 's minimax theorem [45]. 17
C.3. Proofs of Section 4
Note that Proposition 4.1 is consistent with Corollary 3 because the optimal LQG controller corre- sponding to PP" is linear in the past observations.
Proof of Proposition 4.1. By [38, Lemma A.3], the inner problem in (9) admits a maximizer (W*, V*) with XI ž Amin (Xo) as well as W7 = Amin (W)) and Vi m Amin (V)) for all t [T- 1]. Thus, the optimal value of problem (9) and its strong dual (11) does not change if we restrict Gw and Gv to Si and g respectively. We may thus conclude that problem (11) has a maxi- mizer (W*, V*) with V > Amin (V:) > O for all / E [T - 1j. This in turn implies that problem (6) solved by normal distribution F under which the covariance matrix of the observation noise v satisfies V > 0 for every / E E [T-1] I]. As (5) and (6) are strong duals, the optimal solution w* of problem (5) forms a Nash equilibrium with P i.e., u² a best response to P and thus solves the classic LQG problem corresponding to P*. As R₂ > 0 for every t E [T - 1]. this best response u* is unique, and as Vi= for every E [T' - 1], u² is in fact the Kalman filter-based optimal output- feedback strategy corresponding to P (which can be obtained using the techniques highlighted in Appendix $A).
Before proving Proposition 4.2, recall that f(W, V) is called B-smooth for some B > 0f
[V V)- - V.(W", V*I s 3 (IIW- - [[W--'| + IV - V'1 VA, W' 9w, V.V' E gt
where || || denotes the Frobenius norm.
Proof of Proposition 4.2, The function / f(W, V) is concave because the objective function of the inner minimization problem in (11) is linear (and hence concave) in W and V and because concavity is preserved under minimization. To prove that /(W, V) is 3-smooth, we first recall from Proposition 3.3 that it coincides with the optimal value of the inner minimization problem in (10). As 24, Lly. fWW V) can thus be viewed as the optimal value of the classic LQG problem corresponding the normal distribution P determined by the covariance matrices W and V. Hence, f(W, V) coincides with (A.17), where Σ,, for [T- 1]. is a function of (W, V) defined recursively through the Kalman filter equations (6. Note that all inverse matrices in (A.16) are well-defined because any V E G is strictly positive definite. Therefore, Σ, constitutes proper rational function (that is, a ratio of two polyonmials with the polynomial in the denominator being strictly positive) for every & [T-- -r-11 1]. Thus, f(W, V) is infinitely often continuously differentiable on a neighborhood of Siv x Gt.
As f(W, V) is concave and (at least) twice continuously differentiable, it is B-smooth on G x G and only if the largest eigenvalue of the Hessian matrix of -f(W, is bounded above by / throughout Git x gr. Also, the largest eigenvalue of the positive semidefinite Hessian matrix (- V)) coincides with the spectral norm of 2 F(W,VV We may thus set
B =
sup Veat
J(W, 1)
(A.19)
where || ll2 denotes the spectral norm. The supremum in the above maximization problem is finite and attained thanks to Weierstrass' theorem, which applies because f(W V) is twice continuously differentiable and the spectral norm is continuous, while the sets G and 9 are compact by virtue of [38, Lemma A.6]. This observation completes the proof.
D
17

---

D. SDP Reformulation of the Lower Problem (11)
Instead of solving the dual problem (11) with the customized Frank- Wolfe algorithm of Section 4, be reformulated an SDP amenable off-the-shelf solvers. This reformulation obtained by dualizing the inner minimization problem and by exploiting the following preliminary lemma.
Lemma D.1. For any Z 84 and pz 0, the set 9z {z $4 G(Z, 2) S pz) coincides with
{
{
E s4 3Eg E $4 with Tr(Z + 2 2F,) S p
2
z2 Piz
Eg ス 0
Proof of Lemma D.I. By Definition 2, we have
G2 = 1Z e Tr(Z + 2 -
Next, introduce an auxiliary variable E, E s subject to the matrix inequality E² = (2tz21)
By [4, Theorem 1J. this inequality can be recast as E, 3 (2 22+)= ) Hence, we can reformulate the nonlinear matrix inequality in the above representation of Gz Tr(2 + 2 - 2E,) S p A standard Schur complement argument reveals that the inequality E (0) 72ł) is also equivalent
7
2
Eg
0 0.
E
The claim then follows by combining all of these insights.
We are now ready to derive the desired SDP reformulation of problem (11), Proposition D.2. V O, then problem (11) is equivalent the SDP
max Tr(G" (GW) Tr(F(R+ + H QH)
s.t. WE E B (T())) V E ME M, F FE ST Tr(W + W. -2F I Pius Tr(V, + V. - -2F..) S P., At E -1
Emo E ES, Em, E EST, En, ES, Vt E [ -]]
Tr(Wo + No -2Fz)) S Pro
X2 X
Furo
(A.20)
20,
Ezo
In
W W W
Eus In
VVVVVV
Es, 20 20 Ip
20,
We W
Exc
Es,
F (H QTWDT + M/2)T
H HTQGWD + M/2 DWDT +V
20
Wo - Amiu (8o)1,
W411 ł Amin Ww.),,
V
(v)
E - 1].
Here, M denotes the set of all strictly upper block triangular matrices of the form
M1
M1.
MI. M2r
0
M23
TTxTTp
0
MT-1. -.T
where Mr,3 E Duup for every t, E E with <t<sss T.
Proof of Proposition D.2. The proof relies on dualizing the inner minimization problem in (11). Note that strong duality holds because the primal problem is trivially feasible and involves only equality constraints, which implies that any feasible point is in fact a Slater point. In the following we use M e M to denote the Lagrange multiplier of the constraint U e U, which requires all blocks of
18

---

the matrix U above the main diagonal to vanish. The Lagrangian function of the inner minimization problem in (11) can therefore be represented as
Ulq,U, Tr ((DTUT (R (R+ H THHU + + GT +2 QHUDW) + Ir ((U (R+ QH) + (R+ HTQH4 Tr(UMT).
Recall now that R and Q = O, and thus R+ + HH Q1I O. Consequently, L is minimized by qº = 0 for any fixed U and M. In addition, the partial gradient of L with respect U is given by
dc
au 2(R H QHYUDWD" + H QH)UV +2 QGWD" +
Recall also that V E G is strictly positive, which implies that DWD + V > 0 invertible. As we already know that R + H QH > 0 is invertible, as well, C is minimized by
U - (R H T Q11) (HT QGW + M/2) (DWD + VV
for any fixed M. Substituting both q" and U into L yields the dual objective function
((M) f(g", U* M) Tr(G" TT(GTQGW)
Tr ((R+H ((R+H"QH)"(FF (H"QWD" + M/2)(DWD' + VH QGWD" T
The dual of the inner minimization problem in (11) is thus given by MaxMEM (MM) To linearize the dual objective function, we next introduce an auxiliary variable F E STT" subject to the matrix inequality F = (H QGWDI D + M/2)(DWD) D + V-- FEFTTGWDT + M/2)T. By using standard Schur complement reformulation, we can then rewrite the dual problem
Tr(G² (GW) Tr((R HT -H)-1F
s.t.
F F STST
(A.21)
F'
H HTQ M/2 DWD + V
>
LETQGWDT + M/2)T
Next, by replacing the inner problem in (11) with its strong dual (A.21), we can reformulate (11) as max
Tr(G QGW) - Tr((R + (TQH)-FF)
s.t. MEM, PESTI W E E ((((T+)) V SPT F (HT QUTD M G(Xo, 80)² pho G(W,, W.) S ph, G(V,, V.) S ph AI E
H HTQGWD M/2 00 DWD +V
(A.22)
By Proposition 4.1, the inclusion of the constraints Xp Z Amin (8o),, W. 2 Amin(W)) and V. > Amin T for all E [T -1] 1] has no effect on the solution to problem (A.22). In addition, by Lemma 1. each (non-linear) Gelbrich constraint in (A.22) can be reformulated as an equivalent (linear) SDP constraint. Thus, problem (A.22) reduces to (A.20), and the claim follows.
E. Bisection Algorithm for the Linearization Oracle
We now show that the iirection-findingg subproblem (14) can be solved efficiently via bisection. To this end, we first establish that (14) can be reduced the solution of univariate algebraic equation.
Proposition E.1 ([38, Proposition A.4 (ii))). If 2 E 8 Iz E S Iz + and Ps E R++, then
max s.t.
(T L Z)
G(L, 2) S P=
(A.23)
1 Amin [
uniquely solved by L (-1)2(7^1 - [T)) -1, where 7* is the unique solution of 0 - (2, -1- -
(A.24)
in the interval (Amax (Tz), 00).
19

---

In practice, we need to solve the algebraic equation (A.24) numerically. The numerical error in approximating Y should be contained to ensure that L* approximates the exact maximizer of problem (A.23). The next proposition shows that, for any tolerance å (0, 1), d-approximate solution of (A.23) can be computed with an efficient bisection algorithm.
Proposition E.2 (138, Theorem 6.4]). For any fixed Pa E ER++,2 Z E 84+ and I'z E S4, & O, define 92 (Z E $ G(Z, 2) P2, 2 Amin (2)} as the feasible of problem (A.23), and let Z E 2 be any reference covariance matrix. Additionally, let 8 (0, 1) be the desired oracle precision, and define (p()) = +(p2 - [2z (112)>111 -1_ (2,,, Tz)
for any > Amax (Tz). Then, Algorithm A.2 returns in finite time a matrix L 54 with the following properties. (i) Feasibility: 12 E 92 (ū) 5-Suboptimality: (L2 2, z) 2 max LEOt ZT L 2).
Algorithm A.2 Bisection algorithm compute 47
Input: nominal covariance matrix 2 E S4 radius E R++, reference covariance matrix Z E 92
gradient matrix I'z E $ I'y # 0, precision ò E (0, 1), dual objective function ((>) defined in Proposition E.2
1: set λ + Amax (Tz), and let pi an eigenvector for X₁
2: set y +- (1 + (pT ) pp)p and 7 A1 + Tr(2)z//1
3: repeat
4: set 7 +- (7 + 2)2 and L (5) [z) 2 Tz)
5: if de (f)) < 0 then set t else
7 + 3 endif
6: until de (7) > and (L - Z,T'z 2 86(5)
Output: L
In summary, for any Z {Xo, Wo..... .... Wr-1. Vo...., Vr-1), Algorithm A.2 computes a approximate solutions to the direction- finding subproblem (14) with I'z = Vzf(W,
F. Additional Information on Experiments
Generation of Nominal Covariance Matrices. The nominal covariance matrices of the exoge- nous uncertainties are constructed randomly using the following procedure. For each exogenous uncertainty { Wo..... Wf-11 Up. ... -y we denote the dimension of by d and sample matrix Mz E IRdxd from the uniform distribution on the hypercube [o, 1]dxd Next, we define Ez E pd as the orthogonal matrix whose columns represent the orthonormal eigenvectors of the symmetric matrix M2 + M Finally, we set 2 EzAzE, where Az is a diagonal matrix whose main diagonal is sampled uniformly from the interval [1, 2)d The rationale for adopting this cumbersome procedure is to ensure that the covariance matrix Z is positive definite.
Optimality Gap. The optimality gap of the Wolfk-Wolfe algorithm visualized in Figure 1b is calculated as the sum of the surrogate optimality gaps (L&) Z, Z2/WW V)} across all Z E (Xo. Wo ...... Wr-1, Vo...., 1r-1}. For more information on the surrogate optimality gaps see [31].
20


Reviews:
[{"review_id": "DGCch2mjGu", "reviewer": "Reviewer_RqGY", "summary": "- The paper considers a standard LQ setup with uncertainty in the distributions of system noise, observation noise and initial state.\n- Described as Zero-Sum Game between the controller and nature, they show the optimal decisions for both players. Specifically they show that the worst-case distribution is Gaussian and the optimal control law is linear.\n- They provide an numerically efficient approach to solve the distributionally robust LQ control problem based on Frank-Wolfe algorithm.\n-  Simulation experiments are provided to show the computational efficiency of the proposed algorithm compared to MOSEK. ", "questions": "- In general the adaptive linear quadratic control results are provided under assumption of sub-gaussian system noises (Assumption A1, http://proceedings.mlr.press/v19/abbasi-yadkori11a/abbasi-yadkori11a.pdf). Is it possible to extend these results to Sub-Gaussian nominal distribution?", "limitations": "Authors provide some limitations in form of possible extensions and future work in the last section.\n\nSocial Impact: NA", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The results provided in the manuscript make multiple important contribution.\n\n- The optimal control law still remains linear and the worst case distribution is still Gaussian.\n- The proof technique is novel which relies on the \"purified states\" instead of usual dynamic programming approaches, it is interesting to see how this approach can be used in other LQG problems.\n- Simulation results show the computational efficiency of the Frank-Wolfe algorithm over MOSEK.\n", "weaknesses": "- There are no obvious major weakness in the manuscript. Some potential minor weakness are currently mentioned in the form of questions in the comment below. \n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We sincerely appreciate your insightful comment, which has led us to a significant advancement in our work. Specifically, we managed to extend the applicability of our findings to cases where the nominal distribution is an elliptical distribution with finite first- and second-order moments. This extension also involves relaxing the assumption of independence among noise components to uncorrelatedness. Therefore, our results extend to the sub-Gaussian distributions that are elliptical. For more comprehensive information regarding this extension, please refer to the overall response. "}, {"review_id": "BKRf6x3WaV", "reviewer": "Reviewer_P3Ek", "summary": "The paper proposes a distributionally robust version of the output feedback linear quadratic control problem. The goal is to control a system with known partially observed linear dynamics in the face of stochastic disturbances. The stochastic disturbances (both measurement and state disturbances) are drawn from an unknown distribution. This distribution is assumed to belong to some ambiguity set centered at a known gaussian distribution. The disturbances drawn from any distribution in the ambiguity set are assumed to be mean zero, and independent across time, with measurement disturbances independent from the process disturbances. Furthermore, the marginal distribution for the measurement or process disturbance at each time is assumed to be close to the corresponding nominal marginal distribution (as measured by the 2-Wasserstein distance).\n\nUnder this setting, the paper finds the optimal output feedback controller for the worst case distribution of disturbances in the ambiguity set. It is found that the worst case distribution of disturbances belonging to the ambiguity set described above is Gaussian. Therefore, the corresponding optimal controller is a Linear-Quadratic Gaussian controller designed for this worst case distribution. \n\nIt is then shown that the worst case distribution may be determined in a computationally efficient manner. In particular, it is demonstrated that the Frank-Wolfe algorithm converges to the optimal covariance parameters for the worst case distribution. Furthermore, it is shown that each step of the Frank-Wolfe algorithm may be decomposed into simple, easily parallelizable components.\n\n", "questions": "Is there any concrete theoretical connection between distributionally robust linear quadratic control and the other methods for incorporating robustness mentioned above (e.g. mixed H-2/H-inf)? \n\nAre there any practical examples where the proposed method substantially outperforms conventional approaches for incorporating robustness to unknown disturbance distributions? Such an example would make the setting more compelling for practical use. \n\nFrom the experiments, is it possible to detect any clear trends regarding the worst case covariance relative to the nominal covariance of the ambiguity set? E.g. Do we see or expect to see that the worst-case covariances are larger than the central covariances in Loewner order? Combined with the theoretical results, such an observation might justify crude approximate approaches for distributionally robust control design in practice. For example, one could take their estimate for the covariance for the central distribution in the ambiguity set, and simply scale it up by some constant. The resulting covariances could then be used to design a LQG controller.\n", "limitations": "There is no negative societal impact that the authors must address. Limitations to the method that would be worth addressing are mentioned in the weaknesses section. To reiterate, it would be helpful to address:\n-\tWhat the ambiguity set can model, and what it cannot\n-\tIn which settings the proposed method would outperform conventional methods from robust control theory for handling unknown disturbance distributions\n-\tAcknowledging the computational burden of the approach relative to e.g. LQG with a known distribution\n", "rating": 7, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 2, "strengths": "The formulation of the distributionally robust output feedback LQ control problem is novel. It is related to previously published results on distributionally robust state feedback LQ control. As acknowledged by the authors, the output feedback setting brings an extra technical challenge due to the dependence of the optimal state estimator upon the disturbance distribution. \n\nThe problem and corresponding solution are clearly presented, and easy to follow.  The authors highlight the rather surprising result that the worst case disturbance distribution from the prescribed ambiguity set is Gaussian. \n\nGiven the appropriate context, the results here could be a meaningful step in a unification of worst case and stochastic control. \n", "weaknesses": "The contextualization of the setting studied relative to prior work could be improved. In particular, there are a large class of control synthesis approaches for settings where the noise distributions are either not available or are not Gaussian. In particular, consider H-infinity approaches [Zhou et. al, Robust and Optimal Control, 1996], mixed H-2/H-infinity approaches [(Doyle et. al Optimal Control with Mixed H-2/H-infinity Performance Objectives, 1989) (Bernstein and Haddad,  LQG control with an H-infinity performance bound, 1988)], adversarially robust control [Lee et. al, Performance-Robustness Tradeoffs in Adversarially Robust Control and Estimation, 2023], and nonstochastic control [Hazan and Singh, Introduction to Online Nonstochastic Control, 2023]. It would be useful to discuss several of these results and contrast the setting with the distributionally robust setting. \n\nA thorough discussion about the choice for the ambiguity set and what distributions it can model would be beneficial. In particular, considering mean zero disturbances which are independent across time is restrictive. It fails to model e.g. colored noise. \n\nDespite the fact that it was possible to propose a relatively efficient method to optimize for the worst-case covariance of the noise distribution, the required computation time appears to scale quite poorly with the problem horizon, T. Only short horizons, up to T=20, were considered in the simulations, however many control problems have much longer horizons. This appears to limit the practical applicability. \n\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Extension of Theorem 3.5**:  We sincerely appreciate your insightful comment, which has led us to a significant advancement in our work. Specifically, we managed to extend the applicability of our findings to cases where the nominal distribution is an elliptical distribution with finite first- and second-order moments. This extension also involves relaxing the assumption of independence among noise components to uncorrelatedness. Additionally, we can also relax the assumption that the distributions in the ambiguity set have a fixed mean. For more comprehensive information regarding this extension, please refer to the overall response.\n* **Scalability of Algorithm 1**: Thank you very much for pointing us in this direction. Your question helped us realize that our chosen matrix $A$—which was not a convergent matrix ($i.e.$,\nwas not satisfying $\\lim_{t \\to \\infty} (A^t)_{ij} = 0$ for all $i, j \\in [n]$)—was causing numerical instability in the linearization oracle of Algorithm 1 for $T > 50$. To address this issue, we rescaled the $A$ matrix described in the numerical section by a factor of $0.1$, ensuring it is a convergent matrix. Subsequently, we reran our experiment according to the procedures outlined in Section 5, but this time for $T=100$. The results of this modified experiment are shared in Figure 1a of the attached PDF document.\n* **Relation to the relevant literature**:  Indeed, our work is related to the aforementioned literature, which involves minimizing a worst-case objective functional ($e.g.$, transfer function, cost, regret), where the noise perturbations are selected adversarially [Zhou et. al., Hazan and Singh], or considers mixtures of nominal and worst-case perspectives [Doyle et. al., Bernstein and Hadded, Lee et. al.]. In contrast to these approaches, we adopt a distributionally robust approach, evaluating performance based on the worst-case expected cost in view of all noise distributions close to a nominal one. Although we share the same motivation of enhancing control policies' robustness against uncertainty, we are currently unable to identify any rigorous equivalence or theoretical relationship between our approach and the aforementioned methods. We will investigate further and will incorporate a comprehensive discussion of this related literature in the revised paper.\n* **Loewner order of worst-case covariance matrices**: Assume that we fix the optimal affine controller. Then, the expected cost can be reformulated as $\\mathbb E_{\\mathbb{P}}[\\xi^\\top S\\xi]$ for some matrix $S \\succeq 0$, where $\\xi\\sim {\\mathbb{P}}$ encapsulates the uncertainties inherent to the model with  $\\mathbb E_{\\mathbb{P}}[\\xi]=0$ and $\\mathbb E_{\\mathbb{P}}[\\xi \\xi^\\top]=\\Sigma$. Let $\\hat\\Sigma$ denote the nominal covariance matrix. Then, by Theorem 3.5 maximizing the expected cost over a 2-Wasserstein ambiguity set $\\mathcal{W}$ of radius $\\rho$ is equivalent to maximizing $Tr(S \\Sigma)$ over all $\\Sigma\\succeq 0$ satisfying $G(\\Sigma, \\hat \\Sigma)\\leq\\rho^2$. In the following, we will show that the worst-case covariance matrix cannot be smaller than the nominal covariance matrix $\\hat\\Sigma$ with respect to the Loewner order. By applying a linear coordinate transformation, we may assume without loss of generality that $\\hat\\Sigma=I$. Consider now any matrix $\\Sigma\\succeq 0$ in the Gelbrich ball, that is, $G(\\Sigma, I)^2 = Tr(I+\\Sigma-2\\Sigma^{1/2})\\leq \\rho^2$. Using the spectral decomposition $\\Sigma=\\sum_{i=1}^n\\lambda_i v_i v_i^\\top$, where $\\lambda_i\\geq 0$ is the $i$-th eigenvalue and $v_i$ is the corresponding eigenvector (we may assume the eigenvectors are orthonormal), the squared Gelbrich distance can be reformulated as $Tr(I+\\Sigma-2\\Sigma^{1/2}) = \\sum_{i=1}^n (1-\\lambda_i^{1/2})^2$. Note that $\\hat \\Sigma = I$ has eigenvalues 1. In the following, we will show that if $\\Sigma$ has any eigenvalues strictly less than 1, it will not attain the maximum value of $Tr(\\Sigma S)$, that is, there exists $\\Sigma'$ satisfying $G(\\Sigma', I) \\leq \\rho$ with eigenvalues greater or equal than 1. Now let $\\Sigma'=\\sum_{i=1}^n \\max${$1,\\lambda_i$}$v_i v_i^\\top$, which has the same eigenvectors as $\\Sigma$, but all eigenvalues smaller than 1 are set to 1, that is, $\\Sigma'$ has larger eigenvalues than $\\Sigma$ with the same eigenvectors. The matrix $\\Sigma'$ is feasible because $Tr(I+\\Sigma'-2(\\Sigma')^{1/2}) = \\sum_{i=1}^n (1-\\max${$1,\\lambda_i$}$^{1/2})^2 \\leq \\sum_{i=1}^n (1-\\lambda_i^{1/2})^2 = Tr(I+\\Sigma-2\\Sigma^{1/2}) \\leq \\rho^2$. In addition, the expected loss with respect to $\\Sigma'$ is at least as large as the expected loss with respect to $\\Sigma$. Formally, we have $Tr(S\\Sigma') \\geq Tr(S\\Sigma)$ because $S\\succeq 0$ and because $\\Sigma'\\succeq \\Sigma$ by construction. This shows that the worst-case covariance matrix must be non-inferior to the nominal covariance matrix in Loewner order. Similar arguments can be used to show that increasing $\\rho$, increases the worst-case covariance matrix in Loewner order. In Figures 2a and 2b shared in the uploaded PDF file, for $n=m=p=2$, we illustrate the worst-case covariance matrices for 10 different values of $\\rho$ split evenly in the range $[0,1]$ satisfying $\\rho_{x_0} = \\rho_{w_0} = \\ldots= \\rho_{w_{T-1}} = \\rho_{v_0}=\\ldots=\\rho_{v_{T-1}}=\\rho$.  We can see empirically from the plots that indeed the worst-case covariance matrices are inflating as the radius $\\rho$ of the balls increases. Interestingly, Figure 2 indicates that scaling the empirical covariance matrix $\\hat X_0$ might approximate its worst-case in this experiment, but the same scaling wouldn't provide an accurate estimate for the worst-case covariance of $v_0$. "}, {"review_id": "6x0ljlwzse", "reviewer": "Reviewer_BGQp", "summary": "The paper proposes a robust method for controlling LQ systems, where process and observation noise distribution laws are unknown, but noise samples are assumed to be independent, zero mean and have their distributions lying close to a a nominal Gaussian distribution in Wasserstein-2 space.\n\nA pair of relaxed settings is used to prove a strong duality between the minimax\\maximin problems, showing that the worst-case distribution is Gaussian. A numerical method is proposed for computing this distribution for long horizons.\n", "questions": "Some minor comments \\ typos:\n\nIt is somewhat misleading to call $\\mathcal{W}$ a 'Wasserstein ball', e.g. l.33, while it is not even convex, but I guess that's forgivable since this set is defined clearly.\n\nl. 93 dimension should be $p \\times T$. \n\nl.104 it would probably be more didactic to introduce the Wasserstein distance before using it to define  $\\mathcal{W}$.\n", "limitations": "The authors addressed the limitations, but a more detailed discussion about the assumption that the nominal distribution is Gaussian is required.", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": " The paper itself is well-written and gives a thorough overview of the problem.\n\nThe fact that the solution for the minimax problem is given by a Gaussian distribution \\ linear controller, even if not surprising considering the nature of the quadratic cost and $W_2$/Gelbrich distance, is still not clear from the outset and is important enough from a theoretical point of view. Its proof seems correct as well.\n", "weaknesses": "My main concern here, is that the theoretical contribution is limited to one main result, which is confined to a Gaussian-centered ambiguity set, i.e. the real distribution is assumed to be close to a Gaussian one. This calls for an additional study of other nominal distributions (as the Gelbrich distance equals to $W_2$ for other families of distributions [18], and discussion may be confined to linear filters) or mixtures, or at least for a detailed discussion about the practical implications of this choice (beyond that of the 'concluding remarks' in Sec.6) i.e. about whether this hypothesis may (or might not) be practical in real problems.\n\nHowever, even if limited, the contribution is still novel and important enough to warrant acceptance.\n\n", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "* **Extension of Theorem 3.5:** We sincerely appreciate your insightful comment, which has led us to a significant advancement in our work. Specifically, we managed to extend the applicability of our findings to cases where the nominal distribution is an elliptical distribution with finite first- and second-order moments. For more comprehensive information regarding this extension, please refer to the overall response.\n\n* **Wasserstein ball:** Thank you very much for bringing this into our attention. You are correct in asserting that the phrase \"Wasserstein ball\" could lead to confusion due to the non-convex nature of the ambiguity set $\\mathcal W$. The phrase was chosen for its brevity and alignment with customary distributionally robust literature. We will emphasize the non-convex nature of our ambiguity set in the revised paper.\n\n* **Typo:**  Thank you very much for pointing out this typo. We will correct this in the revised manuscript.\n\n* **Definition of Wasserstein distance:** Thank you for bringing this to our attention. We will move the definition of the Wasserstein distance to l.100 in the revised manuscript.\n"}, {"review_id": "4Tbhhoe3ep", "reviewer": "Reviewer_iQxV", "summary": "The paper **Distributionally Robust Linear Quadratic Control** considers the case of finite-time linear quadratic optimal control with process and measurement noises and known time-varying dynamics. The novelty lies in the fact that the distributions of the initial state and noises are unknown but lie in an ambiguity set defined as a 2-Wasserstein ball centered in a known Gaussian distribution and with known radius. The objective is to minimize the expected LQR cost with the distributions chosen adversarially in that ball. The authors call this adversarial problem \"distributionally robust Linear Quadratic Gaussian\". It strict generalizes LQG control, which assumes these distributions are known and Gaussian. \nThe main contribution is twofold. The first part is theoretical, as the authors prove the existence of a Gaussian adversarial distribution, that is, the distributionally robust LQG reduces to classic LQG with unknown Gaussian distributions. A consequence is that the optimal controller is a linear state feedback, like in the classic LQG case. This leads to the second contribution: an algorithm to efficiently estimate the adversarial Gaussian distribution from data. The authors claim that, then, the LQG problem with estimated noise distributions can be solved efficiently by using classical methods.\nThe authors illustrate convergence of their estimation algorithm on a simulated example.", "questions": "My questions are detailed in the above paragraph on weaknesses.  I would appreciate if the authors could respond to these.", "limitations": "The authors have adequately addressed limitations and potential negative societal impact of their work.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "The paper is interesting and the clear exposition makes the reasonings easy to follow. The appendix is well-managed, and I appreciate Appendix A on the solution of classical LQG with a Kalman filter.\nI like the idea of allowing for a whole family of noise distributions rather than assuming a fixed one. The main theoretical result helps mitigate the restrictiveness of the commonly-accepted \"white noise assumption\"; indeed, it shows that allowing for a larger family of noises does not provide any benefits (for distributions close enough to Gaussians). This new problem formulation thus seems relevant. The subsequent algorithm to solve the distributionally robust LQG problem follows naturally and further justifies the interest of the theoretical result. Overall, the reasoning exposed is sound and well-motivated.\nFinally, I appreciate that the authors went the extra mile in Section IV by augmenting the theoretical result with a data-efficient algorithm.", "weaknesses": "The three main weaknesses of the paper are, in my opinion, 1. the lack of thoroughness of the simulation study; 2. the lack of discussion of the effect of hyperparameters; and 3. some arguments are unclear and should be explicited (although I do not question their conclusions). I detail these three points. These points are not critical for acceptance, but I believe the paper would be improved by addressing them.\n\nW1) The simulation study shows convergence of the algorithm on a randomized use case. This is a good sanity check, but I would also appreciate a simulation of the case when the noise distribution is not normal but still within the Wasserstein ball. The theoretical result ensures that Nature's adversarial strategy _is_ normal, but an illustration of the case when Nature is sub-adversarial would be welcome. In particular, is the cost of the learned policy reduced compared to when the noise is normal?\n\nW2) The role of the choice of ball radii is not discussed. While there is obviously no notion of \"optimal\" radius, since it simply corresponds to a degree of robustness, I would like to see the evolution in performance with increasing radii. Intuition tells that performance should decrease as robustness increases, but a confirmation in simulation would be welcome.\n\nW3) Some arguments are unclear to me, although their conclusions seem to be valid. In particular:\n   1. The authors make multiple claims on convexity of sets of probability measures. Examples are lines 112, 151, 189. I understand that these claims are made by considering Borel measures as a subset of the vector space of signed Borel measures with standard addition and scalar multiplication. I believe this superset should be mentioned explicitly at least ones, as it is rather unusual and notions of convexity of a metric space exist (and differ).\n      This is also relevant on line 151, where $\\mathcal{W}$ is claimed to be infinite-dimensional despite not being a vector space.\n   2. With this understanding of convexity, why is the set $\\mathcal{W}$ non convex, as claimed e.g. on lines 112, 151 and 189? As far as I understand, each set $\\mathcal{W}_{z}$ is convex, with $z\\in\\{x_0, w_t, v_t\\}$. Then, $\\mathcal{W}$ should be convex as the cartesian product of these sets. The only way I understand this non-convexity is if $\\mathcal{W}$ is not _equal_ to the cartesian product, but only isometric to this cartesian product by the mapping that computes the marginal distributions. I believe this should be mentioned explicitly, at least in a footnote, as this questions distracted me from more central claims of the paper for a while.\n   3. I am unsure about the claim on line 131 that the controller can compute the fictious states $\\hat x_0,\\dots,\\hat x_t$ from the real observations $y_0, \\dots, y_t$ without knowing the initial state. Indeed, this would imply in particular that the controller can reconstruct the initial state. Since the time origin is arbitrary, any state could be reconstructed. This claim appear unnecessary for the rest of the argument and should be either removed or clarified in my opinion.\n   4. I find the formulation of the first sentence of Proposition 4.1 extremely confusing. In particular, what comes after \"then\" in the first sentence reads as a logical consequence of what precedes whereas it is actually a definition of the symbols $\\mathbb{P}^\\star$ and $V_t^\\star$. I recommend reformulating.\n   5. In Proposition 4.2, I recommend avoiding using the term \"smooth\". While it has a precise meaning in a specific branch of mathematics, it is often only used informally in control. I would prefer the more standard \"infinitely differentiable with $\\beta$-Lipschitz gradient\".", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We sincerely appreciate your insightful comments.\nW1) Thank you for suggesting an investigation into the expected cost under different distributions within the ambiguity set. We note that due to the construction of our distributionally robust control problem, the expected cost of our optimal policy under any distribution (Gaussian or not) in the ambiguity set cannot exceed its worst-case expected cost, $i.e.$, the expected cost under the worst-case distribution. However, it is still valuable to investigate how the expected cost evolves under different distributions, so we conducted a numerical experiment along the lines of your suggestion. \nIn this experiment, we consider the same setup as in Section 5, where the time horizon is set to $T=2$ and the common radius is set to $\\rho = 10$. We generate different distributions that fall within our Wasserstein ambiguity set via a contamination model; specifically, for any $\\varepsilon\\in[0, 1]$, we compute the $\\varepsilon$-contamination distribution ${\\mathbb{P}}^{\\varepsilon}$ as the Gaussian distribution with mean $0$ and covariance matrix $\\Sigma^\\varepsilon=\\varepsilon\\times\\Sigma^\\star+(1-\\varepsilon)\\times\\hat\\Sigma$, where $\\hat \\Sigma,\\Sigma^\\star$ denote the covariances matrices of $\\hat{\\mathbb P}$ and $\\mathbb P^\\star$, respectively. By leveraging the convexity of the squared Gelbrich distance and the equivalence between the Gelbrich and Wasserstein distances in the case of Gaussian distributions, we can show that $\\mathbb{P}^\\varepsilon$ belongs to the Wasserstein ambiguity set for all $\\varepsilon \\in [0, 1]$. We generate Gaussian distributions because checking whether an arbitrary distribution belongs to the Wasserstein ambiguity set is very challenging computationally (in fact, establishing whether an arbitrary discrete distribution falls within ${\\mathcal{W}}$ is tantamount to solving a semi-discrete optimal transport problem, which is known to be #P-Hard [R1, Theorem 2.2]).\n* Figure 1b depicts the expected cost of the robust optimal policy $u^*$ under the distribution $\\mathbb{P}^\\varepsilon$ as a function of $\\varepsilon$. Note that the expected cost increases with $\\varepsilon$ (as the contaminated distribution approaches the worst-case distribution). Additionally, the expected cost is linear with respect to the contamination level $\\varepsilon$, which is aligned with our theoretical results that the objective function of the robust problem admits a linear reformulation in the covariance matrices, as in (9).\n* Figure 1c shows the difference in expected costs between a policy $\\hat{u}$ that is optimal under the nominal case (i.e.,that minimizes the expected cost under $\\hat{\\mathbb P}$) and the robustly optimal policy $u^*$, with both policies evaluated under the contamination distribution $\\mathbb{P}^\\varepsilon$ for different values of $\\varepsilon$. Note that as long as $\\varepsilon\\geq 0.05$, the robustly optimal policy outperforms the nominal one under $\\mathbb{P}^\\varepsilon$, therefore resulting in better performance for the vast majority of contamination levels. In addition, even for $\\varepsilon\\leq 0.05$, the performance of the robust policy is similar to that of the nominal one.\n\n[R1] B. Taskesen, S. Shafieezadeh-Abadeh, and D. Kuhn. \"Semi-discrete optimal transport: Hardness, regularization and numerical solution.\" Math. Prog. 199.1-2 (2023): 1033-1106.\n\nW2) You are correct in noting that the optimal worst-case expected cost is nondecreasing in each of the radii, because increasing any of the radii expands the ambiguity set and relaxes nature's maximization problem. Following your suggestion, we designed a new experiment to quantify the benefits of the robustly optimal policy $u^*$ in comparison to the nominal optimal policy $\\hat u$ (that minimizes the expected cost under the nominal distribution) when both are assessed under (i) the nominal distribution and (ii) the corresponding worst-case distributions, respectively, with respect to different radii. We consider the same setup as in Section 5, with a time horizon $T=2$. We vary the common radius $\\rho$ from 0 to 10 and we estimate the difference between the expected costs of $\\hat{u}$ and the expected costs of $u^*$ under (i) the nominal distribution $\\hat{\\mathbb P}$ and (ii) under their respective worst-case distributions, which depend on the radius $\\rho$. The results, which are shown in Figure 3, illustrate that the gap in expected costs under the worst-case distributions drastically increase as $\\rho$ increases, which indicates that the performance of the nominal optimal policy deteriorates rapidly in comparison to the robust policy in worst-case scenarios. In contrast, the gap in expected costs under $\\hat{\\mathbb P}$ barely changes with the radius $\\rho$, which shows that the robust policy has a robust performance in the nominal scenario, almost matching that of the nominal policy (which is optimal under that scenario).\n\nW3)\n* Indeed, the concept of convexity for sets is meaningful only when these sets are part of a linear (vector) space. In our revised manuscript, we will specify that when referring to the convexity of $\\mathcal{W}$, where all measures are supported on $\\mathbb R^d$, it is implied that $\\mathcal{W}$ is a subset of the linear (vector) space of all signed measures on $\\mathbb R^d$.\n* In our revised manuscript, we will use $\\otimes$ operator instead of $\\times$ and make the distinction clear. \n* In the noise-free system (please see lines 127-128), the fictitious initial state is set to $\\hat{x}_0 = 0$. In addition, the control input $u_t$ can be computed from observations $y_0,\\ldots,y_t$ and does not necessitate knowledge of the true initial state $x_0$. Consequently, the true initial state does not appear in or affect the noise-free system in any manner, and it is not possible to construct the true initial state from this system. We hope our response clarifies your concern.\n* We will clarify the last two points in the revised manuscript."}]

Author rebuttal:
We express our gratitude to the reviewers for their thoughtful comments and questions.

In the reviews, a common concern was raised about our main result in Theorem 3.5. This theorem establishes the optimality of a linear control policy and identifies the worst-case distribution as a Gaussian distribution. However, our proof was specifically formulated under the assumption of a Gaussian nominal distribution.

Thanks to the insightful feedback provided by the reviewers, we have identified the potential for extending Theorem 3.5. Specifically, we can **extend Theorem 3.5** to situations where **the nominal distribution is an elliptical distribution with finite first- and second-order moments** and show that a linear control policy is optimal while the worst-case distribution is an elliptical distribution. In addition, we can **relax the independence assumption** of the noise components **to uncorrelatedness** in this case. 
Finally, we can also **drop the zero mean assumption for the nominal distribution** and **relax the assumption that the distributions in the ambiguity set have a fixed mean**.


The class of elliptical distributions generalizes multivariate Gaussian distributions and spherical distributions [R1] and includes symmetric distributions with light and heavy tails. 
Examples of elliptical distributions include the Laplace, logistic, and $t$-distribution, among others.  
Therefore, even though we cannot extend our results readily to any sub-Gaussian distribution as asked by Reviewer RqGy, we can extend our results to some sub-Gaussian distributions that are also elliptical.

There are two key reasons that allow us to extend Theorem 3.5 to situations where the nominal distribution is elliptical with finite first- and second-order moments. Firstly, when ${\mathbb{P}}$ is an elliptical distribution, the conditional expectation $\mathbb E_{\mathbb{P}}[\xi | \xi' = H\xi ]$ is linear in $\xi'$ [R2, Theorem 4]. This implies that under any fixed elliptical distribution, a linear control policy is optimal. Consequently, this allows us to establish a parallel lower-bound problem to (10), wherein $\mathcal W_{\mathcal N}$ is substituted by a subset of the ambiguity set $\mathcal W$ containing only elliptical distributions sharing the same characteristic function with the nominal one. This observation permits us to focus on linear policies within the inner minimization of the lower-bound problem without loss of generality.
Secondly, the distance between two elliptical distributions, sharing identical characteristic functions, is equivalent to the Gelbrich distance between their corresponding mean vectors and covariance matrices [R3, Theorem 2.4]. This observation further empowers us to substitute the feasible set of the outer maximization with the Gelbrich ambiguity set. The remainder of our proof seamlessly adapts to elliptical distributions, as other arguments do not rely on the Gaussian structure of the nominal distribution. 

Furthermore, we can also extend our numerically efficient method to compute the worst-case distribution, now elliptical, along with the optimal linear control policy, to situations where the nominal distribution is elliptical. Notably, Proposition 4.1's applicability naturally extends the case where the worst-case distribution ${\mathbb{P}}^\star$ becomes an elliptical distribution sharing the same characteristic function with the nominal one. By leveraging our Frank-Wolfe algorithm, we can compute the covariance matrix of this distribution. Subsequently, solving a linear quadratic control problem using ${\mathbb{P}}^\star$ allows us to determine the optimal controller. Remarkably, even in this context where the nominal distribution is elliptical, the separation principle remains applicable [R4], and the recursive equations of the Kalman filter continue to hold for elliptical distributions [R5, R6].

The departure from the fixed mean assumption in the ambiguity set can be inferred from [R7, Theorem 2.7], [R8, Theorem 2.16], and [R7, Theorem 3.5]. However, to maintain focus on the core message of the paper and prevent undue complexity in terms of calculations and notation, we did not include this extension in the submitted version of our manuscript.

We thank the reviewers for recommending us to explore this meaningful direction that strengthens the main result of our paper.

[R1] R.D. Lord, The Use of the Hankel Transform in Statistics I. General Theory and Examples, Biometrika, vol. 41, no. 1/2, 1954, pp. 44–55, JSTOR.


[R2] K.-C. Chu, Estimation and decision for linear systems with elliptical random processes, IEEE Conference on Decision and Control on Adaptive Processes, 1972, pp. 647-651.


[R3] M. Gelbrich, On a formula for the $L^2$-Wasserstein metric between measures on Euclidean and Hilbert spaces, Mathematische Nachrichten, 147 (1990), pp. 185–203.



[R4] Witsenhausen, Hans S., Separation of estimation and control for discrete time systems, Proceedings of the IEEE 59.11 (1971): 1557-1566.

[R5] Basu, A. K., \& Das, J. K. (1994). A Bayesian Approach to Kalman Filter for Elliptically Contoured Distribution and its Application in Time Series Models. Calcutta Statistical Association Bulletin, 44(1–2), 11–28. 

[R6] Girón, F. J., \& Rojano, J. C. (1994). Bayesian Kalman Filtering with Elliptically Contoured Errors. Biometrika, 81(2), 390–395. 


[R7] V. A. Nguyen, S. Shafieezadeh-Abadeh, D. Kuhn, \& P. Mohajerin Esfahani, Bridging Bayesian and minimax mean square error estimation via Wasserstein distributionally robust optimization. Mathematics of Operations Research, 2023, 48(1), 1-37.

[R8] K.-T. Fang, S. Kotz, and K. W. Ng, Symmetric Multivariate and Related Distributions, Chapman \& Hall, 1990.


Meta-review: **Summary:** This paper consider a generalization of the discrete-time, finite-horizon LQG problem, where the noise distributions are unknown and belong to Wasserstein ambiguity sets centered at nominal (Gaussian) distributions. The goal is to minimize a worst-case cost across all distributions in the ambiguity set, including non-Gaussian distributions. The main contribution of this paper are in showing that the optimal control policy is linear (w.r.t. observations), and in providing a numerical approach to compute this optimal policy. 

**Comments:** The review scores are 7, 7, 7, 7. Reviewers are unanimously positive about the technical contributions of this paper. This paper will be of interest to the NeurIPS community, especially for the subareas of control+learning / reinforcement learning. For the final version,  please update your paper based on the comments and the suggestions by the reviewers.

Thanks!

AC
Decision: Accept (spotlight)