{'fraction_non_alphanumeric': 0.0441008869728953, 'fraction_numerical': 0.011768495119171498, 'mean_word_length': 4.358868826340945, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 5, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Active Learning in Gaussian Process State Space Model

We investigate active learning in Gaussian Process state-space models (GPSSM). Our problem is to actively steer the system through latent states by determining its inputs such that the underlying dynamics can be optimally learned by a GPSSM. In order that the most informative inputs are selected, we employ mutual information as our active learning criterion. In particular, we present two approaches for the approximation of mutual information for the GPSSM given latent states. The proposed approaches are evaluated in several physical systems where we actively learn the underlying non-linear dynamics represented by the state-space model.


Introduction
State-space models (SSMs) are a compact representation of dynamical systems as a set of input, output and state variables where the transition could be deterministic or stochastic. SSMs are widely used in practice, where applications range from aerospace-related work (e.g. [16]), medicine (e.g. [1]) to economics and finance (e.g. [56]). The two important aspects of SSMs are learning and control and they have been a fruitful research field in classical engineering (e.g. [18,30]). Learning in SSM is about finding optimal hyperparameters so that the dynamical system is accurately captured whereas control refers to optimising a specific objective function based on inputs and states of a given SSM. The latter is often referred as optimal control. During the last decade, both aspects of SSMs also attracted the attention from the machine learning community (e.g. [19,43]). To capture the uncertainty of system dynamics, Bayesian models are often employed in the context of SSM. For example, a Gaussian Process (GP) prior [44] is mostly used to describe the transition between states, which leads to family of Gaussian Process State-Space Models (GPSSMs) with an impressive amount of progress has been made (e.g. [12,13,20,21,28,29,32,39,53], and more). Learning and inference of GPSSM are an active research topic in machine learning community and this paper focuses on learning of GPSSM.
While different approaches of GPSSM mainly address the modelling issues, the question about data-efficient learning of GPSSM remains largely unanswered. Such questions go naturally to the field of Active Learning (AL), which is the process of strategically selecting and generating new data for supervised learning. The purpose of AL is to learn the proposed model with the least amount of data possible. Though the majority of recent papers in this field focused on classification problems (e.g. [3,48]) and rather few were for regression (e.g. [7,46]), Schreiter et al. [45] applied AL in GP regression and Zimmer et al. [57] extended the scope to learning time-series models. Both papers employed the entropy criterion for exploration, while using a greedy selection based on the maximum predictive variance of the GP. In context of state-space models for predictive control, Buisson et al. [6] proposed an AL approach for learning the transition dynamics via optimising the trajectory. For AL in GPSSM, Capone et al. [8] was the first to introduce the concept here. However, in both cases, Buisson et al. [6] and Capone et al. [8] simplified the settings by assuming the states are observable and measurable. By doing so, the resulting state-space models can be actively learned in the same way as standard GP regression models in a supervised manner.
In contrast to previous work, we present an AL strategy on GPSSM without relying on the assumption that the states are observed. That is, we actively steer the system through the latent state-space by determining its inputs such that the underlying dynamics can be optimally learned by a GPSSM. For the information criterion, we employ an approximated measure of the mutual information. We propose and discuss two different approaches for the approximation of mutual information in the context of unmeasurable, latent state-space. This paper is on learning a stochastic nonlinear dynamical system actively from noisy data via probabilistic SSM. Its main contributions include: -We derive tractable approximate mutual information estimates for GPSSM with latent states based on the state-of-the-art learning scheme and the technique of approximate Gaussian integral. -We propose an AL strategy for GPSSM with latent states based on inputs.
-We conduct experiments for different examples to demonstrate the usefulness of our strategy in physical systems, which is a typical application for SSM.
The reminder of this paper is organised as follow: We give a brief overview of backgrounds in section 2 and introduce GPSSM in section 3. Section 4 shows our active learning strategies and section 5 evaluates our concepts on several test scenarios. A conclusion is given in section 6. Detailed proof of propositions, experimental details and few further remarks are left in supplementary materials.

Background
AL is itself a broad topic and here we refer readers to, for example, Settles [47] and Dasgupta [10] for an in-depth survey of the basic algorithmic and theoretical ideas. There are many different paradigms in AL. For example, Bayesian Active Learning incorporates the Bayesian framework and is often referred to AL with GP models (e.g. [27,54]), because GPs naturally carry the uncertainty measures. Specifically, Houlsby et al. [27] applied AL in GP classification problems, where they defined an acquisition function which estimates the quantity of mutual information between model predictions and model parameters. Other variants include Batch Active Learning and Deep Bayesian Active Learning. While Batch AL refers to making multiple queries in parallel (e.g. [27,31] In dynamics modelling, SSM has been a long-lasting topic with enormous literature across different disciplines (e.g. [18,30,56]). In machine learning, the combination between SSM and GP began arguably from Wang et al. [55], where they learned the latent state of GPSSM via maximum a posteriori. Turner et al.
[53] extended it to learning the transition function. Frigola et al. [21] presented a Bayesian treatment based on particle Markov-Chain-Monte-Carlo and Frigola et al. [20] further introduced a variational inference scheme to overcome the computational complexity of Monte-Carlo based methods. This becomes the state-of-the-art and due to the fact that this approach depends heavily on the approximated function, most recent work introduced their own variational inference scheme with their own advantages and disadvantages (e.g. [12,28,29,39]). For example, Doerr et al. [12] and Ialongo et al. [29] introduced their specific variational inference scheme to account for the dependence between transition function and latent states, which were treated as independent in earlier work.
To the best of our knowledge, AL and GPSSM have not been discussed together until Buisson et al. [6] first brought the topic of AL in learning the transition function using GP. They stated that the dynamical problem, i.e. learning the transition function, is fundamentally different from a static AL in GP problem, because we need to steer the system to certain states through the unknown transition function with a sequence of actions, which is the only component we can control. For each round of exploration, the next input is actively picked so that the model error can be maximally reduced. They proposed their AL strategy using the maximum entropy. Another recent work on AL in GPSSM is by Capone et al. [8]. They proposed taking the most informative data point while employing the mutual information between the most informative point and reference points. The states, in the existing work on AL in GPSSM, are hitherto assumed to be observable such that the model can be trained in a fashion similar to GP regression. Namely, in Capone et al. [8], the fact that states are observable simplify the estimation of the mutual information, resulting in two optimisation problems to solve in parallel.
Our paper differs from the previous work by presenting an AL strategy on GPSSM with states that are latent, resulting in a different estimation of the mutual information as exploration criterion while extending the scope of AL for GPSSM. We leave further quantitative discussion in section 4.

Gaussian Process State-Space Model with Inputs
We consider a discrete-time sequence of T observations y 1:T ≡ {y t } T t=1 , where each observed point y t ∈ Y ⊆ R dy is generated by a corresponding latent variable x t ∈ X ⊆ R dx and the previous step control c t−1 ∈ C ⊆ R dc . The collection of latent states and controls are denoted by x 0:T ≡ {x t } T t=0 and c 0:T ≡ {c t } T t=0 , respectively. These latent variables and controls are assumed to satisfy the Markov property, meaning that any x t+1 can be generated by only conditioning on x t , c t and the transition function f : R dx ×R dc → R dx . To align with previous work, we use Gaussian distribution for both transition and observation density function. We also simplify the mean of the observation function by using linear mapping and tackle multivariate latent states by placing a GP prior on each dimension separately, as in, for instance, Ialongo et al. [29].
GPSSM is defined to be a probabilistic SSM with a GP prior over the transition function f , specified by f ∼ GP(m(·), k(·, ·)), (1) where in equation 1, f is a GP governed by a given mean function m(·) and positive definite covariance function k(·, ·). The initial state p(x 0 ) = N (x 0 |µ 0 , Σ 0 ) in equation 2 is assumed known. C and d are parameters in the linear mapping, whereas Q and R are covariance matrices which capture process and observation noise, respectively. The advantages of linear mapping are that linear mean won't limit the range of systems that can be modelled and this reduces the nonidentifiabilities between transitions and emissions. As a convention from previous work, we denote f t ≡ f (x t−1 , c t−1 ) and equation 3 becomes For brevity, we writex * = (x * , c * ), where * could be a specific index or a collection of indices. Thus, the matrix of covariance functions are denoted by K i:j := (k(x s ,x t )) j s,t=i . For convenience, we write k(x i:j ,x k ) ≡ (k(x i ,x k ), · · · , k(x j ,x k )) as a vector collection of kernel evaluations. It holds that k(x i:j ,x k ) T = k(x k ,x i:j ).
The control inputs c t can be viewed as an augmented latent state. Since c t is not a random variable, we always condition on this quantity. Most previous x0 x1 · · · · · · ct−1 ct y1 · · · yt−1 yt yt+1 papers related to this model were either studied without controls or, if a control is involved, by assuming that this quantity was randomly distributed [37]. We, however, construct an AL strategy based on c t to learn the GPSSM efficiently. Hence, we purposely keep this term in our definitions. Figure 1 shows a graphical model of such a setting. We would like to readdress here that c t is a known quantity where our task is to optimise the value to select. This should not be confused with control in other contexts such as that in RL.
To complete the structure of GPSSM, let us define f 1:T ≡ {f t } T t=1 and we can write the joint density as follows (e.g. chapter 3 in Frigola [19]), p(y 1:T , x 0:T , f 1: In particular, from equation 1, the last factor in equation 6 is derived in the same manner as the posterior in GP regression, given by where and the term f 1:t−1 −m(x 0:t−2 ) is written as
The major challenge is the computation of p(y 1:t ). Throughout this paper, we are going to apply variational inference to train GPSSM, which is a technique based on making assumptions about the posterior over latent variables that leads to a tractable lower bound, often called evidence lower bound (ELBO). There are two reasons for choosing this approach, which is currently state-of-the-art in training GPSSM. First, variational inference is computationally efficient. Second, this approach will lead to an approximation of the posterior, which can be carried to the AL strategy, as we will present further in section 4.
The first step of using variational inference in GPSSM is to introduce M (M T ) inducing points u = u 1: with their corresponding inputs. This is referred as sparse GP technique (see e.g. [22,52]) and the joint density is (8) Then, the ELBO L t to the log marginal likelihood of p(y 1:t ) in equation 6, which is based on KL divergence (see e.g. chapter 8.5 in Cover [9]), is given by log(p(y 1:t )) = L t +KL [q(x 0:t , f 1:t ) p(x 0:t , f 1:t |y 1:t )] for t = 1,· · ·, T , where q(·) is the distribution function to approximate p(·). With further derivations from variational inference methodology, L t is given by Analogously to Ialongo et al. [29], we set which leads to equation 17 of [29] (11) Under sparse GP approximation, we can specify a free Gaussian density on the function values u giving q(u) = N (u|µ u , Σ u ). The only term that requires further specification is q(x i |f i ) and [29] presented q( , wheref i is a free variational parameter. A i , b i and S i are other free variational parameters depending on the choice off i . This generalises a number of previous work. For example, by setting f i = x i in q(·|·) and q(x) = q(x 0 ), we recover the expression from [20] and the optimal q(x) can be solved by calculus of variations. If we setf i = x i , then we have the Gaussian factorised approximation of [39]. In [29], they proposed settingf i = k(x i , u)K −1 uu u and the training turned out to be marginally better if a sufficiently large amount of iterations are allowed.
We choose to setf i = x i throughout this paper, which is the approach by [39]. This is because their approach has shown to have a stable improvement in model accuracy during the training phase. Also, its linear factorising nature, given a fixed and often small number of allowed iterations, is also beneficial for training time. Such setting leads to q( The prediction with GPSSM can be done rather cheaply. Suppose we would like to predict the new observation y based on a new control input c , by defining all training data as D ≡ {y 1:T , x 0:T }, we can show easily [20] that the predictive distribution p(y |c , D) is given by In practice, we sample f (x T , c ), where the latent states x 1:T are already estimated via training. Then, we predict the new observation based on drawing samples from the distribution. This is the approach used by [29].

Active Learning Strategies
The previous section presents the modelling of GPSSM and in this section, we are interested in acquiring a strategy to actively learn the GPSSM by steering the system through the latent state-space. We employ AL as a query strategy to pick the most informative new control input c * t while observing y t , and to learn the model as we explore. Contrary to the two previous work [6,8], we treat the states as latent and unknown and c t is the only controllable variable. For the exploration criterion, we use approximated mutual information. The usage of mutual information for AL has been well motivated by, for example, Krause et al. [34]. They pointed out that the mutual information might lead to a more accurate model than the differential entropy. This quantity was also shown to be the same as minimising the expected uncertainty of the model [15].
First, let us recall from chapter 8.5 of Cover [9] that the mutual information between any two sets of random variables Y and F with joint density p(Y , F ) is defined as and a well-known relationship between mutual information and differential entropy h(·) is given by For steering the system through latent state-space while gathering information for efficiently learning the model, a sensible quantity for exploration is the mutual information between the latest observations y t+1 and latest predicted transition functions f t+1 . This quantity is an extension to the active strategy from Buisson et al.
[6] by using mutual information as a criterion. That is, where c t is placed within the term f t+1 ≡ f (x t , c t ). We useŷ t+1 because this term is provided from predicting GPSSM with a given new control c t . Altogether, c t is the only independent variable. This expression is intractable and there have been enormous efforts in estimating it (e.g. [2, 26,35,40,42]). However, most previous work assumed that we have little or no information on both entropies in the mutual information expression and this comes down to the fact that we do not know the corresponding probability density function (pdf). Their approach includes deriving bounds but extra functions such as critic or variational parameters are often required (e.g. section 2 in [42]). However, should we have some information about this quantity, it is much better to make use of this known information instead. Therefore, we derive our estimate for the mutual information by directly approximating the pdf. The main advantages of our method versus others are: (1) As discussed in the previous paragraph, other methods often require extra functions called critic or variational parameters, or some free parameters to set. Finding the best settings of these extra functions and parameters could be complicated at times. Mcallester et al. [38] also showed that there are statistical limitations in various bounds of mutual information and a better approach is to measure this quantity as a difference in entropies.
(2) Our mathematical derivation is purely via inequalities and numerical approximations, resulting an approximated closed form. Hence, the final expression derived can be easily implemented.
The quantity I(y t ; f t ) is often studied in SSM (e.g. [15]), motivated by the assumption of Markov property. However, if the Markov assumption is not completely valid, as there are (weak) dependencies between the states, the mutual information between all observations y 1:t+1 and all predictions from the transition function f 1:t+1 within a time-horizon might present a better quantity. To the best of our knowledge, we are the first to propose using mutual information between all variables in time, and the AL strategy presented in equation 13 can be viewed as a special case of this general formulation. Thus, another more general AL strategy is Again, the termŷ t+1 is provided from predicting GPSSM. This also leads to the third advantage of our method: (3) To the best of our knowledge, we are the first to estimate the total mutual information by leveraging the inference scheme of the model. This turns out to be computationally efficient as well.
input : Initial T observations y1:T and control inputs c0:T . Initial state x0. output: An optimised GPSSM after taking N exploration steps, the final data set including the explored observations y1:T +N and the corresponding controls c0:T +N .
Train the initial GPSSM for t = T, T + 1, · · · , T + N do 1. Solve c * t via equation (14) (or (13)); 2. Evaluate yt+1 from the system; 3. Update the training data set; 4. Retrain GPSSM to update hyperparameters; 5. t = t + 1; end Algorithm 1: Active Learning in GPSSM on control through approximation of mutual information For clarity, we refer I(ŷ t+1 ; f t+1 ) as latest mutual information (latMI) and I(y 1:t ,ŷ t+1 ; f 1:t+1 ) as total mutual information (totMI). In the remainder of this section, we present how the two quantities can be approximated for GPSSM with latent states. We first present an estimation scheme for the latest mutual information to align with previous work. Then, we proceed to the estimation of our proposed total mutual information. Given the estimates of mutual information, we can formulate the AL strategy for GPSSM, as summarised in algorithm 1.

Computation of latest Mutual Information
To estimate the latest mutual information I(y t ; f t ), we employ the Gaussian approximate integral derived by Girard [24]. Here, for two random variables x and y with p(y|x) = N (y|µ(x), σ 2 (x)) and p(x|u, Σ x ) = N (x|u, Σ x ) for some mean u and variance Σ x , the Gaussian approximation yields where M (u, Σ x ) and V (u, Σ x ) are integral functions to be evaluated (see supplementary materials for details). Based on this approximation, I(y t ; f t ) can be approximately computed as in the following proposition.
Proposition 1. Given the definition of GPSSM from equation 1 -4, as well as the notation of M t and K t defined in equation 7. We define the approximation of the following integrals as and, recursively, for all t = T, T + 1, · · · , T + N , Then, the latest mutual information is approximately Proof (Sketch). Based on the definition of mutual information, the numerator in equation 15 relies on the approximation via moment matching, whereas the denominator in equation 15 can be computed directly. For more details, we refer the interested reader to the supplementary materials.

Computation of total Mutual Information
The ELBO L t given by equation 11 can be utilised to estimate the total mutual information. Since this quantity is estimated via samples, this quantity is sample driven and the following proposition presents our computational approach.
Proposition 2. Given the definition of GPSSM from equation 1 -4, as well as the expression of the ELBO L t from equation 11, for t = T, T + 1, · · · , T + N , if S samples are drawn, the s-th sample estimate (s = 1, · · · , S) of the mutual information between observations y 1:t and the prediction f 1:t , denoted by i s , is bounded by where N s (·) and L t,s are the s-th sample estimate of the normal distribution and ELBO L t , respectively. The total mutual information is then approximated by I(y 1:t ; f 1:t ) ≈ 1 S S s=1 i s .
Proof (Sketch). The proof employs mainly the definition of mutual information. The intractable part is bounded via ELBO, and the other term can be computed directly. Since it is a sample based approach, the bound from a sample becomes an approximation after taking all samples. For more details, we refer the interested reader to the supplementary materials.
Remark 1. In computing equation 16, the sampling nature comes from the specific value of x 0 drawn, and this value will be different between samples. Therefore, the evaluated i s will also be different. The efficiency of this approach heavily relies on the inference of GPSSM and it is possible to use another form of ELBO by defining q(·) differently or even other inference techniques such as a MCMC based approaches. While different approaches will incur their corresponding computational effort, it is better to keep the same approach in training GPSSM and computing total mutual information in order to avoid extra computational cost, as we reuse the posterior with trained hyperparameters in estimating the total mutual information.

Experiments
In this section, our experiments demonstrate how our proposed AL strategies help in learning GPSSM. From our discussion during previous sections, the two latest work from Buisson et al. [6] and Capone et al.
[8] assumed known states, so we cannot compare against them in the setting with unobservable states we focus on. Therefore, the reasonable benchmark we can compare to is the usual setting where models are learned with randomly distributed control inputs (e.g. [37]). We first perform an experiment on a simulated function, then on few non-linear physical problems. 3 The overall goal is to obtain a high accuracy for the GPSSM with as few selected points as possible.

Simulated Function
We first consider a modified kink function based on [29]. Our goal is to learn the dynamics of the system within a range of X = [−3.0, 1.1] for the latent states and C = [0, 1] for the control. We begin with 5 given points and explore 30 steps.
The left diagram in Figure 2 depicts the RMSE of training sessions applying random selection, latest mutual information (latMI) and total mutual information (totMI) strategy, respectively. Each session consists of independent runs. Soon after collecting the first control input, both AL show results not worse than random but latMI is not significant, while totMI is more profound.
To illustrate why we need our AL strategies, especially totMI, we also present a snapshot of the position of the maximum totMI as we explore in the right of Figure 2 Fig. 3. The left diagram shows the results of the pendulum: Mean and standard deviation of RMSE for 30 exploration steps (starting with 5 initial points) are collected from independent trials. latMI is presented in yellow, totMI in red, and random exploration in blue. The right diagram shows the results of the cart-pole: Again, mean and standard deviation of RMSE for 30 exploration steps (starting with 5 initial points) are collected from independent trials. informative control where we can attain maximum totMI. Therefore, optimisation in equation 14 (or 13) is necessary. These diagrams could change as we run different trials, and such evaluation will be more complicated as we increase the dimension of controls or the dynamics become more complex.

Pendulum and Cart-Pole
In these experiments, our aim is to learn the physics of the pendulum and cartpole by actively controlling the input torque and the force of the cart, respectively. We evaluate the state-space model's accuracy via the angular position obtained from the model against the ground-truth, which is calculated via the equation of motions. The description of the physics of both pendulum and cartpole, as well as the experimental settings, are presented in supplementary materials. Figure 3 depicts the RMSE of training sessions applying random selection, latest mutual information and total mutual information strategy for pendulum (left) and cart-pole (right), respectively. Both experiments consist of independent runs. For the pendulum, both totMI and latMI strategy show advantages against the random exploration. However, totMI shows even better results as the exploration proceeds. For cart-pole, all strategies have, in the initial exploration phase, a raise in RMSE, which is subsequently reduced, as the model becomes more stable. latMI is indeed unable to attain noticeable improvement but totMI strategy still shows favourable performance via a consistently lower RMSE.
Since the physics of both systems are common enough that we can actually measure the states, we also compare our work directly with the AL strategy by Capone et al. [8], even though the fairness of such a comparison is debatable. We refer to the supplementary materials for a detailed discussion.

Twin-Rotor Aerodynamical System
As a more realistic setting, we test our AL strategy on a Twin-Rotor Aerodynamical System (TRAS). This is a typical design for control experiments, and its behaviour resembles that of a helicopter. There are two rotors -one horizontal and another vertical -joined by a beam. There is also a counter-weight, which determines a stable equilibrium position. When the system is switched off, the main rotor is lowered. The motion of the system is controlled by the two motor supply voltages to each motor. We demonstrate a simple example of such a system in the left diagram of Figure 4.
Our aim is to learn the physics of the TRAS by actively control the two voltages to the system. We evaluate the model's accuracy via the angular position measured from the model against which is calculated via the equation of motions on both rotors. The physics of the TRAS is based on the description by [4] and we leave the details in supplementary materials. The right diagram of Figure 4 depicts the RMSE of training sessions applying random selection, latMI and totMI strategy, respectively. Each session consists of independent runs. Both latMI and totMI strategies show advantage against random exploration in presence of more sophisticated dynamics. However, GPSSM with totMI reaches an acceptable accuracy with much lower number of data points.

Discussion
We have studied AL in GPSSM by proposing a tractable mutual information estimate, in order to select the most informative control inputs to maximise our model accuracy as we explore. Both latest and total mutual information are proposed. The empirical results show that total mutual information strategy outperforms the other approach and random exploration, which shows that total mutual information should be adopted. There is a promising potential of our approach towards more complex industrial settings. Since we leverage ELBO to derive our estimates, with suitable modifications, we also expect our propositions to be applicable to future work on inference of GPSSM and other probabilistic models with latent space.
Another direction of extension is to look into the use of maximal mutual information in reinforcement learning (RL) community, as such strategy has been used in planning and control. For instance, Ding et al.
[11] applied such strategy in Partially Observable Markov Decision Process (POMDP) framework in order to efficiently learn the RL model. However, there are key differences between the two themes. First, the word control is interpreted very differently. In GPSSM, this refers to an augmented latent state and each input is simply a deterministic real number whereas in POMDP, this is a mapping between an action taken in certain state to maximise the expected return or reward. Second, AL in GPSSM has a very specific problem formulation whereas efficient learning in POMDP is a general framework with a large room of freedom to specify different components, leading to various research topics (e.g. [5,14,25]  Proposition 1. Given the definition of GPSSM from equation 1 -4, as well as the notation of M t and K t defined in equation 7. We define the approximation of the following integrals as and, recursively, for all t = T, T + 1, · · · , T + N , Then, the latest mutual information is approximately Before the proof, we first elaborate the Gaussian approximation technique from Girard [24]. For two random variables x and y with p(y|x) = N (y|µ(x), σ 2 (x)) and p(x|u, Σ x ) = N (x|u, Σ x ) for some mean u and variance Σ x , the Gaussian approximation technique from [24] yields to In particular, if µ(x) is linear (e.g. µ(x) = F x + a) and σ 2 (x) can be represented as a matrix A. The integral reduces to the exact formula of Gaussian integral In the context of GP, M (·, ·) and V (·, ·) can be written based on the expectation of the kernel, which is where β = K −1 y and β = {β i } t i=1 and K refers to the matrix of all covariances of the GP. There are three terms we need to evaluate, denoted by This can be further simplified when the covariance function has a specific form. For example, if the kernel is squared exponential, which can be rewritten as with c = (2π) dx/2 |Λ| 1/2 σ 2 , we can show that Proof. We begin by using the relationship I(y t ; f t ) = h(y t )−h(y t |f t ). Note that ELBO L t , as presented in equation 9, gives a bound of p(y 1:t ) but not a specific p(y t ). There is no clear and simple relationship between one observation and the joint of all observations, so we require an approach which allows us to estimate p(y t ) directly. While there could be different estimates available, we adopt the approximate Gaussian integral from [24] because of reason (1) and (2) stated in section 4. First, we express p(y t ) as p(y t ) = p(y t , x 0:t , f 1:t )dx 0:t df 1:t = p(y t |x 0:t , f 1:t )p(x t |x 0:t−1 , f 1:t )p(f t |x 0:t−1 , f 1:t−1 )p(x 0:t−1 , f 1:t−1 )dx 0:t df 1:t where the third line and terms in the last line follow by the definition of GPSSM in the main text (equation 1-4 and 7). Note that in the forth line, f 1:0 ≡ ∅.
The next step is to integrate the expression in the order of and since M 0 , which depends on x 0 , is in general non-linear, we cannot use the identity of Gaussian integral and this is where we apply the Gaussian approximation by [24]. By interlacing the exact and approximate Gaussian integral, we first have Then, the next integration is an exact Gaussian integral The process continues until the term M t and V t are defined (recall that for M t and K t , these two terms depends on x 0:t and f 1:t ): The last integration with respect to x t , which is an exact Gaussian Integral, is given by Hence, it shows that approximately, p(y t ) ≈ N (y t |CM t + d, R + C(Q + V t )C T ). Since this is a normally distributed expression, the entropy of which is given by It remains to compute h(y t |f t ). This term can be directly computed, as So, the entropy of this distribution is given by Finally, we obtain the approximation of the mutual information between the latest observation and prediction as stated in the proposition, given by Proposition 2. Given the definition of GPSSM from equation 1 -4, as well as the expression of the ELBO L t from equation 11, for t = T, T + 1, · · · , T + N , if S samples are drawn, the s-th sample estimate (s = 1, · · · , S) of the mutual information between observations y 1:t and the prediction f 1:t , denoted by i s , is bounded by where N s (·) and L t,s are the s-th sample estimate of the normal distribution and ELBO L t , respectively. The total mutual information is then approximated by I(y 1:t ; f 1:t ) ≈ 1 S S s=1 i s .
Proof. Again, we first use I(y 1:t ; f 1:t ) = h(y 1:t ) − h(y 1:t |f 1:t ). Since we have the estimate for the probability density function (pdf), a generic Monte-Carlo estimate of the entropy (e.g. , where S is the number of samples drawn andp s (·) is the estimate of the pdf of the random variables. With two estimates of pdf formulas, the total mutual information can be estimated by (log (p s (y 1:t |f 1:t )) − log (p s (y 1:t ))) .
This shows that for the s-th sample, s = 1, · · · , S, the total mutual information estimate can be defined and bounded by i s := log (p s (y 1:t |f 1:t )) − log (p s (y 1:t )) ≤ log (p s (y 1:t |f 1:t )) − L t,s , where the upper bound is based on the inference scheme we have in training the GPSSM, as we have derived an ELBO L t (equation 9) and we use the fact that log (p(y 1:t )) ≥ L t .
The remaining part is to computep s (y 1:t |f 1:t ) and this term can be exactly derived based on the definition of the GPSSM. The outcome is a product of normal random variables given bŷ p s (y 1:t |f 1:t ) = p s (y 1:t , x 0:t |f 1:t )dx 0:t = p s (y t |y 1:t−1 , x 0:t , f 1:t )p s (x t |y 1:t−1 , x 0:t−1 , f 1:t )p s (y 1:t−1 , x 0:t−1 |f 1:t )dx 0:t where the last line is based on using the definition of GPSSM recursively. The remaining derivation is purely mathematical techniques involving Fubini's theorem and Gaussian integral.
With products becomes summations by taking logarithm, we arrive proposition 4.1.
Remark 2. If the dimension of the observation is 1 (i.e. d y = 1), the product of normal random variables can be further simplified. This can be written in a closed form based on Meijer G-function [49]. Recent work showed that it is possible to compute the cumulative distribution function for such distribution as well (e.g. [51]). We postulate that there will be development on the entropy for such distribution from the mathematics community which could further simplify our expression.

A.2 Experimental Details
Here we supplement the experimental section in the main text for further details and results. We also present the equations of motion for the physical problems studied for completeness.
Simulated function The transition function of the modified kink function, based on Ialongo et al. [29], is given by where y t and c t are the target observation and control input, respectively. We set the latent state x t = y t for simplicity.
Pendulum and Cart-Pole -the physics The observed variable here is the angle of the pendulum θ (i.e. y t ≡ θ t ) and the control c is the torque acting on the pendulum. The equations of motion are given by The value of the parameters are g = 9.81 ms −2 , m = 0.1 kg, l = 1 m, b = 0.05 and I = ml 2 .
The cart-pole we used is based on Florian [17]. The observed variables are the angle of the pendulum θ and the position of the cart x. The control c is the + m cb (l cb sin(α v )) 2 + m ts 2 r 2 ts + m ts (l t cos(α v )) 2 + m ms r 2 ms + m ms (l m cos(α v )) 2 .
It appears to have many parameters but all of them are predetermined and the only independent variables are the two voltages applied, namely c v and c h . We try to learn the model by showing that an optimised selection of the voltages is better than randomly selected ones. Table 2 summarises key hyperparameters used for the experiments. For fairness, we use these settings for all random exploration, latest mutual information (latMI) strategy and total mutual information (totMI) strategy.  Running Time Analysis In addition to the predictive performance, we also present the exploration time depending the number of points. The elapsed time is displayed in Figure 5. While it is impossible to beat random exploration because there is no need to optimise for the control, we show that the totMI strategy runs significantly faster than the latMI strategy. This is because we leverage from the inference scheme for the calculation of p(y 1:t ), which saves a significant amount of computational effort. This effect applies to all four experiments we investigate in this paper.  as in GP regression and they resolve this two optimisation problems to be solved in parallel. Therefore, the number of hyperparameters from their settings is rather different from ours. To this end, we keep their initial hyperparameter settings for simplicity. Some settings from [8] are adjusted to be the same as our settings to facilitate a more plausible comparison, including the number of training epochs, taking arithmetic mean for different trials, and the number of points explored before the model is updated.

Hyperparameters settings
We choose Pendulum and Cart-Pole to run the 'comparison'. This is because the true value of the angular velocity in these two problems can be computed. Since we have shown from the previous subsection that in addition to the inferior performance, running time for latMI was prohibitively high in comparison to totMI. Therefore, we only look into totMI strategy. We focus on random exploration, totMI and the approach from [8], denoted by locAL.
Results are summarised in Figure 6. locAL showed a very unstable improvement in RMSE. Specifically, for pendulum, their results appear to be better in average, but not significant. Whereas for cart-pole, the RMSE begins to decay faster. This could be explained by the fact that their approach collected the information of the states after each round of exploration, which allows the model to be stable quicker. However, by the time when all approaches becomes stable, our totMI strategy still demonstrates a stronger performance. Our approach processes a more stable improvement in RMSE and a competitive performance, compared to the approach by [

A.4 Further discussion between GPSSM and Reinforcement Learning
We would also like to address that the strategy of maximising mutual information has also been used in the reinforcement learning (RL) community, especially in planning and control. For instance, Ding et al.
[11] applied such strategy in Partially Observable Markov Decision Process (POMDP) framework in order to efficiently learn the RL model. However, there are key differences between AL in GPSSM and efficient learning in POMDP. First, the word control is interpreted very differently. In GPSSM, this refers to an augmented latent state and each input is simply a real number (or a real-valued vector if we have multiple controls per time step). This quantity is always deterministic because it is a user-input to steer the exploration of the model. On the other hand, in POMDP, this is a mapping between an action taken in certain state to maximise the expected return or reward. Optimising control in this context is essentially finding the optimal policy in RL. Second, AL in GPSSM is a very specific problem formulation and this paper is presented to tackle this whereas efficient learning in POMDP is a general framework with a large room of freedom to specify different components, which leads to various research topics. Afterall, AL for SSM with GPs and RL are two different subfields in machine learning despite GPs are often applied as a tool in RL to help achieve their objective (e.g. [5, 14,25]).