Convergent Policy Optimization for Safe Reinforcement Learning

We study the safe reinforcement learning problem with nonlinear function approximation, where policy optimization is formulated as a constrained optimization problem with both the objective and the constraint being nonconvex functions. For such a problem, we construct a sequence of surrogate convex constrained optimization problems by replacing the nonconvex functions locally with convex quadratic functions obtained from policy gradient estimators. We prove that the solutions to these surrogate problems converge to a stationary point of the original nonconvex problem. Furthermore, to extend our theoretical results, we apply our algorithm to examples of optimal control and multi-agent reinforcement learning with safety constraints.

Reinforcement learning [58] has achieved tremendous success in video games [41, 44, 56, 36, 66] and board games, such as chess and Go [53, 55, 54], in part due to powerful simulators [8, 62]. In contrast, due to physical limitations, real-world applications of reinforcement learning methods often need to take into consideration the safety of the agent [5, 26]. For instance, in expensive robotic and autonomous driving platforms, it is pivotal to avoid damages and collisions [25, 9]. In medical applications, we need to consider the switching cost [7].

A popular model of safe reinforcement learning is the constrained Markov decision process (CMDP), which generalizes the Markov decision process by allowing for inclusion of constraints that model the concept of safety [3]. In a CMDP, the cost is associated with each state and action experienced by the agent, and safety is ensured only if the expected cumulative cost is below a certain threshold. Intuitively, if the agent takes an unsafe action at some state, it will receive a huge cost that punishes risky attempts. Moreover, by considering the cumulative cost, the notion of safety is defined for the whole trajectory enabling us to examine the long-term safety of the agent, instead of focusing on individual state-action pairs. For a CMDP, the goal is to take sequential decisions to achieve the expected cumulative reward under the safety constraint.

Solving a CMDP can be written as a linear program [3], with the number of variables being the same as the size of the state and action spaces. Therefore, such an approach is only feasible for the tabular setting, where we can enumerate all the state-action pairs. For large-scale reinforcement learning problems, where function approximation is applied, both the objective and constraint of the CMDP are nonconvex functions of the policy parameter. One common method for solving CMDP is to formulate an unconstrained saddle-point optimization problem via Lagrangian multipliers and solve it using policy optimization algorithms [18, 60]. Such an approach suffers the following two drawbacks:

First, for each fixed Lagrangian multiplier, the inner minimization problem itself can be viewed as solving a new reinforcement learning problem. From the computational point of view, solving the saddle-point optimization problem requires solving a sequence of MDPs with different reward functions. For a large scale problem, even solving a single MDP requires huge computational resources, making such an approach computationally infeasible.

Second, from a theoretical perspective, the performance of the saddle-point approach hinges on solving the inner problem optimally. Existing theory only provides convergence to a stationary point where the gradient with respect to the policy parameter is zero [28, 37]. Moreover, the objective, as a bivariate function of the Lagrangian multiplier and the policy parameter, is not convex-concave and, therefore, first-order iterative algorithms can be unstable [27].

In contrast, we tackle the nonconvex constrained optimization problem of the CMDP directly. We propose a novel policy optimization algorithm, inspired by [38]. Specifically, in each iteration, we replace both the objective and constraint by quadratic surrogate functions and update the policy parameter by solving the new constrained optimization problem. The two surrogate functions can be viewed as first-order Taylor-expansions of the expected reward and cost functions where the gradients are estimated using policy gradient methods [59]. Additionally, they can be viewed as convex relaxations of the original nonconvex reward and cost functions. In  §4 we show that, as the algorithm proceeds, we obtain a sequence of convex relaxations that gradually converge to a smooth function. More importantly, the sequence of policy parameters converges almost surely to a stationary point of the nonconvex constrained optimization problem.

Related work. Our work is pertinent to the line of research on CMDP [3]. For CMDPs with large state and action spaces, [19] proposed an iterative algorithm based on a novel construction of Lyapunov functions. However, their theory only holds for the tabular setting. Using Lagrangian multipliers, [46, 18, 1, 60] proposed policy gradient [59], actor-critic [33], or trust region policy optimization [51] methods for CMDP or constrained risk-sensitive reinforcement learning [26]. These algorithms either do not have convergence guarantees or are shown to converge to saddle-points of the Lagrangian using two-time-scale stochastic approximations [10]. However, due to the projection on the Lagrangian multiplier, the saddle-point achieved by these approaches might not be the stationary point of the original CMDP problem. In addition, [65] proposed a cross-entropy-based stochastic optimization algorithm, and proved the asymptotic behavior using ordinary differential equations. In contrast, our algorithm and the theoretical analysis focus on the discrete time CMDP. Outside of the CMDP setting, [31, 35] studied safe reinforcement learning with demonstration data, [61] studied the safe exploration problem with different safety constraints, and [4] studied multi-task safe reinforcement learning.

Our contribution. Our contribution is three-fold. First, for the CMDP policy optimization problem where both the objective and constraint function are nonconvex, we propose to optimize a sequence of convex relaxation problems using convex quadratic functions. Solving these surrogate problems yields a sequence of policy parameters that converge almost surely to a stationary point of the original policy optimization problem. Second, to reduce the variance in the gradient estimator that is used to construct the surrogate functions, we propose an online actor-critic algorithm. Finally, as concrete applications, our algorithms are also applied to optimal control (in  §5) and parallel and multi-agent reinforcement learning problems with safety constraints (in supplementary material).

A Markov decision process is denoted by  (S, A, P, γ, r, µ), where Sis the state space, A is the action space, P is the transition probability distribution,  γ ∈ (0, 1)is the discount factor,  r: S × A → R isthe reward function, and  µ ∈ P(S)is the distribution of the initial state  s0 ∈ S, where we denote P(X) as the set of probability distributions over X for any X. A policy is a mapping  π : S → P(A)that specifies the action that an agent will take when it is at state s.

