text:
Distributed Distance-Bounded Network Design Through Distributed Convex Programming

Solving linear programs is often a challenging task in distributed settings. While there are good algorithms for solving packing and covering linear programs in a distributed manner (Kuhn et al.~2006), this is essentially the only class of linear programs for which such an algorithm is known. In this work we provide a distributed algorithm for solving a different class of convex programs which we call"distance-bounded network design convex programs". These can be thought of as relaxations of network design problems in which the connectivity requirement includes a distance constraint (most notably, graph spanners). Our algorithm runs in $O( (D/\epsilon) \log n)$ rounds in the $\mathcal{LOCAL}$ model and finds a $(1+\epsilon)$-approximation to the optimal LP solution for any $0<\epsilon \leq 1$, where $D$ is the largest distance constraint. While solving linear programs in a distributed setting is interesting in its own right, this class of convex programs is particularly important because solving them is often a crucial step when designing approximation algorithms. Hence we almost immediately obtain new and improved distributed approximation algorithms for a variety of network design problems, including Basic $3$- and $4$-Spanner, Directed $k$-Spanner, Lowest Degree $k$-Spanner, and Shallow-Light Steiner Network Design with a spanning demand graph. Our algorithms do not require any"heavy"computation and essentially match the best-known centralized approximation algorithms, while previous approaches which do not use heavy computation give approximations which are worse than the best-known centralized bounds.

Introduction
Distributed network design is a classical type of distributed algorithmic problem, going back at least to the seminal work on distributed MST by Gallager, Humblet, and Spira [15]. By "network design", we mean the class of problems which can be phrased as "given input graph G, find a subgraph H which has some property P , and minimize the cost of H". Clearly different properties P , and different notions of cost, lead to very different problems. One important class of problems are distance-bounded network design problems, where the property P is that certain pairs of vertices are within some distance of each other in H (where distance refers to the shortest-path distance). The most well-known type of distance-bounded network design problems are problems involving graph spanners, in which the distance requirement is that the distance in H for all (or certain) pairs is within a certain factor (known as the stretch) of their original distance in H. But there are many other important versions of distance-bounded network design, such as the bounded diameter problem [12] and the shallow-light Steiner tree/network problems [18].
Many of these problems are NP-hard, so they cannot be solved optimally in polynomial time even in the centralized setting. Thus they have been studied extensively from an approximation algorithms point of view, where we design algorithms which approximate the optimal solution but which run in polynomial time. For many of these problems, a key step in the best-known centralized approximation algorithm is solving a linear programming relaxation of the problem, and then rounding the optimal fractional solution into a feasible integral solution. Interestingly, it is relatively common for the rounding to be "local": if we are in a distributed setting and happen to know the optimal fractional LP solution, then the algorithm used to round this to an integral solution can be accomplished with a tiny amount of extra time (either 0 or a small constant number of rounds). So the bottleneck when trying to make these algorithms distributed is solving the LP, not rounding it.
Solving LPs in distributed settings has received only a small amount of attention, since it unfortunately turns out to be extremely challenging in general. Most notably, Kuhn, Moscibroda, and Wattenhofer [20] gave an efficient distributed algorithm (in the LOCAL model of distributed computation) for packing/covering LPs. Unfortunately, the LPs used for distance-bounded network design are not packing/covering LPs 1 , and hence we are not able to use their techniques. In this paper we show how to solve these LPs (and convex generalizations of them) in the LOCAL model of distributed computation, which almost immediately gives the best-known results for a variety of distance-bounded network design problems.
In particular, for many network design problems (Directed k-Spanner, Basic 3-Spanner, Basic 4-Spanner, Lowest-Degree k-Spanner, Directed Steiner Network with Distance Constraints with spanning demands, and Shallow-Light Steiner Network with spanning demands) we give approximation algorithms which run in O(D log n) rounds (where D is the maximum distance bound) and have the same approximation ratios as in the centralized setting. Previous distributed algorithms for these problems with similar round complexity have either used "heavy" computations (non-polynomial time algorithms) at the nodes (in which case they can often do better than the best computationally-bounded centralized algorithm), or give approximation bounds which are asymptotically worse than the best centralized bounds. See Section 1.2 for more discussion of previous work.

Our Results
We give two main types of results. First, we give a distributed algorithm that (approximately) solves distance-bounded network design convex programs with small round complexity. We then use this result to (almost immediately) get improved distributed approximation algorithms for a variety of network design problems.

Solving convex programs
Stating our main technical result (distributed approximations of distance-bounded network design convex programs) in full generality requires significant technical setup, so we provide an informal description here. See Section 4 for the full definitions and theorem statements (Theorem 9 in particular). But informally, a distance-bounded network design convex program is the following.
We are given a graph G = (V, E), a set S ⊆ V × V , and for each (u, v) ∈ S there is a set of "allowed" u − v paths P u,v . Informally, the integral problem is to find a subgraph H of G so that every (u, v) ∈ S is connected by at least one path from P u,v in H, and the goal is to minimize some notion of "cost". If our notion of "cost" is captured by an objective function g : R |E| ≥0 → R (which is typically linear, but which can be more general convex functions as long as they satisfy a "partitionability" constraint -see Section 4 for the details), then the natural relaxation of this problem is the following convex program, which has a variable x e for every edge and a variable f P for every allowed path.
Informally, the first type of constraint says that an allowed path is included only if all edges in it are included, and the second type of constraint required us to include at least one allowed path for each (u, v) ∈ S. We call this type of convex program a distance-bounded network design convex program. It is clearly not a packing/covering LP due to the first type of constraint, and hence there is no known distributed algorithm to solve this kind of program. However, note that if the maximum length of any allowed path is constant, then there are only a polynomial number of such paths, and hence the size of the convex program is polynomial and so it can be solved in polynomial time in the centralized setting under reasonable assumptions on g (see [16] for details on solving convex programs in polynomial time).
Our main technical result is that we can approximately solve these optimization problems even in a distributed setting. For any path P let ℓ(P ) denote the length of the path (the number of edges in it). Theorem 1. For any constant ǫ > 0, any distance-bounded network design convex program can be solved up to a (1 + ǫ)-approximation in O(D log n) rounds in the LOCAL model, where D = max (u,v)∈S max P ∈Pu,v ℓ(P ). Moreover, if the convex program can be solved in polynomial time in the centralized sequential setting, then the distributed algorithm uses only polynomial-time computations at every node The dependence on ǫ in the above theorem is hidden in the O(·) notation -see Theorem 9 for the full statement.
Our main technique is to use a distributed construction of padded decompositions, a specific type of network decomposition which we explain in detail in Section 3. Padded decompositions have been very useful for metric embeddings and approximation algorithms (e.g., [17,19]), but to the best of our knowledge have not been used before in distributed algorithms (with the exception of [10], which used a special case of them to give a distributed algorithm for the fault-tolerant 2-spanner problem). Very similar decompositions, such as the famous Linial-Saks decomposition [21], have been used extensively in distributed settings, but the guarantees for padded decompositions are somewhat different (and we believe that these decompositions may prove useful in the future when designing distributed approximation algorithms). In Section 3 we give a distributed algorithm in the LOCAL model to construct padded decompositions. These padded decompositions allow us to solve a collection of "local" convex programs with the guarantees that a) most of the demands in S are satisfied in one of the local programs, and b) the solutions of the local convex programs combine into a (possibly infeasible) global solution with cost at most the cost of the global optimum. Then by averaging over O(log n) of these decompositions we get a feasible global solution which is almost optimal.

