Recurrent neural networks for hydrodynamic imaging using a 2D-sensitive artificial lateral line

The lateral line allows fish to perceive a plethora of hydrodynamic phenomena, ranging from movements as small as the water disturbance generated by plankton, to movements as large as river streams (Coombs and van Netten 2005). This phenomenon has been called 'touch at a distance' (Dijkgraaf 1963), and has been associated with behaviours such as schooling, prey detection, rheotaxis, courtship, and station holding (Sutterlin and Waddy 1975, Partridge and Pitcher 1980, Hoekstra and Janssen 1985, Satou et al 1994, Montgomery et al 1997, Coombs and van Netten 2005). Since lateral line perception does not rely on light or sound, it is still effective in environments where these are absent, such as in low-light (e.g. in a dark cave or in blind species) or low-visibility (e.g. murky waters) environments. Fish are capable of detecting sources up to a distance of roughly one fish's body length (Kalmijn 1988, Coombs and Conley 1997, Ćurić-Blake and van Netten 2006).

Lateral lines are comprised of distributed sensors called neuromasts. Two types of neuromasts exist: superficial neuromasts and canal neuromasts. Superficial neuromasts are found on the skin surface all over the fish body, and their number present in one individual ranges from very few up to thousands, depending on the species (Coombs et al 1988). Canal neuromasts are found in fluid-filled canals under the outer skin. The canals themselves can usually be found along the trunk or near the head.

Many different types of ALL sensors exist, although fluid flow is generally measured through deflection of a cantilever structure (Liu et al 2016). These sensors range from Microelectromechanical systems (MEMS) e.g. Asadnia et al (2016) and ionic polymer-metal composites (IMPC) (Abdulsadda and Tan 2013b) to optical guides (Herzog et al 2015) and, recently, all-optical 2D-sensitive deflection sensors (Wolf et al 2018).

The current all-optical sensor also senses fluid flow via a cantilever deflection. Each sensor contains four optical fibres with fibre Bragg grating (FBG) strain sensing elements. These FBGs reflect a specific wavelength of light, depending on periodic variations in the fibre core. When an individual fibre is strained, the periodic variation stretches or contracts, resulting in a measurable wavelength shift (Flockhart et al 2003). The fluid forces acting on the sensor are thus encoded in the sensor deflection, which can be measured as a linear function of change in reflected wavelength peaks (Wolf et al 2018).

We discern between two different types of artificial neural networks: feed-forward networks and recurrent networks. The main difference between the two types of networks is the direction of the information flow.

In feed-forward networks, the input strictly moves forward, sequentially going through all layers until it reaches the output layer, which predicts the location of the source.

Conversely, recurrent neural networks have connections which allow the information to not only strictly go forward, but also allow connecting to the same layer, or even to a previous layer. Thanks to these recurrent connections, these networks can take history of the input into account, and are therefore often used when processing temporal data, such as time series.

2.3.1. State of the art

Artificial neural networks, such as the multilayer perceptron (MLP) used by Abdulsadda and Tan (2011) for locating a dipole source, are capable of learning complex, large, non-linear mappings from high-dimensional data sets.

Using simulated data, Boulogne et al (2017) compared the performance of an MLP, an echo state network (ESN) and the extreme learning machine (ELM) (Huang et al 2004) on localizing a moving source under different levels of signal-to-noise ratio. While the ELM outperformed the other algorithms in low noise conditions, the added complexity of the MLP and the recurrency of the ESN allowed these networks to take past velocity profiles into account, which may have caused this type of network to perform better with higher noise levels. This indicates that, under more noisy conditions (such as when using measured data), recurrent networks may be more suitable compared to feed-forward networks.

Specific to the all-optical sensors used in the current study (figure 1), only the ELM architecture has been applied. In Wolf and van Netten (2019), data from an ALL consisting of four sensors was processed by an ELM with 1500 hidden units. The ELM used a time window of 481 ms as the input, which was down sampled to 27 Hz. The combined system was able to predict the location of a 6 cm diameter sphere moving at 7 cm s⁻¹ with an average error of 3.3 cm in an area of cm.

Zoom In Zoom Out Reset image size

2.3.2. ELM

For our feed-forward neural network architecture, we make use of a variant of the efficient single layer ELM network. This type of neural network has a single tuneable hyper parameter, the size of the hidden layer, which allows for fast optimization while preventing over fitting. The ELM's input-to-hidden weights are randomized and fixed (figure 2(a)). This provides a hidden representation and leaves only the hidden-to-output weights to be learned (Huang et al 2006):

