LORM: a novel reinforcement learning framework for biped gait control

Overall theory

Reinforcement learning aims to generate an agent to provide proper actions according to the observations from the environment. A reinforcement learning structure is composed of an agent to be trained and the environment for training and evaluating. As shown in Fig. 1, the agent repeatedly interacts with the environment in RL training, according to the value function generated by a neural network (Mnih et al., 2013; Mnih et al., 2015) or a simple Q-value table that stores the expected reward for each combination of the state and action. In trajectory T (from the beginning state to the end state of the environment), at the t-th time step, the agent receives observation o_t from the environment, which is based on current state s_t of environment and outputs action a_t accordingly. The environment receives the action and gets into the next state s_t+1 and outputs reward r_t for the action and new observation o_t+1. After some iterations of interaction, the agent will update its weights to get a higher total reward $R_{total}=\sum _{t=1}^n{r}_t$, where n stands for the number of steps in one trajectory. This may be achieved by either getting higher reward in every step or by trying to maintain more steps before the end of the trajectory.

Policy gradient

Policy gradient algorithms are one mainstream family of reinforcement learning algorithms that are capable of complex continuous tasks. Policy gradient algorithms generate optional actions directly without one function or table indicating the reward for different actions. The output of the network is the probability distribution of actions.

(1) $$p(a|s,\theta )=\mathbb{P}[A_t=a|S_t=s,{\theta }_t=\theta ]$$where θ is the current weights of the neural network. The agent will choose the action according to the probability distribution at every step. Thus the probability P(τ|θ) of trajectory τ = (s₁, a₁,…, s_T, a_T) and the total reward can be calculated as:

(2) $$P(\tau |\theta )=p(s_1)\prod _{t=1}^T{p}_{\theta }(a_t|s_t,\theta )p(s_{t+1}|s_t,a_t)$$

(3) $$R(\tau )={\sum }_{t=1}^T{r}_t$$

Given N different trajectories τ¹,τ²,…,τ^N sampled by policy π with weights θ, the expectation of the total reward for one trajectory can be approximated as:

(4) $${\overline{R}}_{\theta }=\sum _{\tau }R(\tau )P(\tau |\theta )\approx {\displaystyle \frac{1}{N}\sum _{n=1}^{N}R({\tau }^n)}$$

The training will optimize the weights θ of policy π to obtain the highest expected reward.

(5) $${\theta }^{\ast }=\underset{\theta }{\arg max}{\overline{R}}_{\theta }$$

During the gradient ascent of backward propagation, the probabilities of actions bringing more rewards will increase while other probabilities are reduced. The gradient can be calculated as:

(6) $$\begin{array}{c}\mathrm{\nabla }{\overline{R}}_{\theta }\approx {\displaystyle \frac{1}{N}\sum _{n=1}^{N}R({\tau }^n)\mathrm{\nabla }\log P({\tau }^n|\theta )}\hfill \\ ={\displaystyle \frac{1}{N}\sum _{n=1}^{N}R({\tau }^n)\sum _{t=1}^{T_n}\mathrm{\nabla }\log p(a_t^n|s_t^n,\theta )}\hfill \\ ={\displaystyle \frac{1}{N}\sum _{n=1}^N\sum _{t=1}^{T_n}R({\tau }^n)\mathrm{\nabla }\log p(a_t^n|s_t^n,\theta )}\hfill \end{array}$$

Advantage function A^θ(s_t, a_t) is introduced replacing R(τ) to further improve the training efficiency:

(7) $$\mathrm{\nabla }{\overline{R}}_{\theta }\approx {\displaystyle \frac{1}{N}\sum _{n=1}^N\sum _{t=1}^{T_n}{A}^{\theta }(s_t,a_t)\mathrm{\nabla }\log p(a_t^n|s_t^n,\theta )}$$

(8) $$A^{\theta }(s_t,a_t)=G_t-b$$