Distributed approximation algorithms for network design
Solving convex programming problems in distributed environments is interesting in its own right, and Theorem 1 is our main technical contribution, but the particular class of convex programs that we can solve are mostly interesting as convex relaxations of interesting combinatorial optimization problems. Many of the problems are NP-hard, but there has been significant work (some quite recent) on designing approximation algorithms for them (see, e.g., [5,6,9,11]). Almost all of these approximations depend on convex relaxations which fall into our class of "distance-bounded network design convex programs". This means that as long as the rounding scheme can be computed locally, we can design distributed versions of these approximation algorithms by using Theorem 1 to solve the appropriate convex relaxation and then using the local rounding scheme.
We are able to use this framework to give distributed approximation algorithms for several problems. Most of them are variations of graph spanners, which were introduced by Peleg and Ullman [25] and Peleg and Schäffer [24], and are defined as follows.
The value k is called the stretch of the spanner.
Before stating our results, we first define the problems. In the Basic k-Spanner problem we are given an undirected graph G and a value k ∈ N. A subgraph H of G is a feasible solution if it is a k-spanner of G, and the objective is to minimize the number of edges in H. For k = 3, 4, the best-known approximation algorithm for this problem isÕ(n 1/3 ) [5,11]. If the input graph G (and the solution H) are directed, then this is the Directed k-Spanner problem, for which the best-known approximation isÕ( √ n) [5]. If the objective is instead to minimize the maximum degree in H then this is the Lowest-Degree k-Spanner problem, for which the best-known approximation isÕ(n (1−1/k) 2 ) [6].
The following theorem contains our results on distributed approximations of graph spanners. Informally, it states that for we can give the same approximations in the the LOCAL model as are possible in the centralized model. We emphasize that our algorithms for these spanner problems both match the best-known centralized approximations and only use polynomial-time computations at each node. There is significant previous work (see Section 1.2) on designing distributed approximation algorithms for these and related problems that has only one of these two properties, but all previous approaches which use only polynomial-time computations necessarily do worse than the best centralized bound (or have much worse round complexity). At a high level, this is because previous approaches (most notably [2]) do not actually use the structure of the centralized algorithm: they only use the efficient centralized approximation as a black box. By going inside the black box and noticing that they all use a similar type of convex relaxation, we can simultaneously get low round complexity, best-known approximation ratios, and efficient local computation.
It turns out that we can use our techniques for an even broader question: Directed Steiner Network with Distance Constraints with a spanning demand graph. In this problem there is a set S ⊆ V × V of demands, and for every demand (u, v) ∈ S there is a length bound L(u, v). The goal is to find a subgraph H so that d H (u, v) ≤ L(u, v) for all (u, v) ∈ S, and the objective is to minimize the number of edges in H. The state of the art centralized bound for this problem is a O(n 3/5+ǫ )-approximation [7], but if we further assume that every vertex u ∈ V is the endpoint of at least one demand in S (which we will refer to as a spanning demand graph) then it is straightforward to see that the centralized algorithm of [5] for Directed k-Spanner can be generalized to give aÕ( √ n)-approximation. Our distributed version of this algorithm also generalizes, giving the following result. Note that Directed k-Spanner and Basic k-Spanner are special cases of this problem, where there is a demand for every edge and the length bound is just k times the original distance. Interestingly, other network design problems which have proved important for distributed systems are also special cases, including the Distance Preserver problem (when , and the Shallow-Light Steiner Network problem (where L(u, v) = D for all (u, v) ∈ S, for some global parameter D). Shallow-Light Steiner Network in particular is a key component in state of the art systems for reliable Internet transport [1], although in that particular application the demand graph is not spanning. Extending our techniques to handle totally general demands by giving a distributed version of [7] is an extremely interesting open question.

Related Work
While distributed solving of convex programs is a natural question, there is little previous work in the LOCAL model. Possibly most related to our results is a line of work on solving positive linear programs (packing and covering LPs). This was introduced by [22], improved by [4], and then essentially optimal upper and lower bounds were given by [20]. Unfortunately, the convex programs we consider are not positive linear programs due to the "capacity" constraints in which some variables appear with positive coefficients while others have negative coefficients.
A special case of our result was proved earlier in [10], who showed how to solve the LP relaxation of Basic 2-Spanner in the LOCAL model in O(log 2 n) rounds (they actually show more than this, by giving a distributed algorithm for the fault-tolerant version of Basic 2-Spanner, but that is not germane to our results). Our techniques are heavily based on [10], which is itself based on the ideas from [20]. In particular, [20] uses a Linial-Saks decomposition [21] to solve "local" versions of the linear program in different parts of the graph, and then combines these appropriately. To make this work for the Basic 2-Spanner LP relaxation, [10] had to use padded decompositions, which can be thought of as a variant of Linial-Saks with slightly different guarantees which, for technical reasons, are more useful for network design LPs. In this paper we extend these techniques further by giving a more general definition of padded decomposition which works for larger distance requirements, showing how to construct them in the LOCAL model, and then showing that the basic "combining" idea from [10] can be extended to handle these more general decompositions and far more general constraints and objective functions.
The major type of combinatorial optimization problem which our techniques allow us to approximate are various versions of graph spanners. There are an enormous number of papers on spanners in both centralized and distributed models, but fewer papers which attempt to find the "best" spanner for the particular given input graphs (most papers on spanners give existential results and algorithms to achieve them, rather than optimization results). These optimization questions (e.g., Basic k-Spanner, Directed k-Spanner, and Lowest-Degree k-Spanner) have been considered quite a bit in the context of centralized approximation algorithms and hardness of approximation [5,6,8,9,11], but almost all of the known centralized results use linear programming relaxations, making them difficult to adapt to distributed settings. Hence there have been only two results on optimization bounds in distributed models: [10] and [2].
Barenboim et al. [2] provided a distributed algorithm using Linial-Saks decompositions that for any integer parameters k, α, gives an O(n 1/α )-approximation for Directed k-Spanner in exp(O(α))+O(k) time. This is an extremely strong approximation bound, and in fact is better than even the best centralized bound. This is possible due to their use of very heavy (exponential time) local computation. Our algorithms, on the other hand, take polynomial time for local computations. Barenboim et al. [2] show that heavy local computations can be removed from their algorithm by using a centralized approximation algorithm for a variant of spanners known as client-server kspanners, and in particular that an f (k)-approximation for client-server k-spanner can be turned into an O(n 1/α f (k))-approximation algorithm running in exp(O(α)) + O(k) rounds in the LOCAL model for minimum k-spanner with only polynomial local computation. So in order to achieve the same asymptotic approximation ratio as the best-known centralized algorithm, the parameter α must be Ω(log n) and hence the running time is polynomial in n, even though k might be a constant. It is essentially known (though not written anywhere) that a variety of other results with slightly different tradeoffs can be achieved through similar uses of Linial-Saks [13] or refinements of Linial-Saks such as [14]. However, since all of these approaches treat the centralized approximation algorithm as a black box, none of them can achieve the same approximation ratio as the centralized algorithm without suffering a much worse (usually polynomial) round complexity that the O(k log n) that we achieve.

Preliminaries and Notation
The distributed setting we will be considering is the LOCAL model [23], in which time passes in synchronous rounds and in each round every node can send an arbitrary message of unbounded size to each of its neighbors in the underlying graph G = (V, E) (as always, we will let n = |V | and m = |E|). We will assume that all nodes know n (or at least know a constant approximation of n). Usually in this model the communication graph is the same as the graph of computational interest; e.g., we will be trying to compute a spanner of the communication graph itself. But for some applications we will want the graph to be directed, in which case we make the standard assumption that communiocation is bidirectional: the graph for which we are trying to compute a convex relaxation / network design problem is directed, but messages can be sent in both directions across a link. In other words, the communication graph is just the undirected version of the given directed graph.
For any pair of nodes u, v ∈ V we define d(u, v) to be the distance between u and v in the com-munication graph (i.e. the length of a shortest path between u and v regardless of edge directions). We define B(u, k) to be an undirected ball of radius k from u in the communication graph. More precisely, If x is a vector then we use x i to denote the i'th component of x. Most of the time our vectors will be indexed by edges in a graph, in which case we will also use the notation (x e ) e∈E .
Given a partition of the vertices V of a graph, we will refer to each part of the partition as a "cluster". For any graph G = (V, E) and set S ⊆ V , we let E(S) denote the set of edges in the subgraph induced by S, i.e., E(S) = {(u, v) ∈ E | u, v ∈ S}. We will frequently need "restrictions" of vectors to induced subgraphs, so for any vector x ∈ R m , we define

Padded decompositions
We will now define and give an algorithm to construct padded decompositions, which are one of the key technical tools that we will use when designing algorithm to solve distance-bounded network design convex programs. In this section all graphs are undirected and all distances are with respect to this undirected graphs (in fact, the definition and our algorithm work more generally for any metric space). Recall that B(u, k) denotes the undirected ball of radius k from node u (in the communication graph).

Definition 2.
Given an undirected graph G, a (k, ǫ)-padded decomposition, where 0 < ǫ ≤ 1, is a probability measure µ over the set of graph partitions (clusterings) that has the following properties: 1) For every P ∈ supp(µ), and every cluster C ∈ P , we have: 2) For every u ∈ V , it holds that Pr(∃C ∈ P | B(u, k) ⊆ C) ≥ 1 − ǫ. That is to say, the probability that all nodes in B(u, k) are in the same cluster is at least 1 − ǫ.
This notion of padded decompositions is standard in metric embeddings and approximation algorithms [17,19], but to the best of our knowledge has not yet been used in distributed algorithms. We first use a centralized algorithm (Algorithm 1) to sample from a (k, ǫ)-padded decomposition, and then describe how it can be implemented in the LOCAL model. Algorithm 1 and its analysis are similar to a partitioning algorithm proposed in [3], which was shown to have a low probability of separating nodes in a close neighborhood.
For any partition P constructed by Algorithm 1, each cluster is clearly C(v) for some v ∈ V . We call this special node v the center of cluster C(v). Later, we will use the center of each cluster for solving locally defined convex programs.
Lemma 4. Algorithm 1 partitions a given undirected graph G = (V, E) into a partition P such that P is sampled from a (k, ǫ)-padded decomposition.
Proof. The first property in Definition 2 is directly implied by the definition of r v for all nodes v ∈ V .
For the second property we consider an arbitrary node u ∈ V , and compute the probability that the ball B(u, k) is not in any of the clusters in P . Consider an arbitrary value 1 ≤ t ≤ n, let v ∈ V be the node such that t = π(v), and let z = z v be the real number sampled by v. Also, for any x, y ∈ V , letd(x, y) = min(d(x, y), r ln n + k). Let us also order the clusters based on their center's position in the permutation, so that C t is the cluster corresponding to t = π(v) (i.e. v is the cluster center of C t ). We define X t to be the event that if B(u, k) is not in the first t − 1 clusters, then it Set the radius r v = min(z v , r ln n + k).
is also not in any of the remaining clusters. We provide a recursive bound on X t based on X t+1 . Then we will get the second property once we show Pr(X 0 ) ≤ ǫ. We need to define the following events: • A t : B(u, k) does not intersect with any of the clusters C 1 , .., C t−1 .
In other words, conditional on the event that B(u, k) is not in any of the first t − 1 clusters, either B(u, k) ⊆ C t , or else one the following two events will occur: M cut t is the event that B(u, k) partially intersects C t , and M ex is the event that B(u, k) does not intersect C t . Now the event X t occurs only when either M cut t occurs or both M ex t and X t+1 occur (i.e. when B(u, k) is not in C t or any of the next clusters). Hence we can write Pr(X t ) ≤ Pr(M cut t ) + Pr(M ex t ) Pr(X t+1 ). Recall that z is independently sampled from the density function p(z v ) = n n−1 e −zv/r r , and thus M cut can be written as follows: Similarly, we can write, We now inductively prove that Pr(X t ) ≤ (2 − t n−1 )( 2k r ). If t < n is the last step, then Pr(X t ) = 0, and thus this bound clearly holds. Assume that the bound is true for X t+1 , we show that then it also holds for X t . We have, Since e −(d(v,u)−k)/r ≥ e −(ln n) ≥ 1/n, we get that Pr(X t ) ≤ 2 − t n−1 2k r . The second property is then implied by the fact that Pr(X 0 ) ≤ 2k r = 2k 2k(1/ǫ) = ǫ.
We will now use an idea similar to the one used in [10] to make Algorithm 1 distributed. In [10] they only considered the special case of k = 1 and ǫ = 1/2, which is why we cannot simply use their result as a black box.

