Machine learning of electro-hydraulic motor dynamics

In this paper we propose an innovative machine learning approach to the hydraulic motor load balancing problem involving intelligent optimisation and neural networks. Two different nonlinear artificial neural network approaches are investigated, and their accuracy is compared to that of a linearised analytical model. The first neural network approach uses a multi-layer perceptron to reproduce the load simulator dynamics. The multi-layer perceptron is trained using the Rprop algorithm. The second approach uses a hybrid scheme featuring an analytical model to represent the main system behaviour, and a multi-layer perceptron to reproduce unmodelled nonlinear terms. Four techniques are tested for the optimisation of the parameters of the analytical model: random search, an evolutionary algorithm, particle swarm optimisation, and the Bees Algorithm. Experimental tests on 4500 real data samples from an electro-hydraulic load simulator rig reveal that the accuracy of the hybrid and the neural network models is comparable, and significantly superior to the accuracy of the analytical model. The results of the optimisation procedures suggest also that the inferior performance of the analytical model is likely due to the non-negligible magnitude of the unmodelled nonlinearities, rather than suboptimal setting of the parameters. Despite its limitations, the analytical linear model performs comparably to the state-of-the-art in the literature, whilst the neural and hybrid approaches compare favourably.


Introduction
A load simulator [1] is a ground system that reproduces the aerodynamic force on the rudder of a flying aircraft. It is a typical example of passive loading servo system, and finds important application in hardware-in-the-loop simulator devices. It allows designers under strict laboratory conditions to foresee and detect problems related to the rudder mechanics and control policy, cutting down costs and development time.
There are different kinds of load simulators: electric, pneumatic, and electro-hydraulic. The latter excels in terms of durability, power to weight ratio, controllability, accuracy, and reliability [2]. Thanks to these advantages, electro-hydraulic load simulators found wide application in the development and testing of systems beyond the aviation and aerospace domain, such as ship fin stabilisers [3], power steering for heavy duty vehicles [4], and robotic manipulators [5,6].
To reproduce aerodynamic flow, the load simulator is typically connected with the actuation system via a rigid shaft. The direct coupling between the two systems produces a strong position disturbance on the load system. This generates an extra torque which needs to be compensated mechanically or via a control system [7][8][9][10][11][12][13].
For control and simulation purposes, it is often useful to identify the dynamics of the load simulator. That is, to model the relationship between the device state, control signal, and torque on the load bearing system (including the extra torque). optimisation (PSO) [32], and the Bees Algorithm (BA) [33]. The latter algorithm showed great promise on various optimisation benchmarks [33][34][35], and its application in this context constitutes a novel contribution to its knowledge. As a term of reference, random guess (RG) was also tested for the optimisation of the parameters of the linearised model.
The non-linear model of the the load simulator dynamics was built using a multi-layer perceptron (MLP) [27] ANN. In this case, the MLP architecture defines the structure of the model, and the connection weights are the parameters to be optimised. The strong points of the MLP approach are the possibility of reproducing arbitrarily complex nonlinear mappings [36], and the existence of well-known inbuilt algorithms to optimise the weights inductively from data points. In this study, a fast weight learning algorithm [37] was used to reduce the ANN training time. The weak point of the MLP approach is the complexity of the ANN, which makes the model non-transparent to the user. That is, the MLP model is a black box.
In the hybrid approach, the analytical model was used to reproduce the main behaviour of the system, and the MLP to fit unmodelled nonlinearities.
The three modelling approaches are presented in Sect. 2, whilst Sect. 3 describes their optimisation. The experimental results are presented in Sect. 4 and discussed in Sect. 5. Section 6 concludes the paper.

Modelling approaches
The overall experimental system is described in Fig. 1. It consists of two subsystems connected by a rigid shaft and a torque sensor: the loading system (on the left in the figure) and load-bearing system (on the right). The overall Fig. 1 Schematic diagram of the experimental rig. The components are: 1 load motor; 2 torque sensor; 3 coupling axis; 4 load bearing motor; 5 angle encoder; 6 load bearing valve; 7 AD/DA card; 8 computer; 9 load valve stiffness of the connecting link includes the stiffness of the torque sensor.
The two subsystems are controlled using closed loop schemes. The angle of the load-bearing system is monitored using an encoder. The signal from the encoder is used by the load bearing controller to achieve servo position control. The load system controller uses angle and torque information to regulate the torque on the load bearing system. The remaining of this section describes the linear and non-linear models of the electro-hydraulic load simulator used in this study.