(9) $$G_t={\sum }_{t^{\mathrm{\prime }}=t}^{T_n}{\gamma }^{t^{\mathrm{\prime }}-t}{r}_{t^{\mathrm{\prime }}}^n$$where G_t is the discounted accumulated reward, γ is the discounted coefficient, and b is the regularizer which will be further discussed in the next subsection.

Actor-critic structure

Actor-critic structure introduces one critic that calculates the value function V^θ(s) of the state s indicating whether the state is good or bad. The value function is used in the calculation of the advantage function. By the prediction of critic, the advantage function of each state can be calculated respectively without the information of the whole trajectory, as shown in Eq. (10). Here r_t + V^θ(s_t+1) acts as the discounted accumulated reward G_t and V^θ(s_t) acts as regularizer b.

(10) $$A^{\theta }(s_t,a_t)=r_t+V^{\theta }(s_{t+1})-V^{\theta }(s_t)$$

The critic is a neural network in deep reinforcement learning, which will be trained together with the agent (called actor) during the training process. The loss for critic can be calculated by temporal-difference (TD) approach:

(11) $$Loss=V^{\theta }-\gamma V^{\theta }(s_{t+1})-r_t$$

Proximal policy gradient

The sampling is time-consuming in the training of reinforcement learning. However, the sampled trajectories expired after the update of policy network π and new trajectories need to be sampled based on the new policy, which resulted in low trajectory utilization and training efficiency. Thus importance sampling is introduced in TRPO (Schulman et al., 2015a) and PPO (Schulman et al., 2017) to use the trajectories repeatedly, improving the training efficiency.

Importance sampling is widely used to estimate the expectation where the sampling of the distribution is difficult to obtain. Given one distribution x∼p(x), the expectation f(x) can be estimated as Eq. (12), where x_i is sampled from p(x).

(12) $$E_{x\sim p}[f(x)]={\displaystyle \frac{1}{n}\sum _{i}f(x_i)}$$

However, when the sample of x_i is difficult to access, importance sampling can be used by introducing another distribution x∼q(x) with a relatively small difference with p(x). The expectation f(x) with distribution p(x) can be calculated as Eq. (13), where x_i is sampled from q(x). Based on importance sampling, PPO uses one neural network to interact with the environment collecting samplings and another neural network to update weights, which is called an off-policy algorithm. The training efficiency and sample utilization are largely improved. In addition, the difference of distribution p(x) and q(x) should be small to ensure the accuracy.

(13) $$E_{x\sim p}[f(x)]=E_{x\sim q}[f(x){\displaystyle \frac{p(x)}{q(x)}}]={\displaystyle \frac{1}{n}\sum _{i}f(x_i){\displaystyle \frac{p(x_i)}{q(x_i)}}}$$

PPO maintains two actors π_θ(Pi) and π_θ′ (OldPi) with the same network structure. The actor π_θ is used to update weights while the π_θ′ is used to interact with the environment sampling the training data. The gradient for optimization is calculated based on importance sampling:

(14) $$\mathrm{\nabla }{\overline{R}}_{\theta }=E_{\tau \sim p_{\theta }(\tau )}[R(\tau )\mathrm{\nabla }\log p_{\theta }(\tau )]=E_{\tau \sim p_{{\theta }^{\mathrm{\prime }}}(\tau )}\left[{\displaystyle \frac{p_{\theta }(\tau )}{p_{{\theta }^{\mathrm{\prime }}}(\tau )}R(\tau )\mathrm{\nabla }\log p_{\theta }(\tau )}\right]$$

PPO maintains one critic neural network for state evaluation. Introducing advantage function to replace R(τ), the gradient is further calculated as:

(15) $$\begin{array}{c}\mathrm{\nabla }{\overline{R}}_{\theta }=E_{(s_t,a_t)\sim {\pi }_{\theta }}[A^{\theta }(s_t,a_t)\mathrm{\nabla }\log p_{\theta }(a_t^n|s_t^n)]\\ =E_{(s_t,a_t)\sim {\pi }_{{\theta }^{\mathrm{\prime }}}}\left[{\displaystyle \frac{P_{\theta }(s_t,a_t)}{P_{{\theta }^{\mathrm{\prime }}}(s_t,a_t)}{A}^{\theta }(s_t,a_t)\mathrm{\nabla }\log p_{\theta }(a_t^n|s_t^n)}\right]\\ \approx E_{(s_t,a_t)\sim {\pi }_{{\theta }^{\mathrm{\prime }}}}\left[{\displaystyle \frac{P_{\theta }(s_t,a_t)}{P_{{\theta }^{\mathrm{\prime }}}(s_t,a_t)}{A}^{{\theta }^{\mathrm{\prime }}}(s_t,a_t)\mathrm{\nabla }\log p_{\theta }(a_t^n|s_t^n)}\right]\\ =E_{(s_t,a_t)\sim {\pi }_{{\theta }^{\mathrm{\prime }}}}\left[{\displaystyle \frac{p_{\theta }(a_t|s_t)}{p_{{\theta }^{\mathrm{\prime }}}(a_t|s_t)}{\displaystyle \frac{p_{\theta }(s_t)}{p_{{\theta }^{\mathrm{\prime }}}(s_t)}{A}^{{\theta }^{\mathrm{\prime }}}(s_t,a_t)\mathrm{\nabla }\log p_{\theta }(a_t^n|s_t^n)}}\right]\\ \approx E_{(s_t,a_t)\sim {\pi }_{{\theta }^{\mathrm{\prime }}}}\left[{\displaystyle \frac{p_{\theta }(a_t|s_t)}{p_{{\theta }^{\mathrm{\prime }}}(a_t|s_t)}{A}^{{\theta }^{\mathrm{\prime }}}(s_t,a_t)\mathrm{\nabla }\log p_{\theta }(a_t^n|s_t^n)}\right]\end{array}$$

The objective function for gradient ascent can be derived from ∇R_θ:

(16) $$J^{{\theta }^{\mathrm{\prime }}}(\theta )=E_{(s_t,a_t)\sim {\pi }_{{\theta }^{\mathrm{\prime }}}}\left[{\displaystyle \frac{p_{\theta }(a_t|s_t)}{p_{{\theta }^{\mathrm{\prime }}}(a_t|s_t)}{A}^{{\theta }^{\mathrm{\prime }}}(s_t,a_t)}\right]$$

The clip coefficient ε is furthermore introduced to limit the update step, which can keep the difference between Pi and OldPi small enough for importance sampling. The final objective function for gradient ascent is calculated as Eq. (17), where k denotes the number of current iteration. The overall structure of PPO is shown in Fig. 2. In this work, PPO from OpenAI Baselines (Dhariwal et al., 2017) is implemented as the training algorithm with hyperparameters as shown in Table 1. The neural network structure of Pi and OldPi is shown in Fig. 3, and is composed of one input layer, one output layer, and two hidden layers. The algorithm also introduced generalized advantage estimation (GAE) (Schulman et al., 2015b) in the calculation of advantage function, improving the training performance.

(17) $$J_{PPO}^{{\theta }^k}(\theta )=\sum _{(s_t,a_t)}min({\displaystyle \frac{p_{\theta }(a_t|s_t)}{p_{{\theta }^k}(a_t|s_t)}{A}^{{\theta }^k}(s_t,a_t),clip\left({\displaystyle \frac{p_{\theta }(a_t|s_t)}{p_{{\theta }^k}(a_t|s_t)},1-\epsilon ,1+\epsilon )A^{{\theta }^k}(s_t,a_t)}\right)}$$

Task reward

Task reward, including velocity reward and deviation penalty, was used to prompt the agent to interact properly with the environment and complete the task. The velocity reward encouraged the robot to walk at the expected velocity stably. Two different velocity rewards were designed according to the tasks. For the task of walking as fast as possible, the velocity reward was calculated as shown in Eq. (20).