Lemma 5.
There is an algorithm in the LOCAL model that runs in O( k ǫ ln n) rounds and samples from a (k, ǫ)-padded decomposition (so every node knows the cluster that it is in).
Proof. Without loss of generality, we assume that all nodes have unique IDs 2 . The sequence of IDs in ascending order will determine the permutation π used in Algorithm 1, i.e. if ID u < ID v then π(u) < π(v). The algorithm proceeds as follows until all nodes have been assigned to a cluster: each node u ∈ V chooses a radius r u based on the distribution defined in Algorithm 1. Then every u ∈ V simultaneously sends a message containing ID u to all nodes in B(u, r u ). After receiving all the messages, each node chooses the node with the smallest ID as the cluster center. Then Lemma 4 implies that the clusters satisfy the properties of a (k, ǫ)-padded decomposition. Since the radius that each node chooses is O((k/ǫ) log n), and each node only communicates with nodes within its radius, the running time in the LOCAL mode is O((k/ǫ) log n).

Distributed distance bounded network design convex programming
In this section we prove Theorem 1, giving an algorithm similar to [10] which can almost optimally solve distance-bounded network design convex programs. We first make all definitions formal in Section 4.1, and in particular define formally the class of objective functions where our results hold. Then in Section 4.2 we give a distributed algorithm which solved these programs up to arbitrarily small error.

Distance bounded network design convex programs
We will first describe a general class of objective functions that our algorithm applies to. For a graph G = (V, E) and a set S ⊆ V , we let E(S) denote the set of edges in the subgraph of G induced by S. Recall that for a vector x ∈ R m (where m = |E|), we define x S = (x S e ) e∈E ∈ R m to be the vector where x S e = 0 if e ∈ E(S) and x S e = x e if e ∈ E(S).
Definition 3. Given a graph G = (V, E), a function g : R m → R is convex partitionable with respect to G if g is a non-decreasing 3 and convex function with the following property: for all partitions σ = {σ 1 , ..., σ ℓ } of nodes in V , there exists a non-decreasing function h σ : R ℓ → R, s.t. g(x) = h σ (g(x σ 1 ), g(x σ 2 ), ..., g(x σ ℓ )) for all x = (x e ) e∈E where x e = 0 for any edge e with endpoints in different clusters of σ (equality does not need to hold for vectors x with nonzero values on edges between clusters).
Convex partitionable functions for graphs are a natural class of functions for distributed computing purposes. Moreover, this class includes many types of objective functions that are of interest in network design problems, including p-norms and linear functions. For example, if the function g is the p-norm with p ∈ Z ≥0 , then it is easy to verify that by setting the function h σ to also be the p-norm for any partition σ of V , the conditions of Definition 3 are satisfied. Note that an unweighted sum is the 1-norm, and the max function is the infinity norm, and hence they will also satisfy the conditions of Definition 3. Similarly, in case of linear functions, it is easy to see that conditions of Definition 3 are satisfied by setting h σ to be the unweighted sum.
There are also other, less trivial examples. For example, it is not hard to show the p-norm of the degree vector (rather than just the edge vector) is also convex partitionable with respect to G. An important special case of this is the ∞-norm of the degree vector, i.e., the maximum degree. Since we use this objective in some of our applications (i.e., for the Lowest-Degree k-Spanner problem), we give a short proof of this case in Lemma 6. For an integral vector x ∈ R m we can write g(x) = max v∈V deg(v). By generalizing this notation to all x ∈ R m , we can define fractional node degrees as deg(v) = u:(v,u)∈E x (v,u) 4 . u) ). Then we can set h σ (y) = max i∈[ℓ] (y i ), y ∈ R ℓ , where y i is the i-th coordinate of y. Let σ(v) ∈ σ be the cluster that node v belongs to. For all

Lemma 6. Given a graph
It is also easy to see that the function h σ is convex and non-decreasing. Hence h σ satisfies the conditions in Definition 3. Now that this class of functions has been defined, we can formally define the class of distancebounded network design convex programs. Definition 4. Let S ⊆ V × V be a set of pairs in the graph G = (V, E), and for any pair (u, v) ∈ S let P u,v be a set of paths from u to v, which we sometimes call the set of "allowed" paths. Let g be a non-decreasing convex-partitionable function of x = (x e ) e∈E with g( 0) = 0. Then we call a convex program of the following form a distance bounded network design CP : As we will see in Section 5, many network design problems use linear (or convex) programming relaxations that satisfy the conditions of Definition 4. A key parameter of such a program is the length of the longest allowed path D = max (u,v)∈S max p∈Pu,v ℓ(p) (where ℓ(p) is the length of path p).

Distributed Algorithm
In order to solve these convex programs in a distributed manner, we will first use padded decompositions to form a local problem using a simple distributed algorithm. Let P be a partition sampled from a (k, ǫ)-padded decomposition (in particular, obtained by Lemma 5), where 0 < ǫ ≤ 1. Recall that for each cluster C ∈ P , E(C) = {(u, v) ∈ E | u, v ∈ C}. We define G(C) to be the subgraph induced by C.