Policy gradient method. Let  {πθ : S → P(A)}be a parametrized policy class, where  θ ∈ Θis the parameter defined on a compact set  Θ. This parameterization transfers the original infinite dimensional policy class to a finite dimensional vector space and enables gradient based methods to be used to maximize (1). For example, the most popular Gaussian policy can be written as π(·|s, θ) = N�µ(s, θ), σ(s, θ)�, where the state dependent mean  µ(s, θ)and standard deviation σ(s, θ)can be further parameterized as  µ(s, θ) = θ⊤µ · x(s) and σ(s, θ) = exp�θ⊤σ · x(s)�with x(s)being a state feature vector. The goal of an agent is to maximize the expected cumulative reward


where  s0 ∼ µ, and for all  t ≥ 0, we have st+1 ∼ P(· | st, at) and at ∼ π(· | st). Given a policy  π(θ),we define the state- and action-value functions of  πθ, respectively, as


The policy gradient method updates the parameter  θthrough gradient ascent


where �∇θR(θk)is a stochastic estimate of the gradient  ∇θR(θk)at k-th iteration. Policy gradient method, as well as its variants (e.g. policy gradient with baseline [58], neural policy gradient [64, 39, 16]) is widely used in reinforcement learning. The gradient  ∇θR(θ)can be estimated according to the policy gradient theorem [59],


Actor-critic method. To further reduce the variance of the policy gradient method, we could estimate both the policy parameter and value function simultaneously. This kind of method is called actor-critic algorithm [33], which is widely used in reinforcement learning. Specifically, in the value function evaluation (critic) step we estimate the action-value function  Qθ(s, a)using, for example, the temporal difference method TD(0) [20]. The policy parameter update (actor) step is implemented as before by the Monte-Carlo method according to the policy gradient theorem (3) with the action-value Qθ(s, a)replaced by the estimated value in the policy evaluation step.

Constrained MDP. In this work, we consider an MDP problem with an additional constraint on the model parameter  θ. Specifically, when taking action at some state we incur some cost value. The constraint is such that the expected cumulative cost cannot exceed some pre-defined constant. A constrained Markov decision process (CMDP) is denoted by  (S, A, P, γ, r, d, µ), whered: S × A → Ris the cost function and the other parameters are as before. The goal of an the agent in CMDP is to solve the following constrained problem


where  D0is a fixed constant. We consider only one constraint  D(θ) ≤ D0, noting that it is straightforward to generalize to multiple constraints. Throughout this paper, we assume that both the reward and cost value functions are bounded:��r(st, at)�� ≤ rmaxand��d(st, at)�� ≤ dmax. Also, the parameter space  Θis assumed to be compact.

In this section, we develop an algorithm to solve the optimization problem (4). Note that both the objective function and the constraint in (4) are nonconvex and involve expectation without closedform expression. As a constrained problem, a straightforward approach to solve (4) is to define the following Lagrangian function


and solve the dual problem


However, this problem is a nonconvex minimax problem and, therefore, is hard to solve and establish theoretical guarantees for solutions [2]. Another approach to solve (4) is to replace  J(θ)and  D(θ)by surrogate functions with nice properties. For example, one can iteratively construct local quadratic approximations that are strongly convex [52], or are an upper bound for the original function [57]. However, an immediate problem of this naive approach is that, even if the original problem (4) is feasible, the convex relaxation problem need not be. Also, these methods only deal with deterministic and/or convex constraints.

In this work, we propose an iterative algorithm that approximately solves (4) by constructing a sequence of convex relaxations, inspired by [38]. Our method is able to handle the possible infeasible situation due to the convex relaxation as mentioned above, and handle stochastic and nonconvex constraint. Since we do not have access to  J(θ)or  D(θ), we first define the sample negative cumulative reward and cost functions as


Given  θ, J∗(θ) and D∗(θ)are the sample negative cumulative reward and cost value of a realization (i.e., a trajectory) following policy  πθ. Note that both  J∗(θ)and  D∗(θ)are stochastic due to the randomness in the policy, state transition distribution, etc. With some abuse of notation, we use J∗(θ) and D∗(θ)to denote both a function of  θand a value obtained by the realization of a trajectory. Clearly we have  J(θ) = E�J∗(θ)�and D(θ) = E�D∗(θ)�.

We start from some (possibly infeasible)  θ0. Let  θkdenote the estimate of the policy parameter in the k-th iteration. As mentioned above, we do not have access to the expected cumulative reward J(θ). Instead we sample a trajectory following the current policy  πθkand obtain a realization of the negative cumulative reward value and the gradient of it as  J∗(θk)and  ∇θJ∗(θk), respectively. The cumulative reward value is obtained by Monte-Carlo estimation, and the gradient is also obtained by Monte-Carlo estimation according to the policy gradient theorem in (3). We provide more details on the realization step later in this section. Similarly, we use the same procedure for the cost function and obtain realizations  D∗(θk) and ∇θD∗(θk).

We approximate  J(θ) and D(θ) at θkby the quadratic surrogate functions


where  τ > 0is any fixed constant. In each iteration, we solve the optimization problem


where we define


with the initial value  J(0)(θ) = D(0)(θ) = 0. Here  ρkis the weight parameter to be specified later. According to the definition (5) and (6), problem (7) is a convex quadratically constrained quadratic program (QCQP). Therefore, it can be efficiently solved by, for example, the interior point method. However, as mentioned before, even if the original problem (4) is feasible, the convex relaxation problem (7) could be infeasible. In this case, we instead solve the following feasibility problem


