Introducing Neuromodulation in Deep Neural Networks to Learn Adaptive Behaviours

TL;DR

This paper introduces a neuromodulation-inspired architecture for deep neural networks that enhances their ability to learn and adapt to changing environments, inspired by biological cellular neuromodulation mechanisms.

Contribution

The paper proposes a novel neuromodulation-based neural network architecture specifically designed for adaptive behavior learning in dynamic environments.

Findings

Neuromodulation improves task adaptation in navigation benchmarks.

The approach outperforms state-of-the-art methods in meta-reinforcement learning.

Neuromodulation enhances the flexibility and robustness of artificial agents.

Abstract

Animals excel at adapting their intentions, attention, and actions to the environment, making them remarkably efficient at interacting with a rich, unpredictable and ever-changing external world, a property that intelligent machines currently lack. Such an adaptation property relies heavily on cellular neuromodulation, the biological mechanism that dynamically controls intrinsic properties of neurons and their response to external stimuli in a context-dependent manner. In this paper, we take inspiration from cellular neuromodulation to construct a new deep neural network architecture that is specifically designed to learn adaptive behaviours. The network adaptation capabilities are tested on navigation benchmarks in a meta-reinforcement learning context and compared with state-of-the-art approaches. Results show that neuromodulation is capable of adapting an agent to different tasks and that neuromodulation-based approaches provide a promising way of improving adaptation of artificial systems.

Tables

Table 1. Table 1 : Value of the hyper-parameters that are kept constant for every benchmark in this paper.

$B$	$50$
$\lambda$	$0.98$
$\gamma$	$0.998$
${\beta }_0$	$1$
${\beta }_{min}$	$1/30$
${\beta }_{max}$	$30$
$d_{targ}$	$0.003$
$a_{l{r}_0}$	$2\ast {10}^{-4}$
${\omega }_1$	$0.9$
${\omega }_2$	$0.999$
$ϵ$	${10}^{-8}$
$e_{actor}$	$20$
$crb$	$2$
$cmb$	$25$
$c_{lr}$	$6\ast {10}^{-3}$
$T$	$200$
$e_{critic}$	$10$
$c_{v_0},c_{z_0},a_{v_0},a_{z_0}$	$0$
$\eta$	$50$

Equations

$$\mathcal{X}=[-3.0,3.0{]}^2$$

$$\mathcal{A}=\mathbb{R}$$

$$\mathbf{\alpha }[\mathbf{i}]\sim \mathbb{U}[-1.0,1.0],\forall i\in [1,2]$$