Linearised analytical model
This section briefly outlines the linear model that approximates the relationship at time k between the current plant output T g (k) (output torque of the load-on motor), its previous values, and the current an past values of the two plant inputs: d (k) (turn angle of the motor drive shaft in radians), and u m (k) (load input from the controller). The derivation of the model is fully detailed in "Appendix 1".
The linearised model consists of the following equation: Equation (1) includes the following 14 independent variables : the torques at time steps k − 1, … , k − 4 (4 variables), load inputs u m at time steps k, … , k − 4 (5 variables), and angles d at time steps k, … , k − 4 (5 variables). The contribution of each of the 14 independent variables is weighted by a number of constant parameters that are detailed in "Appendix 1", and depend on the 9 primary parameters listed in Table 1. The optimal value of these 9 parameters is not known. The boundaries showed in the table were set based on expertise and the physical properties of the system.

Multi-layer perceptron model
The analytical model described in Sect. 2.1 approximates linearly the non-linear relationship between the output torque and plant inputs. Given the difficulty of representing such relationship, the linearised model offers a reasonable trade-off between modelling accuracy and mathematical complexity. An alternative approach to circumvent the difficulty of identifying the load simulator dynamics is to employ a non-linear neural network model.
Artificial neural networks (ANNs) [27] are widely used learning systems composed of an input layer of units (neurons), a variable number of hidden layers, and an output layer. The number of units per layer and their connectivity defines the ANN topology.
The multi-layer perceptron (MLP) [38] is arguably the best known ANN model. It is composed of three or four layers of units. Each element is fully connected to all the units in the immediately previous and following layer, with no inter-layer connections between units. Each connection is associated to a weight, which allows modulating the signal passed from one unit to the other. The plant inputs are collected by the first (input) layer and processed in feedforward manner, one layer after the other (feedforward model). The input layer is composed of as many units as the plant inputs, and acts as a buffer. The units of the hidden layer(s) perform a non-linear mapping of the weighted sum of their inputs (i.e. the outputs of the previous layer). Depending on the implementation, the units of the last (output) layer perform a linear or non-linear mapping of the weighted sum of their inputs (i.e. the outputs of the last hidden layer). The output of the units of the last layer constitutes the ANN response to a given input pattern. In the specific case of the hydraulic motor plant, there is only one output unit giving the motor output torque. The transfer function of each unit, which determine the unit response at a weighted input, is the hyperbolic tangent function (tanh). Given a sufficient number of non-linear hidden units, an MLP is able to approximate any continuous mapping to an arbitrary degree of precision [39]. The ANN response is fitted to the input-output data patterns by adjusting the weights of the unit connections. The supervised learning approach changes (trains) the weights based on the error between the ANN and the actual plant output. In the dataset, the values of the plant output have been divided by a factor 200. This way, all the actual output values present in the dataset are in the [− 1, 1] range, and they can be compared with the MLP output.
The MLP is fed to the same 14 inputs used by the analytical model (Eq. 1), and is required to approximate the the torque T g (k) at time step k.

Hybrid model
A common divide-and-conquer strategy is to combine the two models described in Sects. 2.1 and 2.2 into one hybrid model (HM) [40,41]. In detail, the analytical linearised model described in Eq. (1) is first optimised. In absence of noise, the error (t) of the optimised analytical model consists of the unmodelled non-linear dynamics of the system. That is, for each pattern p, (p) is equal to the difference between the desired value T g (p) and the analytical model output T g (p): Rearranging (2), the desired output T g (p) can be expressed as: An MLP is then trained to reproduce the analytical model error ( MLP(p) = (p) , where MLP(p) is the MLP output). Once the MLP has learned to reproduce the error, the plant dynamics are obtained substituting the MLP output in (3).