In particular, we relax the infeasible constraint and find  θkas the solution that gives the minimum relaxation. Due to the specific form in (6),  D(k)(θ)is decomposable into quadratic forms of each component of  θ, with no terms involving  θi · θj. Therefore, the solution to problem (9) can be written in a closed form. Given  θkfrom either (7) or (9), we update  θk by


where  ηkis the learning rate to be specified later. Note that although we consider only one constraint in the algorithm, both the algorithm and the theoretical result in Section 4 can be directly generalized to multiple constraints setting. The whole procedure is summarized in Algorithm 1.


Obtaining realizations  J∗(θk) and ∇θJ∗(θk).We detail how to obtain realizations  J∗(θk)and ∇θJ∗(θk)corresponding to the lines 3 and 4 in Algorithm 1. The realizations of  D∗(θk)and ∇θD∗(θk)can be obtained similarly.

First, we discuss finite horizon setting, where we can sample the full trajectory according to the policy  πθ. In particular, for any  θk, we use the policy  πθkto sample a trajectory and obtain  J∗(θk)by Monte-Carlo method. The gradient  ∇θJ(θ)can be estimated by the policy gradient theorem [59],


Again we can sample a trajectory and obtain the policy gradient realization  ∇θJ∗(θk)by Monte-Carlo method.

In infinite horizon setting, we cannot sample the infinite length trajectory. In this case, we utilize the truncation method introduced in [48], which truncates the trajectory at some stage T and scales the undiscounted cumulative reward to obtain an unbiased estimation. Intuitively, if the discount factor γis close to 0, then the future reward would be discounted heavily and, therefore, we can obtain an accurate estimate with a relatively small number of stages. On the other hand, if  γis close to 1, then the future reward is more important compared to the small  γcase and we have to sample a long trajectory. Taking this intuition into consideration, we define T to be a geometric random variable with parameter  1 − γ: Pr(T = t) = (1 − γ)γt. Then, we simulate the trajectory until stage T and use the estimator  Jtruncate(θ) = −(1 − γ) · �Tt=0 r(st, at), which is an unbiased estimator of the expected negative cumulative reward  J(θ), as proved in proposition 5 in [43]. We can apply the same truncation procedure to estimate the policy gradient  ∇θJ(θ).

Variance reduction. Using the naive sampling method described above, we may suffer from high variance problem. To reduce the variance, we can modify the above procedure in the following ways. First, instead of sampling only one trajectory in each iteration, a more practical and stable way is to sample several trajectories and take average to obtain the realizations. As another approach, we can subtract a baseline function from the action-value function  Qθ(s, a)in the policy gradient estimation step (11) to reduce the variance without changing the expectation. A popular choice of the baseline function is the state-value function  V θ(s)as defined in (2). In this way, we can replace  Qθ(s, a)in (11) by the advantage function  Aθ(s, a) defined as


This modification corresponds to the standard REINFORCE with Baseline algorithm [58] and can significantly reduce the variance of policy gradient.

Actor-critic method. Finally, we can use an actor-critic update to improve the performance further. In this case, since we need unbiased estimators for both the gradient and the reward value in (5) and (6) in online fashion, we modify our original problem (4) to average reward setting as



Let  V Jθ (s)and  V Dθ (s)denote the value and cost functions corresponding to (2). We use possibly nonlinear approximation with parameter w for the value function:  V Jw (s) and vfor the cost function: V Dv (s). In the critic step, we update w and v by TD(0) with step size  βw and βv; in the actor step, we solve our proposed convex relaxation problem to update  θ. The actor-critic procedure is summarized in Algorithm 2. Here J and D are estimators of  J(θk)and  D(θk). Both of J and D, and the TD error  δJ, δD can be initialized as 0.

The usage of the actor-critic method helps reduce variance by using a value function instead of Monte-Carlo sampling. Specifically, in Algorithm 1 we need to obtain a sample trajectory and calculate  J∗(θ) and ∇θJ∗(θ)by Monte-Carlo sampling. This step has a high variance since we need to sample a potentially long trajectory and sum up a lot of random rewards. In contrast, in Algorithm 2, this step is replaced by a value function  V Jw (s), which reduces the variance.

In this section, we show almost sure convergence of the iterates obtained by our algorithm to a stationary point. We start by stating some mild assumptions on the original problem (4) and the choice of some parameters in Algorithm 1.

Assumption 1 The choice of  {ηk}and  {ρk}satisfy  limk→∞�k ηk = ∞, limk→∞�k ρk = ∞and  limk→∞�k η2k + ρ2k < ∞. Furthermore, we have  limk→∞ ηk/ρk = 0 and ηkis decreasing.

Assumption 2 For any realization,  J∗(θ) and D∗(θ)are continuously differentiable as functions of θ. Moreover,  J∗(θ), D∗(θ), and their derivatives are uniformly Lipschitz continuous.

Assumption 1 allows us to specify the learning rates. A practical choice would be  ηk = k−c1and ρk = k−c2with  0.5 < c2 < c1 < 1. This assumption is standard for gradient-based algorithms. Assumption 2 is also standard and is known to hold for a number of models. It ensures that the reward and cost functions are sufficiently regular. In fact, it can be relaxed such that each realization is Lipschitz (not uniformly), and the event that we keep generating realizations with monotonically increasing Lipschitz constant is an event with probability 0. See condition iv) in [67] and the discussion thereafter. Also, see [45] for sufficient conditions such that both the expected cumulative reward function and the gradient of it are Lipschitz.

The following Assumption 3 is useful only when we initialize with an infeasible point. We first state it here and we will discuss this assumption after the statement of the main theorem.

Assumption 3 Suppose  (θS, αS)is a stationary point of the optimization problem


We have that  θSis a feasible point of the original problem (4), i.e.  D(θS) ≤ D0.

We are now ready to state the main theorem.