$$\mathbf{\alpha }[3]\sim \mathbb{U}[-\pi ,\pi [$$

$$u_0^k\sim \mathbb{U}[-1.5\ast \pi ,1.5\ast \pi ]\forall k\in \{x,y\}$$

$$x_0=[\mathbf{\alpha }[1]-a_0^x,\mathbf{\alpha }[2]-a_0^y]$$

$$(x,y)\in target\iff \sqrt{(x-\mathbf{\alpha }[1]{)}^2+(y-\mathbf{\alpha }[2]{)}^2}\le 0.4.$$

$$u_{t+1}^x=u_t^x+0.25\ast (\sin (a_t)+\sin (\mathbf{\alpha }[3]+n_t))$$

$$u_{t+1}^y=u_t^x+0.25\ast (\cos (a_t)+\cos (\mathbf{\alpha }[3]+n_t)).$$

$$u_{t+1}^k\leftarrow \{\begin{array}{ll}{u}_{t+1}^k-4& \text{if}{u}_{t+1}^k>2\\ u_{t+1}^k+4& \text{if}{u}_{t+1}^k<-2\\ u_{t+1}^k& \text{otherwise.}\end{array}$$

$$\rho (a_t,x_t,x_{t+1})=\{\begin{array}{ll}100& \text{if}(u_t^x,u_t^y)\in target\\ -2& \text{otherwise.}\end{array}$$

$$\mathcal{X}=[-2.5,2.5{]}^4$$

$$\mathcal{A}=\mathbb{R}$$

$$\mathbf{\alpha }[\mathbf{i}]\sim \mathbb{U}[-1.0,1.0],\forall i\in [1,2,3,4]$$

$$\mathbf{\alpha }[5]\sim \mathbb{U}\{-1,1\}$$

$$u_0^k\sim \mathbb{U}[-1.5,1.5]\forall k\in \{x,y\}.$$

$$x_0=[{\alpha }^{x_1}-u_0^x,{\alpha }^{y_1}-u_0^y,{\alpha }^{x_2}-u_0^x,{\alpha }^{y_2}-u_0^y].$$

$$\sqrt{(x-{\alpha }^{x_i}{)}^2+(y-{\alpha }^{y_i}{)}^2}\le 0.4.$$

$$u_{t+1}^x=u_t^x+\sin (a_t\ast \pi )\ast 0.25$$

$$u_{t+1}^y=u_t^y+\cos (a_t\ast \pi )\ast 0.25.$$

$$u_{t+1}^k\leftarrow \{\begin{array}{ll}{u}_{t+1}^k-4.0& \text{if}{u}_{t+1}^k>2\\ u_{t+1}^k+4.0& \text{if}{u}_{t+1}^k<-2\\ u_{t+1}^k& \text{otherwise.}\end{array}$$

$$\rho (a_t,x_t,x_{t+1})=\{\begin{array}{ll}100\ast \mathbf{\alpha }[5]& \text{if}(u_t^x,u_t^y)\in targe{t}_1\wedge (u_t^x,u_t^y)\notin targe{t}_2\\ -50\ast \mathbf{\alpha }[5]& \text{if}(u_t^x,u_t^y)\in targe{t}_2\wedge (u_t^x,u_t^y)\notin targe{t}_1\\ 0& \text{if}(u_t^x,u_t^y)\notin targe{t}_1\wedge (u_t^x,u_t^y)\notin targe{t}_2.\end{array}$$

$$\rho (a_t,x_t,x_{t+1})=\{\begin{array}{ll}100\ast \mathbf{\alpha }[5]& \text{if}\sqrt{(u_t^x-p^{x_1}{)}^2+(u_t^y-p^{y_1}{)}^2}\le \sqrt{(u_t^x-p^{x_2}{)}^2+(u_t^y-p^{y_2}{)}^2}\\ -50\ast \mathbf{\alpha }[5]& \text{otherwise.}\end{array}$$

$$R_{\mathcal{M}}^{{\pi }_{\theta }}=\underset{T\to \infty }{\lim }\sum _{t=0}^T{\gamma }^t{r}_t$$

$$\underset{a_t\sim {\pi }_{\theta }\forall t}{\underset{x_0\sim P_{x_0}(\cdot )}{\mathbb{E}}}{R}_{\mathcal{M}}^{{\pi }_{\theta }}$$

$${\pi }_{\mathcal{M}}^{\ast }\in \underset{\pi \in \Pi }{\mathrm{argmax}}\underset{a_t\sim {\pi }_{\theta }\forall t}{\underset{x_0\sim P_{x_0}(\cdot )}{\mathbb{E}}}{R}_{\mathcal{M}}^{\pi }$$

$$J_{\mathcal{M}}^{{\pi }_{\theta }}(h_t)=\sum _{j=t}^{\infty }{\gamma }^{j-t}{\rho }_{\mathcal{M}}(x_j,a_j\sim {\pi }_{\theta }(h_j),x_{j+1})$$

$$V^{\pi }=\underset{\mathcal{M}\sim \mathcal{D}}{\underset{a_t\sim {\pi }_{\theta }\forall t}{\underset{x_0\sim P_{x_0}(\cdot )}{\mathbb{E}}}}.$$

$$T{D}_i=r_i+\gamma \ast c_{{\psi }_k}(h_{i+1})-c_{{\psi }_k}(h_i),\forall i\in [0,\dots ,L]$$

$$GA{E}_j^i=\sum _{t=j}^L(\gamma \ast \lambda {)}^{t-j}\ast T{D}_i^j,\forall j\in [1,\dots ,E],i\in [0,\dots ,L^{\prime }]$$

Full text

1 Introduction

The field of machine learning has seen tremendous progress made during the past decade, predominantly owing to the improvement of deep neural network (DNN) algorithms. DNNs are networks of artificial neurons whose interconnections are tuned to reach a specific goal through the use of an optimization algorithm, mimicking the role of synaptic plasticity in biological learning. This approach has led to the emergence of highly efficient algorithms that are capable of learning and solving complex problems. Despite these tremendous successes, it remains difficult to learn models that generalise or adapt themselves efficiently to new, unforeseen problems based on past experiences. This calls for the development of novel architectures specifically designed to enhance adaptation capabilities of current DNNs.

In biological nervous systems, adaptation capabilities have long been linked to neuromodulation, a biological mechanism that acts in concert with synaptic plasticity to tune neural network functional properties. In particular, cellular neuromodulation provides the ability to continuously tune neuron input/output behaviour to shape their response to external stimuli in different contexts, generally in response to an external signal carried by biochemicals called neuromodulators [1, 2]. Neuromodulation regulates many critical nervous system properties that cannot be achieved solely through synaptic plasticity [3, 4]. It has been shown as being critical to the adaptive control of continuous behaviours, such as in motor control, among others [3, 4]. In this paper, we introduce a new neural architecture specifically designed for DNNs and inspired from cellular neuromodulation, which we call NMN, standing for “Neuro-Modulated Network”.

At its core, the NMN architecture comprises two neural networks: a main network and a neuromodulatory network. The main network is a feed-forward DNN composed of neurons equipped with a parametric activation function whose parameters are the targets of neuromodulation. It allows the main network to be adapted to new unforeseen problems. The neuromodulatory network, on the other hand, dynamically controls the neuronal properties of the main network via the parameters of its activation functions. Both networks have different inputs: the neuromodulatory network processes feedback and contextual data whereas the main network is in charge of processing other inputs.

Our proposed architecture can be related to previous works on different aspects. In [5], the authors take inspiration from Hebbian plasticity to build networks with plastic weights, allowing them to tune their weights dynamically. In [6] the same authors extend their work by learning a neuromodulatory signal that dictates which and when connections should be plastic. Our architecture is also related to hypernetworks [7], in which a network’s weights are computed through another network. Finally, other recent works focused on learning fixed activation functions [8, 9].

2 NMN architecture

The NMN architecture revolves around the neuromodulatory interaction between the neuromodulatory and main networks. We mimic biological cellular neuromodulation [10] in a DNN by assigning the neuromodulatory network the task to tune the slope and bias of the main network activation functions.

Let $\sigma(x):\mathbb{R}\rightarrow\mathbb{R}$ denote any activation function and its neuromodulatory capable version $\sigma_{\text{NMN}}(x,\mathbf{z};\mathbf{w}_{s},\mathbf{w}_{b})=\sigma\left(\mathbf{z}^{T}(x\mathbf{w}_{s}+\mathbf{w}_{b})\right)$ where $\mathbf{z}\in\mathbb{R}^{k}$ is a neuromodulatory signal and $\mathbf{w}_{s},\mathbf{w}_{b}\in\mathbb{R}^{k}$ are two parameter vectors of the activation function, respectively governing a scale factor and an offset. In this work, we propose to replace all the main network neuron activation functions with their neuromodulatory capable counterparts. The neuromodulatory signal $\mathbf{z}$, where size $k$ is a free parameter, is shared for all these neurons and computed by the neuromodulatory network as $\mathbf{z}=f(\mathbf{c})$, where $\mathbf{c}$ is a vector representing contextual and feedback inputs. The function $f$ can be any DNN taking as input such vector $\mathbf{c}$. For instance, $\mathbf{c}$ may have a dynamic size (e.g. more information about the current task becomes available as time passes), in which case $f$ could be parameterised as a recurrent neural network (RNN) or a conditional neural process [11], enabling refinement of the neuromodulatory signal as more data becomes available. The complete NMN architecture and the change made to the activation functions are depicted in Figure 1.

Notably, the number of newly introduced parameters scales linearly with the number of neurons in the main network whereas it would scale linearly with the number of connections between neurons if the neuromodulatory network was affecting connection weights, as seen for instance in the context of hypernetworks [7]. Therefore our approach can be extended more easily to very large networks.

3 Experiments

3.1 Setting

A good biologically motivated framework to which the NMN can be applied and evaluated is meta-reinforcement learning (meta-RL), as defined in [12]. In contrast with classical reinforcement learning (RL), which is formalised as the interaction between an agent and an environment defined as a Markov decision process (MDP), the meta-RL setting resides in the sub-division of an MDP as a distribution $\mathcal{D}$ over simpler MDPs. Let $t$ denote the discrete time, $\mathbf{x}_{t}$ the state of the MDP at time $t$, $\mathbf{a}_{t}$ the action taken at time $t$ and $r_{t}$ the reward obtained at the subsequent time-step. At the beginning of a new episode $i$, a new element is drawn from $\mathcal{D}$ to define an MDP, referred to as $\mathcal{M}$, with which the meta-RL agent interacts for $T\in\mathbb{N}$ time-steps afterwards. The only information that the agent collects on $\mathcal{M}$ is through observing the states crossed and the rewards obtained at each time-step. We denote by $\mathbf{h}_{t}=\left[\mathbf{x}_{0},\mathbf{a}_{0},r_{0},\mathbf{x}_{1},\ldots,\mathbf{a}_{t-1},r_{t-1},\mathbf{x}_{t}\right]$ the history of the interaction with $\mathcal{M}$ up to time-step $t$. As in [12], the goal of the meta-learning agent is to maximise the expected value of the discounted sum of rewards it can obtain over all the time-steps and episodes.

3.2 Training

In [12], the authors tackle this meta-RL framework by using an advantage actor-critic (A2C) algorithm. This algorithm revolves around two distinct parametric functions: the actor and the critic. The actor represents the policy used to interact with the MDPs, while the critic is a function that rates the performance of the agent policy. All actor-critic algorithms follow an iterative procedure that consists of the three following steps.

Use the policy to interact with the environment and gather data. 2. 2.

Update the actor parameters using the critic ratings. 3. 3.

Update the critic parameters to better approximate a value function.

In [12], the authors chose to model the actor and the critic with RNNs, taking $\mathbf{h}_{t}$ as the input. In this work, we propose comparing the NMN architecture to standard RNN by modelling both the actor and the critic with NMN. To this end, we define the feedback and contextual inputs $\mathbf{c}$ (i.e. the neuromodulatory network inputs) as $\mathbf{h}_{t}\setminus\mathbf{x}_{t}$ while the main network input is defined as $\mathbf{x}_{t}$. Note that $\mathbf{h}_{t}$ grows as the agent interacts with $\mathcal{M}$, motivating the usage of a RNN as neuromodulatory network. A graphical comparison between both architectures is shown on Figure 2.

To be as similar as possible to the neuronal model proposed by [10], the main network is a fully-connected neural network built using saturated rectified linear unit (sReLU) activation functions $\sigma(x)=\min(1,\max(-1,x))$, except for the final layer (also neuromodulated), for which $\sigma(x)=x$. In Section 4, we also report results obtained with sigmoidal activation functions which are often appreciably inferior to those obtained with sReLUs, further encouraging their use.

We built our models such that both standard RNN and NMN architectures have the same number of recurrent layers/units and a relative difference between the numbers of parameters that is lower than $2\%$. Both models are trained using an A2C algorithm with generalized advantage estimation [13] and proximal policy updates [14]. Finally, no parameter is shared between the actor and the critic. We motivate this choice by noting that the neuromodulatory signal might need to be different for the actor and the critic. For completeness and reproducibility, we provide a formal description of the algorithms used as supplementary material (5.2). This material aims mainly to describe and discuss standard RL algorithms in the context of meta-RL and, to a lesser extent, it aims to provide full implementation details. We also provide the exact neural architectures used for each benchmark as supplementary material (5.3).

3.3 Benchmarks description

We carried out our experiments on three custom benchmarks: a simple toy problem and two navigation problems with sparse rewards. These benchmarks were built to evaluate our architecture in environments with continuous action spaces. For conciseness and clarity, we only provide a mathematical definition of the first benchmark. The two other benchmarks are briefly textually depicted and further details are available as supplementary material (see 5.1). Figures 3, 4 and 5 are a graphical representation of each of the benchmarks.

Benchmark 1.

We define the first benchmark (made of a 1-D state space and action space) through a random variable $\alpha$, informative enough to distinguish all different MDPs in $\mathcal{D}$. With this definition, $\alpha$ represents the current task and drawing $\alpha$ at the beginning of each episode amounts to sampling a new task in $\mathcal{D}$. At each time-step, the agent observes a biased version $x_{t}=p_{t}+\alpha$ of the exact position of a target $p_{t}$ belonging to the interval $[-5-\alpha,5-\alpha]$, with $\alpha\sim\mathbb{U}[-10,10]$. The agent outputs an action $a_{t}\in\left[-20,20\right]$ and receives a reward $r_{t}$ which is equal to $10$ if $|a_{t}-p_{t}|<1$ and $-|a_{t}-p_{t}|$ otherwise. In case of positive reward, $p_{t+1}$ is re-sampled uniformly in its domain, else $p_{t+1}=p_{t}$. This benchmark is represented on Figure 3.

Benchmark 2.

The second benchmark consists of navigating towards a target in a 2-D space with noisy movements. Similarly to the first benchmark, we can distinguish all different MDPs in $\mathcal{D}$ through a three-dimensional random vector of variables $\boldsymbol{\alpha}$. The target is placed at $(\boldsymbol{\alpha}[1],\boldsymbol{\alpha}[2])$ in the 2-D space. At each time-step, the agent observes its relative position to the target and outputs the direction of a move vector $\mathbf{m}_{t}$. A perturbation vector $\mathbf{w}_{t}$ is then sampled uniformly in a cone, whose main direction $\boldsymbol{\alpha}[3]\sim\mathbb{U}[-\pi,\pi[$, together with the target’s position, define the current task in $\mathcal{D}$. Finally the agent is moved following $\mathbf{m}_{t}+\mathbf{w}_{t}$ and receives a reward ($r_{t}=-0.2$). If the agent reaches the target, it instead receives a high reward ($r_{t}=100$) and is moved to a position sampled uniformly in the 2-D space. This benchmark is represented on Fig 4

Benchmark 3.

The third benchmark also involves navigating in a 2-D space, but which contains two targets. As for the two previous benchmarks, we distinguish all different MDPs in $\mathcal{D}$ through a five-dimensional random vector of variables $\boldsymbol{\alpha}$. The targets are placed at positions $(\boldsymbol{\alpha}[1],\boldsymbol{\alpha}[2])$ and $(\boldsymbol{\alpha}[3],\boldsymbol{\alpha}[4])$. At each time-step, the agent observes its relative position to the two targets and is moved along a direction given by its action. One target, defined by the task in $\mathcal{D}$ through $\boldsymbol{\alpha}[5]$, is attributed a positive reward ($100$) and the other a negative reward ($-50$). In other words, $\boldsymbol{\alpha}[5]$ is a Bernoulli variable that determines which target is attributed the positive reward and which is attributed the negative one. As for benchmark 2, once the agent reaches a target, it receives the corresponding reward and is moved to a position sampled uniformly in the 2-D space. This benchmark is represented on Figure 5.

4 Results

Learning.

From a learning perspective, a comparison of the sum of rewards obtained per episode by NMNs and RNNs on the three benchmarks is shown in Figure 6. Results show that, on average, NMNs learn faster (with respect to the number of episodes) and converge towards better policies than RNNs (i.e., higher rewards for the last episodes). It is worth mentioning that, NMNs show very stable results, with small variances over different random seeds, as opposed to RNNs. To put the performance of the NMN in perspective, we note that an optimal Bayesian policy would achieve an expected sum of rewards of $4679$ on benchmark 1 (see 5.5 for proof) whereas NMNs reach, after $20000$ episodes, an expected sum of rewards of $4534$. For this simple benchmark, NMNs manage to learn near-optimal Bayesian policies.

Adaptation.

From an adaptation perspective, Figure 7 shows the temporal evolution of the neuromodulatory signal $\mathbf{z}$ (part A), of the scale factor (for each neuron of a hidden layer, part B) and of the rewards (part C) obtained with respect to $\alpha$ for $1000$ episodes played on benchmark $1$. For small values of $t$, the agent has little information on the current task, leading to a non-optimal behaviour (as it can be seen from the low rewards). Of greatest interest, the signal $\mathbf{z}$ for the first time-steps exhibits little dependence on $\alpha$, highlighting the agent uncertainty on the current task and translating to noisy scale factors. Said otherwise, for small $t$, the agent learned to play a (nearly) task-independent strategy. As time passes, the agent gathers further information about the current task and approaches a near-optimal policy. This is reflected in the convergence of $\mathbf{z}$ (and thus scale factors) with a clear dependency on $\alpha$ and also in wider-spread values of $\mathbf{z}$. For a large value of $t$, $\mathbf{z}$ holding constant between time-steps shows that the neuromodulatory signal is almost state-independent and serves only for adaptation. We note that the value of $\mathbf{z}$ in each of its dimensions varies continuously with $\alpha$, meaning that for two similar tasks, the signal will converge towards similar values. Finally, it is interesting to look at the neurons scale factor variation with respect to $\alpha$ (B). Indeed, for some neurons, one can see that the scale factors vary between negative and positive values, effectively inverting the slope of the activation function. Furthermore, it is interesting to see that some neurons are inactive (scale factor almost equal to [math], leading to a constant activation function) for some values of $\alpha$.

For benchmark 2, let us first note that $\mathbf{z}$ seems to code exclusively for $\boldsymbol{\alpha}[3]$. Indeed, $\mathbf{z}$ converges slowly with time with respect to $\boldsymbol{\alpha}[3]$, whatever the value of $\boldsymbol{\alpha}[1]$ and $\boldsymbol{\alpha}[2]$ (Figure 8). This, could potentially be explained by the fact that one does not need the values of $\boldsymbol{\alpha}[1]$ and $\boldsymbol{\alpha}[2]$ to compute an optimal move. The graphs on Figure 8 are projected on the dimension $\boldsymbol{\alpha}[3]$, allowing the same analysis as for benchmark 1.

The results obtained for benchmark 2 (Figure 8) show similar characteristics. Indeed, despite the agent receiving only noisy information on $\boldsymbol{\alpha}[3]$ at each time-step (as perturbation vectors are sampled uniformly in a cone centered on $\boldsymbol{\alpha}[3]$), $\mathbf{z}$ quasi-converges slowly with time (part A). The value of $\mathbf{z}$ in each of its dimensions also varies continuously with $\boldsymbol{\alpha}[3]$ (as for the first benchmark) resulting also in continuous scale factors variations. This is clearly highlighted at time-step $100$ on Figure 8 where the scale factors of some neurons appear highly asymmetric, but with smooth variations with respect to $\boldsymbol{\alpha}[3]$. Finally, let us highlight that for this benchmark, the agent continues to adapt even when it is already performing well. Indeed, one can see that after $40$ time-steps the agent is already achieving good results (part C), even though $\mathbf{z}$ has not yet converged (part A), which is due to the stochasticity of the environment. Indeed, the agent only receives noised information on $\alpha$ and thus after $40$ time-steps it has gathered sufficient information to act well on the environment, but insufficient information to deduce a near-exact value of $\boldsymbol{\alpha}[3]$. This shows that the agent can perform well, even while it is still gathering relevant information on the current task.

It is harder to interpret the neuromodulatory signal for benchmark 3. In fact, for that benchmark, we show that the signal seems to code not only for the task in $\mathcal{D}$ but also for the state of the agent in some sense. As $\boldsymbol{\alpha}$ is five-dimensional, it would be very difficult to look at its impact on $\mathbf{z}$ as a whole. Rather, we fix the position of the two references in the 2-D space and look at the behaviour of $\mathbf{z}$ with respect to $\boldsymbol{\alpha}[5]$. In Figure 9 adaptation is clearly visible in the rewards obtained by the agent (part C) with very few negative rewards after $30$ time-steps. We note that for later time-steps, $\mathbf{z}$ tends to partially converge (A) and :

• •

some dimensions of $\mathbf{z}$ are constant with respect to $\boldsymbol{\alpha}[5]$, indicating that they might be coding for features related to $\boldsymbol{\alpha}[1,2,3,4]$.

• •

Some other dimensions are well correlated to $\boldsymbol{\alpha}[5]$, for which similar observations than for the two other benchmarks can be made. For example, one can see that some neurons have a very different scale factors for the two possible different values of $\boldsymbol{\alpha}[5]$ (B).

• •

The remaining dimensions do not converge at all, implying that these are not related to $\boldsymbol{\alpha}$, but rather to the state of the agent.

These results suggest that in this case, the neuromodulation network is used to code more complex information than simply that required to differentiate tasks, making $\mathbf{z}$ harder to interpret. Despite $\mathbf{z}$ not converging on some of its dimensions, we stress that freezing $\mathbf{z}$ after adaptation will not strongly degrade the agent’s performance. That is, the features coded in $\mathbf{z}$ that do not depend on $\boldsymbol{\alpha}$ are not critical to the performance of the agent. To illustrate this, we will analyse the behaviour of the agent within an episode when freezing and unfreezing the neuromodulation signal and when changing task. This behaviour is shown on Figure 10, for which:

1. (a)

Shows the behaviour of the agent when $\mathbf{z}$ is locked to its initial value. This plot thus shows the initial ”exploration” strategy used by the agent; that is, the strategy played by the agent when it has not gathered any information on the current task. 2. (b)

Shows the behaviour of the agent after unlocking $\mathbf{z}$, that is when the agent is able to adapt freely to the current task by updating $\mathbf{z}$ at each time-step. 3. (c)

Shows the behaviour of the agent when locking $\mathbf{z}$ at a random time-step after adaptation. $\mathbf{z}$ is thus fixed at a value which fits well the current task. As one can see, the agent continues to navigate towards the correct target. The performance is however a slightly degraded as the agent seems to lose some capacity to avoid the wrong target. This further suggests that, in this benchmark (as opposed to the two others), the neuromodulation signal does not only code for the current task but also for the current state, in some sense, that is hard to interpret. 4. (d)

Shows the same behaviour as in (c) as $\mathbf{z}$ is still locked to the same value, but the references are now switched. As there is no adaptation without updating $\mathbf{z}$; the agent is now always moving towards to wrong target. 5. (e)

Shows the behaviour of the agent when unlocking $\mathbf{z}$ once again. As one can see, the agent is now able to adapt correctly by updating $\mathbf{z}$ at each time-step, and thus it navigates towards the correct target once again.

Robustness study.

Even though results are quite promising for the NMN, it is interesting to see how it holds up with another type of activation function as well as analysing its robustness to different main networks’ architectures.

Sigmoid activation functions.

Figure 11 shows the comparison between having sigmoids as the main network’s activation function instead of sReLUs. As one can see, sigmoid activation functions lead to worse or equivalent results to sReLUs, be it for RNNs or NMNs. In particular, the NMN architecture seems more robust to the change of activation function as opposed to RNNs, as the difference between sReLUS and sigmoids is often far inferior for NMNs than RNNs (especially for benchmark 2).

Architecture impact.

Figure 12 shows the learning curve, on benchmark 1, for different main network architectures ([math], $1$ and $4$ hidden layers in the main network respectively). As one can see, RNNs can, in fact, reach NMNs’ performances for a given architecture (no hidden layer in this case), but seem relatively dependant on the architecture. On the contrary, NMNs seem surprisingly consistent with respect to the number of hidden layers composing the main network.

5 Conclusions

In this work, we adopt a high-level view of a nervous system mechanism called cellular neuromodulation to improve artificial neural networks adaptive capabilities. The results obtained on three meta-RL benchmark problems showed that this new architecture was able to perform better than classical RNN. The work reported in this paper could be extended along several lines.

First, it would make sense to explore other types of machine-learning problems where adaptation is required. Supervised meta-learning would be an interesting track to follow.

Second, research work could also be carried out to further improve the NMN introduced here. For instance, one could introduce new types of parametric activation functions which are not linear, or even spiking neurons. It would also be of interest to look at sharing activation function parameters per layer.

Furthermore, it would be a logical progression to tackle other benchmarks to see if the observations made here hold true. More generally, analysing the neuromodulatory signal further in depth (and its impact on activation functions) with respect to different more complex tasks would be worth-while.

Finally, let us emphasize that even if the results obtained by our NMN are good and also rather robust with respect to a large choice of parameters, further research is certainly still needed to better characterise NMNs performances.

Supporting information

5.1 Detailed description of benchmark 2 and 3

Before defining the three benchmark problems, let us remind that for each benchmark, the MDPs that belong to the support of $\mathcal{D}$, which generates the different tasks (see Section 3.1), have transition probabilities and reward functions that differ only according to the value of a scalar $\alpha$. Drawing an MDP according to $\mathcal{D}$ will amount for all the benchmark problems to draw a value of $\alpha$ according to a probability distribution $P_{\alpha}(\cdot)$ and to determine the transition function and the reward function that correspond to this value. Let us also denote by $\mathcal{X}$ and $\mathcal{A}$ the state and action spaces respectively.

5.1.1 Benchmark 2

State space and action space:

[TABLE]

[TABLE]

Probability distribution of $\alpha$:

[TABLE]

[TABLE]

where $\mathbb{U}[a,b]$ stands for a uniform distribution between $a$ and $b$.

Initial state distribution:

The initial state $x_{0}$ is drawn through $2$ auxiliary random variables that represent the $x$ and $y$ initial coordinates of the agent and are denoted $u^{x}_{0},u^{y}_{0}$. At the beginning of an episode, those variables are drawn as follows:

[TABLE]

From those four auxiliary variables, we define $x_{0}$ as:

[TABLE]

The distribution $P_{x_{0}}(\cdot)$ is thus fully given by the distributions over the auxiliary variables.

Transition function:

Fist, let $target$ be the set of points $(x,y)\in\mathbb{R}^{2}$ such that

[TABLE]

When taking action $a_{t}$ in state $x_{t}$ drawing the state $x_{t+1}$ from the transition function amounts to first compute $u^{x}_{t+1}$ and $u^{y}_{t+1}$ according to the following procedure:

If $(u^{x}_{t},u^{y}_{t})\in target$ then $u^{k}_{t+1}\sim\mathbb{U}[-1.5,1.5]~{}\forall k\in\{x,y\}\quad.$ 2. 2.

If the preceding condition is not met, an auxiliary variable $n_{t}\sim\mathbb{U}[\frac{-\pi}{4},\frac{\pi}{4}]$ is drawn to compute $u^{x}_{t+1}$ and $u^{y}_{t+1}$ through the following sub-procedure:

1. (a)

Step one:

[TABLE]

[TABLE]

One can see that taking an action $a_{t}$ moves the agent in a direction which is the vectoral sum of the intended move $\mathbf{m}_{t}$ of direction $a_{t}$ and of a perturbation vector $\mathbf{p}_{t}$ of direction $\alpha+n_{t}$ sampled through the distribution over $n_{t}$. 2. (b)

Step two: In the case where the coordinates computed by step one lay outside $S[-2;2]^{2}$, they are corrected so as to model the fact that when the agent reaches an edge of the 2D space, it is moved to the opposite edge from which it continues its move. More specifically, $\forall k\in\{x,y\}$:

[TABLE]

Once $u^{x}_{t+1}$ and $u^{y}_{t+1}$ have been computed, $x_{t+1}$ is set equal to $[\boldsymbol{\alpha[1]}-u^{x}_{t+1},\boldsymbol{\alpha[2]}-u^{y}_{t+1}]$.

Reward function:

The reward function can be expressed as follows:

[TABLE]

5.1.2 Benchmark 3

State space and action space:

[TABLE]

[TABLE]

Probability distribution of $\alpha$:

[TABLE]

[TABLE]

Note that $\boldsymbol{\alpha}[1,2,3,4]$ define the 2-D positions of two targets. For clarity, we will refer to these values respectively by $\alpha^{x_{1}}$,$\alpha^{y_{1}}$,$\alpha^{x_{2}}$ and $\alpha^{y_{2}}$.

Initial state distribution:

The initial state $x_{0}$ is drawn through two auxiliary random variables that represent the $x$ and $y$ initial coordinates of the agent and are denoted $u^{x}_{0},u^{y}_{0}$. At the beginning of an episode, those variables are drawn as follows:

[TABLE]

From those six auxiliary variables, we define $x_{0}$ as:

[TABLE]

Transition function:

For all $i\in\{1,2\}$ let $target_{i}$ be the set of points $(x,y)\in\mathbb{R}^{2}$ such that

[TABLE]

. When taking action $a_{t}$ in state $x_{t}$, drawing the state $x_{t+1}$ from the transition function amounts to first compute $u^{x}_{t+1}$ and $u^{y}_{t+1}$ according to the following procedure:

If $\exists i\in\{1,2\}:(u^{x}_{t},u^{y}_{t})\in target_{i}$, which means that the agent is in one of the two targets, then $u^{k}_{t+1}\sim\mathbb{U}[-1.5,1.5]~{}\forall k\in\{x,y\}$ 2. 2.

If the preceding condition is not met, $u^{x}_{t+1}$ and $u^{y}_{t+1}$ are computed by the following sub-procedure:

1. (a)

Step one:

[TABLE]

[TABLE]

This step moves the agent in the direction it has chosen. 2. (b)

Step two: In the case where the coordinates computed by step one lay outside $[-2;2]^{2}$, they are corrected so as to model the fact that when the agent reaches an edge of the 2D space, it is moved to the opposite edge from which it continues its move. More specifically, $\forall k\in\{x,y\}$:

[TABLE]

Once $u^{x}_{t+1}$ and $u^{y}_{t+1}$ have been computed, $x_{t+1}$ is set equal to $[\alpha^{x_{1}}-u^{x}_{t+1},\alpha^{y_{1}}-u^{y}_{t+1},\alpha^{x_{2}}-u^{x}_{t+1},\alpha^{y_{2}}-u^{y}_{t+1}]$.

Reward function:

In the case where $(u^{x}_{t},u^{y}_{t})$ either belongs to only $target_{1}$, only $target_{2}$ or none of them, the reward function can be expressed as follows:

[TABLE]

In the case where $(u^{x}_{t},u^{y}_{t})$ belongs to both $target_{1}$ and $target_{2}$, that is $(u^{x}_{t},u^{y}_{t})\in target_{1}\wedge(u^{x}_{t},u^{y}_{t})\in target_{2}$, the reward function can be expressed as follows:

[TABLE]

That is, we consider that the agent belongs to the target to which it is closer to the centre.

5.2 Advantage actor-critic with generalized advantage estimation

In our meta-RL setting, both the actor and the critic are parametric functions that are defined on the trajectories’ histories. With $\theta\in\Theta$ and $\psi\in\Psi$ the parameters of the actor and critic ($\Theta$ and $\Psi$ are the actor and critic parameters spaces), respectively, we define $\pi_{\theta}$ and $c_{\psi}$ as the policy and critic functions. Let $\pi_{\theta_{k}}$ and $c_{\psi_{k}}$ be the models for the policy and the critic after $k$ updates of the parameters $\theta$ and $\psi$, respectively. To update from $\theta_{k}$ to $\theta_{k+1}$ and $\psi_{k}$ to $\psi_{k+1}$, the actor-critic algorithm uses the policy $\pi_{\theta_{k}}$ to select actions during $B$ MDPs drawn sequentially from $\mathcal{D}$, where $B\in\mathbb{N}_{0}$ is a parameter of the actor-critic approach. This interaction between the actor-critic algorithm and the meta-RL problem is presented in a tabular version in Algorithm 1 (5.2.1).

Using the $L\in\mathbb{N}_{0}$ first elements of each trajectory generated from the interaction with the $B$ MDPs and the values of $\theta_{k}$ and $\psi_{k}$, the algorithm computes $\theta_{k+1}$ and $\psi_{k+1}$. To this end, the algorithm exploits the set $[h_{B*k,L},\ldots,h_{B*(k+1)-1,L}]$, which we denote as $H_{k}$. Note that we use a replay buffer for updating $\psi$, thus for this update we also use several previous sets $H_{k-1}$, $H_{k-2}$, etc… A tabular version of the algorithm that details how MDPs are drawn and played, as well as how the set $H_{k}$ is built, is presented in Algorithm 2 of Appendix 5.2.1. Let $R^{\pi_{\theta}}_{\mathcal{M}}$ denote the sum of discounted rewards obtained when playing policy $\pi_{\theta}$ on task $\mathcal{M}$. That is,

[TABLE]

where $r_{t}$ are the rewards gathered at each time-step. To have a properly performing actor-critic algorithm, the value chosen for $L$ has to be chosen sufficiently large to produce an accurate estimation of the returns $R^{\pi_{\theta_{k}}}_{\mathcal{M}_{i}}~{}\forall i\in[B*k,\ldots,B*(k+1)-1]$ obtained by the policy $\pi_{\theta_{k}}$.

When used in a classical RL setting, an AC algorithm should interact with its environment to find the value of $\theta$ that leads to high values of the expected return given a probability distribution over the initial states. This expected return is written as:

[TABLE]

where $\mathcal{M}$ denotes the Markov Decision Process with which the AC algorithm interacts. When working well, actor critic algorithms produce a sequence of policies $\pi_{\theta_{1}}$, $\pi_{\theta_{2}}$, $\pi_{\theta_{3}}$, $\ldots$ whose expected returns increase as the iterative process evolves and eventually reaches values close to those obtained by $\pi_{\theta^{*}_{\mathcal{M}}}$ with $\theta^{*}_{\mathcal{M}}=\underset{\theta\in\Theta}{\arg\max}\underset{a_{t}\sim\pi_{\theta}\forall t}{\underset{x_{0}\sim P_{x_{0}}(\cdot)}{\mathbb{E}}}R^{\pi_{\theta}}_{\mathcal{M}}$, which, if $\pi^{\theta}$ is flexible enough, are themselves close to those obtained by an optimal policy $\pi^{*}_{\mathcal{M}}$ defined as:

[TABLE]

where $\Pi$ is the set of all admissible policies.

Let $h_{t}=\{x_{0},a_{0},r_{0},\ldots,x_{t}\}$ be a trajectory generated by policy $\pi_{\theta}$ on $\mathcal{M}$ and let $J_{\mathcal{M}}^{\pi_{\theta}}(h_{t})$ be the expected sum of discounted rewards that can be obtained while starting from $h_{t}$ and playing the policy $\pi_{\theta}$ in this environment, that is:

[TABLE]

where $\rho_{\mathcal{M}}(x_{j},a_{j},x_{j+1})$ is the reward function of task $\mathcal{M}$. In a classical RL setting, and again for an efficient AC algorithm, the value of the critic for $h_{t}$, $c_{\psi}(h_{t})$, also converges to $J^{\pi_{\theta^{*}_{\mathcal{M}}}}_{\mathcal{M}}(h_{t})$. We also note that in such a setting, the critic is updated at iteration $k+1$ in a direction that provides a better approximation of $J^{\pi_{\theta_{k}}}_{\mathcal{M}}(\cdot)$. Now, let us go back to our meta-RL problem and let $V^{\pi}$ denote the expected sum of returns that policy ${\pi}$ can obtain on this problem:

[TABLE]

Let $\theta^{*}\in\underset{\theta\in\Theta}{\arg\max}V^{\pi_{\theta}}$. When interacting with our meta-RL problem, a performant AC algorithm should, in principle, converge towards a policy $\pi_{\hat{\theta}^{*}}$, leading to a value of $V^{\pi_{\hat{\theta}^{*}}}$ close to $V^{\pi_{\theta}^{*}}$ that is itself close to $\underset{\pi\in\Pi}{\max}V^{\pi}$. A policy $\pi^{*}$ such that $\pi^{*}\in\underset{\pi\in\Pi}{\arg\max}V^{\pi}$ is called a Bayes optimal policy in a Bayesian RL setting where the distribution $\mathcal{D}$ is assumed to be known. If we are working with policies that are, indeed, able to quickly adapt to the environment, we may also expect that the policy $\pi_{\hat{\theta}^{*}}$ learned by the algorithm is such that, when applied on an $\mathcal{M}$ belonging to the support of $\mathcal{D}$, it leads to a value of $J^{\pi_{\hat{\theta}^{*}}}_{\mathcal{M}}(h_{t})$ close to $\underset{\pi\in\Pi}{\max}J^{\pi}_{\mathcal{M}}(h_{t})$ as $t$ increases. In other words, once the agent has gathered enough information to adapt to the current MDP, it should start behaving (almost) optimally. This is the essence of meta-RL.

We may also expect that, in such case, the value of the critic for $h_{t}$ when the budget is exhausted closely estimates the expected value of the future discounted rewards that can be obtained when using policy $\pi^{\hat{\theta}^{*}}$ and after having already observed a trajectory $h_{t}$. Therefore, we may also expect that once the episode budget is exhausted, $c_{\psi}(h_{t})$:

will be close to $\underset{\mathcal{M}\sim\mathcal{D}}{\mathbb{E}}J^{\pi_{\hat{\theta}^{*}}}_{\mathcal{M}}(h_{t})\simeq\underset{\mathcal{M}\sim\mathcal{D}}{\mathbb{E}}\underset{\pi\in\Pi}{\max}J^{\pi}_{\mathcal{M}}(h_{t})$ if $h_{t}=\{x_{0}\}$; 2. 2.

will, as $t$ increases, tend to get closer to $\underset{\pi\in\Pi}{\max}J^{\pi}_{\mathcal{M}}(h_{t})\simeq J^{\pi_{\hat{\theta}^{*}}}_{\mathcal{M}}(h_{t})$ where $\mathcal{M}$ can be any environment belonging to the support of $\mathcal{D}$ used to generate $h_{t}$.

Existing actor-critic algorithms mainly differ from each other by the way the actor and critic are updated. While in early actor-critic algorithms the critic was directly used to compute the direction of update for the actor’s parameters (see for example the REINFORCE policy updates [15]), now it is more common to use an advantage function. This function represents the advantage in terms of return of selecting specific actions given a trajectory history (or simply a state when AC algorithms are used in a standard setting) over selecting them following the policy used to generate the trajectories. Here, we use generalised advantage estimations (GAE), as introduced in [13]. More recently, it has been shown that avoiding too large policy changes between updates can greatly improve learning ([16], [14]). Therefore, while in classical AC algorithms the function used to update the actor aims at representing directly the gradient of the actor’s return with respect to its parameters, we rather update the actor’s parameters $\theta$ by minimising a loss function that represents a surrogate objective. We have selected as surrogate function one that is similar to the one introduced in [14] with an additional loss term that proved to improve (albeit slightly) the performances of PPO in all cases.

As our actor and critic are modelled by differentiable functions, they are both updated through gradient descent. We now proceed to explain the losses use to compute the gradient for both the actor and the critic.

Actor update

First, we define the temporal error difference term for any two consecutive time-steps of any trajectory:

[TABLE]

where $\psi_{k}$ denotes the critic’s parameters for playing the given trajectory. This temporal difference term represents, in some sense, the (immediate) advantage obtained, after having played action $a_{j}$ over what was expected by the critic. If $c_{\psi_{k}}(\cdot)$ was the true estimate of $J^{\pi_{\theta_{k}}}_{\mathcal{M}}(\cdot)$ and if the policy played was $\pi_{\theta_{k}}$, the expected value of these temporal differences would be equal to zero. We now define the GAE’s terms that will be used later in our loss functions:

[TABLE]

where $\lambda\in[0,1]$ is a discount factor used for computing GAEs, $TD^{j}_{i}$ is the value of $TD_{i}$ for trajectory $j$, and where $L^{\prime}$ is another hyper-parameter of the algorithm, chosen in combination with $L$ in order to have a value of $GAE_{i,j}$ that accurately approximates $\sum_{t=j}^{\infty}(\gamma*\lambda)^{k-j}*TD^{j}_{i}$ $\forall i,j$. Note that the value chosen for $L^{\prime}$ also has to be sufficiently large to provide the loss function with a sufficient number of GAE terms. These GAE terms, introduced in [13], represent the exponential average of the discounted future advantages observed. Thanks to the fact that GAE terms can catch the accumulated advantages of a sequence of actions rather than of a single action, as it is the case with the temporal difference terms, they can better represent the advantage of the new policy played by the AC algorithm over the old one (in terms of future discounted rewards).

In the loss function, we will actually not use the advantage terms as defined by Equation 5, but normalised versions in order to have advantages that remain in a similar range regardless of rewards magnitude. Thanks to this normalisation, the policy learning rate does not have to be tuned according to the loss magnitude. However, this normalisation does not mask actions that have led to higher or lower returns than average. We normalize as follows $\forall k\in[1,\ldots,E]$ with $E$ the number of actor and critic updates that have been carried:

[TABLE]

[TABLE]

[TABLE]

where \mathop{\mathchoice{{\vphantom{\leavevmode\hbox{\set@color\displaystyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\displaystyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\textstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\textstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=8.00009pt}\cr\scriptstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptscriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=6.00006pt}\cr\scriptscriptstyle\sum\cr}}}}\displaylimits is the symbol we use to represent the average sum operator (i.e. \mathop{\mathchoice{{\vphantom{\leavevmode\hbox{\set@color\displaystyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\displaystyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\textstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\textstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=8.00009pt}\cr\scriptstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptscriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=6.00006pt}\cr\scriptscriptstyle\sum\cr}}}}\displaylimits_{x=1}^{m}f(x)=\sum_{x=1}^{m}\frac{f(x)}{m}). To define the loss functions used to compute $\theta_{k+1}$ and $\psi_{k+1}$, only the GAE terms corresponding to time-steps $[0,\ldots,L^{\prime}]$ of episodes $[B*k,\ldots,B*(k+1)-1]$ are computed. A tabular version of the algorithm used to compute these terms is given in Algorithm 3 of Appendix 5.2.1 111Although not explicitly written in the text for clarity, we use a normalisation technique when computing discounted sums for the AC algorithm update. In fact, when carrying an update of the AC algorithm, if rewards appear in discounted sums, they are multiplied by $(1-\gamma)$. This has for effect that the discounted sum values remain of the same magnitude regardless of $\gamma$. The implications of this normalization are two-fold. (i) The critic does not directly approximate $J^{\pi}_{\mathcal{M}}(\cdot)$ but rather $(1-\gamma)*J^{\pi}_{\mathcal{M}}(\cdot)$. (ii) Second, for the temporal differences to remain coherent with this normalisation, $r_{i}$ must also be multiplied by $(1-\gamma)$ when computing $TD_{i}$. Those two small changes are included in Algorithm 3..

Once advantages have been computed, the values of $\theta_{k+1}$ are computed using updates that are strongly related to PPO updates with a Kullback Leibler (KL) divergence implementation [14]. The loss used in PPO updates is composed of two terms: a classic policy gradient term and a penalisation term. Let us now present a standard policy gradient loss, note that from now on, we will refer to the value of $a_{t}$, $x_{t}$ and $h_{t}$ at episode $i$ by $a^{i}_{t}$, $x^{i}_{t}$ and $h^{i}_{t}$ respectively.

[TABLE]

where $\mathcal{B}_{k}$ is the set of all pairs $(i,t)$ for which $i,t\in([B*k,\ldots,B*(k+1)-1]*[0,\ldots,L^{\prime}])$, that is, the set containing the first $L^{\prime}$ time-steps of the $B$ trajectories played for iteration $k$ of the actor-critic algorithm.

One can easily become intuitive about Equation 6 as, given an history $h^{i}_{t}$, minimising this loss function tends to increase the probability of the policy taking actions leading to positive advantages (i.e. $GAE^{i}_{t}>0$) and decreases its probability to take actions leading to negative advantages (i.e. $GAE^{i}_{t}<0$). It has been found that to obtain good performances with this above-written loss function, it was important to have a policy that does not change too rapidly from one iteration to the other. Before explaining how this can be achieved, let us first give an explanation on why it may be important to have slow updates of the policy. Let us go back to the loss function given by Equation 6. Minimising this loss function will give a value for $\theta_{k+1}$ that will lead to higher probabilities of selecting actions corresponding to high values of the advantages $GAE^{i}_{t}$. A potential problem is that these advantages are not really related to the advantages of the would-be new policy $\pi_{\theta_{k+1}}$ over $\pi_{\theta_{k}}$ but are instead related to the advantages of policy $\pi_{\theta_{k}}$ over $\pi_{\theta_{k-1}}$. Indeed, the advantages $GAE^{i}_{t}$ are computed using the value function $c_{\psi_{k}}$, whose parameters have been updated from $\psi_{k-1}$ in order to better approximate the sum of discounted rewards obtained during the episodes $[B*(k-1),\ldots,B*k-1]$. It clearly appears that $\psi_{k}$ has, in fact, been updated to approximate discounted rewards obtained through the policy $\pi_{\theta_{k-1}}$ (used to play episodes for update $k-1$). A solution to this problem is to constraint the minimisation to reach a policy $\pi_{\theta_{k+1}}$ that does not stand too far from $\pi_{\theta_{k}}$. We may reasonably suppose that the advantage function used in (6) still correctly reflects the real advantage function of $\pi_{\theta_{k+1}}$ over $\pi_{\theta_{k}}$. To achieve this, we add a penalisation term $\mathcal{P}(\theta)$ to the loss function. In the PPO approach, the penalisation term is $\mathcal{P}_{ppo}(\theta)=\beta_{k}*d(\theta)$, where:

• i)

$\beta_{k}$ is an adaptive weight

• ii)

d(\theta)=\mathop{\mathchoice{{\vphantom{\leavevmode\hbox{\set@color\displaystyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\displaystyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\textstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=11.00012pt}\cr\textstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=8.00009pt}\cr\scriptstyle\sum\cr}}}{{\vphantom{\leavevmode\hbox{\set@color\scriptscriptstyle\sum}}\ooalign{\smash{\vrule height=1.0pt,depth=6.00006pt}\cr\scriptscriptstyle\sum\cr}}}}\displaylimits_{[i,t]\in\mathcal{B}_{k}}{}[KL(\pi_{\theta_{k}(.|h_{i,t})},\pi_{\theta(.|h_{i,t})})], where $KL$ is the Kullback-Leibler divergence, detailed later on. This term penalises policies that are too different from $\pi_{\theta_{k}}$.

We note that the $\beta_{k}$ dynamical updates use a hyper-parameter $d_{targ}\in\mathbb{N}_{0}$ called the divergence target. The update is done through the following procedure (note that, unlike updates of $\beta$ proposed in [14], we constrain $\beta$ to remain in the range $[\beta_{min},\beta_{max}]$; we explain later why):

[TABLE]

With this update strategy, the penalisation term will tend to evolve in a way such that the KL divergence between two successive policies does not tend to go beyond $d_{targ}$ without having to add an explicit constraint on $d$, as was the case in Trust Region Policy Optimization (TRPO) updates [16], which is more cumbersome to implement.

As suggested in [17], adding another penalisation term (squared hinge loss) to $\mathcal{P}_{PPO}$ to further penalise the KL divergence, in cases where it surpasses $2*d_{targ}$, improved algorithm performance. The final expression of the penalisation term is:

[TABLE]

where $\delta$ is a hyper-parameter that weights the third loss term. The loss function $\mathcal{L}_{policy}$ that we minimise as a surrogate objective becomes:

[TABLE]

We now detail how to compute the KL divergence. First, let us stress that we have chosen to work with multi-variate Gaussian policies for the actor. This choice is particularly well suited for MDPs with continuous action spaces. The approximation architecture of the actor will therefore not directly output an action, but the means and standard deviations of an m-dimensional multi-variate Gaussian from which the actor’s policy can be defined in a straightforward way. For each dimension, we bound the multi-variate Gaussian to the support, $\mathcal{U}$, by playing the action that is clipped to the bounds of $\mathcal{U}$ whenever the multi-variate Gaussian is sampled outside of $\mathcal{U}$. In the remaining of this paper, we will sometimes abusively use the terms ”output of the actor at time $t$ of episode $i$” to refer to the means vector $\mu_{i,t}^{\theta_{k}}$ and the standard deviations vector $\sigma_{i,t}^{\theta_{k}}$ that the actor uses to define its probabilistic policy at time-step $t$ of episode $i$. Note that we have chosen to work with a diagonal covariance matrix for the multi-variate Gaussian distribution. Its diagonal elements correspond to those of the vector $\sigma_{i,t}^{\theta_{k}}$. We can then compute the KL divergence in each pair $[i,t]$ following the well-established formula:

[TABLE]

where $\Sigma_{\theta_{k},i,t},\Sigma_{\theta,i,t}$ are the diagonal covariance matrices of the two multi-variate Gaussian distributions $\pi_{\theta_{k}}(\cdot|h^{i}_{t}),\pi_{\theta}(\cdot|h^{i}_{t})$ that can be derived from $\sigma_{i,t}^{\theta_{k}}$ and $\sigma_{i,t}^{\theta}$. The loss function $\mathcal{L}_{vanilla}$ can be expressed as a function of $\Sigma_{\theta_{k},i,t}$, $\Sigma_{\theta,i,t}$, $\mu^{\theta}_{i,t}$ and $\mu^{\theta_{k}}_{i,t}$ when working with a multi-variate Gaussian. To this end, we use the log-likelihood function $\ln{(\pi_{\theta}(a_{i,t}|h^{i}_{t}))}$, which gives the log-likelihood of having taken action $a^{i}_{t}$ given a trajectory history $h^{i}_{t}$. In the case of a multi-variate Gaussian, $\ln{(\pi_{\theta}(a^{i}_{t}|h^{i}_{t}))}$ is defined as:

[TABLE]

where $m$ is the dimension of the action space and where $|\Sigma_{\theta,i,t}|$ represents the determinant of the matrix. From this definition, one can rewrite $\mathcal{L}_{vanilla}$ as:

[TABLE]

By merging equation (11), (10) and equation (8), one gets a loss $\mathcal{L}_{policy}$ that depends only on $\Sigma_{\theta_{k},i,t}$, $\Sigma_{\theta,i,t}$, $\mu^{\theta}_{i,t}$ and $\mu^{\theta_{k}}_{i,t}$.

Critic update

The critic is updated at iteration $k$ in a way to better approximate the expected return obtained when following the policy $\pi_{\theta_{k}}$, starting from a given trajectory history. To this end, we use a mean-square error loss as a surrogate objective for optimizing $\psi$. First, we define $\hat{R}^{i}_{j}=\sum_{k=j}^{L}\gamma^{k-j}*r^{i}_{j}$ $\forall i,j\in[B*k,\ldots,B*(k+1)-1],[0,\ldots,L]$. From the definition of $\hat{R}^{i}_{j}$ we express the loss as:

[TABLE]

where (i) $crb\in\mathbb{N}_{0}$ is a hyper-parameter; (ii) $\mathcal{B}_{k-crb}$ is the set of all pairs $(i,t)$ for which $i,j\in([B*(k-crb),\ldots,B*(k+1)-1]*[0,\ldots,L^{\prime}])$. The set $\mathcal{B}_{k-crb}$ used in (12) contains all the pairs from the current trajectory batch and from the $crb$ previous trajectory batches. We call this a replay buffer whose length is controlled by $crb$ which stands for ”criticreplaybuffer”. Minimising $\mathcal{L}_{critic}$ does not lead to updates such that $c_{\psi}$ directly approximates the average expected return of the policy $\pi_{\theta_{k}}$. Rather, the updates are such that $c_{\psi}$ directly approximates the average expected return obtained by the last $crb+1$ policies played. We found out that using a replay buffer for the critic smoothed the critic’s updates and improved algorithm performances.

Note that the loss (12) is only computed on the $L^{\prime}<<L$ first time-steps of each episode, as was the case for the actor. The reason behind this choice is simple. The value function $c_{\psi_{k}}$ should approximate $R^{i}_{j}=\sum_{t=j}^{+\infty}\gamma^{t-j}*r^{i}_{t}$ for every $h^{i}_{j}$, where $R^{i}_{j}$ the infinite sum of discounted rewards that are attainable when ”starting” from $h^{i}_{j}$. However, this approximation can become less accurate when $j$ becomes close to $L$ since we can only guarantee $\hat{R}^{i}_{j}$ to stand in the interval: $[R^{i}_{j}-\frac{\gamma^{L-j}}{1-\gamma}R_{max},R^{i}_{j}-\frac{\gamma^{L-j}}{1-\gamma}R_{min}]$. Hence this choice of $L^{\prime}$.

Gradients computation and update

The full procedure is available as Tabular versions in 5.2.1. As a summary, we note that both the actor and critic are updated using the Adam procedure ([18]) and back-propagation through time ([19]). The main difference between both updates is that the actor is updated following a full-batch gradient descent paradigm, while the critic is updated following a mini-batch gradient descent paradigm.

5.2.1 Tabular version

5.3 Architecture details

For conciseness, let us denote by $f_{n}$ a hidden layer of $n$ neurons with activation functions $f$, by $\rightarrow$ a connection between two fully-connected layers and by $\multimap()$ a neuromodulatory connection (as described in Section 2).

Benchmark 1.

The architectures used for this benchmark were as follows:

• •

RNN : $GRU_{50}\rightarrow ReLU_{20}\rightarrow ReLU_{10}\rightarrow I_{1}$

• •

NMN : $GRU_{50}\rightarrow ReLU_{20}\multimap(SReLU_{10}\rightarrow I_{1})$

Benchmark 2 and 3.

The architectures used for benchmark 2 and 3 were the same and as follows:

• •

RNN : $GRU_{100}\rightarrow GRU_{75}\rightarrow ReLU_{45}\rightarrow ReLU_{30}\rightarrow ReLU_{10}\rightarrow I_{1}$

• •

NMN : $GRU_{100}\rightarrow GRU_{75}\rightarrow ReLU_{45}\multimap(ReLU_{30}\rightarrow ReLU_{10}\rightarrow I_{1})$

5.4 Hyper-parameter values

5.5 Bayes optimal policy for benchmark 1

A Bayes optimal policy is a policy that maximises the expected sum of rewards it obtains when playing an MDP drawn from a known distribution $\mathcal{D}$. That is, a Bayes optimal policy $\pi^{*}_{bayes}$ belongs to the following set:

[TABLE]

with $P_{\mathcal{M}}$ being the state-transition function of the MDP $\mathcal{M}$ and $R^{\pi}_{\mathcal{M}}$ the discounted sum of reward obtained when playing policy $\pi$ on $\mathcal{M}$.

In the first benchmark, the MDPs only differ by a bias, which we denote $\alpha$. Drawing an MDP according to $\mathcal{D}$ amounts to draw a value of $\alpha$ according to a uniform distribution of $\alpha$ over $[-\alpha_{max},\alpha_{max}]$, denoted by $\mathbb{U}_{\alpha}$, and to determine the transition function and the reward function that correspond to this value. Therefore, we can write the previous equation as:

[TABLE]

with $\mathcal{M}(\alpha)$ being a function giving as output the MDP corresponding to $\alpha$ and $\Pi$ the set of all possible policies.

We now prove the following theorem.

Theorem 1

The policy that selects:

at time-step $t=0$ the action $a_{0}=x_{0}+\frac{\gamma*(\alpha_{max}+4.5)}{1+\gamma}$ 2. 2.

at time-step $t=1$

1. a)

if $r_{0}=10$, the action $a_{1}=x_{1}+a_{0}-x_{0}$ 2. b)

else if $|r_{0}|>\alpha_{max}-(a_{0}-x_{0})\quad\wedge\quad a_{0}-x_{0}>0$, the action $a_{1}=a_{0}+r_{0}$ 3. c)

else if $|r_{0}|>\alpha_{max}-(x_{0}-a_{0})\quad\wedge\quad a_{0}-x_{0}<0$, the action $a_{1}=a_{0}-r_{0}$ 4. d)