where denotes the desired (teacher) outputs.

Zoom In Zoom Out Reset image size

Using the Moore–Penrose generalized inverse , the network is able to find optimal weights in a single learning step. An added benefit to the Moore–Penrose generalized inverse is that it produces the smallest norm solution, which also aids avoiding over fitting. This smallest norm, similar to L2-regularization (Goodfellow et al 2016), penalizes large weights and thus penalizes (over-)dependence on single input nodes.

2.3.3. OS-ELM

One consequence of the single learning step, calculating the pseudo-inverse, is that this step becomes computationally intensive with larger data sets. Eventually, it is no longer computationally viable to perform the pseudo-inverse in one step due to hardware, i.e. working memory, constraints. The online sequential ELM (OS-ELM), which functions similarly to the ELM, can deal with larger data sets by implementing a sequential version of the least squares solution algorithm for the inverse operation.

Liang et al (2006) employ a sequential implementation by first calculating output weights based on an initial, sufficiently large, input batch. This initial batch produces a hidden activation or representation ; consequently, its optimal output weights are calculated using the left pseudo-inverse via:

This produces an initial weights estimate , optimized for this initial subset. Subsequently, the network is given new input samples in small batches, and the output weight vector is updated based on the discrepancy between the target values of the next subset and the estimate based on the previously found weights :

This second step repeats until all training samples have been processed. In the special case where the initial batch size equals the number of training examples, the OS-ELM implements the original ELM (Liang et al 2006).

2.3.4. LSTM

Data is gathered in a 1200 800 260 mm (w l h) water tank. A sphere with a diameter of 6 cm is moved horizontally in this tank by an adapted XY plotter (Makeblock) at a velocity of 13 cm s⁻¹ to create hydrodynamic stimuli. The plotter is controlled with an Arduino-compatible mainboard, which reports its location on set intervals.

The sensor array consists of eight sensors with an inter-sensor distance of 64 mm. One of the eight sensors malfunctioned (see figure 3(a)), and its data is therefore excluded from further analysis, functionally making it a seven-sensor ALL.

Zoom In Zoom Out Reset image size

Each sensor is protected by a perforated 32 mm plastic tube (figure 1) as a failsafe to avoid impacts with the object. The perforation pattern is chosen to be uniform with a ring of eight equidistant 9 mm diameter holes, sandwiched between two rings of eight equidistant 3 mm holes. This allows the bulk of the flow field to enact the sensor, while still providing structural integrity for protection. These perforated tubes may slightly affect the sensor sensitivity (Wolf et al 2018) and introduce a dampening effect; lowering the resonance peak magnitude and frequency, and reducing overall sensitivity.

The measured sensor deflections are obtained via an optical sensing interrogator (Micron Optics si255, Atlanta, USA) at a sampling frequency of 250 Hz. Each sensor is calibrated via a calibration procedure described by Wolf et al (2018). Given that the sensors were calibrated for fluid flow speed outside of their protective tubes, we opt to use the measured deflections as inputs for the neural networks. The expected measured velocity profiles are low frequency (<1 Hz) stimuli, we may therefore assume that the relation between the deflections and fluid velocity is close to linear (Wolf et al 2018).

Both the XY plotter mainbord and the optical sensing interrogator are controlled via a Matlab R2016a script. This allows synchronizing the sensor data with the stimulus location, providing accurate teacher labels for the neural networks.

Object movement is constrained to the (2D) xy-plane; its depth is kept constant during all experiments. The sphere is moved along 256 unique linear paths in a region of interest in front of the array, see figure 3(b). For each path, the sphere is moved back and forth 6 times between the left and right bound, creating a total of 3072 runs.

In addition to the recorded sensor signal, the plotter also provides timestamps for movement initiation, ramp-up, ramp-down, and movement end. These timestamps are used to synchronize and extract the parts of the signal where the velocity is constant; i.e. in between ramp-up and ramp-down. We combine the timestamps with the coordinates given to the plotter to automatically generate location labels for these segments with a constant velocity.

After each plotter movement, the plotter pauses for 20 s in order to minimize the effect of the previous movement on the hydrodynamic situation during subsequent movements.

3.2.1. Preprocessing