Optimisation of models
The experimental rig shown in Fig. 1 was used to generate three sets (series) of data, each containing 1500 samples recorded over a 3 s span (2 ms sampling frequency). Each data sample consists of two input features and the corresponding model output, respectively the load input from the controller at time step t − 1 ( u m (t − 1) ), the load bearing position (angle) at time step t − 1 ( d (t − 1) ), and the output torque of the load-on motor at time step t ( T g (t) ). The three sets differ by the frequency of the load input u m generated by the controller, which is a sinusoidal signal of respectively 1, 7 and 9 Hz. The measured output torque is always within the range T g ∈ (−150, 150) N m , except for one outlier at T g = 194 N m.
The test set is composed of 1500 data patterns. It was formed concatenating the first 500 instances of each data series (1, 7 and 9 Hz). The training set contains the remaining 3000 data patterns (1000 patterns from each of the three sets). The difficulty of the data series is uneven, as some high-frequency disturbances are present in a small section of the training set.
The three models predict the output torque using information ( d , u m , T g (t) ) from the current and the previous four time steps (Eq. 1). For this reason, the first five time steps of each data set are used only as inputs to the model. That is, for each data set, 995 and 495 data instances are used respectively for training and testing the models.

Optimisation of the analytical model
For this task, three different global optimisation techniques were tested: EA, BA and PSO. These techniques use different nature-inspired population-based metaheuristics. All three algorithms use the same representation scheme, fitness evaluation function, initialisation procedure, and are parameterized to sample an equal number of solutions in the search space.
Namely, candidate solutions are encoded as real-valued numerical strings composed of the 9 parameters listed in Table 1. The fitness F of a candidate solution is calculated as the Root Mean Square (RMS) error of the model, that is: where N is the number of data samples (995 for training, 495 for final testing), T g (i) is the desired (plant) output for data pattern i, and M(i) is the model output for data pattern i. The initialisation procedure randomly draws the parameter values of the candidate solutions with uniform probability from the range of values defined in Table 1.
In order to obtain the same sampling of the search space, the EA, PSO, and BA used the same number of individuals and were run for the same number of learning cycles. Preliminary experiments showed fast convergence of the three optimisation procedures. For this reason, each run was limited to 100 generations, and the population size set to 100 individuals. These figures correspond to a total of 104 candidate solutions tested per optimisation procedure.
The above conditions ensure that the three algorithms search the same solution landscape (encoding and fitness function), start from comparable states (initialisation procedure), and are given the same search opportunities (equal sampling opportunities). To define a baseline of performance for the optimisation algorithms, RG was also implemented. RG uses the same encoding scheme, fitness evaluation function, and initialisation procedure as the EA, BA, and PSO algorithms, and is allowed the same sampling opportunities. RG randomly picks candidate solutions from the search space with uniform probability, and retains the best solution encountered. A summary of the parameters used in the optimization algorithms is present in Table 2.

Evolutionary algorithm
The EA [31] used in this study employed the Roulette Wheel [31] selection procedure to allocate the reproduction opportunities proportionally to the fitness of the solutions, and generational replacement with elitism [31] to renew the population. Genetic mutations were simulated by adding some 'noise' to each gene. This noise was randomly sampled from a normal distribution centred in zero. That is: where og i (j) is the ith gene of offspring j, pg i (j) is the ith gene of parent j, and N(0, j ) is a random number sampled from a normal probability distribution of mean 0 and variance j . The variance j is individual specific; it is initialised at 1 and undergoes adaptation (via mutation) like the other genes. The parameter K i rescales the magnitude of the mutation to the range of parameter i according to the following formula: where par max j and par min j are respectively the upper and lower bound of the ith parameter (see Table 1).
Cycles of selection and mutation are repeated until a given number of iterations has elapsed. Due to the lack of crossover and the adaptation of the mutation width, the EA employed in this study is close to the Evolutionary Programming paradigm [31]

Particle swarm optimisation
A standard Particle swarm optimiser [32] was employed in this study. At each iteration, the position (x) and velocity (v) of each particle i were updated according to the following equations: where , and are fixed parameters, 1 and 2 are random numbers drawn with uniform probability in the interval [0,1], PB is the best solution so far encountered by the particle, NB is the best solution so far encountered by the social neighbours of the particle, and t is the current iteration of the optimisation procedure. In this work the whole population was considered connected. That is, the number of social neighbours of an individual corresponded with the population size (100 individuals). The standard values = 1, = 2, = 2 were used.