and otherwise the action $a_{1}=a_{0}+r_{0}+1$ 3. 3.

for the remaining time-steps:

1. a)

if $r_{0}=10$, the action $a_{t}=x_{t}+a_{0}-x_{0}$ 2. b)

else if $r_{1}=10$, the action $a_{t}=x_{t}+a_{1}-x_{1}$ 3. c)

and otherwise the action $a_{t}=x_{t}+i_{t}$ where $i_{t}$ is the unique element of the set $\{a_{0}-x_{0}+r_{0};a_{0}-x_{0}-r_{0}\}\cap\{a_{1}-x_{1}+r_{1};a_{1}-x_{1}-r_{1}\}$

is Bayes optimal for benchmark 1 defined in SubSection 3.3.

**Proof ** Let us denote by $\pi^{*}_{theorem1}$ the policy described in this theorem. To prove this theorem, we first prove that in the set of all possible policies $\Pi$ there are no policy $\pi$ which leads to a higher value of

[TABLE]

than $\pi^{*}_{theorem1}$. Or equivalently:

[TABLE]

Afterwards, we prove that the policy $\pi^{*}_{theorem1}$, generates for each time-step $t\geq 2$ a reward equal to $R_{max}$ which is the maximum reward achievable, or written alternatively as:

[TABLE]