Theorem 4 Suppose the Assumptions 1 and 2 are satisfied with small enough initial step size  η0. Suppose also that, either  θ0is a feasible point, or Assumption 3 is satisfied. If there is a subsequence {θkj}of  {θk}that converges to some �θ, then there exist uniformly continuous functions �J(θ)and �D(θ)satisfying


Furthermore, suppose there exists  θ such that �D(θ) < D0(i.e. the Slater’s condition holds), then �θ isa stationary point of the original problem (4) almost surely.

The proof of Theorem 4 is provided in the supplementary material.

Note that Assumption 3 is not necessary if we start from a feasible point, or we reach a feasible point in the iterates, which could be viewed as an initializer. Assumption 3 makes sure that the iterates in Algorithm 1 keep making progress without getting stuck at any infeasible stationary point. A similar condition is assumed in [38] for an infeasible initializer. If it turns out that  θ0is infeasible and Assumption 3 is violated, then the convergent point may be an infeasible stationary point of (12). In practice, if we can find a feasible point of the original problem, then we proceed with that point. Alternatively, we could generate multiple initializers and obtain iterates for all of them. As long as there is a feasible point in one of the iterates, we can view this feasible point as the initializer and Theorem 4 follows without Assumption 3. In our later experiments, for every single replicate, we could reach a feasible point, and therefore Assumption 3 is not necessary.

Our algorithm does not guarantee safe exploration during the training phase. Ensuring safety during learning is a more challenging problem. Sometimes even finding a feasible point is not straightforward, otherwise Assumption 3 is not necessary.

Our proposed algorithm is inspired by [38]. Compared to [38] which deals with an optimization problem, solving the safe reinforcement learning problem is more challenging. We need to verify that the Lipschitz condition is satisfied, and also the policy gradient has to be estimated (instead of directly evaluated as in a standard optimization problem). The usage of the Actor-Critic algorithm reduces the variance of the sampling, which is unique to Reinforcement learning.

We apply our algorithm to the linear-quadratic regulator (LQR), which is one of the most fundamental problems in control theory. In the LQR setting, the state dynamic equation is linear, the cost function is quadratic, and the optimal control theory tells us that the optimal control for LQR is a linear function of the state [23, 6]. LQR can be viewed as an MDP problem and it has attracted a lot of attention in the reinforcement learning literature [12, 13, 21, 47].

We consider the infinite-horizon, discrete-time LQR problem. Denote  xtas the state variable and  utas the control variable. The state transition and the control sequence are given by


where  vt and wtrepresent possible Gaussian white noise, and the initial state is given by  x0. The goalis to find the control parameter matrix F such that the expected total cost is minimized. The usual cost function of LQR corresponds to the negative reward in our setting and we impose an additional quadratic constraint on the system. The overall optimization problem is given by


where  Q1, Q2, R1,and  R2are positive definite matrices. Note that even thought the matrices are positive definite, both the objective function J and the constraint D are nonconvex with respect to the parameter F. Furthermore, with the additional constraint, the optimal control sequence may no longer be linear in the state  xt. Nevertheless, in this work, we still consider linear control given by (13) and the goal is to find the best linear control for this constrained LQR problem. We assume that the choice of A, B are such that the optimal cost is finite.

Random initial state. We first consider the setting where the initial state  x0 ∼ Dfollows a random distribution D, while both the state transition and the control sequence (13) are deterministic (i.e. vt = wt = 0). In this random initial state setting, [24] showed that without the constraint, the policy gradient method converges efficiently to the global optima in polynomial time. In the constrained case, we can explicitly write down the objective and constraint function, since the only randomness comes from the initial state. Therefore, we have the state dynamic  xt+1 = (A − BF)xtand the objective function has the following expression ([24], Lemma 1)


where  PFis the solution to the following equation


The gradient is given by


We start from some  F0and apply our Algorithm 1 to solve the constrained LQR problem. In iteration k, with the current estimator denoted by  Fk, we first obtain an estimator of  PFkby starting from Q1and iteratively applying the recursion  PFk ← Q1 + F ⊤k R1Fk + (A − BFk)⊤PFk(A − BFk)until convergence. Next, we sample an  x∗0from the distribution D and follow a similar recursion given by (17) to obtain an estimate of  SFk. Plugging the sample  x∗0and the estimates of  PFkand SFkinto (14) and (16), we obtain the sample reward value  J∗(Fk)and  ∇F J∗(Fk), respectively. With these two values, we follow (5) and (8) and obtain  J(k)(F). We apply the same procedure to the cost function D(F) with  Q1, R1replaced by  Q2, R2to obtain  D(k)(F). Finally we solve the optimization problem (7) (or (9) if (7) is infeasible) and obtain  Fk+1 by (10).

Random state transition and control. We then consider the setting where both  vtand  wtare independent standard Gaussian white noise. In this case, the state dynamic can be written as xt+1 = (A − BF)xt + ϵt where ϵt ∼ N(0, I + BB⊤). Let PFbe defined as in (15) and  SF be thesolution to the following Lyapunov equation


The objective function has the following expression ([68], Proposition 3.1)


and the gradient is given by


Although in this setting it is straightforward to calculate the expectation in a closed form, we keep the current expectation form to be in line with our algorithm. Moreover, when the error distribution is more complicated or unknown, we can no longer calculate the closed form expression and have to sample in each iteration. With the formulas given by (18) and (19), we again apply our Algorithm 1. We sample  x ∼ N(0, SF )in each iteration and solve the optimization problem (7) or (9). The whole procedure is similar to the random initial state case described above.

Other applications. Our algorithm can also be applied to constrained parallel MDP and constrained multi-agent MDP problem. Due to the space limit, we relegate them to supplementary material.


Figure 1: An experiment on constrained LQR problem. The iterate starts from an infeasible point and then becomes feasible and eventually converges.


Table 1: Comparison of our method with Lagrangian method