Bees Algorithm
The standard version [33] of the Bees Algorithm was used in this study. The experiments were performed setting the number of elite sites ne = 2 , the number of best sites nb = 6 , and allocating nre = 20 and nrb = 10 foragers to perform local exploitative search respectively in the elite and best sites. The parameter for site abandonment was set to slim = 10.

Optimisation of neural network model
The MLP was trained using the Rprop [37] algorithm. Rprop is a fast gradient descent technique that uses only the sign of the gradient, ignoring its module. Rprop increases the step size of the adjustment of a given weight if the partial error derivative along that weight keeps the same sign from the previous iteration, and reduces it otherwise. The model inputs ( d and u m ) were normalised using the mean-variance procedure: where v and v are respectively the mean and standard deviation of the input v = d , u m .
The MLP structure was experimentally optimised, varying the number of hidden layers (1 or 2) and units per hidden layer (from 1 to 100). The weight decay parameter of the Rprop algorithm was also experimentally optimised. The following procedure was used to find the best ANN parameters. The training set (3000 patterns) was split into one training subset containing 2000 patterns, and one validation subset containing the remaining 1000 patterns. Different instances of MLP structures were trained on the training subset, and their performance was evaluated on the validation subset. For each parameter setting, the training procedure was repeated 20 times, and the results averaged. At the end, the parameter setting that provided the best results on the validation set was retained as the 'optimal' one. If the difference in the results obtained by two different settings was not statistically significant (Mann-Whitney test), the setting that minimised the MLP structure was preferred.
Using the above procedure, the MLP structure of the non-linear model was optimised to 2 hidden layers of 77 units each. The weight decay parameter was set to 5 × 10 −4 .
Once the MLP structure had been optimised, 20 different instances of the ANN were trained on the whole training set, and the final performance evaluated on the test set. Since the ANN was optimised and trained using only data from the training set, the results of the final evaluation on the test set are guaranteed to be unbiased.
Due to the relatively fast convergence shown in the preliminary experiments, the training procedure was stopped after 500 epochs.

Optimisation of the hybrid model
The hybrid model used the best analytical linear model obtained in Sect. 3.1. An MLP was optimised and trained to compensate (i.e. reproduce) the error (t) (Eq. 2) of the analytical model.
The procedure described in Sect. 3.2 was used to optimise and train the ANN. The optimal MLP structure for the hybrid model was found to be one hidden layer of 80 hidden units, with no weight decay in the Rprop procedure.
A summary of the parameters used in the MLP and the Hybrid Model can be found in Table 3.

Experimental results
The results presented in Tables 4 and 5 show the five number summary of the RMS error (5) attained by the different models on the training and test set. For each optimisation procedure, the tables summarise the results of 20 independent runs. The learning curves are plotted in Fig. 2. Pairwise two-sided Mann-Whitney tests were performed to assess the statistical significance of the differences among the results obtained. The significance level was set to 5%, and the p values are shown in Tables 6 and 7.
Overall, all algorithms gave very consistent results in the optimisation of the analytical model. The BA attained the smallest RMS error on the test set, whilst the PSO gave the most consistent performances. The pairwise Mann-Whitney tests (Tables 6, 7) revealed that the performance of the BA is significantly better than the performance of the other algorithms. In all cases, the RMS error on the test set is smaller than the RMS error on the training set. This result is due to the fact that the analytical model was not able to reproduce the high-frequency disturbances present in a small section of the training set.
The MLP achieved significantly higher accuracy results than the analytical model. The measures of the RMS error obtained on the training and test sets are comparable (Fig. 3), indicating that the ANN model was able to reproduce the high-frequency disturbances in the training set.
According to the results of pairwise Mann-Whitney tests (Tables 6, 7), the modelling accuracies obtained using the Hybrid Model are statistically comparable to those obtained by the MLP model.