By merging (14) and (15), we have that

[TABLE]

which proves the theorem.

$\triangleright$ Part 1.

Let us now prove inequality (14). The first thing to notice is that for a policy to maximise expression (13), it only needs to satisfy two conditions for all $x_{0}$. The first one: to select an action $a_{1}$, which knowing the value of $(x_{0},a_{0},r_{0},x_{1})$, maximises the expected value of $r_{1}$. We denote by $V_{1}(x_{0},a_{0},r_{0},x_{1})$ the maximum expected value of $r_{1}$ that can be obtained knowing the value of $(x_{0},a_{0},r_{0},x_{1})$. The second one: to select an action $a_{0}$ knowing the value of $x_{0}$ that maximises the expected value of the sum $r_{0}+\gamma V_{1}(x_{0},a_{0},r_{0},x_{1})$. We now show that the policy $\pi_{theorem1}$ satisfies these two conditions.

Let us start with the first condition that we check by analysing four cases, which correspond to the four cases a), b), c), d) of policy $\pi_{theorem1}$ for time step $t=1$.

1. a)

If $r_{0}=10$, the maximum reward that can be obtained, we are in a context where $a_{0}$ belongs to the target interval. It is easy to see that, by playing $a_{1}=x_{1}+a_{0}-x_{0}$, we will obtain $r_{1}$ equal to 10. This shows that in case a) for time step $t=1$, $\pi_{theorem1}$ maximises this expected value of $r_{1}$. 2. b)