Lemma 7.
For each cluster C sampled from a (k, ǫ)-padded decomposition, there is a distributed algorithm running in O( k ǫ log n) rounds so that every cluster center knows G(C). Proof. The first property of (k, ǫ)-padded decompositions implies that for all nodes u ∈ C, we have where v is the center of cluster C. Each node u ∈ C that determines v as the center of the cluster it belongs to, will send the information of its incident edges to v. Since there is no bound on the size of the messages being forwarded, this can be done in O((k/ǫ) log n) time.
Let CP(G) be a distance bounded network design CP defined on graph G = (V, E). We will define local convex programs based on a partition P of G that is sampled from a (D, λ)-padded decomposition. The value of 0 < λ ≤ 1 will be set later based on the parameters of our distributed algorithm. For each C ∈ P , let CP(C) be CP(G) defined on G(C), but where only demands corresponding to any pair (u, v) ∈ S in which B(u, D) is fully contained in C are included. We denote the set of these demands by N (C), more precisely, The objective will then be to minimize g(x) = g (x e ) e∈E(C) . In other words CP (C) is defined as follows: There is a technical subtlety about computing the function g on each cluster, which is the fact that a solution x C , f C of CP(C) is only defined on G(C). While in practice x C is a vector defined only on edges in E(C), in our analysis we will assume that x C is a vector in R m (where |E| = m) and x C e = 0 for all e ∈ E(C) (this is possible for all the functions we care about). It will become clear that this distinction is also needed when we try to compare the local solutions to the global solutions. The following lemma is similar to Lemma 3.8 in [10], and we show that it holds for our modified definition of local convex programs and for generalized objective functions that satisfy Definition 3.
Lemma 8. Let x * , f * be an optimal solution of CP (G) and let x * C = (x * e ) e∈E(C) . For each cluster C ∈ P , let x C ,f C be an optimal solution of CP(C). Then g(x C ) ≤ g(x * C ).
Proof. We argue that the vector x * C , f * C , where x * e C = x * e for all e ∈ E(C) and f * p C = f * p for all p ∈ P u,v , is a feasible solution to CP(C). By definition of N (C) we have that for any (u, v) ∈ N (C) all paths in P u,v also appear in G(C), and therefore x * C , f * C satisfies both capacity and flow constraints of CP(C) for pairs (u, v) ∈ E(C) since they were satisfied in CP(G). Since we assumed that x C ,f C is an optimal solution of CP(C), this implies that g(x C ) ≤ g(x * C ).
We now provide in Algorithm 2 a distributed algorithm for solving CP(G), by having cluster centers solve CP(C) of their cluster using a sequential algorithm in each iteration, and then averaging over the solutions for each edge. We assume that all nodes know the values of D and ǫ. Let C u,i denote the cluster that node u belongs to in the i-th iteration, and let is the fractional CP value for e = (u, v), and p ∈ P u,v .
Since the objective is a function of edge vectors, what we mean by having a distributed solution to a distance bounded network design CP is that each node u will know the value x e for all the edges e incident to u. It is not hard to see that the algorithm could be modified so that every node u can also know the flow value f p for each path p.
Algorithm 2: Distributed algorithm for approximating distance bounded network design CPs.
2 Sample from (D, λ)-padded decompositions t times by Lemma 5, and let P i be the partition obtained in the i'th run. 3 For each cluster C ∈ P i , the center of cluster C computes G(C) (see Lemma 7). 4 The center of each cluster C ∈ P i solves CP (C) and sends the solution x C,i , f C,i to all nodes u ∈ C. 5 for e = (u, v) ∈ E do 6 Let I u,v = {i | ∃C ∈ P i : u, v ∈ C}. // these are the iterations in which both endpoints are in same cluster Theorem 9. Algorithm 2 takes O((D/ǫ) log n) rounds to terminate, and it will compute a solution of cost (1 + ǫ)CP * to a bounded distance network design CP (Definition 4) with high probability, where CP * is the optimal solution and 0 < ǫ ≤ 1. Moreover, if the convex program can be solved sequentially to arbitrary precision in polynomial time, then all of the node computations are also polynomial time.
Proof. Correctness: We first show that with high probability the valuesx e , e ∈ E form a feasible solution. Here we only need to show that a feasible solution for the flow values exist, and do not require nodes to compute these values. Let I u = {i : ∃C ∈ P i , B(u, D) ⊆ C}, i.e. I u is the set of iterations in which B(u, D) is contained in a cluster, and let I u,v be the set of iterations in which both u and v are in the same cluster. Since we need to implement Algorithm 2 in a distributed manner, we use I u,v in our implementation, while the analysis is based on I u . We can do so since by definition we have I u ⊆ I u,v , for any (u, v) ∈ E.
For any p ∈ P u,v , we set the flow values to bef p = 1 |Iu| i∈Iu f C u,i ,i p . We will show that this gives a feasible flow. First we argue that enough flow is being sent. For all (u, v) ∈ S, we have, We have used the fact that for each i ∈ I u the solution corresponding to the CP of the cluster containing u satisfies the constraint that p∈Pu,v:e∈p f C u,i ,i p ≥ 1, because for each such i we know that (u, v) ∈ N (C).
Next, we will argue that the capacity constraints are also satisfied. The second property of (D, λ)-padded decompositions implies that Pr(i ∈ I u ) By linearity of expectations we have E[|I u |] ≥ t(1 − λ). Since each sampling is performed independently, by Chernoff bound for δ = ǫ/2, we get, Hence by a union bound on all nodes we have that with high probability |I u | > t/(1 + ǫ). Therefore, for all (u, v) ∈ S, e ∈ E, we have (w.h.p.), Upper bound: We will now show that the upper bound holds. Let x * , f * be an optimal solution to CP(G). We havex e = min(1, 1+ǫ t i∈Ie x C u,i ,i e ), and for each e = (u, v) and 1 ≤ i ≤ t, if i ∈ I e , andx i e = 0 otherwise. Note that 0 < (1 + ǫ)/t < 1, and since g is a convex function and g( 0) = 0, by Jensen's inequality we have g 1+ǫ t x ≤ 1+ǫ t g(x). Then for x = (x e ) e∈E we can write: In the final inequality, since g is convex, we used Jensen's inequality to take the sum out of the function. It is now enough to show that in each iteration i, it holds g((x i e ) e∈E ) ≤ g(x * ). Let x i = ((x i e ) e∈E ), and let P i = {C 1 , C 2 , ..., C ℓ } be the partition of V . Since g is a convex partitionable function w.r.t. G, there exists a nondecreasing and convex function h : R m → R for which we can write g(x i ) = h P i (g(x i,C 1 ), g(x i,C 2 ), ..., g(x i,C ℓ )), sincex i e = 0 by definition for edges which go between clusters (for simplicity we are denotingx i C j byx i,C j ).
Recall that x * C is the vector in which x * C e = x * e for all e ∈ E(C) and x * C e = 0 otherwise. By Lemma 8 we get that for all C ∈ P i , g(x i,C ) ≤ g(x * C ). Now we consider a vectorx, defined by settingx e = x * e for all edge e with both endpoints in the same cluster, andx e = 0 otherwise. Since we assumed h P i to be nondecreasing, we get, For the last inequality we have used the fact that g is non-decreasing, and that for all e ∈ E, x e ≤ x * e (since eitherx e = x * e orx e = 0). By plugging this into the above inequalities, we will get , which implies the claim that Algorithm 2 gives a (1 + ǫ)approximation to the optimal solution.
Time Complexity: The decomposition step and sending the information within a cluster takes O((D/ǫ) log n) rounds since the diameter of each cluster is O((D/λ) log n) = O((D/ǫ) log n). Since each decomposition is independent, we can do all of them in parallel, so steps 1-4 of the algorithm only take O((D/ǫ) log n) rounds in total. Clearly the rest of the algorithm can be done in a constant number of rounds. Hence in total w.h.p. the algorithm will take O((D/ǫ) log n) rounds.

Distributed Approximation Algorithms for Network Design
In this section, we will focus on several network design problems which can be approximated by first solving a convex relaxation using Algorithm 2 and then locally rounding the solution. For that purpose, we will describe how each problem has a distance bounded network design CP relaxation (Definition 4), and will then show that existing rounding schemes are local.