In order to determine the optimal network hyperparameters and preprocessing methods, we perform a parameter sweep for both. We first determine a baseline parameter set which gives a reasonable performance, through trial and error. Then, optimal values for all parameters were found by modifying one variable at a time while keeping others constant, which identifies the effect of that single variable on the performance of the complete pipeline. No grid-wise search is performed due to the excessive time it would take to test each combination.

We further limit the preprocessing parameter sweep to exclude effective time windows that exceed one second, for which we have several reasons. First, longer time windows effectively reduce the amount of training examples, making it harder for generalization to occur. Secondly, longer window sizes may include past data that contains less information and might also be irrelevant for the current location, making it harder to learn temporal patterns from the data. Finally, for time-sensitive tasks, having a smaller window enables a faster response from the neural networks.

3.3.1. LSTM

The LSTM network was implemented in Python using Keras (Chollet et al 2015), a deep learning framework running on top of Tensorflow (Abadi et al 2015). Our LSTM implementation used the Adagrad optimizer with a clipnorm of 1 (Goodfellow et al 2016), which limits calculated gradients such that they cannot exceed a maximum norm of 1, which helps in combatting overfitting.

Data is fed to the LSTM network in batches of windows. Each complete motion in the data set is split up into windows; the default size of the windows is 30 samples, but this value is varied during the parameter sweep. The default window thus has a size of , for the x- and y - deflections of each of the 7 functioning sensors. Each window is accompanied by a single teacher label corresponding to the ground truth location of the source at the end of that window, creating a many-to-one mapping.

During the parameter sweep, we evaluate the LSTM performance through k-fold cross-validation, with k  =  5; this results in a train/test/validation split of 80/10/10%. The order of the 3072 recorded motions is shuffled and each complete motion is assigned to one of the sets.

Based on the parameter sweep (table 1), the optimal LSTM network architecture has 100 hidden nodes, a learning rate () of 0.05, a learning rate decay of , no dropout, and uses the hyperbolic tangent activation function (act. fun.). The signal processing pipeline for the LSTM is described by a window size (win sz.) of 15 samples, a stride of 2, a down sampling rate (ds. rate) of 16 (giving the a sampling frequency of 15.625 Hz) and uses, z-score per sensor normalization (norm.).

Table 1. Parameter sweep for the LSTM and signal processing. Indicated here are the default parameters and the variation options for each parameter. From these options, the optimal choices are highlighted in bold.

Parameter	Default	Options
Nodes	50	10, 20, 50, 100, 200
	0.05	0.1, 0.05, , ,
Decay		0, , , ,
Dropout	0.20	0, 0.1, 0.2, 0.35, 0.5
Act. fun.	ReLU	ReLU, sigmoid, tanh
win sz.	30	1, 2, 4, 8, 15, 30
Stride	1	1, 2, 4, 8, 15, 30
ds. rate	8	4, 8, 16
Norm.	4 (z-score)	1, 2, 3, 4, 5

3.3.2. OS-ELM

For the OS-ELM, we similarly divide complete motions from the data set into separate windows, but we also flatten this window to a 1D vector. For the initial training phase, we create a batch of windows, where N is the number of hidden units. For the sequential training steps, the batch size was set to 64 windows.

Since there is no validation set needed for the OS-ELM training phase, we apply an 80/20% train/test split. Again, with five-fold cross validation, we shuffled the order of the data set and assigned each fold to their respective train or test set.

Table 2 shows the results for the OS-ELM parameter sweep. Based on the parameter sweep, our final ELM network architecture had 30 000 hidden nodes, and uses the ReLU() activation function. The processing pipeline includes a window size of 15 samples, a stride of 1, a down sampling factor of 16 (giving a sampling frequency of 15.625 Hz), and uses z-score normalization.

Table 2. Parameter sweep for the OS-ELM and and signal processing. Indicated here are the default parameters and the variation options for each parameter. From these options, the optimal choices are highlighted in bold.

Parameter	Default	Options
Nodes	5000	{1, 2, 5, 10, 20, 30}
act. fun.	sigmoid	ReLU, sigmoid, tanh
win sz.	30	1, 2, 4, 8, 15, 30
Stride	1	1, 2, 4, 8, 15, 30
ds. rate	8	4, 8, 16
Norm.	4 (z-score)	1, 2, 3, 4, 5