If $|r_{0}|>\alpha_{max}-(a_{0}-x_{0})\quad\wedge\quad a_{0}-x_{0}>0$ and $r_{0}\neq 10$ it is easy to see that the value of $\alpha$ to which the MDP corresponds can be inferred from $(x_{0},a_{0},r_{0})$ and that the action $a_{1}=a_{0}+r_{0}$ will fall in the middle of the target interval, leading to a reward of 10. Hence, in this case also, the policy $\pi_{theorem1}$ maximises the expected value of $r_{1}$. 3. c)

If $|r_{0}|>\alpha_{max}-(x_{0}-a_{0})\quad\wedge\quad a_{0}-x_{0}<0$ and $r_{0}\neq 10$, we are also in a context where the value of $\alpha$ can be inferred directly from $(x_{0},a_{0},r_{0})$ and the action $a_{1}=a_{0}-r_{0}$ targets the centre of the target interval, leading to a reward of $10$. Here again, $\pi_{theorem1}$ maximises the expected value of $r_{1}$. 4. d)

When none of the three previous conditions is satisfied, $a$ is not satisfied and so $x_{1}=x_{0}$, we need to consider two cases: $(a_{0}-x_{0})\geq 0$ and $(a_{0}-x_{0})<0$. Let us first start with $(a_{0}-x_{0})\geq 0$. In such a context, $\alpha\in\{a_{0}-x_{0}+r_{0};a_{0}-x_{0}-r_{0}\}=\{a_{0}-x_{0}-|a_{0}-x_{0}-\alpha|,a_{0}-x_{0}+|a_{0}-x_{0}-\alpha|\}$ and where:

$P(\alpha=a_{0}-x_{0}-|a_{0}-x_{0}-\alpha||x_{0},a_{0},r_{0},x_{1})=0.5$ 2. 2)

$P(\alpha=a_{0}-x_{0}+|a_{0}-x_{0}-\alpha||x_{0},a_{0},r_{0},x_{1})=0.5\quad.$

Let us now determine the action $a_{1}$ that maximises $\hat{r}_{1}$, the expected value of $r_{1}$ according to $P(\alpha|x_{0},a_{0},r_{0},x_{1})$. Five cases, represented on Figure 13, have to be considered:

$a_{1}<a_{0}-|a_{0}-x_{0}-\alpha|-1$. Here $\hat{r}_{1}=a_{1}-a_{0}$ and the maximum of $\hat{r}_{1}$ is equal to $-|a_{0}-x_{0}-\alpha|-1$. 2. 2)

$a_{1}\in[a_{0}-|a_{0}-x_{0}-\alpha|-1,a_{0}-|a_{0}-x_{0}-\alpha|+1]$. Here we have $\hat{r}_{1}=\frac{1}{2}(10+a_{0}-|a_{0}-x_{0}-\alpha|-a_{1})$ whose maximum over the interval is $5.5-|a_{0}-x_{0}-\alpha|$ which is reached for $a_{1}=a_{0}+|a_{0}-x_{0}-\alpha|-1$. 3. 3)