Discussion
Three different approaches were tested to reproduce the dynamics of an electro-hydraulic load simulator. These approaches fit the model parameters to a set of experimental data points, measured at three different frequencies of the load input. The MLP-based non-linear model outperformed the baseline linearised analytical model. The RMS error of the non-linear model on the test data corresponds to about 0.5% of the output variable range, and is about 3.5 times smaller than the error of the linearised model. This result confirms that the non-linear terms in the load simulator  The four algorithms (RG, EA, PSO, and BA) used to optimise the baseline analytical model gave similar results: the difference between the most (BA) and least (RG) accurate model obtained is extremely small, amounting to about 3% of the RMS error magnitude. This result indicates that the current error is due to the limits of the linearised approach, rather than inadequacies of the optimisation procedures tested.
The RMS error of the linearised model is higher on the training set of examples, where some short high frequency disturbances are included. The results obtained using the MLP model on the training and test set are comparable. This suggests that the effect of the non-linear terms on the overall system dynamics is larger at high frequencies of the load input signal. In one further test, the training and test set were swapped. That is, the ANN was trained on the test set (not containing the high frequency disturbances) and its performance was validated on the training set (containing the high frequency disturbances). In this case (results not shown), there was a sharp drop in performance between the training accuracy and the test accuracy. The results of this test indicate that the high frequency disturbances contain non-linearities that are not present at low frequency.
Combining the analytical and MLP approach did not bring benefits in terms of modelling accuracy. The difference between the RMS error attained by the hybrid and MLP model is not statistically significant. Overall, the experimental results indicate that the non-linear MLP model is preferable in terms of simplicity and accuracy.
Comparing the proposed approaches with the existing literature is difficult because of the scarcity of quantitative results published, and the different data sets used. Perhaps the most similar tests were carried out by Zhang and Dong [21], who reported a tracking error within 5% of the loadon motor output torque. However, the present study is the only one where the model is optimised and tested at different load input frequencies. Conservatively, it can be said that the accuracy of the linearised model presented in this study is comparable to the accuracy of the stateof-the-art in the literature. The MLP-based model and the HM appear to obtain error magnitudes one order of magnitude smaller than the state-of-the-art in the literature.

Conclusions
The accuracy of analytical linear models of electrohydraulic load simulators is constrained by the extent of unmodelled nonlinearities, and the difficulty of optimising the numerous parameters that characterise them. A small number of nonlinear models has been reported in the literature, mainly based on ANN approaches. These ANNs were part of larger control systems, and whilst the performance of the whole systems was described, no information was given on the actual accuracy of the models.
This study used experimental data describing the dynamics of an electro-hydraulic motor, and fitted to these data the parameters of three kinds of models:  analytical, ANN, and hybrid. The study made the following contributions: • the upper limit of the accuracy of the analytical linear model was estimated using four intelligent parameter optimisation methods. • the difference in accuracy between the linear and nonlinear models was quantified and shown to be significant • it was shown a pure ANN model can be trained with accuracy comparable to a hybrid analytical + ANN model To the best of the authors' knowledge, this study is to date the most exhaustive attempt at identifying the dynamics of an electro-hydraulic load simulator, involving different kinds of linear and non-linear models and parameter optimisation algorithms. Despite its limitations, the accuracy of the optimised analytical model is comparable with the state-of-the-art linear models in the literature. The MLP and HM outperform the best models reported in the literature. Given the advantages of electro-hydraulic load simulators in terms of reliability, durability, and performances, the results of this study are deemed significant for the understanding and simulation of such systems. Further work should extend the experimental tests to higher frequencies of the load input, and more complex (e.g. recurrent [27]) ANN structures.

Funding
The authors received no financial support for the research, authorship, and/or publication of this article.

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.

Linearised analytical model: full description
The load system is driven by a servo valve, which is used to control the output torque. The flow equation of the servo valve of the load system is:  The output torque equation is: where G (N m/rad) is the rigidity of the connecting link; d (rad) is the turn angle of the motor drive shaft; The spool displacement of the servo valve and the input voltage satisfy the following relational expression: where u m is the load input from the controller, and K sv and K p are the correlation coefficients of the servo valve.
Combining the Laplace transforms of equations Eqs. (11) to (15), the following equation is obtained: In order to be able to identify the model using the discrete dataset of experimental measurements, Eq. (19) needs to be transformed into discrete form. This is achieved by performing a bilinear transformation of the following equation: where T is the sampling interval and s is the complex variable. The final discrete equation is obtained: where (