Nov 24, 2016 - Suppose that someone currently writing a paper at his/her office desk wants to get out of the building; the office is located in the th...

3 downloads 26 Views 961KB Size

No documents

Multiscale Inverse Reinforcement Learning using Diffusion Wavelets Jung-Su Ha and Han-Lim Choi Dept. of Aerospace Engineering, KAIST {wjdtn1404, hanlimc}@kaist.ac.kr

Abstract This work presents a multiscale framework to solve an inverse reinforcement learning (IRL) problem for continuous-time/state stochastic systems. We take advantage of a diffusion wavelet representation of the associated Markov chain to abstract the state space. This not only allows for effectively handling the large (and geometrically complex) decision space but also provides more interpretable representations of the demonstrated state trajectories and also of the resulting policy of IRL. In the proposed framework, the problem is divided into the global and local IRL, where the global approximation of the optimal value functions are obtained using coarse features and the local details are quantified using fine local features. An illustrative numerical example on robot path control in a complex environment is presented to verify the proposed method.

1

Introduction

In this paper, we address an inverse reinforcement learning (IRL) problem (or often equivalently called an inverse optimal control problem) for robots operated in a complex environment over a long time horizon, where its forward problem is a continuous time/continuous stochastic optimal control problem called linearly solvable optimal control (LSOC). The objective of an IRL problem is to recover the value and cost functions as well as the optimal policy when experts’ demonstrations are given. A general method to solve IRL problem for standard MDP often involves the procedure of solving the corresponding forward problem in every iteration [7], while a method for LSOC does not [2]. The IRL solution method for LSOC is formulated as a convex optimization problem, where its gradient and Hessian are obtained analytically. However, the optimization tends to be intractable as the size of the problem increases and moreover, in real situation, a demonstration data set may be not sufficient to represent whole state space. Finding a sparse structure of the problem and representing the problem with few meaningful bases are essential for obtaining the solution efficiently. To address a large-scale problems effectively, it is conceivable that the hallmark of human intelligence could provide some insights - in particular, this work notes multiscale and hierarchical structure of human decision making. Suppose that someone currently writing a paper at his/her office desk wants to get out of the building; the office is located in the third floor of the building and there is one set of staircases and an elevator. Then, what would this person’s control policy look like? This person would not try to figure out what he/she should do for all possible situations he/she could face like the standard value function-based approach; instead, he/she would figure out which building gate he/she would use, whether to take the elevator or the stairs, which door he would exit from the room (if there are more than one), etc. A detailed plans such as “which particular start he/she should put their left foot,” would be determined later in the process of executing a piece of overall plan, for example, “go downstairs using the staircase.” It should be noted that this human-like decision takes advantage of the underlying (multiscale) hierarchical structure of state space; in the above example, detailed aspects such a particular certain sequence of stairs to put on is abstracted by just a single notion of staircase. 30th Conference on Neural Information Processing Systems (NIPS 2016), Barcelona, Spain.

Robot dynamics & Environment geometry

Time & State Discretization

Hierarchical Abstraction

Sampling and MC Construction

𝑢𝑢∗ 𝑡𝑡 , ∀𝑡𝑡 ∈ [𝑡𝑡0 , 𝑡𝑡0 + 𝜏𝜏]