$a_{1}\in[a_{0}-|a_{0}-x_{0}-\alpha|+1,a_{0}+|a_{0}-x_{0}-\alpha|-1]$. In this case $\hat{r}_{1}=-|a_{0}-x_{0}-\alpha|$ and is independent from $a_{1}$. 4. 4)

$a_{1}\in[a_{0}+|a_{0}-x_{0}-\alpha|-1,a_{0}+|a_{0}-x_{0}-\alpha|+1]$. The expected reward is $\hat{r}_{1}=\frac{1}{2}(10+a_{0}-|a_{0}-x_{0}-\alpha|-a_{1})$ whose maximum over the interval is $5.5-|a_{0}-x_{0}-\alpha|$ which is reached for $a_{1}=a_{0}+|a_{0}-x_{0}-\alpha|+1$. 5. 5)

$a_{1}>a_{0}+|a_{0}-x_{0}-\alpha|+1$. In this case the expected reward is $\hat{r}_{1}=a_{0}-a_{1}$ and the maximum of $\hat{r}_{1}$ is equal to $-|a_{0}-x_{0}-\alpha|-1$.

From 1), 2), 3), 4) and 5) one can see that, given the conditions considered here, an optimal policy can either play $a_{1}=a_{0}+|a_{0}-x_{0}-\alpha|-1$ or $a_{1}=a_{0}-|a_{0}-x_{0}-\alpha|+1$. In the following we will fix $a_{1}$ to $a_{0}+|a_{0}-x_{0}-\alpha|+1$ when $a_{0}-x_{0}\geq 0$. Let us also observe that the expected value of $r_{1}$ is equal to $5.5-|a_{0}-x_{0}-\alpha|$. Up to now in this item d), we have only considered the case where $(a_{0}-x_{0})>0$. When $(a_{0}-x_{0})\leq 0$, using the same reasoning we reach the exact same expression for the optimal action to be played and for the maximum expected return of $r_{1}$. This is due to the symmetry that exists between both cases. Since $\pi_{theorem1}$ plays the action $a_{1}=a_{0}+r_{0}+1=a_{0}-|a_{0}-x_{0}-\alpha|+1$ in the case d) at time step $1$, it is straightforward to conclude that, in this case, it also plays an action that maximises the expected value of $r_{1}$.