Directed k-Spanner
Dinitz and Krauthgamer [9] introduced a linear programming relaxation for Directed k-Spanner which is just a distance-bounded network design CP with demands pairs S = E, allowed paths P u,v which are the directed paths from u to v of length at most k, and objective function g(x) = e∈E x e . They showed that this LP can be solved in polynomial time (approximately if k is non-constant). We will denote this LP by LP (G). Clearly, LP(G) is a distance bounded network design CP with D = k. Hence, Theorem 9 implies that we can use Algorithm 2 to approximately solve this LP in O(k log n) rounds in the LOCAL model.
We now provide in Algorithm 3 a distributed rounding scheme that gives an O(n 1/2 log n)approximation for Directed k-Spanner. This algorithm matches the best centralized approximation ratio known [5], and is just the obvious distributed version of the algorithm proposed in [5]. The difference is that here we truncate the shortest-path trees at depth k (as opposed to full shortest-path trees), and nodes choose whether to become a tree root independently (rather than chosen without replacement as in [5].
Add e to E ′ with probability min(n 1/2 · ln n · x e , 1). 4 for v ∈ V do // Random tree sampling 5 Choose p uniformly at random from [0, 1]. 6 if p < 3 ln n √ n then 7 T in v ← shortest path in-arborescence rooted at v truncated at depth k. 8 T out v ← shortest path out-arborescence rooted at v truncated at depth k. 9 Output E ′ ∪ (∪ v∈V (T in v ∪ T out v )). // Each node knows its portion of the output.
The following lemma is essentially from [5], with the proof requiring only slight technical changes due to the slightly different algorithms. We sketch it for completeness. Lemma 10. Given a directed graph G, LP(G) as defined, and a fractional solution LP * to LP (G), the output of Algorithm 3 has size O(n 1/2 · (n + LP * ) log n).
Proof. Let N s,t be the subgraph of G induced by the nodes on paths in P s,t . Edge e ∈ E is called a thick edge if |N s,t | ≥ n 1/2 , and otherwise it is called a thin edge. The set E ′ in Algorithm 3 satisfies the spanner property for all thin edges (as argued in [5]), and the random tree sampling phase satisfies the spanner property for the thick edges. Each thick edge (s, t) is spanned if at least one node in N s,t performs the random tree sampling. This probability is at least 1 − (1 − 3 ln n n 1/2 ) n 1/2 ≥ 1 − 1/n 3 . Then a union bound on all the edges (of size at most O(n 2 )) implies that w.h.p. all thick edges are spanned. We now argue that the output is an O(n 1/2 log n)-approximation algorithm: at most O(n 1/2 log n) arborescences are chosen with high probability (each arborscence has O(n) edges), and we argued that |E ′ | = O(n 1/2 log n · LP * ). Hence, the overall size of the output is O(n 1/2 log n · (n + LP * )).
It is also easy to see that this algorithm can be implemented in the LOCAL model. Proof. Each node v in G has received the fractional solutions x e corresponding to all edges e ∈ E incident to v. The randomized rounding step can be performed locally: the node with the smaller ID flips a coin, and exchanges the coin flip result with its corresponding neighbors. In order to form T in i and T out i , v performs a distributed BFS algorithms by forming a shortest path tree while keeping track of the distance from v. When the distance counter reaches k, the tree construction terminates.
We now immediately get our main result for Directed k-Spanner. Proof. We first run Algorithm 2 to solve LP(G) up to a constant factor (by setting ǫ = 1/2), which takes time O(k log n) with high probability (Theorem 9). Since each cluster center can solve the local LP in polynomial time, all computations are polynomial time. We then use Algorithm 3 to round the fractional solutions of LP(G), which takes O(k) time. Since the size of a k-spanner is at least Ω(n), Algorithm 3 then outputs an O(n 1/2 ln n)-approximation to the minimum (Lemma 10).

Basic 3-Spanner and Basic 4-Spanner
If the input graph is undirected then stronger approximations are possible. In particular, for stretch 3 and 4, there areÕ(n 1/3 )-approximations due to [5] (for stretch 3) and [11] (for stretch 4). Without going into details, both of these algorithms use the same LP relaxation as in Directed k-Spanner, but round the LP differently. So in order to give distributed versions of these algorithms, we only need to modify Algorithm 3 to use the appropriate rounding algorithm (and change some of the other parameters in the shortest-path arborescence sampling). Fortunately, both of these algorithms use rounding schemes which are highly local. Informally, rather than sample each edge independently with probability proportional to the (inflated) fractional value as in Algorithm 3, these algorithms sample a value independently at each vertex and then include an edge if a particular function of the values of the two endpoints (different in each of the algorithms) passes some threshold. Clearly this is a very local rounding algorithm: once we have solved the LP relaxation using Theorem 9, each node can draw its random value and then spend one more round to exchange a message with each of its neighbors to find out their values, and thus determine which of the edges have been included by the rounding. Thus the total running time is dominated by the time needed to solve the LP, which in these cases is O(log n) using Theorem 9.

Lowest-Degree k-Spanner
We now turn our attention to Lowest-Degree k-Spanner: Given a graph G = (V, E) and a value k, we want to find a k-spanner that minimizes the maximum degree. We will use the relaxation and rounding scheme proposed by Chlamtáč and Dinitz [6]. The linear programming relaxation used in [6] is very similar to the Directed and Basic k-spanner LP relaxation described earlier, with the difference being that a new variable λ is added to represent the maximum degree, and so the objective is to minimize λ and constraints are added to force λ to upper bound the maximum fractional degree. Proof. It is easy to see that the LP relaxation proposed in [6] can be written as a distance bounded network design CP where the objective is max v∈V deg(v) = max v∈V u:{v,u}∈E x {v,u} (we do not need to use their extra variable λ, since we can instead directly write the objective). Lemma 6 implies that this function is convex partitionable w.r.t. G, and hence the Lowest-Degree k-Spanner problem can be approximately solved (to within a constant factor) by using Algorithm 2 with ǫ = 1/2. Next, we use the following rounding scheme proposed in [6]: each edge e ∈ E is included in the spanner with probability x 1/k e . It is clear that this can be done in a constant number of rounds, and hence the overall algorithm takes O(k log n) rounds (by Theorem 9) in the LOCAL model. In [6], it was shown that this leads to aÕ(∆ (1−1/k) 2 )-approximation solution of the problem.

Directed Steiner Network with Distance Constraints
It is well-known that the centralized rounding of [5] for Directed k-Spanner is more general than is actually stated in their paper. In particular, the randomized rounding for "thin" edges gives the same guarantee even when each demand has a possibly different distance constraint. This fact was used, e.g., in [7] in their algorithms for Distance Preserver, Pairwise k-Spanner, and Directed Steiner Network with Distance Constraints. The difficulty in extending the algorithm of [5] is not in the LP rounding, but rather because the arborescence sampling technique used to handle thick edges in [5] (and in our Algorithm 3) assumes that n is a lower bound on the optimal cost. This assumption is true for Directed k-Spanner, but false for variants where there might be a tiny number of demands. However, it is easy to see that if we assume the demand graph is spanning (i.e., assume that every node is an endpoint of at least one demand) then the optimal solution must have at least n/2 edges, and hence we can again just use [5] to get aÕ( √ n)approximation for Directed Steiner Network with Distance Constraints as long as the demand graph is spanning. While this is in the centralized setting, since our algorithm for Directed k-Spanner is just a lightly modified distributed version of [5] (the only difficulty in the distributed setting is solving the LP, which is why that is the main technical contribution of this paper), we can easily modify it to give the same approximation for Directed Steiner Network with Distance Constraints with spanning demand graphs. The only change is that we use D = max (u,v)∈S L(u, v) instead of k when solving the linear programming relaxation (using Theorem 9) and when truncating the shortest-path arborescences that we sample (note that we have to assume that D is global knowledge, which is reasonable for spanner problems and for Shallow-Light Steiner Network but may be less reasonable for other special cases of Directed Steiner Network with Distance Constraints). This implies Theorem 3, and all of the interesting special cases (Shallow-Light Steiner Network, Distance Preserver, Pairwise k-Spanner, etc.) which it includes.

Conclusion
In this paper we presented a distributed algorithm for solving distance-bounded network design convex programs in the LOCAL model of distributed computation. This is one of the few classes of convex programs (along with positive linear programs and a few other special cases) for which we now have distributed algorithms. This class is particularly interesting since many state-of-the-art approximation algorithms for distance-bounded network design problems work by rounding one of these convex relaxations. So if we can solve the relaxation in a distributed fashion then it is often straightforward to give a distributed approximation algorithm which simply solves the relaxation using our new distributed algorithm and then does the appropriate rounding. Using this framework, we provide distributed approximation algorithms for a variety of problems (e.g., for Directed k-Spanner, Lowest-Degree k-Spanner, Basic 3-Spanner, and Basic 4-Spanner) which use only polynomial-time computations at each node, achieve approximation ratios asymptotically equal to the centralized algorithm, and have very low round complexity. Previous approaches to these problems either use exponential-time computations at local nodes, or if they use only polynomial-time computations require either significantly larger round complexity or give asymptotically worse approximations.


extfieldsofstudy:
Mathematics
Computer Science

text:
Markov Decision Processes with Continuous Side Information

We consider a reinforcement learning (RL) setting in which the agent interacts with a sequence of episodic MDPs. At the start of each episode the agent has access to some side-information or context that determines the dynamics of the MDP for that episode. Our setting is motivated by applications in healthcare where baseline measurements of a patient at the start of a treatment episode form the context that may provide information about how the patient might respond to treatment decisions. We propose algorithms for learning in such Contextual Markov Decision Processes (CMDPs) under an assumption that the unobserved MDP parameters vary smoothly with the observed context. We also give lower and upper PAC bounds under the smoothness assumption. Because our lower bound has an exponential dependence on the dimension, we consider a tractable linear setting where the context is used to create linear combinations of a finite set of MDPs. For the linear setting, we give a PAC learning algorithm based on KWIK learning techniques.

Introduction
Consider a basic sequential decision making problem in healthcare, namely that of learning a treatment policy for patients to optimize some health outcome of interest. One could model the interaction with every patient as a Markov Decision Process (MDP). In precision or personalized medicine, we want the treatment to be personalized to every patient. At the same time, the amount of data available on any given patient may not be enough to personalize well. This means that modeling each patient via a different MDP will result in severely suboptimal treatment policies. The other extreme of pooling all patients' data results in more data but most of it will perhaps not be relevant to the patient we currently want to treat. We therefore face a trade-off between having a large amount of shared data to learn a single policy, and, finding the most relevant policy for each patient. A similar trade-off occurs in other applications involving humans as the agents environment such as online tutoring and web advertising.
A key observation is that in many personalized decision making scenarios, we have some side information available about individuals which might help us in designing personalized policies and also help us pool the interaction data across the right subsets of individuals. Examples of such data include laboratory data or medical history of patients in healthcare, user profiles or history logs in web advertising, and student profiles or historical scores in online tutoring. Access to such side information should let us learn better policies even with a limited amount of interaction with individual users. We refer to this side-information as contexts and adopt an augmented model called Contextual Markov Decision Process (CMDP) proposed by Hallak et al. (2015). We assume that contexts are fully observed and available before the interaction starts for each new MDP. 1 In this paper we study the sample complexity of learning in CMDPs in the worst case. We consider two concrete settings of learning in a CMDP with continuous contexts. In the first setting, the individual MDPs vary in an arbitrary but smooth manner with the contexts, and we propose our Cover-Rmax algorithm in Section 3 with PAC bounds. The innate hardness of learning in this general case is captured by our lower bound construction in Section 3.1. To show that it is possible to achieve significantly better sample complexity in more structured CMDPs, we consider another setting where contexts are used to create linear combinations of a finite set of fixed but unknown MDPs. We use the KWIK framework to devise the KWIK LR-Rmax algorithm in Section 4.1 and also provide a PAC upper bound for the algorithm.

Problem setup and notation
We start with basic definitions and notations for MDPs, and then introduce the contextual case.
Definition 2.1 (Markov Decision Processes). A Markov Decision Process (MDP) is defined as a tuple (S, A, p(·|·, ·), r(·, ·), µ) where S is the state space and A is the action space; p(s |s, a) defines the transition probability function for a tuple (s, a, s ) where s, s ∈ S, a ∈ A; and µ defines the initial state distribution for the MDP.
We consider the case of fixed horizon (denoted H) episodic MDPs. We use π h (s) to denote a policy's action for state s at timestep h. For each episode, an initial state is observed according to the distribution µ and afterwards, for 0 ≤ h < H, the agent chooses an action a h = π h (s h ) according to a (non-stationary) policy. There is a reward r h and then a next state s h+1 according to the reward and the transition functions. For policy π we define its value as follows: An optimal policy π * is one that achieves the largest possible value (called optimal value and denoted V * M ). Next we define the contextual model which is similar to the definition given by Hallak et al. (2015): Definition 2.2 (Contextual MDP). A contextual Markov Decision Process (CMDP) is defined as a tuple (C, S, A, M) where C is the context space, S is the state space, and A is the action space. M is a function which maps a context c ∈ C to MDP parameters M(c) = {p c (·|·, ·), r c (·, ·), µ c (·)}.
We denote the MDP for a context c by M c . We make the simplifying assumption that the initial state distribution is independent of the context and is same for all. We assume throughout the paper that the rewards are bounded between 0 and 1. We denote |S|, |A| by S, A respectively. We also assume that the context space is bounded, and for any c ∈ C the 1. Observe context c t ∈ C.
2. Choose a policy π t (based on c t and previous episodes).
3. Experience an episode in M ct using π t .
We do not make any distributional assumptions over the context sequence. Instead, we allow the sequence to be chosen in an arbitrary and potentially adversarial manner. A natural criteria for judging the efficiency of the algorithm is to look at the number of episodes where it performs sub-optimally. The main aim of the PAC analysis is to bound the number of episodes where we have V πt M c t < V * M c t − , i.e., the value of the algorithm's policy is not -optimal (Dann and Brunskill, 2015). Although, we do give PAC bounds for the Cover-Rmax algorithm given below, the reader should make note that, we have not made explicit attempts to achieve the tightest possible result. We use the Rmax (Brafman and Tennen holtz, 2002) algorithm as the base of our construction to handle exploration-exploitation because of its simplicity. Our approach can also be combined with the other PAC algorithms (Strehl and Littman, 2008;Dann and Brunskill, 2015) for improved dependence on S, A and H.

Cover-Rmax
In this section, we present the Cover-Rmax algorithm and provide a PAC bound for it under a smoothness assumption. The key motivation for our contextual setting is that sharing information among different contexts is helpful. Therefore, it is natural to assume that the MDPs corresponding to similar contexts will themselves be similar. This can be formalized by the following smoothness assumption: Definition 3.1 (Smoothness). Given a CMDP (C, S, A, M), and a distance metric over the context space φ(·, ·), if for any two contexts c 1 , c 2 ∈ C, we have the following constraints: |r c 1 (s, a) − r c 2 (s, a)| ≤ L r φ(c 1 , c 2 ) then, we call this a smooth CMDP with smoothness parameters L p and L r .
We assume that the distance metric and the constants L p and L r are known. This smoothness assumption allows us to use a minimally tweaked version of Rmax (Brafman and Tennenholtz, 2002) and provide an analysis for smooth CMDPs similar to existing literature on MDPs(Kearns and Singh, 2002;Strehl et al., 2009;Strehl and Littman, 2008). If we know the transition dynamics and the expected reward functions for each state-action pair in a finite MDP, we can easily compute the optimal policy. The idea of Rmax is to distinguish the state-action pairs as known or unknown: a state-action pair is known if we have visited it for enough number of times, so that the empirical estimates of reward and transition probabilities are near-accurate due to sufficient data. A state s becomes known when all for all actions a the pairs (s, a) become known. Rmax then constructs an auxiliary MDP which encourages optimistic behaviour by assigning maximum reward (hence the name Rmax) to the remaining unknown states. When we act according to the optimal policy in the auxiliary MDP, one of the following must happen: 1) we exploit the information available and achieve near-optimal value, or, 2) we visit unknown states and accumulate more information efficiently.
Formally, for a set of known states K, we define an (approximate) induced MDPM K in the following manner. Let n(s, a) and n(s, a, s ) denote the number of observations of state-action pair (s, a) and transitions (s, a, s ) respectively. Also, let R(s, a) denote the total reward obtained from state-action pair (s, a). For each s ∈ K, define the values (2) For each s / ∈ K, define the values as pM K (s |s, a) = I{s = s} and rM K (s, a) = 1. We use the certainty equivalent policy computed for this induced MDP and perform balanced wandering (Kearns and Singh, 2002) for unknown states. Balanced wandering ensures that all actions are tried equally and fairly for unknown states. Assigning maximum reward to the unknown states pushes the agent to visit these states and provides the necessary exploration impetus. The generic template of Rmax is given in Algorithm 1.
For the contextual case, we have an infinite number of such MDPs. The idea behind our algorithm is that, one can group close enough contexts and treat them as a single MDP. Utilizing the boundedness of the context space C, we create a cover of C with finitely many balls B r (o i ) of radius r centered at o i ∈ R d . By tuning the radius r, we can control the bias introduced by ignoring the differences among the MDPs in the same ball. Doing so allows us to pool together the data from all MDPs in a ball, so that we avoid the difficulty of infinite MDPs and instead only deal with finitely many of them. The size of the cover, i.e., the number of balls can be measured by the notion of covering numbers (see e.g., Wu (2016)), defined as N (C, r) = min{|Y| : C ⊆ ∪ y∈Y B r (y)}.
The resulting algorithm, Cover-Rmax, is obtained by using the subroutines in Algorithm 2, and we state its sample complexity guarantee in Theorem 3.2.
Algorithm 1: Rmax Template for CMDP 1 Initialize(S, A, C, , δ); 2 for each episode t = 1, 2, · · · do 3 Receive context c t ∈ C; 4 Set K, M K using P redict(c t , s, a) for all (s, a). π ← π * M K ; Choose a h : (s h , a h ) is unknown; Theorem 3.2 (PAC bound for Cover-Rmax). For any input values 0 < , δ ≤ 1 and a CMDP with smoothness parameters L p and L r , with probability at least 1 − δ, the Cover-Rmax algorithm produces a sequence of policies {π t } which yield at most non--optimal episodes, where N = N (C, r 0 ) and r 0 = min( 8HLp , 8Lr ).
Proof sketch. We first of all carefully adapt the analysis of Rmax by Kakade (2003) to get the PAC bound for an episodic MDP. Let m be the number of visits to a state-action pair after which the model's estimatep(·|s, a) for p(·|s, a) has an 1 error of at most /4H and reward estimater(s, a) has an absolute error of at most /4. We can show that: Lemma 3.3. Let M be an MDP with the fixed horizon H. Ifπ is the optimal policy for M K as computed by Rmax, then for any starting state s 0 , with probability at least 1 − 2δ, we have Vπ M ≥ V * M − 2 for all but O( mSA ln 1 δ ) episodes. Now instead of learning the model for each contextual MDP separately, we combine the data within each ball. Therefore, we have to take care of two things: choose the radius r 0 for a fine enough cover and a value of m which is the number of visits after which a state becomes known for a ball. For satisfying the conditions of Lemma 3.3 for all MDPs within a ball, we need the radius r of the cover to be r ≤ r 0 = min 8HLp , 8Lr and the value of . Using Lemma 3.3, we obtain an upper bound on number of nonepisodes in a single ball as O H 2 SA 3 S + ln SA δ ln 1 δ with probability at least 1 − δ. Setting the individual failure probability to be δ/N (C, r 0 ) and using the union bound, we get the stated PAC bound. 2 We observe that the PAC bound has linear dependence on the covering number of the context space. In case of a d-dimensional Euclidean metric space, the covering number would be of the order O( 1 r d ). However, we show in Section 3.1, that, the dependence would be at least linear, and hence, indicate the difficulty of optimally learning in such cases.