Transition data from demonstrations, {x,x'}

Diffusion Wavelets

Local IRL

Compute Control Receding Horizon Control

Adaptive Basis Addition

Global IRL

Recursive solution of Compressed IRL

For every 𝜏𝜏 = ℎ × 𝑘𝑘𝑅𝑅𝑅𝑅𝑅𝑅

For every ℎ × 𝑘𝑘𝐿𝐿𝐿𝐿

If demonstrations for new task are given

Figure 1: Proposed framework Same intuition can be applied to the inverse inference problem. By utilizing it, the key contribution of this paper is to present a systematic framework for solving IRL, as depicted in Fig. 1. The framework consists of five phases: (i) In the discretization phase, Markov chain associated to the robot dynamics is constructed by sampling a set finite set of states in state-space. (ii) In the abstraction phase, the hierarchical bases structure is obtained using the diffusion wavelet method. (iii) In global IRL phase, an IRL problem constructed only using the coarse bases (or on “abstract-state”) is solved; this IRL is much more tractable to handle than using the original bases set. (iv) In the local planning phase, the finer bases located in focused regions where the demonstration data visit frequently are sought for and detailed solution associated with these focused regions are computed. (vi) In the control phase, a continuous control sequence is computed and applied to the robot in a receding horizon fashion. Rest of the paper is primarily focused on elaborating the details of this framework, followed by numerical example for validation of the method.

2 2.1

Inverse Reinforcement Learning Linearly-solvable optimal control

Consider a stochastic dynamics whose deterministic drift term is affine in control input: dx = f (x)dt + G(x)(udt + σdw) dx

(1)

du

where x ∈ X ⊂ R , u ∈ R and w are a state, control vector and an du -dimensional Brownian motion process, respectively. Let functions q : X → R and πc : X → Rdu be an instantaneous state cost rate and a control policy, respectively and define an instantaneous cost rate as l(x, u) := q(x) + 2σ1 2 uT u. Then the cost functional is given as: Z tf 1 πc Jconti (x) = lim E l(x(t), πc (x(t))dt . (2) tf →∞ tf 0 The problem with the cost function (2) and dynamics (1) is called the infinite horizon average cost stochastic optimal control (SOC) problem or is referred as linearly-solvable optimal control (LSOC) problem since its solution is obtained from the linear partial differential equation [6]. If the time-axis is discretized by a time step h, the transition probability of one step without/with any control input is defined as: x[k + 1] ∼ p(·|x[k]), x[k + 1] ∼ π(·|x[k]),

(3)

which are called the passive and controlled dynamics, respectively. The passive and controlled dynamics are approximated as N (y; µ(h), Σ(h)), where N is a Gaussian distribution with a mean µ(h) and covariance Σ(h). For small h, the Kullback-Leibler divergence between two distributions is approximated as DKL (π(·|x)||p(·|x)) = 2σh2 u0 u. Therefore, the cost functional (2) is written in the discrete time setting: "K # X 1 π J (x) = lim E hq(x[k]) + DKL (π(·|x[k])||p(·|x[k])) . (4) K→∞ K k=0

2

Moreover, the state space can be discretized by sampling a set of states X = {xn } [3]. Transition probability matrix for passive dynamics P , where Pnm means a transition probability from xn to xm , is approximated via Gaussian distribution as: N (xm : µ(h), Σ(h)) , 0 m0 N (xm : µ(h), Σ(h)

Pnm = P

(5)

where µ(h) and Σ(h) are computed by integrating moment dynamics of linearized SDE for t ∈ [0, h]: ˙ µ(t) ˙ = Aµ(t) + c, Σ(t) = AΣ(t) + Σ(t)A0 + BB 0 (6) df from µ(0) = xn , Σ(0) = 0, where A = dx , B = σG(xn ) and c = f (xn ) − Axn . One x=xn can truncate tails of Gaussian distribution to make P sparse. With a set of discrete states, X, the state-space as well as time-axis discretized version of SOC is formulated as the Markov decision process (MDP) of which cost function is given by (4). This type of MDP is called the linearly-solvable MDP [6] and its solution is known to converge to SOC solution as |X| → ∞ and h → 0. Define the optimal cost-to-go value c := minπ J π (x), the differential cost-to-go function v(x), the desirability function z(x) = exp(−v(x)), and the linear operator G[z](x) = P (differential) 0 0 p(x |x)z(x ). Then z(x) satisfies the following linear Bellman equation: 0 x exp(−c)z(x) = exp(−hq(x))G[z](x),

(7)

and the optimal policy is obtained analytically: π ∗ (x0 |x) =

p(x0 |x)z(x0 ) . G[z](x)

(8)

For more details of problem formulation and discrete approximation method, we would refer the reader to [3] and references therein. 2.2

Inverse reinforcement learning for LSOC problem

While the objective of (forward) SOC problem is to find the optimal control policy for the given system and cost function, the objective of inverse reinforcement learning (IRL) problem is to recover the value and cost functions as well as the optimal policy when experts’ demonstrations are given [2, 7]. Suppose a dataset of transitions {xn , x0n }n=1,··· ,N is obtained from the optimal policy (8). Let v be a vector representation of value function v(·). Then, the negative log-likelihood of the dataset is given as: L[v] = aT v + bT log(P exp(−v)), (9) 0 where each component of a and b represent visitation counts of xn and xn , respectively. Since L is convex and its gradient and Hessian are computed analytically, it can be minimized by applying iterative second order convex optimization methods. Once the value function v is obtained, the cost function q(·) and the optimal policy π ∗ (·|·) can be recovered directly from (7) and (8), respectively. In real situation, however, a and b are sparse since the dataset is not sufficient, so it is impossible to compute v(·) over whole state space. Also, the optimization procedure in (9) gets intractable for problems having the large size of the state space. To represent the problem efficiently, a linear value function approximation can be used: ˆ = Φw, v (10) where each column of Φ represents a feature (or basis) and w is weight. Note that L[w] is also convex. In this work, we obtain the hierarchical structure of feature sets which is naturally induced from the passive (diffusion) dynamics of the system and utilize those set of features to solve IRL problem efficiently.

3 3.1

Multiscale Inverse Reinforcement Learning Multiscale feature extraction: Diffusion Wavelets

From now on, we consider T = P 0 for notion simplicity; then Tnm represents a transition probability from xm to xn . The Markov chain, T , obtained by discretizing a diffusion process ((1) with u = 0) is known to have some interesting properties: local, smoothing and contractive [1]. From any initial 3

point, δm , the state (numerically) transitions to only a few its neighbors (i.e., T δm has a small support) and T j δm is a smooth probability distribution. Also since ||T ||2 ≤ 1, a dimension of a subspace, Vj , which is -spanned by {T j δm }xm ∈X monotonically decreases as j increases and V0 ⊇ V1 ⊇ · · · ⊇ Vj ⊇ · · · ; especially for an irreducible Markov chain, dim(Vj ) → 1 as j increases and a limit of Vj corresponds to the stationary distribution of the Markov chain. Let Wj be an orthogonal complement of Vj+1 into Vj , i.e., Vj = Vj+1 ⊕ Wj and suppose the orthonormal bases Φj and Ψj span Vj and Wj , respectively. By using aforementioned properties of T , Diffusion wavelets constructs a hierarchical structure of a set of well-localized bases Φj and Ψj called scaling and wavelet functions, respectively, in order that the subspace spanned by feature set 2 j−1 j [Φj ]Φ0 is j-close to the subspace spanned by {T 1+2+2 +···+2 δm = T 2 −1 δm }xm ∈X . Roughly speaking, Φj and Ψj represent smooth bump function and oscillatory function, respectively. We omit the procedure of Diffusion wavelets algorithm because of the space limitation and would refer the readers to [1, 3] for more details. Let [B]Φj be a set of vectors B represented on a basis Φj , where the columns of [B]Φj are the coordinates of the vectors B in the coordinates Φj . A set of features at level j can be written in the original coordinate (or can be unpacked) as: Φj = [Φj ]Φ0 = [Φj−1 ]Φ0 [Φj ]Φj−1 = [Φ1 ]Φ0 · · · [Φj−1 ]Φj−2 [Φj ]Φj−1 ,

(11)

which is represented as a |X| × |Xj | matrix. Note that each column of [Φj ]Φ0 can be viewed as an “abstract-state" of the original Markov chain. At the scale j, there are only |Xj | meaningful combinations of states and each combination, [Φj ]Φ0 , represents “abstract-state”. 3.2

IRL with hierarchical multiscale feature sets

Rather than solving original |X|-dimensional optimization problem, we can treat the lowerdimensional coarsened problem. Suppose a set of “abstract-state" at level j, Φj , is utilized as a set of features, which means the problem is viewed in a lower resolution with 2j time scale. Then, ˆ j = Φj wj , and the optimization v is approximated as a linear combination of this feature set as v problem (9) is also written as: L[wj ] = aT Φj wj + bT log(P exp(−Φj wj )).

(12)

The compressed problem (12) is much more tractable than the original problem (9) if |Xj | << |X|. Note that due to its localization property, the features are naturally interpretable (see Fig. 2 (a)–(d)); thus user can choose the appropriate level, considering trade-off between the size of the problem and the solution quality. Also, the hierarchical structure of diffusion wavelet tree can be utilized to solve the problem more efficiently; the solution of jth level, wj , can provide an initial guess to ˜ j−1 = [Φj ]Φj−1 wj (j − 1)th level problem; that is, the optimization at (j − 1)th level starts from w ˆ j is not sharply changed through the scale j, this initial guess by unpacking the jth level solution. If v would be near the optimum of (j − 1)th level problem and the optimization procedure would rapidly converge to the minimum. Suppose we solve lth level problem and have the approximate solution vl by using features Φl . Then, we can achieve more exact solution by considering additional features from wavelet functions, Ψ1:(l−1) which are orthogonal to Φl (note that V0 = Vl ⊕ Wl−1 ⊕ Wl−2 ⊕ · · · ⊕ W0 ). The wavelet bases are also built as being well-localized. In this work, we utilize the intuition that the important region where optimal policy frequently visits should have highest resolution [3, 2]. The wavelet functions are evaluated as: s = bT |Ψ1:(l−1) |, (13) where the score s ∈ R|X|−|Xl | represents how each feature overlap with visitation counts of x0n . By adding features with high scores and solving the corresponding problem, the value, cost functions and policy will have higher resolution in important region. Also, user can easily choose the additional features for their objective since the wavelet features are also interpretable. Finally, the continuous-time optimal policy can be extracted from the multiscale quantification of the optimal value structure for the interval τ = h × kRHC in a receding horizon control fashion. The control can be computed by matching the 1st order moment of original SDE (1) and the optimal policy (8) for the MDP: u∗ (t) = −σ 2 GT eA

T

(τ −t)

4

Σ(τ )−1 (µ(τ ) − ynew ),

(14)

5

5

5

4

4

4

4

3

3

3

3

2

2

2

2

1

1

1

0

-1

-1

-1

-1

-2

-2

-2

-2

-3

-3

-3

-3

-4

-4

-4

-5 5 -5

-5 5 -5

-4

-3

-2

-1

0

1

2

3

4

-4

x

-3

-2

-1

0

1

2

3

4

-4

-4

-3

-2

x

(a)

4.2

-3

-2

-1

2

3

4

5

5 3.75 3.7

3

1

3.65

1

3.6

0

3.55

-1

3.5

-2

3.45

3.8

3.8

0

y 3.7

2 3.9

1

-1

1

4

2

0

0

4

3.9

y

0.15

-2

-4

4.1

2

0.2

-5 -5

x

3

0.3

0.25

5

(d)

4

2

-1

4

4

3

0

3

4.1

3

1

2

5

4.2

4

0.35

1

(c)

5

0.4

4

0

x

(b)

5

-1

y

-5 -5

y

1

0

y

0

y

0

y

y

5

3.7

-1

-2

3.6 -2

3.6

-3

3.5 -3

3.5 -3

-4

3.4 -4

3.4 -4

0.1

0.05

-4 -5 -5

0

0

5

-5 -5

0

(e) q(x)

-5 -5

5

0

x

5

0.4

4

0.35

3

5

0.2

y

y

0

0.15

-2

0.9

2

0.8

1

0.7 0.6

0 -1

0.4

-2 -3 0.05

-4 -5 -5

0

0

5

ˆ8 (h) v 0.8

2000

0.6

1500 0.4 1000 0.2

500

0.3

0.1

-3

5

2500 1

3

0.5

-1

0

x

1.1

4

# of bases

0.25

1

3.35

-5 -5

x

0.3

2

5

ˆ3 (g) v

(f) v(x)

x

3.4

RMS error

-3

0.2

-4 -5 -5

ˆ3 (i) q x

0.1

0 0

2

4

6

8

10

12

14

0 16

Level 0

5

ˆ8 (j) q x

(k)

Figure 2: Scaling functions at (a - b) level 3 and (c - d) level 8. (f - h) Exact and approximate value functions. (e, i - j) Exact and approximate cost functions. (k) The number of features, RMS error of the solution at each level. P where ynew = y∈X yP r(x[kRHC ] = y|x[0] = xcur , π ˆ ∗ ) denotes the expected state after τ when the system follows the (approximate) optimal policy (8) from xcur .

4

Numerical Example

We consider a simple two-dimensional stochastic single integrator in the fractal-like environment. The environment consists of 5 groups of rooms where one group is made of 5 square rooms as shown in Fig. 2; one can observe that the environment has 2 level self-similarity, which makes the problem have a multiscale nature. The dynamics is given by f (x) = 0, G(x) = I2 ; that is, the position of a robot in the configuration space, x ∈ X , is controlled by the velocity input, u ∈ R2 while being disturbed by white noise. We set h = 0.1 and σ = 1. In order to discretize the state space, 100 samples are obtained from each room, therefore there are 2500 discrete state in total. The transition data set is obtained from true occupancy measure induced by the optimal policy, which is equivalent to using infinite number of samples. Fig. 2 (a) - (d) shows multiscale features: some scaling functions in the diffusion wavelet tree at level 3 and 8. It is seen that at level 3 and 8, where roughly 8h and 256h are considered as 1-step, scaling functions roughly represent each small room and one group of 5-rooms, respectively; it has no meaning to make a distinction within a room or a group of 5 rooms at those levels. Fig. 2 (g, i) and (h, j) depict the approximated value and cost functions at level 3 and 8, respectively and Fig. 2 (k) shows the number of features and the RMS error of value functions at each scale. At level 3, only 421 basis functions out of 3000 are used, but the value and cost functions are recovered quite exactly; at level 8, where 25 bases are used, even though the solution contains some error, it interprets the information of optimal policy and the preference of cost function between the groups of 5 rooms. We omit the computational time at each scale because of the space limitation but we observe that the computational time decreases as the number of feature increases; It is obvious that there is a trade-off between the solution quality and the computational cost and the appropriate level can be chosen by observing the features at that level. 5

Acknowledgments This work was supported by Agency for Defense Development (under contract #UD150047JD).

References [1] Ronald R Coifman and Mauro Maggioni. Diffusion wavelets. Applied and Computational Harmonic Analysis, 21(1):53–94, 2006. [2] Krishnamurthy Dvijotham and Emanuel Todorov. Inverse optimal control with linearly-solvable mdps. In Proceedings of the 27th International Conference on Machine Learning, 2010. [3] Jung-Su Ha and Han-Lim Choi. Multiscale abstraction, planning and control using diffusion wavelets for stochastic optimal control problems. In arXiv preprint, 2016. [4] Mauro Maggioni and Sridhar Mahadevan. Fast direct policy evaluation using multiscale analysis of markov diffusion processes. In Proceedings of the 23rd international conference on Machine learning, pages 601–608. ACM, 2006. [5] Sridhar Mahadevan and Mauro Maggioni. Value function approximation with diffusion wavelets and laplacian eigenfunctions. In Advances in neural information processing systems, pages 843–850, 2005. [6] Emanuel Todorov. Efficient computation of optimal actions. Proceedings of the national academy of sciences, 106(28):11478–11483, 2009. [7] Brian D Ziebart, Andrew L Maas, J Andrew Bagnell, and Anind K Dey. Maximum entropy inverse reinforcement learning. In AAAI, pages 1433–1438, 2008.

6