We verify the effectiveness of the proposed algorithm through experiments. We focus on the LQR setting with a random initial state as discussed in Section 5. In this experiment we set  x ∈ R15and u ∈ R8. The initial state distribution is uniform on the unit cube:  x0 ∼ D = Uniform�[−1, 1]15�. Each element of A and B is sampled independently from the standard normal distribution and scaled such that the eigenvalues of A are within the range  (−1, 1). We initialize  F0as an all-zero matrix, and the choice of the constraint function and the value  D0are such that (1) the constrained problem is feasible; (2) the solution of the unconstrained problem does not satisfy the constraint, i.e., the problem is not trivial; (3) the initial value  F0is not feasible. The learning rates are set as  ηk = 23k−3/4and ρk = 23k−2/3. The conservative choice of step size is to avoid the situation where an eigenvalue of A − BFruns out of the range  (−1, 1), and so the system is stable. 5

Figure 1(a) and 1(b) show the constraint and objective value in each iteration, respectively. The red horizontal line in Figure 1(a) is for  D0, while the horizontal line in Figure 1(b) is for the unconstrained minimum objective value. We can see from Figure 1(a) that we start from an infeasible point, and the problem becomes feasible after about 100 iterations. The objective value is in general decreasing after becoming feasible, but never lower than the unconstrained minimum, as shown in Figure 1(b).

Comparison with the Lagrangian method. We compare our proposed method with the usual Lagrangian method. For the Lagrangian method, we follow the algorithm proposed in [18] for safe reinforcement learning, which iteratively applies gradient descent on the parameter F and gradient ascent on the Lagrangian multiplier  λfor the Lagrangian function until convergence.

Table 1 reports the comparison results with mean and standard deviation based on 50 replicates. In the second and third columns, we compare the minimum objective value and the number of iterations to achieve it. We also consider an approximate version, where we are satisfied with the result if the objective value exceeds less than 0.02% of the minimum value. The fourth and fifth columns show the comparison results for this approximate version. We can see that both methods achieve similar minimum objective values, but ours requires less number of policy updates, for both minimum and approximate minimum version.

[1] Joshua Achiam, David Held, Aviv Tamar, and Pieter Abbeel. Constrained policy optimization. In International Conference on Machine Learning, pages 22–31, 2017.

[2] Leonard Adolphs. Non convex-concave saddle point optimization. Master’s thesis, ETH Zurich, 2018.

[3] Eitan Altman. Constrained Markov decision processes, volume 7. CRC Press, 1999.

[4] Haitham Bou Ammar, Rasul Tutunov, and Eric Eaton. Safe policy search for lifelong reinforcement learning with sublinear regret. In International Conference on Machine Learning, pages 2361–2369, 2015.

[5] Dario Amodei, Chris Olah, Jacob Steinhardt, Paul Christiano, John Schulman, and Dan Man´e. Concrete problems in ai safety. arXiv preprint arXiv:1606.06565, 2016.

[6] Brian DO Anderson and John B Moore. Optimal control: linear quadratic methods. Courier Corporation, 2007.

[7] Yu Bai, Tengyang Xie, Nan Jiang, and Yu-Xiang Wang. Provably efficient q-learning with low switching cost. arXiv preprint arXiv:1905.12849, 2019.

[8] Marc G Bellemare, Yavar Naddaf, Joel Veness, and Michael Bowling. The arcade learning environment: An evaluation platform for general agents. Journal of Artificial Intelligence Research, 47:253–279, 2013.

[9] Felix Berkenkamp, Matteo Turchetta, Angela Schoellig, and Andreas Krause. Safe model-based reinforcement learning with stability guarantees. In Advances in Neural Information Processing Systems, pages 908–918, 2017.

[10] Vivek S Borkar. Stochastic approximation with two time scales. Systems & Control Letters, 29(5):291–294, 1997.

[11] Craig Boutilier. Planning, learning and coordination in multiagent decision processes. In Proceedings of the 6th conference on Theoretical aspects of rationality and knowledge, pages 195–210. Morgan Kaufmann Publishers Inc., 1996.

[12] Steven J Bradtke. Reinforcement learning applied to linear quadratic regulation. In Advances in neural information processing systems, pages 295–302, 1993.

[13] Steven J Bradtke, B Erik Ydstie, and Andrew G Barto. Adaptive linear quadratic control using policy iteration. In Proceedings of the American control conference, volume 3, pages 3475–3475. Citeseer, 1994.

[14] Leo Breiman. Random forests. Machine learning, 45(1):5–32, 2001.

[15] Lucian Busoniu, Robert Babuska, and Bart De Schutter. A comprehensive survey of multiagent reinforcement learning. IEEE Transactions on Systems, Man, And Cybernetics-Part C: Applications and Reviews, 38 (2), 2008, 2008.

[16] Qi Cai, Zhuoran Yang, Jason D Lee, and Zhaoran Wang. Neural temporal-difference learning converges to global optima. arXiv preprint arXiv:1905.10027, 2019.

[17] Tianyi Chen, Kaiqing Zhang, Georgios B Giannakis, and Tamer Bas¸ar. Communication-efficient distributed reinforcement learning. arXiv preprint arXiv:1812.03239, 2018.

[18] Yinlam Chow, Mohammad Ghavamzadeh, Lucas Janson, and Marco Pavone. Risk-constrained reinforcement learning with percentile risk criteria. Journal of Machine Learning Research, 18(167):1–167, 2017.

[19] Yinlam Chow, Ofir Nachum, Edgar Duenez-Guzman, and Mohammad Ghavamzadeh. A lyapunov-based approach to safe reinforcement learning. arXiv preprint arXiv:1805.07708, 2018.