Lower Bound
We prove a lower bound on the number of sub-optimal episodes for any learning algorithm in a smooth CMDP which shows that a linear dependence on the covering number of the context space is unavoidable. As far as we know, there is no existing way of constructing PAC lower bounds for continuous state spaces with smoothness, so we cannot simply augment the state representation to include context information. Instead, we prove our own lower bound in Theorem 3.4 which builds upon the work of Dann and Brunskill (2015) on lower bounds for episodic finite MDPs and of Slivkins (2014) on lower bounds for contextual bandits.
Theorem 3.4 (Lower bound for smooth CMDP). There exists constants δ 0 , 0 , such that for every δ ∈ (0, δ 0 ) and ∈ (0, 0 ), any algorithm that satisfies a PAC guarantee for ( , δ) and computes a sequence of deterministic policies for each context, there is a hard CMDP (C, S, A, M) with smoothness constant L p = 1, such that where B is the number of sub-optimal episodes and 1 = 1280H e 4 (H−2) .
Proof. The overall idea is to embed multiple MDP learning problems in a CMDP, such that the agent has to learn the optimal policy in each MDP separately and cannot generalize across them. We show that the maximum number of problems that can be embedded scales  (Dann and Brunskill, 2015). The initial state 0 moves to a uniform distribution over states 1 to n regardless of the action, and states +/− are absorbing with 1 and 0 rewards respectively. States 0 to n have 0 reward for all actions. Each state i ∈ [n] essentially acts as a hard bandit instance, whose A actions move to + and − randomly. Action a 0 satisfies p(+|i, a 0 ) = 1 2 + 2 and there is at most one other action a i with p(+|i, a i ) = 1 2 + . Any other action a j satisfies p(+|i, a j ) = 1 2 .
with the covering number, and the result follows by incorporating known PAC lower bound for episodic MDPs. We start with the lower bound for learning in episodic MDPs. See Figure 1 and its caption for details. The construction is due to Dann and Brunskill (2015) and we adapt their lower bound statement to our setting in Theorem 3.5.
Theorem 3.5 ( Lower bound for episodic MDPs (Dann and Brunskill, 2015) ). There exists constants δ 0 , 0 , such that for every δ ∈ (0, δ 0 ) and ∈ (0, 0 ), any algorithm that satisfies a PAC guarantee for ( , δ) and computes a sequence of deterministic policies, there is a hard instance M hard so that E[B] = Ω SA 2 , where B is the number of sub-optimal episodes. The constants can be chosen as Now we discuss how to populate the context space with these hard MDPs. Note in Figure 1 that, the agent does not know which action is the most rewarding (a i ), and the adversary can choose i to be any element of [A] (which is essentially choosing an instance from a family of MDPs). In our scenario, we would like to allow the adversary to choose the MDP independently for each individual packing point to yield a lower bound linear in the packing number. However, this is not always possible due to the smoothness assumption, as committing to an MDP at one point may restrict the adversary's choices in another point.
To deal with this difficulty, we note that any pair of hard MDPs differ from each other by O( ) in transition distributions. Therefore, we construct a packing of C with radius r = 8 , defined as a set of points Z such that any two points in Z are at least r away from 3. The lower bound here differs from that in the original paper by H 2 , because our value is normalized (see Eq. (1)), whereas they allow the magnitude of value to grow with H. each other. The maximum size of such Z is known as the packing number : which is related to the covering number as N (C, r) ≤ D(C, r). The radius r is chosen to be O( ) so that arbitrary choices of hard MDP instances at different packing points always satisfy the smoothness assumption (recall that L p = 1). Once we fix the MDPs for all c ∈ Z, we specify the MDP for c ∈ C \ Z as follows: for state i and action a, Essentially, as we move away from a packing point, the transition to +/− become more uniform. We can show that: Claim 3.6. The CMDP defined above is satisfies Definition 3.1 with constant L p = 1. 4 We choose the context sequence given as input to be repetitions of an arbitrary permutation of Z. By construction, the learning at different points in Z are independent, so the lower bound is simply the lower bound for learning a single MDP (Theorem 3.5) multiplied by the cardinality of Z (the packing number). Using the well known relation that N (C, r) ≤ D(C, r), we have the desired lower bound. We refer the reader to the appendix for proof of Claim 3.6 and a more detailed analysis.