(20) $$r_{sim}=50\times min(x_{pos}-{\hat{x}}_{pos},L_{max})$$where x_pos denotes the position in the forward direction at the current time-step while ${\hat{x}}_{pos}$ denotes the position at the last time-step. L_max is introduced as the ceiling of the reward to guarantee stability. The velocity reward for the task of tracking specific velocity is calculated as shown in Eq. (21a):

(21a) $$r_{sim}=exp[-{10}^5\times (x_{pos}-{\hat{x}}_{pos}-x_t{)}^2]$$

(21b) $$x_t=v_t\times t$$where x_pos and ${\hat{x}}_{pos}$ have the same meaning as Eq. (20). x_t is the expected forward distance in one time step, which can be converted from target velocity v_t by Eq. (21b), where t denotes the length of one time step.

The deviation penalty is implemented to guarantee the accuracy of direction. When the robot deviates from the target direction, the total reward will decrease. The deviation penalty is defined as follows:

(22) $$Los{s}_y=\{\begin{array}{c}0,-20\le pitch\le 20,\hfill \\ -1,else.\hfill \end{array}$$

The deviation is not significant when the pitch of the robot is in the range from −20 to 20, thus the penalty is zero. When the absolute value exceeds this range, the penalty is set to −1 to command the robot to adjust the direction. The total environment reward is the sum of velocity reward and deviation penalty:

(23) $$R_T=r_{sim}+Los{s}_y$$

Imitation reward

Imitation reward R_I encourages the robot to imitate the reference motion in each phase of the gait cycle. In this case, the policy is constrained by the prior information obtained from the reference, which significantly improves the training efficiency. Imitation is especially important at the beginning of the training process where little knowledge is accumulated and the reward is almost zero in every training epoch. The imitation reward is composed of imitation rewards for joints, the center of mass, and Euler angles:

(24) $$R_I=w_I^t{r}_I^t+w_I^c{r}_I^c+w_I^o{r}_I^o$$where w^t_I, w^c_I and w^o_I are separate weights of the three parts of the reward. Based on sufficient experiments, weights are set as follows: w^t_I = 0.6, w^c_I = 0.3 , w^o_I = 0.1.

The reward for joints is calculated as shown in Eq. (25), where ${\hat{q}}_t^i$ and qⁱ_t respectively stand for the angles of the i-th joint of the reference and the simulation environment, w_it is the weight indicating the importance of the i-th joint. Different joints contribute differently to the walk. For instance, the joints of ankles are more important than joints of shoulders. Thus the weights of left and right shoulder are set to 0.05 while other joints of legs are 0.09, and add up to 1.

(25) $$r_I^t=exp[-{\sum }_i{w}_{it}(||{\hat{q}}_t^i-q_t^i|{|}^2)]$$

The reward for the center of mass is calculated as shown in Eq. (26), where ${\hat{q}}_c^j$ and q^j_c denotes the reference and simulation position of the center of mass at the current time-step.

(26) $$r_I^c=exp[-{\sum }_j(||{\hat{q}}_c^j-q_c^j|{|}^2)]$$

The reward for the Euler angle is calculated as shown in Eq. (27), where ${\hat{q}}_o^j$ and q^j_o denotes the reference and simulation Euler angle at the current time-step.

(27) $$r_I^o=exp[-{\sum }_k(||{\hat{q}}_o^k-q_o^k|{|}^2)]$$

Random state initialization

The environment state is initialized at the beginning of every episode. Traditional RL training is based on fixed state initialization (FSI), resetting the environment to a fixed beginning state. In such methods, agent will learn the policy serially. For example, in the walking task, the robot will learn to stand stable at the beginning, then learn the next sub-task of walking. However, at the beginning of the training, the robot falls directly in most episodes, thus the latter sub-tasks cannot be trained, resulting in a low training efficiency. To make full use of the reference motion and improve the training efficiency, the random state initialization (RSI) is introduced into training. This method has been widely adopted (Peng et al., 2018; Nair et al., 2018). At the beginning of one episode, the environment is initialized as one random phase of the reference gait cycle. In this way, the agent learns earlier and later sub-tasks simultaneously with a much higher efficiency.

Symmetrization of actions and states

The pattern of gait is defined as the cycling switch of the swing leg and instance leg between the two legs. In the first half of the cycle, the left leg is the swing leg and in the second half of the cycle the right leg is the swing leg. The two halves are symmetrical, which provides useful prior knowledge for the training. The symmetrization of reference and the symmetrization of state and action are proposed to improve the training efficiency and performance.

Symmetrization of reference: The gait cycle of reference motion is not symmetric, decreasing the performance of the trained agent, as shown in Fig. 5, where the curves of AnkleR and AnkleL are not perfectly symmetrical. To achieve absolute symmetry, the phases 9–17 are derived based on the phases 0–8, leveraging the symmetry of the two halves of the gait cycle. The resulting gait cycle has a perfect mirror symmetry in corresponding phase pairs: 0–9, 1–10, and so on.

Symmetrization of state and action: The state space is composed of the phase information (1D), the simulation state (35D) and the reference state (35D). The robot only needs to learn the first half of the gait cycle due to the symmetry of the gait cycle, then the second half can be generated according to symmetry. We propose the compression of state space and action space to accelerate the training process based on this innovation as shown in Fig. 6. At phases 0–8, the states are unchanged and are sent to the RL agent directly, while states at phases 9–17 are processed into the symmetry state as follows: For the simulation state, mirror symmetry to the 20-dimension joint angles was obtained replacing the original ones. Other dimensions remained unchanged. Then the processed symmetry state was sent to the RL agent. The 20-dimension joint angles of the reference state is also replaced by the symmetric ones. Thus, all states indicating joint angles are symmetric. For states at symmetric phase pairs, such as phases 0 and 9, 1 and 10, the output action generated by the RL agent will be highly similar. At phase 0–8, the output of the RL agent will be sent to the adder directly. At phase 9–17, the output of the RL agent will be recovered by symmetrization firstly. Then the symmetry output will be sent into the adder. The final target angles for the joints are the sum of the output action and the reference angles for both the first and second halves of gait cycle.

Noise in training

Gaussian noise was introduced into the output of the RL agent in the training process. During training, the RL agent was trained to overcome uncertain noise. Thus the training efficiency slightly decreased while the robustness in various environments and resistance to the external force of the trained model were significantly improved.

Walking on plains

The position and orientation during the walking procedure are shown in Fig. 13, where the task was to walk as fast as possible. The robot achieved an average speed of 0.488 m/s, with 22.48 m forward distance in 1,440 time steps. The velocity performance achieved a 5.8 time improvement compared with the reference motion generated by the traditional controller. The absolute deviation was 0.31 m and the relative deviation was 1.4%. The relative deviation decreased 87.3% compared with the original reference motion, which is 11.0%. In addition, the max Euler angle deviation from the expected direction (euler_y in Fig. 13B) was 16.2°, however, the robot corrected the direction automatically and simultaneously to maintain accuracy while the reference motion kept the deviated direction with no adjustment. The comparison of velocity performance is shown in Table 6. The compared five works include: the official document of Darwin-op, Li, Li & Cui, 2016, Gil, Calvo & Sossa (2019), and two different algorithms in Xi & Chen (2020), in which the target platform is the same as or similar with our platform. Our work and the first two compared works are tested on the platform Darwin-op, while the last three works on the platform NAO. For Xi & Chen (2020), the highest velocity among various cases is selected as the absolute velocity on the plain. Normalization of the velocity is introduced to ensure the precision of the comparison, which divides the absolute velocity by the height of the controlled robot. The normalized velocity of our work is 2 times compared with the rated max velocity, and approximately 10 times compared with other algorithms. The result shows that our work fully explored the potential of the robot, getting better performance compared with traditional algorithms and other reinforcement learning algorithms combined with dynamics models. To our knowledge, we have achieved the highest velocity on the platform Darwin-op.

For the task of tracking specific speed, the comparison of the actual speed and target speed is shown in Fig. 14, where the equivalent target speed (forward distance in one time step) is x_t = 0.01 m in (a) and x_t = 0.005 m in (b). The tracking accuracy achieved above 95% in both conditions.

Walking on slopes

The sequence of walking uphill is shown in Fig. 15, where the slope was 0.1 rad. The task was tracking specific velocity, rather than walking as fast as possible, to ensure the stability of walking. The robot kept the velocity of 0.295 m/s, with an accuracy of 94.0% compared with the target velocity of 0.313 m/s. The directional error was adjusted simultaneously, with a max value of 25.4°. The sequence of walking downhill is shown in Fig. 16, where the slope was 0.1 rad and the task was tracking specific velocity. The robot maintained a velocity of 0.278 m/s, with an accuracy of 88.8% compared with the target velocity of 0.313 m/s. The slight drop in accuracy was caused by the high velocity when walking downhill, which results in greater difficulty to keep balance. The velocity performance outperforms the rated max velocity on plain though the velocity is limited to around 0.313 m/s during the training. The comparison of velocity performance on slopes is shown in Table 7, where the compared works come from Xi & Chen (2020). The velocity performance is more than 5 times compared with the reference work. The robot maintained the directional error within 18° in the tested 1,440 time steps, which is better than the original reference motion. In addition, the robot adjusted its pose before losing balance to maintain a stable walking gait, which is the benefit of RL based controller.

Uneven terrain

To further evaluate the robustness of the proposed method, an environment with the uneven ground was constructed in Webots as shown in Fig. 17. The walking sequence is shown in Fig. 18. For the task of walking as fast as possible, the robot achieved a max velocity of 0.453 m/s, with a max position deviation of 0.234 m, and a max direction deviation of 17.5°. The comparison of velocity performance is shown in Table 8. The three compared works come from Liu, Ning & Chen, 2018, Yi et al. (2016), and Morisawa et al. (2012). The velocity performance is largely improved compared with the reference works. For the task of tracking specific velocity, the robot maintained the velocity of 0.299 m/s, with an accuracy of 95.5% compared with the target velocity 0.313 m/s. The max position deviation was 0.385 m and the max directional deviation was 21.4°. The trained robot was able to adapt well to the uneven ground and complete both tasks with high quality.

Series of slopes

A series of slopes combining plain-uphill-plain-downhill-plain was designed to further validate the adaptability of the proposed method. The walking sequence is shown in Fig. 19. The robot automatically adjusted to the gait policy to adapt different slopes and perform better. For instance, the angles of AnkleR and AnkleL during walking are shown in Fig. 20. The angle of AnkleR is higher during walking uphill and lower downhill. Correspondingly, the angle of AnkleL is lower during walking uphill and higher during the downhill. Both two ankles lift higher when the robot is walking uphill and lift lower when downhill.

External force

In the validation of resistance to external force, an environment providing external force for eight time steps in every 200 time steps was constructed. The robot was able to walk stably with the external force up to 10N from forward and back, or up to 6N from the left of right. The Euler angle of the robot is shown in Fig. 21 and the linear velocity is shown in Fig. 22. The robot made adjustments to recover the sudden change of the Euler angle and linear velocity caused by a random external force. The robot achieved a velocity of 0.363 m/s, which is also improved significantly from the original reference motion. The comparison of performance with other works is shown in Table 9. The compared works come from Liu, Ning & Chen (2018) and Smaldone et al. (2019). Our work has a resistance to stronger external force from both x-axis and y-axis with a better velocity performance.

Variety of capable environments

Our method is capable in the widest range of environments including: plains, slopes, uneven terrains and walking with random disturbance of external force. The comparison of capable environments is shown in Table 10, where the compared works are those mentioned in the previous comparisons.