[20] Christoph Dann, Gerhard Neumann, and Jan Peters. Policy evaluation with temporal differences: A survey and comparison. The Journal of Machine Learning Research, 15(1):809–883, 2014.

[21] Sarah Dean, Horia Mania, Nikolai Matni, Benjamin Recht, and Stephen Tu. On the sample complexity of the linear quadratic regulator. arXiv preprint arXiv:1710.01688, 2017.

[22] Nelson Dunford and Jacob T Schwartz. Linear operators part I: general theory, volume 7. Interscience publishers New York, 1958.

[23] Lawrence C Evans. An introduction to mathematical optimal control theory. Lecture Notes, University of California, Department of Mathematics, Berkeley, 2005.

[24] Maryam Fazel, Rong Ge, Sham Kakade, and Mehran Mesbahi. Global convergence of policy gradient methods for the linear quadratic regulator. In International Conference on Machine Learning, pages 1466–1475, 2018.

[25] Jaime F Fisac, Anayo K Akametalu, Melanie N Zeilinger, Shahab Kaynama, Jeremy Gillula, and Claire J Tomlin. A general safety framework for learning-based control in uncertain robotic systems. IEEE Transactions on Automatic Control, 2018.

[26] Javier Garcıa and Fernando Fern´andez. A comprehensive survey on safe reinforcement learning. Journal of Machine Learning Research, 16(1):1437–1480, 2015.

[27] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio. Generative adversarial nets. In Advances in neural information processing systems, pages 2672–2680, 2014.

[28] Ivo Grondman, Lucian Busoniu, Gabriel AD Lopes, and Robert Babuska. A survey of actor-critic reinforcement learning: Standard and natural policy gradients. IEEE Transactions on Systems, Man, and Cybernetics, Part C (Applications and Reviews), 42(6):1291–1307, 2012.

[29] Jingyu He, Saar Yalov, and P Richard Hahn. XBART: Accelerated Bayesian additive regression trees. In The 22nd International Conference on Artificial Intelligence and Statistics, pages 1130–1138, 2019.

[30] Jingyu He, Saar Yalov, Jared Murray, and P Richard Hahn. Stochastic tree ensembles for regularized supervised learning. Technical report, 2019.

[31] Jessie Huang, Fa Wu, Doina Precup, and Yang Cai. Learning safe policies with expert guidance. arXiv preprint arXiv:1805.08313, 2018.

[32] John L Kelley. General topology. Courier Dover Publications, 2017.

[33] Vijay R Konda and John N Tsitsiklis. Actor-critic algorithms. In Advances in neural information processing systems, pages 1008–1014, 2000.

[34] R Matthew Kretchmar. Parallel reinforcement learning. In The 6th World Conference on Systemics, Cybernetics, and Informatics. Citeseer, 2002.

[35] Jonathan Lacotte, Yinlam Chow, Mohammad Ghavamzadeh, and Marco Pavone. Risk-sensitive generative adversarial imitation learning. arXiv preprint arXiv:1808.04468, 2018.

[36] Dennis Lee, Haoran Tang, Jeffrey O Zhang, Huazhe Xu, Trevor Darrell, and Pieter Abbeel. Modular architecture for starcraft ii with deep reinforcement learning. In Fourteenth Artificial Intelligence and Interactive Digital Entertainment Conference, 2018.

[37] Yuxi Li. Deep reinforcement learning: An overview. arXiv preprint arXiv:1701.07274, 2017.

[38] An Liu, Vincent Lau, and Borna Kananian. Stochastic successive convex approximation for non-convex constrained stochastic optimization. arXiv preprint arXiv:1801.08266, 2018.

[39] Boyi Liu, Qi Cai, Zhuoran Yang, and Zhaoran Wang. Neural proximal/trust region policy optimization attains globally optimal policy. arXiv preprint arXiv:1906.10306, 2019.

[40] Volodymyr Mnih, Adria Puigdomenech Badia, Mehdi Mirza, Alex Graves, Timothy Lillicrap, Tim Harley, David Silver, and Koray Kavukcuoglu. Asynchronous methods for deep reinforcement learning. In International conference on machine learning, pages 1928–1937, 2016.

[41] Volodymyr Mnih, Koray Kavukcuoglu, David Silver, Andrei A Rusu, Joel Veness, Marc G Bellemare, Alex Graves, Martin Riedmiller, Andreas K Fidjeland, Georg Ostrovski, et al. Human-level control through deep reinforcement learning. Nature, 518(7540):529, 2015.

[42] Arun Nair, Praveen Srinivasan, Sam Blackwell, Cagdas Alcicek, Rory Fearon, Alessandro De Maria, Vedavyas Panneershelvam, Mustafa Suleyman, Charles Beattie, Stig Petersen, et al. Massively parallel methods for deep reinforcement learning. arXiv preprint arXiv:1507.04296, 2015.

[43] Santiago Paternain. Stochastic Control Foundations of Autonomous Behavior. PhD thesis, University of Pennsylvania, 2018.

[44] Peng Peng, Ying Wen, Yaodong Yang, Quan Yuan, Zhenkun Tang, Haitao Long, and Jun Wang. Multiagent bidirectionally-coordinated nets: Emergence of human-level coordination in learning to play starcraft combat games. arXiv preprint arXiv:1703.10069, 2017.

[45] Matteo Pirotta, Marcello Restelli, and Luca Bascetta. Policy gradient in lipschitz markov decision processes. Machine Learning, 100(2-3):255–283, 2015.

[46] LA Prashanth and Mohammad Ghavamzadeh. Variance-constrained actor-critic algorithms for discounted and average reward mdps. Machine Learning, 105(3):367–417, 2016.

[47] Benjamin Recht. A tour of reinforcement learning: The view from continuous control. Annual Review of Control, Robotics, and Autonomous Systems, 2018.