Contextual Linear Combination of MDPs
From the previous section, it is clear that for a contextual MDP with just smoothness assumptions, exponential dependence on context dimension is unavoidable. Further, the computational requirements of our Cover-Rmax algorithm scales with the covering number of the context space. As such, in this section, we focus on a more structured assumption about the mapping from context space to MDPs and show that we can achieve substantially improved sample efficiency.
The specific assumption we make in this section is that the model parameters of an individual MDP M c is the linear combination of the parameters of d base MDPs, i.e., . . .
We use P (s, a, s ) and R(s, a) as shorthand for the d × 1 vectors that concatenate the parameters from different base MDPs for the same s, a (and s ). The parameters of the base MDPs (p i and r i ) are unknown and need to be recovered from data by the learning agent, and the combination coefficients are directly available which is the context vector c itself. This assumption can be motivated in an application scenario where the user/patient responds according to her characteristic distribution over d possible behavioural patterns. A mathematical difficulty here is that for an arbitrary context vector c ∈ R d , p c (·|·, ·) is not always a valid transition function and may violate non-negativity and normalization constraints. Therefore, we require that c ∈ ∆ d−1 , that is, c stays in the probability simplex so that p c (·|·, ·) is always valid. 5

KWIK LR-Rmax
We first explain how to estimate the model parameters in this linear setting, and then discuss how to perform exploration properly.
Model estimation Recall that in Section 3, we treat the MDPs whose contexts fall in a small ball as a single MDP, and estimate its parameters using data from the local context ball. In this section, however, we have a global structure due to our parametric assumption (d base MDPs that are shared across all contexts). This implies that data obtained at a context may be useful for learning the MDP parameters at another context that is far away, and to avoid the exponential dependence on d we need to leverage this structure and generalize globally across the entire context space.
Due to the linear combination setup, we use linear regression to replace the estimation procedure in Equation 2: in an episode with context c, when we observe the state-action pair (s, a), a next-state s next will be drawn from p c (·|s, a). 6 Therefore, the indicator of whether s next is equal to s forms an unbiased estimate of p c (s |s, a), i.e., whenever we observe a transition tuple (s, a, s next ) under context c, and their relationship is governed by a linear prediction rule with P (s, a, s ) being the coefficients. Hence, to estimate P (s, a, s ) from data, we can simply collect the feature-label pairs that correspond to this particular (s, a, s ) tuple, and run linear regression to recover the coefficients. The case for reward function is similar, hence, not discussed. If the data is abundant (i.e., (s, a) is observed many times) and exploratory (i.e., the design matrix that consists of the c vectors for (s, a) is well-conditioned), we can expect to recover P (s, a, s ) accurately. But how to guarantee these conditions? Since the context is chosen adversarially, the design matrix can indeed be ill-conditioned.
Observe, however, when the matrix is ill-conditioned and new contexts lie in the subspace spanned by previously observed contexts, we can make accurate predictions despite the inability to recover the model parameters. An online linear regression (LR) procedure will take care of this issue, and we choose KWIK LR (Walsh et al., 2009) as such a procedure.
The original KWIK LR deals with scalar labels, which can be used to decide whether the estimate of p c (s |s, a) is sufficiently accurate (known). A (s, a) pair then becomes known if (s, a, s ) is known for all s . This approach, however, generally leads to a loose analysis, 5. ∆n is the n-simplex {x ∈ R n+1 : n+1 i=1 xi = 1, x ≥ 0}. 6. Here we use snext to denote the random variable, and s to denote a possible realization.
because there is no need to predict p c (s |s, a) for each individual s accurately: if the estimate of p c (·|s, a) is close to the true distribution under L 1 error, the (s, a) pair can already be considered as known. We extend the KWIK LR analysis to handle vector-valued outputs, and provide tighter error bounds by treating p c (·|s, a) as a whole. Below we introduce our extended version of KWIK LR, and explain how to incorporate the knownness information in Rmax skeleton to perform efficient exploration.
Identifying known (s, a) with KWIK LR The KWIK LR-Rmax algorithm we propose for the linear setting still uses Rmax template (Algorithm 1) for exploration: in every episode we build the induced MDP M K , and act greedily according to its optimal policy with balanced wandering. The major difference from Cover-Rmax lies in how the set of known states K are identified and how M K is constructed, which we explain below (see pseudocode in Algorithm 3).
At a high level, the algorithm works in the following way: when constructing M K , we query the KWIK procedure for estimatesp c (·|s, a) andr c (s, a) for every pair (s, a) using P redict(c, s, a). The KWIK procedure either returns ⊥ (don't know), or returns estimates that are guaranteed to be accurate. If ⊥ is returned, then we consider (s, a) as unknown and associate s with R max reward for exploration. Such optimistic exploration ensures significant probability of observing (s, a) pairs on which we have predicted ⊥. If we do observe such pairs in an episode, we call U pdate with feature-label pairs formed via Equation 5 to make progress on estimating parameters for unknown state-action pairs.
Next we walk through the pseudocode and explain how P redict and U pdate work in detail. Then we prove an upper bound on the number of updates that can happen (i.e., the if condition holds on Line 11), which forms the basis of our analysis of KWIK LR-Rmax.
In Algorithm 3, we initialize matrices Q and W for each (s, a) using Initialize(·) and update them over time. Let C t (s, a) be the design matrix at episode t, where each row is a context c τ such that (s, a) was observed in episode τ < t. By matrix inverse rules, we can verify that the update rule on Line 10 essentially yields Q t (s, a) = (I + C t C t ) −1 , where Q t (s, a) is the value of Q(s, a) in episode t. This is the inverse of the (unnormalized and regularized) empirical covariance matrix, which plays a central role in linear regression analysis. The matrix W accumulates the outer product between the feature vector (context) c and the one-hot vector label y = ({I[s next = s ]} ∀s ∈S ) . It is then obvious that Q t (s, a)W t (s, a) is the linear regression estimate of P (s, a) using the data up to episode t.
When a new input vector c t comes, we check whether Q(s, a)c t 2 is below a predetermined threshold α S (Line 5). Recall that Q(s, a) is the inverse covariance matrix, so a small Q(s, a)c t 2 implies that the estimate Q t (s, a)W t (s, a) is close to P (s, a) along the direction of c t , so we predict p ct (·|s, a) = c t P (s, a) ≈ c t Q(s, a)W (s, a); otherwise we return ⊥. The KWIK subroutine for rewards is similar hence omitted. To ensure that the estimated transition probability is valid, we project the estimated vector onto ∆ S−1 , which can be done efficiently using existing techniques (Duchi et al., 2008).
Below we state the KWIK bound for learning the transition function; the KWIK bound for learning rewards is much smaller hence omitted here. We use the KWIK bound for scalar linear regression from Walsh et al. (2009) and the property of multinomial samples to get our KWIK bound.
Proof sketch. (See full proof in the appendix.) We provide a direct reduction to KWIK bound for learning scalar values. The key idea is to notice that for any vector v ∈ R S : So conceptually we can view Algorithm 3 as running 2 S scalar linear regression simultaneously, each of which projects the vector label to a scalar by a fixed linear transformation f . We require every scalar regressor to have ( , δ/2 S ) KWIK guarantee, and the 1 error guarantee for the vector label follows from union bound.  Proof. When the KWIK subroutine (Algorithm 3) makes non-"⊥" predictionsp c (s, a, s ), we require that p c (·|s, a) − p c (·|s, a) 1 ≤ /8H.
After projection onto ∆ S−1 , we have: Further, the update to the matrices Q and W happen only when an unknown state action pair (s, a) is visited and the KWIK subroutine still predicts ⊥ (Line 10). The KWIK bound states that after a fixed number of updates to an unknown (s, a) pair, the parameters will always be known with desired accuracy. The number of updates m can be obtained by setting the desired accuracy in transitions to /8H and failure probability as δ/SA in Theorem 4.1: We now use Lemma 3.3 where instead of updating counts for number of visits, we look at the number of updates for unknown (s, a) pairs. On applying a union bound over all state action pairs and using Lemma 3.3, it is easy to see that the sub-optimal episodes are bounded by O mSA ln 1 δ with probability at least 1 − δ. The bound in Theorem 4.2 is obtained by substituting the value of m.
We see that for this contextual MDP, the linear structure helps us in avoiding the exponential dependence in context dimension d. The combined dependence on S and d is now O(max{d 4 S, d 2 S 3 }).

Related work
Transfer in RL with latent contexts The general definition of CMDPs captures the problem of transfer in RL and multi-task RL. See Taylor and Stone (2009) and Lazaric (2011) for surveys of empirical results. Recent papers have also advanced the theoretical understanding of transfer in RL. For instance, Brunskill and Li (2013) and Hallak et al. (2015) analyzed the sample complexity of CMDPs where each MDP is an element of a finite and small set of MDPs, and the MDP label is treated as the latent (i.e., unseen) context. Mahmud et al. (2013) consider the problem of transferring the optimal policies of a large set of known MDPs to a new MDP. The commonality of the above papers is that the MDP label (i.e., the context) is not observed. Hence, their methods have to initially explore in every new MDP to identify its label, which requires the episode length to be substantially longer than the planning horizon. This can be a problematic assumption in our motivating scenarios, where we interact with a patient / user / student for a limited period of time and the data in a single episode (whose length H is the planning horizon) is not enough for identifying the underlying MDP. In contrast to prior work, we propose to leverage observable context information to perform more direct transfer from previous MDPs, and our algorithm works with arbitrary episode length H.
RL with side information Our work leverages the available side-information for each MDP, which is inspired by the use of contexts in contextual bandits (Langford and Zhang, 2008;Li et al., 2010). The use of such side information can also be found in RL literature: Ammar et al. (2014) developed a multi-task policy gradient method where the context is used for transferring knowledge between tasks; Killian et al. (2016) used parametric forms of MDPs to develop models for personalized medicine policies for HIV treatment.
RL in metric space For smooth CMDPs (Section 3), we pool observations across similar contexts and reduce the problem to learning policies for finitely many MDPs. An alternative approach is to consider an infinite MDP whose state representation is augmented by the context, and apply PAC-MDP methods for metric state spaces (e.g., C-PACE proposed by Pazis and Parr (2013)). However, doing so might increase the sample and computational complexity unnecessarily, because we no longer leverage the structure that a particular component of the (augmented) state, namely the context, remains the same in an episode. Concretely, the augmenting approach needs to perform planning in the augmented MDP over states and contexts, which makes its computational/storage requirement worse than our solution: we only perform planning in MDPs defined on S, whose computational characteristics have no dependence on the context space. In addition, we allow the context sequence to be chosen in an adversarial manner. This corresponds to adversarially chosen initial states in MDPs, which is usually not handled by PAC-MDP methods.
KWIK learning of linear hypothesis classes Our linear combination setting (Section 4) provides an instance where parametric assumptions can lead to substantially improved PAC bounds. We build upon the KWIK-Rmax learning framework developed in previous work (Li et al., 2008;Szita and Szepesvári, 2011) and use KWIK linear regression as a sub-routine. For the resulting KWIK LR-Rmax algorithm, its sample complexity bound inherently depends on the KWIK bound for linear regression. It is well known that even for linear hypothesis classes, the KWIK bound is exponential in input dimension in the agnostic case (Szita and Szepesvári, 2011). Therefore, the success of the algorithm relies on the validity of the modelling assumption.
Abbasi-Yadkori and Neu (2014) studied a problem similar to our linear combination setting, and proposed a no-regret algorithm by combining UCRL2 (Jaksch et al., 2010) with confidence set techniques from stochastic linear optimization literature (Dani et al., 2008;Filippi et al., 2010). Our work takes an independent and very different approach, and we provide a PAC guarantee which is not directly comparable to regret bound. Still, we observe that our dependence on A is optimal for PAC whereas theirs is not ( √ A is optimal for bandit regret analysis and they have A); on the other hand, their dependence on T (the number of rounds) is optimal, and our dependence on 1/ , its counterpart in PAC analysis, is suboptimal. It is an interesting future direction to combine the algorithmic ideas from both papers to improve the guarantees.

Conclusion
In this paper, we present a general setting of using side information for learning near-optimal policies in a large and potentially infinite number of MDPs. The proposed Cover-Rmax algorithm is a model-based PAC-exploration algorithm for the case where MDPs vary smoothly with respect to the observed side information. Our lower bound construction indicates the necessary exponential dependence of any PAC algorithm on the context dimension in a smooth CMDP. We also consider another instance with a parametric assumption, and using a KWIK linear regression procedure, present the KWIK LR-Rmax algorithm for efficient exploration in linear combination of MDPs. Our PAC analysis shows a significant improvement with this structural assumption.
The use of context based modelling of multiple tasks has rich application possibilities in personalized recommendations, healthcare treatment policies, and tutoring systems. We believe that our setting can possibly be extended to cover the large space of multi-task RL quite well with finite/infinite number of MDPs, observed/latent contexts, and deterministic/noisy mapping between context and environment. We hope our work spurs further research along these directions. We now need to compute the required resolution of the cover and the number of transitions m which will guarantee the approximation for the value functions as required in the previous lemma. The following result is the key result: Lemma A.5 (Cover approximation). For a given CMDP and a finite cover, i.e., C = ∪ and |r c 1 (s, a) − r c 2 (s, a)| ≤ /8 if we visit every state-action pair m = 128(S ln 2+ln SA δ )H 2 2 times in a ball B i summing observations over all c ∈ B i , then, for any policy π and with probability at least 1 − 2δ, the approximate MDPM i corresponding to B i computed using empirical averages will satisfy Proof. For each visit to a state action pair (s, a), we observe a transition to some s ∈ S for context c t ∈ B i in t th visit with probability p ct (s, a). , the error becomes /4H. One can easily verify using similar arguments that, the error in rewards for any context c ∈ B i is less than /4. By using the simulation lemma A.2, we get the desired result.

B.2 Lower bound for smooth CMDP
In the lower bound construction, we made a key argument that the defined contextual mapping with specific instances, requires any agent to learn the models of a significant number of MDPs separately to get decent generalization. The construction populates a set of packing points in the context space with hard MDPs and argues that these instances are independent of each other from the algorithm's perspective. To formalize this statement, let Z be the 8 -packing as before. The adversary makes the choices of the instances at each context c ∈ Z, as follows: Select an MDP from the family of hard instances described in Figure 1 where the optimal action from each state in {1, . . . , n} is chosen randomly and independently from the other assignments. The parameter deciding the difference in optimality of actions in Figure 1 is taken as 160H e 4 (H−2) which is required by the construction of Theorem 3.5.
We denote these instances by the set I and an individual instance by I z where z = {z 1 , z 2 , . . . , z |Z| } is the set of these random choices made by the adversary. By construction, we have a uniform distribution Γ ≡ Γ 1 × Γ 2 × . . . Γ |Z| over these possible set of instances I z . From claim 5.3, any assignment of optimal actions to these packing points would define a valid smooth contextual MDP. Further, the independent choice of the optimal actions makes MDPs at each packing point at least as difficult as learning a single MDP. Formally, let the sequence of transitions and rewards observed by the learning agent for all packing points be T ≡ {τ 1 , τ 2 , . . . , τ |Z| }. Due to the independence between individual instances, we can see that: where P Γ [τ i ] denotes the distribution of trajectories τ i . Thus, observing trajectories from MDP instances at other packing points does not let the algorithm deduce anything about the nature of the instance chosen for one point. With respect to this distribution, learning in the contextual MDP is equivalent to or worse than simulating a single MDP learning algorithm at each of these packing points. For any given contextual MDP algorithm Alg, we have: where Alg* is an optimal single MDP learning algorithm. The expectation is with respect to the distribution over the instances I z and the algorithm's randomness. From Theorem 3.5, we can lower bound the expectation on the right hand side of the inequality by Ω SA 2 . The total number of mistakes is lower bounded as: Setting |Z| = D(C, 1 ) ≤ N (C, 1 ) gives the stated lower bound with 1 = 8 .
Appendix C. Proof of the Theorem 4.1 In this section, we present the proof of our KWIK bound for learning transition probabilities. Our proof uses a reduction technique that reduces the vector-valued label setting to the scalar setting, and combines the KWIK bound for scalar labels given by Walsh et al. (2009). For our purpose, β = 2 as |η f | ≤ 2 and M = √ d as θ f 2 = P (s, a)f 2 ≤ √ d. Now set δ = δ 2 S in Theorem C.1. In the KWIK linear regression algorithm, the known status for a context c is checked in the same manner as done in Line 5 in Algorithm 3. Therefore