Now that the first condition for $\pi_{theorem1}$ to maximise expression (13) has been proved, let us turn our attention to the second one. To this end, we will compute for each $x_{0}\in\mathcal{X}$, the action $a_{0}\in\mathcal{A}$ that maximises:

[TABLE]

and show that this action coincide with the action taken by $\pi_{theorem1}$ for time step $t=0$. First let us observe that for this optimisation problem, one can reduce the search space $\mathcal{A}$ to $[x_{0}-\alpha_{max}+1,x_{0}+\alpha_{max}-1]\subset\mathcal{A}$. Indeed, an action $a_{0}$ that does not belong to this latter interval would not give more information about $\alpha$ than playing $a_{0}=x_{0}-\alpha_{max}+1$ or $x_{0}+\alpha_{max}-1$ and lead to a worse expected $r_{0}$. This reduction of the search space will be exploited in the developments that follow.

However, we should first remember that $\mathbb{U}_{\alpha}=\mathbb{U}[-\alpha_{max},\alpha_{max}]$ and that the function $V_{1}(x_{0},a_{0},r_{0},x_{1})$ can be written as follows:

if $r_{0}=10$, $V_{1}$ is equal to $R_{max}=10$ 2. 2.

else if $|r_{0}|>\alpha_{max}-(a_{0}-x_{0})\quad\wedge\quad a_{0}-x_{0}>0$ and $r_{0}\neq 10$, then $V_{1}$ is equal to $R_{max}=10$ 3. 3.

else if $|r_{0}|>\alpha_{max}-(x_{0}-a_{0})\quad\wedge\quad a_{0}-x_{0}<0$ and $r_{0}\neq 10$, then $V_{1}$ is equal to $R_{max}=10$ 4. 4.

and otherwise $V_{1}$ is equal to $5.5-|a_{0}-x_{0}-\alpha|$.

We note that the value of $V_{1}(x_{0},a_{0},r_{0},x_{1})$ does not depend on the state $x_{1}$, which allows us to rewrite expression (16) as follows:

[TABLE]

and since the expectation is a linear operator:

[TABLE]

Let us now focus on the second term of this sum:

[TABLE]

We note that when $a_{0}-x_{0}\geq 0$ the function $V_{1}$ can be rewritten under the following form:

if $\alpha\in[-\alpha_{max},2*(a_{0}-x_{0})-\alpha_{max}[$, $V_{1}$ is equal to 10 2. 2.

else if $\alpha\in[2*(a_{0}-x_{0})-\alpha_{max},a_{0}-x_{0}-1]$, $v_{1}$ is equal to $5.5+\alpha-(a_{0}-x_{0})$ 3. 3.

else if $\alpha\in[a_{0}-x_{0}-1,a_{0}-x_{0}+1]$, $V_{1}$ is equal to $10$ 4. 4.

else if $\alpha\in]a_{0}-x_{0}+1,\alpha_{max}]$, $V_{1}$ is equal to $5.5-\alpha+(a_{0}-x_{0})$.

From here, we can compute the value of expression (19) when $a_{0}-x_{0}\geq 0$. We note that due to the symmetry that exists between the case $a_{0}-x_{0}\geq 0$ and $a_{0}-x_{0}\leq 0$, expression (19) will have the same value for both cases. Since we have:

[TABLE]

where $p_{\alpha}$ is the probability density function of $\alpha$, we can write:

[TABLE]

And thus, by computing the integrals, we have:

[TABLE]

Let us now analyse the first term of the sum in equation (18), namely $\underset{\alpha\sim\mathbb{U}_{\alpha}}{\mathbb{E}}(r_{0})$.

We have that:

[TABLE]

which can be rewritten as:

[TABLE]

Due to the reduction of the search space, we can assume that $a_{0}$ belongs to $[x_{0}-\alpha_{max}+1,x_{0}+\alpha_{max}-1]$, we can write:

[TABLE]

Given that $R_{max}=10$, we have:

[TABLE]

and therefore:

[TABLE]

To find the action $a_{0}$ that maximises (16), one can differentiate (18) with respect to $a_{0}$:

[TABLE]

This derivative has a single zero value equal to:

[TABLE]

It can be easily checked that it corresponds to a maximum of expression (16) and since it also belongs to the reduced search space $[x_{0}-\alpha_{max}+1,x_{0}+\alpha_{max}-1]$, it is indeed the solution to our optimisation problem. Since $\pi_{theorem1}$ plays this action at time $t=0$, Part 1 of this proof is now fully completed.

$\triangleright$ Part 2.

Let us now prove that the policy $\pi^{*}_{theorem1}$ generates for every $t\geq 2$ rewards equal to $R_{max}=10$. We will analyse three different cases, corresponding to the three cases a), b) and c) of policy $\pi_{theorem1}$ for time step $t\geq 2$.

1. a)

If $r_{0}=10$, we are in a context where $a_{0}$ belong to the target interval. It is straightforward to see that, by playing $a_{t}=x_{t}+a_{0}-x_{0}$, the action played by $\pi_{theorem1}$ in this case, we will get a reward $r_{t}$ equal to 10. 2. b)

If $r_{1}=10$ and $r_{0}\neq 10$, one can easily see that playing action $a_{t}=x_{t}+a_{1}-x_{1}$, the action played by $\pi_{theorem1}$, will always generate rewards equal to $10$. 3. c)

If $r_{0}\neq 10$ and $r_{1}\neq 10$, it is possible to deduce from the first action $a_{0}$ that the MDP played corresponds necessarily to one of these two values for $\alpha$: $\{a_{0}-x_{0}+r_{0};a_{0}-x_{0}-r_{0}\}$. Similarly, from the second action played, one knows that $\alpha$ must also stand in $\{a_{1}-x_{1}+r_{1};a_{1}-x_{1}-r_{1}\}$. It can be proved that because $a_{0}\neq a_{1}$ (a property of our policy $\pi_{theorem1}$), the two sets have only one element in common. Indeed if these two sets had all their elements in common, either this pair of equalities would be valid:

[TABLE]

or this pair of equalities would be valid:

[TABLE]

By summing member by member the two equations of the first pair, we have:

[TABLE]

Taking into account that $x_{0}=x_{1}$ because none of the two actions yielded a positive reward, it implies that $a_{0}=a_{1}$, which results in a contradiction. It can be shown in a similar way that another contradiction appears with the second pair. As a result the intersection of these two sets is unique and equal to $\alpha$. From here, it is straightforward to see that in this case c), the policy $\pi_{theorem1}$ will always generate rewards equal to $R_{max}$.

From Theorem 1, one can easily prove the following theorem.

Theorem 2

The value of expected return of a Bayes optimal policy for benchmark 1 is equal to $\frac{3*\gamma^{2}*(\alpha_{max}+4.5)^{2}}{2*\alpha_{max}*(1+\gamma)}+\frac{21+\alpha_{max}^{2}+\gamma*(10+11*\alpha_{max}-\alpha_{max}^{2})}{2*\alpha_{max}}+\frac{\gamma^{2}}{1-\gamma}*10$.

**Proof ** The expected return of a Bayes optimal policy can be written as follows:

[TABLE]

From the proof of Theorem 1, it is easy to see that:

$\underset{\underset{{x_{\cdot}}\sim P_{\mathcal{M}}\mbox{\tiny{(\raisebox{1.29167pt}{.},\raisebox{1.29167pt}{.})}}}{\underset{{a_{\cdot}}\sim\pi^{*}_{bayes}\mbox{\tiny{(\raisebox{1.29167pt}{.})}}}{\underset{{x_{0}}\sim P_{x_{0}}}{\mathcal{M}\sim\mathcal{D}}}}}{\mathbb{E}}\sum_{t=0}^{1}\gamma^{t}*r_{t}=\frac{3*\gamma^{2}*(\alpha_{max}+4.5)^{2}}{2*\alpha_{max}*(1+\gamma)}+\frac{21+\alpha_{max}^{2}+\gamma*(10+11*\alpha_{max}-\alpha_{max}^{2})}{2*\alpha_{max}}$ 2. 2.

$\underset{\underset{{x_{\cdot}}\sim P_{\mathcal{M}}\mbox{\tiny{(\raisebox{1.29167pt}{.},\raisebox{1.29167pt}{.})}}}{\underset{{a_{\cdot}}\sim\pi^{*}_{bayes}\mbox{\tiny{(\raisebox{1.29167pt}{.})}}}{\underset{{x_{0}}\sim P_{x_{0}}}{\mathcal{M}\sim\mathcal{D}}}}}{\mathbb{E}}\sum_{t=2}^{\infty}\gamma^{t}*r_{t}=\frac{\gamma^{2}}{1-\gamma}10$

which proves Theorem 2.