[48] Chang-han Rhee and Peter W Glynn. Unbiased estimation with square root convergence for sde models. Operations Research, 63(5):1026–1043, 2015.

[49] Andrzej Ruszczy´nski. Feasible direction methods for stochastic programming problems. Mathematical Programming, 19(1):220–229, 1980.

[50] Joris Scharpff, Diederik M Roijers, Frans A Oliehoek, Matthijs TJ Spaan, and Mathijs Michiel de Weerdt. Solving transition-independent multi-agent mdps with sparse interactions. In AAAI, pages 3174–3180, 2016.

[51] John Schulman, Sergey Levine, Pieter Abbeel, Michael Jordan, and Philipp Moritz. Trust region policy optimization. In International Conference on Machine Learning, pages 1889–1897, 2015.

[52] Gesualdo Scutari, Francisco Facchinei, Peiran Song, Daniel P Palomar, and Jong-Shi Pang. Decomposition by partial linearization: Parallel optimization of multi-agent systems. IEEE Transactions on Signal Processing, 62(3):641–656, 2013.

[53] David Silver, Aja Huang, Chris J Maddison, Arthur Guez, Laurent Sifre, George Van Den Driessche, Julian Schrittwieser, Ioannis Antonoglou, Veda Panneershelvam, Marc Lanctot, et al. Mastering the game of go with deep neural networks and tree search. nature, 529(7587):484, 2016.

[54] David Silver, Thomas Hubert, Julian Schrittwieser, Ioannis Antonoglou, Matthew Lai, Arthur Guez, Marc Lanctot, Laurent Sifre, Dharshan Kumaran, Thore Graepel, et al. A general reinforcement learning algorithm that masters chess, shogi, and go through self-play. Science, 362(6419):1140–1144, 2018.

[55] David Silver, Julian Schrittwieser, Karen Simonyan, Ioannis Antonoglou, Aja Huang, Arthur Guez, Thomas Hubert, Lucas Baker, Matthew Lai, Adrian Bolton, et al. Mastering the game of go without human knowledge. Nature, 550(7676):354, 2017.

[56] Peng Sun, Xinghai Sun, Lei Han, Jiechao Xiong, Qing Wang, Bo Li, Yang Zheng, Ji Liu, Yongsheng Liu, Han Liu, et al. Tstarbots: Defeating the cheating level builtin ai in starcraft ii in the full game. arXiv preprint arXiv:1809.07193, 2018.

[57] Ying Sun, Prabhu Babu, and Daniel P Palomar. Majorization-minimization algorithms in signal processing, communications, and machine learning. IEEE Transactions on Signal Processing, 65(3):794–816, 2016.

[58] Richard S Sutton and Andrew G Barto. Reinforcement learning: An introduction. MIT press, 2018.

[59] Richard S Sutton, David A McAllester, Satinder P Singh, and Yishay Mansour. Policy gradient methods for reinforcement learning with function approximation. In Advances in neural information processing systems, pages 1057–1063, 2000.

[60] Chen Tessler, Daniel J Mankowitz, and Shie Mannor. Reward constrained policy optimization. arXiv preprint arXiv:1805.11074, 2018.

[61] Matteo Turchetta, Felix Berkenkamp, and Andreas Krause. Safe exploration in finite markov decision processes with gaussian processes. In Advances in Neural Information Processing Systems, pages 4312–4320, 2016.

[62] Oriol Vinyals, Timo Ewalds, Sergey Bartunov, Petko Georgiev, Alexander Sasha Vezhnevets, Michelle Yeo, Alireza Makhzani, Heinrich K¨uttler, John Agapiou, Julian Schrittwieser, et al. Starcraft ii: A new challenge for reinforcement learning. arXiv preprint arXiv:1708.04782, 2017.

[63] Hoi-To Wai, Zhuoran Yang, Zhaoran Wang, and Mingyi Hong. Multi-agent reinforcement learning via double averaging primal-dual optimization. arXiv preprint arXiv:1806.00877, 2018.

[64] Lingxiao Wang, Qi Cai, Zhuoran Yang, and Zhaoran Wang. Neural policy gradient methods: Global optimality and rates of convergence. arXiv preprint arXiv:1909.01150, 2019.

[65] Min Wen and Ufuk Topcu. Constrained cross-entropy method for safe reinforcement learning. In Advances in Neural Information Processing Systems, pages 7461–7471, 2018.

[66] Sijia Xu, Hongyu Kuang, Zhi Zhuang, Renjie Hu, Yang Liu, and Huyang Sun. Macro action selection with deep reinforcement learning in starcraft. arXiv preprint arXiv:1812.00336, 2018.

[67] Yang Yang, Gesualdo Scutari, Daniel P Palomar, and Marius Pesavento. A parallel decomposition method for nonconvex stochastic multi-agent optimization problems. IEEE Transactions on Signal Processing, 64(11):2949–2964, 2016.

[68] Zhuoran Yang, Yongxin Chen, Mingyi Hong, and Zhaoran Wang. On the global convergence of actor-critic: A case for linear quadratic regulator with ergodic cost. arXiv preprint arXiv:1907.06246, 2019.

[69] Kaiqing Zhang, Zhuoran Yang, Han Liu, Tong Zhang, and Tamer Bas¸ar. Finite-sample analyses for fully decentralized multi-agent reinforcement learning. arXiv preprint arXiv:1812.02783, 2018.

A.1 Constrained Parallel Markov Decision Process

We consider the parallel MDP problem [34, 42, 17] where we have a single-agent MDP task and N workers, where each worker acts as an individual agent and aims to solve the same MDP problem. In the parallel MDP setting, each agent is characterized by a tuple  (S, A, P, γ, ri, di, µi), where each agent has the same but individual state space, action space, transition probability distribution, and the discount factor. However, the reward function, cost function, and the distribution of the initial state s0 ∈ Scould be different for each agent, but satisfy  E[ri(s, a)] = r(s, a), E[di(s, a)] = d(s, a), andE[µi(s, a)] = µ(s, a). Each agent i generates its own trajectory  {si0, ai0, si1, ai1, . . . }and collects its own reward/cost value  {ri0, di0, ri1, di1, . . . }.