The constant speed section of each motion lasts between 4.3 s for horizontal and 4.5 s for diagonal motion. This results in a data set of input examples for the OS-ELM. For the LSTM framework, given the stride of 2 rather than 1, the data set is roughly half that size.

Although the number of hidden nodes for the OS-ELM may seem high, the number of input examples (windows) is roughly 7 times higher. This makes it unlikely for overfitting to occur.

We create a separate data set in order to systematically determine the influence of noise. In particular, we are interested in the localization robustness of the two types of neural networks with respect to different noise levels. This data set consists of 60 repetitions of both backwards and forward motion along two unique paths parallel to the array, at 6.15 cm and 8.15 cm distance, respectively, producing 240 recorded motions. These distances are indicated via red dashes in figure 3(a).

For augmenting this data set, we employ a controlled method of adding noise to our data; we want to stay as close to the real environment as possible, while still being able to quantify how much noise is added. We therefore first record a large data set without any stimulus or plotter activity, from which we sample noise specific to each sensor. We augment the data set by combining recordings with the following nine noise augmentation amplification levels: {0, 0.1, 0.25, 0.5, 1, 1.5, 2, 3, 4}, producing a final augmented set of 2160 motions. This range is selected on a pseudo-logarithmic scale to allow testing a large range of noise levels.

We keep the previously determined optimal hyper parameters for both type of neural networks and thus do not re-optimize for this augmented data set. For determining the performance of each network, we employ five-fold proportional stratified sampling; each of the five folds has an equal amount of randomly picked complete motions from each of the nine augmented sets. The training data is shuffled before using these folds for training the networks.

In the first series of experiments, we compared the performance of a feed-forward neural network (OS-ELM) to the performance of a recurrent neural network (LSTM) on the task of localizing a source using measured artificial lateral line (ALL) data.

Both networks were optimized with respect to their hyper parameters and preprocessing methods using five-fold cross validation, as described in the Methods section.

Figure 4 shows an example of a single run, where both networks predict the location of an object moving parallel to the array. This figure shows that the LSTM network output is less erratic and closer to the target path compared to the ELM prediction.

Zoom In Zoom Out Reset image size

When we combine and interpolate between all locations in the region of interest, we can show a spatial mapping of the localisation error. Figure 5 shows the interpolated error magnitude at each location. From these spatial plots, we observe that the localisation error further away from the array increases for both neural networks.

Zoom In Zoom Out Reset image size

The actual distribution of localisation errors, regardless of relative location, is shown in figure 6(a). Here, we binned the localisation error for both types of networks, with a limit of 15 cm. The average and standard deviation of the localisation errors are cm and cm for the LSTM and OS-ELM, respectively. While the LSTM errors are concentrated around the average, the OS-ELM has a large tail in the distribution (figure 6(a)), with 4.2% of the errors exceeding 15 cm (not shown).

Zoom In Zoom Out Reset image size

These results show that the LSTM network significantly outperforms the OS-ELM network on this 2D source localisation task using measured velocity profiles from a 2D-sensitive array.

To further determine the suitability of both neural networks for other circumstances or other types of artificial lateral lines, we also investigated the effect of noise with respect to localisation.

In our performance analysis, we differentiate between noise augmentation levels within the test set. Figure 6(b) and table 3 both provide an insight in the effect of added noise to the overall performance. As the noise level increases, so do the LSTM and OS-ELM network's average localization errors, although the LSTM seems relatively less affected by higher noise levels.

Table 3. Error per noise level for LSTM and OS-ELM networks after stratified five-fold cross-validation. Shown here are the mean Euclidean error and standard deviation. Finally, the average error over all noise levels is shown.

Noise level	LSTM error ( cm)	OS-ELM error ( cm)
0.00
0.10
0.25
0.50
1.00
1.50
2.00
3.00
4.00
Average

Figure 6(b) shows box plots of OS-ELM and LSTM performance, separated per noise level. For all noise levels considered, the median LSTM error is significantly lower than the OS-ELM error.

Table 3 shows the overall and per-noise-level error for both networks. On the lowest noise level (no noise), the LSTM achieved an average error of 0.55 cm, whereas the ELM had an average error of 1.75 cm. For the highest noise level (4), these errors increase to 1.10 cm and 5.64 cm for the LSTM and ELM networks, respectively. Again, the error distribution of the OS-ELM shows a long tail towards higher errors, making it less reliable than the LSTM in this comparison.