The hope is that by solving the single-agent problem using N agents in parallel, the algorithm could be more stable and converge much faster [40]. Intuitively, each agent i may have a different initial state and will explore different parts of the state space due to the randomness in the state transition distribution and the policy. It also helps to reduce the correlation between agents’ behaviors. As a result, by running multiple agents in parallel, we are more likely to visit different parts of the environment and get the experience of the reward/cost function values more efficiently. This mimics the strategy used in tree-based supervised learning algorithms [14, 29, 30].

Following the settings in [17], we have N agents (i.e., N workers) and one central controller in the system. The global parameter is denoted by  θ, and we consider the constrained parallel MDP problem where the goal is to solve the following optimization problem:


During the estimation step, the controller broadcasts the current parameter  θkto each agent and each agent samples its own trajectory and obtains estimators for function value/gradient of the reward/cost function. Next, each agent uploads its estimators to the central controller and the central controller takes the average over these estimators, and then follow our proposed algorithm to solve for the QCQP problem and update the parameter to  θk+1. This process continues until convergence.

A.2 Constrained Multi-agent Markov Decision Process

A natural extension of the (single-agent) MDP is to consider a model with N agents termed multi-agent Markov decision process (MMDP). Recently this kind of problem has been attracting more and more attention. See [15] for a comprehensive survey. Most of the work on multi-agent MDP problems consider the setting where the agents share the same global state space, but each with their own collection of actions and rewards [11, 63, 69]. In each stage of the system, each agent observes the global state and chooses its own action individually. As a result, each agent receives its reward and the state evolves according to the joint transition distribution. An MMDP problem can be fully collaborative where all the agents have the same goal, or fully competitive where the problem consists of two agents with an opposite goal, or the mix of the two.

Here we consider a slightly different setting where each agent has its own state space. The only connection between the agents is that the global reward is a function of the overall states and actions. Furthermore, each agent has its own constraint which depends on its own state and action only. This problem is known as Transition-Independent Multi-agent MDP and is considered in [50]. Specifically, each agent’s task is characterized by a tuple  (Si, Ai, P i, γ, di, µi)with each component defined as usual. Note that  P i : Si × Ai → P(Si)and  di : Si × Ai → Rare functions of each agent’s state and action only and do not depend on other agents. Denote  S = Πi∈N Siand  A = Πi∈N Aias the joint state space and action space. The global reward function is given by  r: S × A → Rthat depends on the joint state and action. Here we consider the fully collaborative setting where all the agents have the same goal. Under this setting, the policy set of each agent is parameterized as {πiθi : Si → P(Ai)}and we denote  θ = [θ1, . . . , θN]as the overall parameters and  πθas the overall policy. In the following, we use N = {1, 2, . . . , N} to denote the N agents. Denote  ait as the action chosen by agent  i at stage t and at = Πi∈N ait as the joint action chosen by all the agents. The goal of this constrained MMDP is to solve the following problem


Inspired by the parallel implementation ([38], Section V), our algorithm applies naturally to constrained MMDP problem with some modifications. This modified procedure can also be viewed as a distributed version of the original algorithm. The overall problem (20) can be viewed as a large “single-agent” problem where the constraints are decomposable into N parts. In this case, instead of solving a large QCQP problem in each iteration, each agent could solve its own QCQP problem in a distributed manner which is much more efficient. As before, we denote the sample negative reward and cost function as


Note that the constraint function is naturally decomposable into N parts. We also “manually” split the objective function into N parts, so that each agent could solve its own QCQP problem in a distributed manner. As before, we define


With this surrogate functions, each agent then solves its own convex relaxation problem


or, alternatively, solves for the feasibility problem if (21) is infeasible


This step can be implemented in a distributed manner for each agent and is more efficient than solving the overall problem with the overall parameter  θ. Finally, the update rule for each agent i is as usual


This process continues until convergence.

According to the choice of the surrogate function (5) and Assumption 2, it is straightforward to verify that the function  J(k)(θ)defined in (8) is uniformly strongly convex in  θfor each iteration t.

Moreover, both  J(k)(θ) and ∇θJ(k)(θ)are Lipschitz continuous functions.

From Lemma 1 in [49] we have


for some constant  L0and the error term  ϵthat goes to 0 as  k1, k2go to infinity. This shows that the function sequence  J(kj)(θ)is equicontinuous. Since  Θis compact and the discounted cumulative reward function is bounded by  rmax/(1 − γ), we can apply Arzela-Ascoli theorem [22, 32] to prove existence of �J(θ)that converges uniformly. Moreover, since we apply the same operations on the constraint function  D(θ)as to the reward function  J(θ)in Algorithm 1, the above properties also hold for  D(θ).

The rest of the proof follows in a similar way as the proof of Theorem 1 in [38]. Under Assumptions 1 - 3, the technical conditions in [38] are satisfied by the choice of the surrogate functions (5) and (6). According to Lemma 2 in [38], with probability one we have


This shows that, although in some of the iterations the convex relaxation problem (7) is infeasible, and we have to solve the alternative problem (9), the iterates  {θk}converge to the feasible region of the original problem (4) with probability one. Furthermore, with probability one, the convergent point �θis the optimal solution to the following problem


The KKT conditions for (22) together with the Slater condition show that the KKT conditions of the original problem (4) are also satisfied at �θ. This shows that �θis a stationary point of the original problem almost surely.

Designed for Accessibility and to further Open Science