1. Introduction

The human nervous system is an efficient and powerful ‘biological computer’. It consumes 20 Watts of power, $<1\%$ of the power consumption of the average kettle and can process an enormous stream of information in real-time [1–3]. Key to the vast capability of the human nervous system lies is how and where such information is processed. An electronic computer processes information at a centralised processing unit (von Neumann architecture) with a separate memory device. The human nervous system is, however, based on a distributed processing architecture in which each neuron “remembers and processes” a stream of (electro-chemical) signals as it travels to the brain for “inferring”. Inspired by the distributed computing architecture of the human nervous system, several neuromorphic computing architectures have been developed [4–7]. The neuromorphic computing system aims to realise the distributed computing architecture of the human nervous system as hardware, enabling a hardware-level integration of Artificial Intelligence. Implementation of neuromorphic computing as a photonic system could potentially lead to a computing system which is 6-8 orders of magnitude faster and 4-8 orders of magnitude more efficient than digital electronics [6].

The photonic reservoir computer (PhRC) is a relatively new member of the neuromorphic computing family implemented as a physical photonic/optical system. Implementing a reservoir computer in photonics leads to a design for very high speed processing exploiting the availability of ultrawide bandwidth optical communication commponents. It offers a fundamentally different approach to other neuromorphic photonic systems, in terms of both architecture and how the system is trained (or optimised) [4–8]. PhRC only trains/optimises the Read-Out layer. The kernel of the PhRC is not trained/optimised, thus retaining the characteristic of the kernel as a random and semi-chaotic system. Other neuromorphic systems require both kernel and read-out training/optimisation [4–7]. Consequently, the PhRC offers a more straightforward implementation and scalability than other neuromorphic photonic architectures [4–8]. Several implementations of PhRC have been reported, for example as a swirl configuration of semiconductor amplifiers [9,10], or as a photonic crystal cavity [11–13]. In this paper, we consider an electro-optical reservoir computer (EORC) implementation with a single delay line and a Mach-Zehnder modulator (MZM) as the non-linear neuron activation function [14,15].

Humans’ powerful computing capacity has mainly been used to interpret sensory information, $~99.99\%$ of the information processing capacity is for sensory processing while $<0.01\%$ is for attention, decision-making, motion and language [16]. Motivated by such a fact, in our previous work [12,13,17], we have numerically demonstrated a new methodology to process ultrafast high-throughput molecular vibrational signals (spectroscopic) using a chaotic photonic crystal billiard system as its kernel. The present paper reports an alternative method which is fully-bio-inspired, meaning that the sensor system is biomimetic of the human eyes and its output will be processed through a neuromorphic EORC system. Numerical simulation successfully demonstrates that both chemical discrimination and concentration determination tasks can be achieved.

The paper is organised as follows: Section 2 describes the system as a whole and its composition, namely the biomimetic sensor and the EORC and the modelling approach used. Section 3 discusses the results obtained. The first part considers the dynamics of the EORC system, which is modelled using a realistic model considering the filtering effects observed in a practical setup. The second part demonstrates the application of such a system to perform the chemical discrimination and concentration determination tasks. The impact of EORC parameters on the accuracy is investigated and discussed. Some concluding remarks are then provided in Section 4.

2. Biomimetic sensing by neuromorphic signal processing

Figure 1(a) schematically shows the apparatus considered in this paper. It is comprised of a biomimetic sensor serving as the sensing component, an electro-optic reservoir computer (EORC) serving as signal processing component and a read-out component serving as an inferring component of the EORC output. Figure 1(b) schematically depicts the implementation of Fig. 1(a) using a single delay-line configuration and a Mach-Zehnder modulator (MZM) as the nonlinear neuron activation kernel. The biomimetic sensor supplies a three-channel reading to the arbitrary waveform generator which performs time-multiplexing on the sensing information for input to the EORC system.

 

Fig. 1. Schematic illustration of the system as a whole; (b) Implementation of (a) in a single physical node configuration with electro-optic with Mach-Zehnder modulator (MZM) as nonlinear neuron activation kernel; The multi-channels output of the biomimetic sensor has been representatively shown in a Red-Green-Blue (RGB) colouring scheme. (c) Signal flow diagram translating the flow of signal in (b); (d) Schematic of the architecture of algorithmic reservoir computer.

Download Full Size | PDF

2.1 Biomimetic sensor system

This sub-section describes the operating principle and the modelling procedures of the sensing component used in this work. It is inspired by the biomimetic design which previously has been demonstrated successfully for chemical [18–20] and biological sample [21] detection. The biomimetic sensor mimics the operation of the human eyes, which are capable of discriminating over 100,000 shades of (visible) colour over an environmental background despite being equipped with only three main types of cone cells [22,23]. Key to this ability is the spectral sensitivity of these cone cell types which are spectrally broad and overlapping [24]. Research on colour sensitivity of mantis shrimp (Haptosquilla trispinosa), which has 12 different types of photosensitive cell covering different light wavelengths and polarisations, however, found that higher number of unique photoreceptors does not mean an improvement in colour discimination [25]; Mantis’s spectral discrimination $\Delta \lambda$ is in the region of 15 to 25 nm in comparison to humans’ which is in the region of 1 to 8 nm [25,26] over the visible light region. In [18–21], a trichromat, three-colour, biomimetic sensor has been constructed in the mid-infrared range of light and a multi-variate statistics tool based on a patented comparative-discrimination spectral detection (CDSD) technique [20,27] has been used to perform the vectorial signal processing. In the present work, we have adopted a biomimetic sensor with a similar trichromat architecture, but the EORC is used for signal processing. The EORC will be described in sub-section 2.2.

As mentioned, the biomimetic sensor used in this paper comprised of three channels whose spectral responses are different, see Fig. 2. Figure 2(a) schematically shows the implementation of the biomimetic sensor constructed in [18–20]. The photodetectors are equipped with a filter at different central wavelength, in Fig. 2(a) the filters have been representatively shown in Red-Green-Blue (RGB) layer in front of the detectors. In the present work, we have specifically considered that each channel has a Gaussian spectral response peaking at different wavelengths, i.e., 3.37, 3.45 and 3.53 $\mu$m and with a half-power bandwidth (HPBW) of 60 nm. Such sensors are readily available, for example as dielectric filtered photodetectors [28]. The specific parameter values of the sensor have been chosen as such they are broad and spectrally overlapped [18–21], see Fig. 2(b). The choice of peak spectral response around the 3 $\mathrm {\mu }$m of the mid-infrared range is motivated by the richness of carbon-hydrogen stretching fingerprints in this region [29]. Note that in previous works [18–21], the sensor has also operated at similar wavelengths of interest.

 

Fig. 2. Schematic design of the three-channel biomimetic sensor; (b) the absorption spectra $\alpha$ and the filter response $\mathcal {D}_m(f)$ of the photodetector (PD).

Download Full Size | PDF

For demonstration, in this work methanol, ethanol and isopropanol and their mixture will be considered as chemicals to be detected. The complex-valued frequency-dependent refractive indices of the pure substances have been taken from [30] whereas the refractive indices of mixture are calculated using the Lorentz-Lorenz rule [31–33] by,

2.2 Electro-optical reservoir computer (EORC)

In the previous sub-section, the operating principle and model description of the biomimetic sensor was described. In this sub-Section, its physical implementation using electro-optical (EO) modulator system and its numerical model are described. The EORC is a physical implementation of the reservoir computer (RC) computational framework which is a generalisation of a traditional artificial neural-network (ANN) based on the recurrent neural-network (RNN) [7,37–39]. For detail on the framework of the RC, the reader is referred to [39]. Figure 1(d) shows a schematic of the architecture of an RC scheme. An RC system consists of three main parts, namely an input, a kernel and a Read-Out layer. The origin of the name ‘reservoir computer’ comes from the architecture of the RC which features a ‘kernel’ that is a reservoir containing a large number of neurons.

Figure 1(b) schematically shows the EORC system considered which is a combination of optical (in blue) and (electronic) radio frequency (RF) (in black) systems. The EORC is a feedback system in which the signal is delayed by a time $\tau$ due to the presence of the long delay line. The delayed feedback signal is combined with the present-time input signal to modulate the MZM completing the feedback loop. The following sub-sections describe the implementation of each part of the RC as an EORC system.

2.2.1 Input node by signal time multiplexing

The role of the input node is to distribute the input signal $u(t)$ to the reservoir kernel. In the algorithmic reservoir computer case, one introduces a numerical input weight $\mathbf {W}^{\mathrm {in}}\in \mathbb {R}^{N_x\times N_u}$ which is a matrix of size, $N_x\times N_u$, where $N_x$ and $N_u$ are the number of neuron nodes and the number of input signals, respectively. For the translation system as an EORC, that is considered here, neuron nodes are created artificially by means of time multiplexing and masking during the information encoding process; such artificial temporal nodes have also been referred to as virtual neuron nodes [14]. Practically this process done by the arbitrary waveform generator (AWG), see Fig. 1(b).

The signal time-multiplexing is achieved by performing a sample and hold operation, that is the measurement from each channel is serialised for a duration of $\tau$. For the three-channel output biomimetic sensor considered in this paper, i.e., $\mathbf {P}(n)=[P_1;P_2;P_3 ]^{\mathrm {T}}$, the serialised signal $j(t)$ is illustrated in Fig. 3(a). This time-multiplexing process allows us to encode the sensor information as a time-series signal.

 

Fig. 3. (a) Serialised information signal $j(t)$ which is based on the outputs of the biomimetic sensor $P_1$, $P_2$ and $P_3$; (b) Illustration of the mask signal $m(t)=m(t+\tau )$; (c) Input signal generated upon application of mask on serialised information signal; (d) Illustration of the time-demultiplexing process to generate virtual nodes.

Download Full Size | PDF

Subsequently, a periodic mask signal $m(t)=m(t+\tau )$ is generated whose values are randomly chosen from some random generator with mean of 1. This mask is piece-wise constant over the period of $\theta =\tau /N_x$. Here, the number of mask points is equal to the number of neuron nodes $N_x$, because, as mentioned, the mask serves as the serialised input weight of the information signal with the virtual node. In practice, one can increase the number of masks $N_x$ up to the sampling speed of the digitising oscilloscope, without any penalty on EORC system’s latency. If a further number of masks is necessary, the length of delay time $\tau$ can be extended which consequently introduces additional latency into the EORC scheme. However, as it will be demonstrated in Section 3 an EORC scheme operating with $\tau =1$ $\mu$s and $N_x=50$ leads to a moderate mask temporal length $\theta =0.02$ $\mu$s which can easily be resolved by a laboratory standard oscilloscope with typical 1 GSa/s performance. Fig. 3(b) illustrates the mask signal for $N_x=10$. The input signal is obtained by application of mask signal to the time-multiplexed signal, i.e., $u(t)=m(t)j(t)$, see Fig. 3(c). The subsequent process applied upon input signal $u(t)$ by the ‘neuron’ in the reservoir is described in the following.

2.2.2 Kernel by a single delay-line non-linear node

The input signal $u(t)$, here simulated by the Arbitrary Waveform Generator, will be combined with the delayed feedback signal $x(t-\tau )$ and transformed via a non-linear modulation process occurring in the neuron kernel. Throughout this paper, we consider a neuron kernel implementation based on the single physical node delay line with a Mach-Zehnder modulator (MZM) as the non-linear neuron activation function, see Fig. 1(b). Implementation of a RC based on this system and its variant have been reported previously in the context of all-optical computing and communication applications [14,40]. Initially, a similar setup to Fig. 1(a) was proposed by Ikeda [41–43] to demonstrate the existence of multiple stability state and chaos in a non-linear optical ring resonator system, mathematically, modelled by the Ikeda map given by,

Equation  (5) is the well-known non-linear transmittance of a MZM due to two-travelling wave interference process [44]. Parameter $I_0$ is the intensity of continuous-wave laser which is being modulated by the driving RF signal $v(t)$. Parameter $V_\pi ^{\mathrm {RF}}$ is the half-wave voltage of the MZM i.e. the voltage which would phase shift the light by $\pi$. Parameter $\Phi _0$ is the initial phase difference between the two-travelling wave in the MZM. Commercial MZMs are usually fabricated with $\Phi _0\approx -\pi /4$, and so operates near the (almost) linear region of Eq. (5); this operation is commonly referred to as the ‘Quad’ point. For numerical simplicity, here, the half-wave voltage is taken as $V_\pi ^{\mathrm {RF}}=\pi /2$.

The original Ikeda map Eq. (4) has been adapted to model a realistic experimental setup [15,45]. Specifically, low-pass and bandpass filter effects are included to describe the inclusion of the photodetector and RF signal amplifier, respectively [15,45]. This leads to the following delay differential equations (DDEs),

The DDEs (Eqs. 6,7) have been represented as a signal flow diagram in Fig. 1(c). Figure 1(c) translates the mathematical operation performed at each stage of the signal’s flow and uses the notation that will be used throughout this paper. Parameter $\gamma$ has been added to represent the effective gain of the system. Operator $h(t)$ represents the time-domain operation of the filter effect, implementing DDE in Eq. (6) or Eq. (7). In this paper, we have implemented a time-stepping algorithm, performing the integration of Eqs. (6) or (7) up to second-order accuracy, based on the $\mathcal {Z}$-bilinear transformation technique. This technique has previously used in conjunction with the transmission-line method to model non-linear photonic materials described by a Duffing model [48,49]. Specific implementation of the filters Eqs.  (6) and (7) is provided in Appendix 5.1. Other methods [43,50,51] could also be used.

The activation state signal $x(t)$, as the observable, is monitored and demultiplexed for further weighting by the Read-out layer. The subsequent sub-section describes the time-demultiplexing process, and the weight applied to the activation state signal.

2.2.3 Read-Out layer

To allow weighting by the Read-Out layer, the activation state signal $x(t)$ is sampled and demultiplexed. Figure 3(d) schematically illustrates the time-demultiplexing process. The neuron activation state signal $x(t)$ is sampled with a sampling-time of $\theta$ and spun with a winding period of $\tau$. The demultiplexed activation signal is partitioned as $\mathbf {x}(\tilde {t})=[x_1 (\tilde {t}),x_2 (\tilde {t} ),\ldots,x_k (\tilde {t})]^{\mathrm {T}}$ where $k=1,2,3,\ldots,N_x$ is the index of the virtual nodes and each entries of matrix $\mathbf {x}$ are given by $x_k (\tilde {t}\tau +k\theta )=x(t)$ with $\tilde {t}=0,1,\ldots$ being the winding index.

In the following step, the demultiplexed signal is weighted by the read-out layer,

3. Results and discussion

In Section 2, the sensing component of the system in the form of a biomimetic sensor with 3 output channels, and the architecture of the EORC system and its model which is used to process the sensor signal have been described. This Section starts with results and discussion of the dynamics of the single delay line feedback system (sub-sections 3.1) before describing in sub-sections 3.2 and 3.3 the application of such a system as an EORC for chemical sensing detection.

3.1 Dynamics of electro-optic reservoir computer under limited bandwidth operation

A delay feedback system, including the Ikeda map [43], is known to exhibits chaotic behaviour. This sub-section will discuss the influence of the filtering effects that would be present in a practical EORC system, due to the inclusion of photodetector and RF signal amplifier. For this purpose, we consider the case that there is no information signal to be processed, i.e., $u(t)=0$.

Figure 4(a) shows the bifurcation diagram of the original Ikeda map (Eq. (4)). The bifurcation diagram shows all possible values for neuron activation signal $x(t)$ at the steady state as a function of laser intensity parameter $I_0$, also see Fig. 1(c) for notations used. Specific parameters used for Fig. 4(a) are $\tau =1$ $\mu$s, $\gamma =1$ and $\Phi _0=-0.85\pi$ similar to that use in [45, Fig. 3]. Such a delay $\tau =$ 1 $\mu$s can, in practice, be achieved by using a $\approx$200 m spool of a standard SMF-28 fibre. Signals $x(t)$ for intensities of $I_1=1.5$, $I_2=1.9$ and $I_3=2.1$ are shown in Fig. 4(c,d,e); these specific intensity points are marked in Fig. 4(a) for reference. When operating at $I_1$, the neuron signal is single valued around 1.5, while for the case of $I_2$ and $I_3$ the signal $x(t)$ are periodic and multivalued, i.e., two and four stationary states, respectively. Figure 4(a) further shows that while at low intensity the steady state of neuron activation signal $x(t)$ is single valued, it splits as $I_0$ increases; the splitting point is called the bifurcation point. Furthermore, Fig. 4(a) shows that by increasing the parameter $I_0$ the steady state solution splits further and begins to behave chaotically, also see Fig. 4(f) for the case of $I_5=2.3$ and also, noting the bistable response at $I_0=1.16$. For detail on the bifurcation and bistablilty phenomena, and the route to chaos of the Ikeda map system, see [45].

 

Fig. 4. (a,b) Bifurcation diagram of Ikeda map (4) and (6), respectively. (c,d,e,f,g,h,i) The temporal signal $x(t)$ at the steady-state for different laser intensity parameters of $I_1=1.5$, $I_2=1.9$, $I_3=2.1$, $I_4=2.14$, $I_5=2.3$, $I_6=2.38$ and $I_7=2.95$ for Ikeda map in Eq. (4). Subfigure (j) is an enlarged view of (n). (k,l,m) Temporal signal $x(t)$ at the steady-state for $I_4$, $I_6$ and $I_7$ for low-pass feedback Ikeda system in Eq. (6). (n) Transient state of $x(t)$ for Ikeda map in Eq. (4) in blue and low-pass feedback Ikeda system in Eq. (6) in red. For convenience, the intensity parameters considered are marked in (a,b).

Download Full Size | PDF

The bifurcation diagram for the low-pass feedback (Eq. (6)) case is shown in Fig. 4(b). Figure 4(b) used the same parameters as for Fig. 4(a) in addition to the high-frequency cut-off (-3dB) parameter of $f_H=10$ GHz; such a high cut-off is currently possible due to the availability of ultrafast photodetectors, e.g., InGaAs detectors [46]. It is noted that in the absence of any filtering effect, the original Ikeda mapin Eq. (4) operates at infinite bandwidth and as such the response of the system is instantaneous. In contrast, the DDE (Eq. (6)), i.e., low-pass feedback system, operates with a limited bandwidth. Note that DDE (Eq. (6)) is a causal system whose response is dependent on a past event. Visual comparison of Fig. 4(a) and Fig. 4(b) shows that the low-pass filter effect has small impact for low intensity $I_0<2.1$. The temporal signals $x(t)$ for intensity parameters of $I_1$, $I_2$ and $I_3$ are also overlapped with that of the original Ikeda map in Fig. 4(c,d,e). Deviation from the original Ikeda map in the low-pass feedback case is noted at a high-intensity whose value depends on the cut-off frequency $f_H$; for $f_H=10$ GHz in Fig. 4(b), this occurs at $I_4\approx 2.14$. Comparing Fig. 4(g) and 4(k), there is a noticeable presence of ripples at the stationary point in the low-pass feedback case but not in the infinite bandwidth case. This ripple originates from the causal response of the low-pass filter which is not instantaneous resulting a dispersion effect of the signal. This rationale is confirmed by analysing Fig. 4(n,j). Figure 4(j) is an enlargement of Fig. 4(n) which shows the transient of $x(t)$ when operated at $I_4$. It shows that, in the low-pass feedback system, signal $x(t)$ exhibits a typical dispersive ripple at signal’s transition whose impact accumulates as $x(t)$ evolves over time.

The bifurcation diagram for the case of the bandpass feedback system is shown in Fig. 5(a); the temporal signal $x(t)$ for different intensity parameters $I_i$, $i=1,2,3,4,5,6$ are shown in Fig. 5(b-g). The bandpass filter parameters are $f_L=10$ kHz and $f_H=10$ GHz, for the low and high-frequency cut-off, see Appendix 5.1. In addition, the MZM is operated at the Quad operational condition, i.e., $\Phi _0=-\pi /4$. At low intensity $I_0$, signal $x(t)$ reaches a steady state value of $x(t\rightarrow \infty )=I_0/2$. This is because bandpass filter removes the zero-frequency (DC) component of signal $v(t)$, as such $v(t\rightarrow \infty )=0$ leads to $x(t\rightarrow \infty )=I_0 \cos ^{2} (\pi /4)$. Visual inspection of Fig. 5(b–d) shows that the time required for the system to reach steady state becomes longer as the intensity increases. Then, at $I_0=1$, the neuron activation signal begins to oscillate. Figure 5(e) which depicts the signal $x(t)$, at $I_0=1.2$, oscillates at the steady state and is exhibits the skewed waveform typically found in strong non-linear media due to the presence of high-order harmonics [53].

 

Fig. 5. Bifurcation diagram of the bandpass filter feedback of Eq. (7). Blue bullets points are the ‘fixed points’ solution of (7), see Appendix 5.2. The line $x=I_0$ is included for reference; (b,c,d,e,f,g) Temporal signal $x(t)-I_0/2$ at different intensity parameter of $I_1=0.5$, $I_2=0.7$, $I_3=0.9$, $I_4=1.2$, $I_5=2$ and $I_6=3$. These intensity parameters also are marked in (a) for reference.

Download Full Size | PDF

Furthermore, visual inspection of Fig. 5(e) shows that the signal $x(t)$ oscillates at longer period $\approx$ 25 $\mu$s than the round-trip period $\tau =1$ $\mu$s. To investigate the origin of this low-frequency oscillation, the spectrum of $x(t)$ is displayed in Fig. 6(a) for different low-frequency cut-offs of the bandpass filter, namely $f_L=10$, 20 and 30 kHz. The resulting temporal signals $x(t)$ are shown in Fig. 6(b-d), respectively. Figure 6(a) displays different characteristic behaviour in low- and high- frequency regions. The high-frequency region is characterised by the resonance of the feedback system whose fundamental frequency occurs at $1/\tau =1$ MHz and harmonics of which are produced from the nonlinear process of the MZM. Note that these harmonics are similar regardless the different $f_L$ parameters.

 

Fig. 6. Impact of low-frequency cut-off $f_L$ on the dynamic of bandpass feedback system beyond the bifurcation point. (a) Spectra of $x(t)$ for different cases of $f_L=10$, 20 and 30 kHz. The temporal evolution of $x(t)$ is shown in (b,c,d), respectively. Note that in (b,c,d), $x(t)$ is shown from $t=0$ to steady state. In all cases the intensity parameter is $I_0=1.2$.

Download Full Size | PDF

The characteristic behaviour in the low-frequency region is, however, significantly different for the different cases of $f_L$. Figure 6(a) shows that the fundamental frequencies of the low-frequency region are increasing, i.e., 38, 55, 67 kHz for increasing $f_L=10$, 20 and 30 kHz, confirming the origin of the low-frequency oscillation as the interplay of the $f_L$ parameter and the nonlinear feedback effect. Analytical investigations on similar bandpass feedback systems [54,55] have also noted the existence of this novel interplay effect, which is attributed to the stability of the zero ‘fixed point’. As such, the fixed points of (7) are shown in Fig. 5(a) as bullets points; Appendix 5.2 describes the fixed points and the bifurcation point of Eq. (7) for Quad operation. For $\gamma =1$, the bifurcation point, which occurs at $I_0=1$, marks the beginning of such oscillation in $x(t)$. Further increase of intensity leads to the generation of more harmonics which manifest themselves as ripples, as shown in Fig. 5(f). Similar to the case shown in Fig. 4(a) and 4(b), further increase of intensity leads to chaotic dynamics of the system. Although not shown in Fig. 5(a), we note the existence of a metastable regime at $I_0\approx 5$ after the chaotic regime.

3.2 Sample discrimination task

In sub-section 3.1, the impact of filter (low-pass and bandpass) on the dynamics of the EORC was discussed. In this sub-section, the EORC is trained and tasked to discriminate three chemicals, namely methanol, ethanol and isopropanol.

To achieve such a task, the input signal has been generated following the scheme described in Section 2.2. Moreover, the teacher signal $\mathbf {y}_\mathrm {teacher}$ has been defined using one-hot-encoding, as

Optimum operation of the EORC, which leads to the minimum normalised mean square error (NMSE) $\varepsilon$, occurs around the bifurcation point $I_0=1/\gamma$. It is because that around the bifurcation point the EORC is sensitive to perturbation as such is suitable for sensing/detection application. To show the relation between the intensity and gain parameter for the optimum operation of the EORC, Fig. 7(a) depicts the error $\varepsilon$ on a logarithmic scale as a function of intensity $I_0$ and effective gain $\gamma$ parameters for $N_x=30$, which is calculated by

 

Fig. 7. Impact of EORC parameters $I_0$ and $\gamma$ to the accuracy of the trichotomy task. (a) The NMSE $\varepsilon$, in logarithmic scale, is shown as the colour as function of $I_0$ and $\gamma$. The Bifurcation points are included for reference. (b) Error as function of intensity parameter $I_0$ for gain parameter $\gamma =5$; (c) Impact of number of nodes $N_x$ to the NMSE $\varepsilon$.

Download Full Size | PDF

3.3 Chemical concentration determination task

In sub-section 3.2, the EORC has been trained and tasked to detecting pure methanol, ethanol and isopropanol. In this sub-section, the EORC is trained and tasked to determine the concentration of methanol, ethanol and isopropanol mixture. To achieve such a task, the teacher signal in Eq. (9) is modified to

Here, the EORC is only trained to recognise the concentration of few possible combinations. Figure 8(b) shows the training points used which are spaced by 20$\%$ within the ternary diagram. This specific task highlights the key capability of the artificial neural network system to perform data-interpolation and regression tasks.

 

Fig. 8. The EORC application for concentration determination task. (a) RGB pseudo-colouring of the output of the biomimetic sensor $\mathbf {P}(n)=[P_1;P_2:P_3]^{\mathrm {T}}$ for different concentration of chemical mixture shown as a ternary diagram. Note that for visual contrast $\mathbf {P}(n)=[P_1;P_2;P_3]^{\mathrm {T}}$ has been normalised, see sub-Section 2.1. (b) Training samples of the EORC at discrete concentration combinations. (c) RGB pseudo-colouring of the steady state output signal $\mathbf {y}(t)$ for different concentration of chemical mixtures.

Download Full Size | PDF

To find the optimum operation of EORC, the error $\varepsilon$ is shown as a function of intensity $I_0$ and gain $\gamma$ parameters in Fig. 8. It is similar to Fig. 8(a), with the difference being that it is based on the training at points shown in Fig. 8(b). From Fig. 8(b), the optimum operational parameters of the EORC interestingly are observed to occur at the low-intensity and high-gain $\gamma$ parameter along the bifurcation line $\gamma =1/I_0$. As such, for this task the following parameters of the EORC will be used, $I_0=0.197$, $\gamma =5$ and $N_x=50$.

Ternary diagram, Fig. 8(d), shows the steady state of output signal $\mathbf {y}(\tilde {t})$ displayed using a RGB pseudo-colouring scheme. It is important to note that no post-processing was applied to generate Fig. 8(d). By visual inspection, one can observe three distinct regions in Fig. 8(d) clearly indicating the presence of three underlying chemical specimens. Furthermore, one can see the change of shades of the colour-encoding as the concentration of the underlying chemicals are varied. In contrast with Fig. 8(a), there is almost no colour contrast for the case of small methanol concentration, indicating that discriminating ethanol and isopropanol mixture is a challenging task. This is perhaps understandable since, as, as seen in Fig. 2(b), the absorption spectra $\alpha$ of pure ethanol and isopropanol are similar.

4. Conclusions

This paper proposed a novel fully bio-inspired apparatus to monitor and determine the concentration of chemical mixture. It implements a human eyes-like sensor and neuromorphic optical system based on electro-optical implementation of reservoir computer (EORC) with the single physical node delay line and Mach-Zehnder nonlinear element for information processing. Here, a delay differential equation including a realistic bandpass filtering effect was used to demonstrate this sensing application. The sensing apparatus proposed in this paper – novel combination of biomimetic sensor and the EORC – offers significant advantages compares to conventional spectroscopic method, e.g., Fourier transform infrared or Raman spectroscopy, in terms of the speed of response and the ability to be trained. Two kinds of sensing application were simulated successfully, namely, to discriminate known (trained) samples of pure chemicals and to determine the concentration of mixture chemical. The latter application highlights the true capabilities of EORC as member of artificial neural network family to perform interpolation and regression task which cannot be achieved in a standard spectroscopic analysis. Experimental demonstration is currently in progress and soon will be reported. Time-stepping algorithms based on $\mathcal {Z}$-bilinear transform technique have been developed to numerically solved the governing equations. Further optimisation that was not explored here is the optimisation of the biomimetic sensor, for example the number of channels and the wavelength of the filter, we believe such an optimisation will further increase the accuracy of the whole system. Finally, we believe the sensing apparatus proposed in this work will enable the advancement of human-inspired artificial visual sensory systems, for example in robotics, which will be high-throughput and require fast hardware-level integration of machine learning.

A. Appendices

A.1 Digital filter model of Eqs. (6) and (7) by $\mathcal {Z}$-bilinear transform technique

Low-pass digital filter: Performing Laplace transform on the DDE for Eq. (6),

(12)$$\mathcal{L}[h(t)] = \frac{v(s)}{r(s)} = H(s)= \frac{1}{1+\delta s}.$$

Performing $\mathcal {Z}$-bilinear transformation $s\leftarrow \frac {2}{\Delta T} \left (\frac {1-z^{-1}}{1+z^{-1}}\right )$ on Eq. (12), it can be shown

(13)$$v(z)=K_1 r(z)+z^{{-}1} \{K_1 K_2 v(z)+K_1 r(z)\},$$

where

(14)$$K_1=\left(\frac{2\delta}{\Delta T}+1\right)^{{-}1}; \qquad K_2=\frac{2\delta}{\Delta T}-1,$$

and $\Delta T$ denote the numerical accuracy of the model, throughout this paper $\Delta T=\tau /N_x$, with $N_x$ denote the number of virtual nodes.

Band-pass digital filter: Performing Laplace transform on the DDE in Eq. (7),

(15)$$\mathcal{L}[h(t)] = \frac{v(s)}{r(s)} = H(s)= \frac{\frac{2m}{\omega_0} s}{1+\frac{2m}{\omega_0} s + \frac{1}{\omega_0^{2}}s^{2}}.$$

Performing $\mathcal {Z}$-bilinear transformation $s\leftarrow \frac {2}{\Delta T} \left (\frac {1-z^{-1}}{1+z^{-1}}\right )$ on Eq. (15), it can be shown

(16)$$v(z)=K_1 K_2 r(z) + z^{{-}1} \left\{ K_2 K_3 v(z) + z^{{-}1} \{ K_2K_4 v(z) - K_1K_2 r(z) \}\right\},$$

where

(17)$$\begin{aligned}&K_1 = \frac{4m}{\omega_0 \Delta T} ;{\kern 7pt} K_2 = \left(1+ \frac{4m}{\omega_0 \Delta T} + \frac{4}{\omega_0^{2} \Delta T^{2}} \right)^{{-}1}; \\& K_3={-}2 + \frac{8}{\omega_0^{2} \Delta T^{2}} ; {\kern 7pt}K_4 ={-}1+ \frac{4m}{\omega_0 \Delta T} - \frac{4}{\omega_0^{2} \Delta T^{2}},\end{aligned}$$

and $\Delta T$ denote the numerical accuracy of the model, throughout this paper $\Delta T=\tau /N_x$, with $N_x$ denote the number of virtual nodes.

Depending on the value of parameter $m$, Eq. (15) can model either a narrow or broad band-pass filter [45]. With $m<1$, Eq. (15) has complex eigenvalues exhibiting a high-Q factor oscillator, as such operates as a narrow band-pass filter centred at $\omega _0$. For $m>1$, the eigenvalues of Eq. (15) are real and negative valued, such system is a strongly damped oscillator. The transfer function (Eq. (15)) with $m>1$ exhibits a broad band-pass filter effect with low-frequency cut-off (-3 dB), $f_L=\omega _0/(4\pi m)$ and the high-frequency cut-off (-3 dB) $f_H=(m \omega _0)/\pi$, from which filter parameters $m=\frac {1}{2} \sqrt {f_H/f_L}$ and $\omega _0=2\pi \sqrt {f_L f_H}$. The digital filter circuit of Eq. (13) for the low-pass and Eq. (16) for the bandpass filter is depicted in Fig. 9(a and b), respectively.

 

Fig. 9. Digital filter design for the low-pass filter (a) and the bandpass (b) filter used to model DDEs in Eqs. (6) and (7), respectively.

Download Full Size | PDF

A.2 Stationary and bifurcation point of Eq. (7) at the quad operational condition

Taking an adiabatic limit at the steady state reduces Eq. (7) to

(18)$$x(t+\tau)= v(t+\tau)+I_0 \cos^{2} \Phi_0 = I_0 \cos^{2} \left(\gamma v(t)+ \Phi_0 \right).$$

For specific ‘Quad’ operational point at positive slope of Eq. (5), one can calculate the ‘fixed point’ $\bar {v}(I_0,\Phi _0=-\pi /4)$ by solving,

(19)$$\bar{v}=I_0 \cos^{2}\left( \gamma \bar{v} -\frac{\pi}{4}\right) - \frac{I_0}{2}.$$

In the presence of band-pass filter, the steady state of signal $v(t)$ will has a mean value of zeros. Thus, one can expand the non-linear term in its series form, from which the stationary point of Eq. (18) can be approximated by,

(20)$$\bar{v}=I_0\left(\frac{1}{2} +\gamma \bar{v} -\frac{2}{3} \gamma^{3} \bar{v}^{3}\right) -\frac{I_0}{2} \rightarrow \bar{v} =\left\{0;\pm \frac{3}{2\gamma} \sqrt{\gamma - \frac{1}{I_0}} \right\}.$$

Substituting Eq. (20) back to Eq. (18),

(21)$$\bar{x} =\left\{\frac{I_0}{2};\frac{I_0}{2} \pm \frac{3}{2\gamma} \sqrt{\gamma - \frac{1}{I_0}} \right\}.$$

Specifically, the bifurcation point occurs at,

(22)$$\gamma = \frac{1}{I_0}.$$

The analytical solution for the general case of Eq. (7) can be found in [54].

Funding

Engineering and Physical Sciences Research Council (EP/P013708/1, EP/T010762/1, EP/V048937/1).

Acknowledgment

Thanks to TM Benson for proof-reading the manuscript.

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. W. B. Levy and V. G. Calvert, “Computation in the human cerebral cortex uses less than 0.2 watts yet this great expense is optimal when considering communication costs,” BioRxiv (2020).

2. C. Wilson, “Your brain is like 100 billion mini-computers all working together,” New Sci. 18 (2018).

3. L. Beaulieu-Laroche, E. H. Toloza, M.-S. van der Goes, M. Lafourcade, D. Barnagian, Z. M. Williams, E. N. Eskandar, M. P. Frosch, S. S. Cash, and M. T. Harnett, “Enhanced dendritic compartmentalization in human cortical neurons,” Cell 175(3), 643–651.e14 (2018). [CrossRef]  

4. P. R. Prucnal, B. J. Shastri, and M. C. Teich, Neuromorphic Photonics (CRC Press, 2017).

5. B. J. Shastri, A. N. Tait, T. F. de Lima, M. A. Nahmias, H.-T. Peng, and P. R. Prucnal, “Principles of neuromorphic photonics,” arXiv preprint arXiv:1801.00016 (2017).

6. T. F. De Lima, B. J. Shastri, A. N. Tait, M. A. Nahmias, and P. R. Prucnal, “Progress in neuromorphic photonics,” Nanophotonics 6(3), 577–599 (2017). [CrossRef]  

7. G. Tanaka, T. Yamane, J. B. Héroux, R. Nakane, N. Kanazawa, S. Takeda, H. Numata, D. Nakano, and A. Hirose, “Recent advances in physical reservoir computing: a review,” Neural Networks 115, 100–123 (2019). [CrossRef]  

8. D. V. Christensen, R. Dittmann, B. Linares-Barranco, A. Sebastian, M. Le Gallo, A. Redaelli, S. Slesazeck, T. Mikolajick, S. Spiga, S. Menzel, I. Valov, G. Milano, C. Ricciardi, S.-J. Liang, F. Miao, M. Lanza, T. J. Quill, S. T. Keene, A. Salleo, J. Grollier, D. Markovic, A. Mizrahi, P. Yao, J. J. Yang, G. Indiveri, J. P. Strachan, S. Datta, E. Vianello, A. Valentian, J. Feldmann, X. Li, W. H. Pernice, H. Bhaskaran, S. Furber, E. Neftci, F. Scherr, W. Maass, S. Ramaswamy, J. Tapson, P. Panda, Y. Kim, G. Tanaka, S. Thorpe, C. Bartolozzi, T. A. Cleland, C. Posch, S.-C. Liu, G. Panuccio, M. Mahmud, A. N. Mazumder, M. Hosseini, T. Mohsenin, E. Donati, S. Tolu, R. Galeazzi, M. E. Christensen, S. Holm, D. Ielmini, and N. Pryds, “2022 roadmap on neuromorphic computing and engineering,”, Neuromorphic Computing and Engineering, accepted manuscript (2022).

9. K. Vandoorne, J. Dambre, D. Verstraeten, B. Schrauwen, and P. Bienstman, “Parallel reservoir computing using optical amplifiers,” IEEE Trans. Neural Netw. 22(9), 1469–1481 (2011). [CrossRef]  

10. K. Vandoorne, W. Dierckx, B. Schrauwen, D. Verstraeten, R. Baets, P. Bienstman, and J. V. Campenhout, “Toward optical signal processing using photonic reservoir computing,” Opt. Express 16(15), 11182 (2008). [CrossRef]  

11. F. Laporte, A. Katumba, J. Dambre, and P. Bienstman, “Numerical demonstration of neuromorphic computing with photonic crystal cavities,” Opt. Express 26(7), 7955 (2018). [CrossRef]  

12. S. Phang, P. Bienstman, A. B. Seddon, C. J. Mellor, and T. M. Benson, “A chaotic microresonator structure for an optical implementation of an artificial neural network,” in Integrated Optics: Devices, Materials, and Technologies XXIV, S. M. García-Blanco and P. Cheben, eds. (SPIE, 2020).

13. S. Phang, P. D. Sewell, A. Vukovic, and T. M. Benson, “The optical reservoir computer: a new approach to a programmable integrated optics system based on an artificial neural network,” in Integrated Optics Volume 2: Characterization, Devices and Applications (Institution of Engineering and Technology, 2020), pp. 361–380.

14. L. Appeltant, M. C. Soriano, G. Van der Sande, J. Danckaert, S. Massar, J. Dambre, B. Schrauwen, C. R. Mirasso, and I. Fischer, “Information processing using a single dynamical node as complex system,” Nat. Commun. 2(1), 468 (2011). [CrossRef]  

15. Y. Paquot, F. Duport, A. Smerieri, J. Dambre, B. Schrauwen, M. Haelterman, and S. Massar, “Optoelectronic reservoir computing,” Sci. Rep. 2(1), 287 (2012). [CrossRef]  

16. T. Wu, A. J. Dufford, M.-A. Mackie, L. J. Egan, and J. Fan, “The capacity of cognitive control estimated from a perceptual decision making task,” Sci. Rep. 6(1), 34025 (2016). [CrossRef]  

17. S. Phang, D. Furniss, C. Mellor, G. Roelkens, A. B. Seddon, P. Bientsman, and T. M. Benson, “Neuromorphic sensing via temporal signal signature processed by photonic reservoir computer,” in Optical Biopsy XIX: Toward Real-Time Spectroscopic Imaging and Diagnosis, R. R. Alfano, S. G. Demos, and A. B. Seddon, eds. (SPIE, 2021).

18. K. J. Major, K. J. Ewing, M. K. Poutous, J. S. Sanghera, and I. D. Aggarwal, “Filter-based chemical sensors for hazardous materials,” in Detection and Sensing of Mines, Explosive Objects, and Obscured Targets XIX, S. S. Bishop and J. C. Isaacs, eds. (SPIE, 2014).

19. K. J. Major, M. K. Poutous, K. F. Dunnill, P. C. Deguzman, J. S. Sanghera, I. D. Aggarwal, and K. J. Ewing, “Biomimetic optical-filter detection system for discrimination of infrared chemical signatures,” Anal. Chem. 88(23), 11491–11497 (2016). [CrossRef]  

20. K. J. Major, M. K. Poutous, I. D. Aggarwal, J. S. Sanghera, and K. J. Ewing, “Analytical procedure to assess the performance characteristics of a non-spectroscopic infrared optical sensor for discrimination of chemical vapors,” Appl. Opt. 57(30), 8903 (2018). [CrossRef]  

21. K. J. Ewing, K. Major, and J. S. Sanghera, “A biomimetic optical approach to skin cancer detection,” in Optical Fibers and Sensors for Medical Diagnostics and Treatment Applications XX, I. Gannot, ed. (SPIE, 2020).

22. M. Flinkman, H. Laamanen, P. Vahimaa, and M. Hauta-Kasari, “Number of colors generated by smooth nonfluorescent reflectance spectra,” J. Opt. Soc. Am. A 29(12), 2566 (2012). [CrossRef]  

23. J. M. M. Linhares, P. D. Pinto, and S. M. C. Nascimento, “The number of discernible colors in natural scenes,” J. Opt. Soc. Am. A 25(12), 2918 (2008). [CrossRef]  

24. J. L. Schnapf, T. W. Kraft, and D. A. Baylor, “Spectral sensitivity of human cone photoreceptors,” Nature 325(6103), 439–441 (1987). [CrossRef]  

25. H. H. Thoen, M. J. How, T.-H. Chiou, and J. Marshall, “A different form of color vision in mantis shrimp,” Science 343(6169), 411–413 (2014). [CrossRef]  

26. R. L. D. Valois and G. H. Jacobs, “Primate color vision,” Science 162(3853), 533–540 (1968). [CrossRef]  

27. M. K. Poutous, I. D. Aggarwal, K. J. Major, J. S. Sanghera, and K. Ewing, “Comparative discrimination spectral detection system and method for the identification of chemicals with overlapping spectral signatures,” US Patent 9, 857, 295 (2018).

28. Infrarec, “FPI Detectors,” https://www.infratec.co.uk/downloads/en/sensor-division/flyer/infratec-fpi-b-en-mail.pdf. Accessed: 2021-11-25.

29. F. Borondics, M. Jossent, C. Sandt, L. Lavoute, D. Gaponov, A. Hideur, P. Dumas, and S. Février, “Supercontinuum-based fourier transform infrared spectromicroscopy,” Optica 5(4), 378 (2018). [CrossRef]  

30. T. L. Myers, R. G. Tonkyn, T. O. Danby, M. S. Taubman, B. E. Bernacki, J. C. Birnbaum, S. W. Sharpe, and T. J. Johnson, “Accurate measurement of the optical constants n and k for a series of 57 inorganic and organic liquids for optical modeling and detection,” Appl. Spectrosc. 72(4), 535–550 (2018). [CrossRef]  

31. H. A. Lorentz, The Theory of Electrons and its Applications to the Phenomena of Light and Radiant Heat, vol. 29 (GE Stechert & Company, 1916).

32. N. An, B. Zhuang, M. Li, Y. Lu, and Z.-G. Wang, “Combined theoretical and experimental study of refractive indices of water–acetonitrile–salt systems,” J. Phys. Chem. B 119(33), 10701–10709 (2015). [CrossRef]  

33. T. G. Mayerhöfer, S. Pahlow, and J. Popp, “The Bouguer-Beer Lambert law: shining light on the obscure,” ChemPhysChem 21(18), 2029–2046 (2020). [CrossRef]  

34. L. Sojka, Z. Tang, D. Jayasuriya, M. Shen, D. Furniss, E. Barney, T. M. Benson, A. B. Seddon, and S. Sujecki, “Ultra-broadband mid-infrared emission from a pr3+/dy3+ co-doped selenide-chalcogenide glass fiber spectrally shaped by varying the pumping arrangement [invited],” Opt. Mater. Express 9(5), 2291 (2019). [CrossRef]  

35. C. R. Petersen, U. Møller, I. Kubat, B. Zhou, S. Dupont, J. Ramsay, T. Benson, S. Sujecki, N. Abdel-Moneim, Z. Tang, D. Furniss, A. Seddon, and O. Bang, “Mid-infrared supercontinuum covering the 1.4–13.3 μm molecular fingerprint region using ultra-high na chalcogenide step-index fibre,” Nat. Photonics 8(11), 830–834 (2014). [CrossRef]  

36. T. Sylvestre, E. Genier, A. N. Ghosh, P. Bowen, G. Genty, J. Troles, A. Mussot, A. C. Peacock, M. Klimczak, A. M. Heidt, J. C. Travers, O. Bang, and J. M. Dudley, “Recent advances in supercontinuum generation in specialty optical fibers [invited],” J. Opt. Soc. Am. B 38(12), F90 (2021). [CrossRef]  

37. N. Ganesh, “Rebooting neuromorphic design-a complexity engineering approach,” in 2020 International Conference on Rebooting Computing (ICRC) (IEEE, 2020), pp. 80–89.

38. G. V. der Sande, D. Brunner, and M. C. Soriano, “Advances in photonic reservoir computing,” Nanophotonics 6(3), 561–576 (2017). [CrossRef]  

39. M. Lukoševičius and H. Jaeger, “Reservoir computing approaches to recurrent neural network training,” Comput. Sci. Rev. 3(3), 127–149 (2009). [CrossRef]  

40. A. Argyris, J. Bueno, and I. Fischer, “Photonic machine learning implementation for signal recovery in optical communications,” Sci. Rep. 8(1), 8487 (2018). [CrossRef]  

41. K. Ikeda, “Multiple-valued stationary state and its instability of the transmitted light by a ring cavity system,” Opt. Commun. 30(2), 257–261 (1979). [CrossRef]  

42. K. Ikeda, H. Daido, and O. Akimoto, “Optical turbulence: Chaotic behavior of transmitted light from a ring cavity,” Phys. Rev. Lett. 45(9), 709–712 (1980). [CrossRef]  

43. J. C. Sprott and J. C. Sprott, Chaos and Time-series Analysis (Oxford University Press, 2001).

44. B. E. Saleh and M. C. Teich, Fundamentals of Photonics (John Wiley & Sons, 2019).

45. L. Larger, “Complexity in electro-optic delay dynamics: modelling, design and applications,” Philos. Trans. R. Soc., A 371(1999), 20120464 (2013). [CrossRef]  

46. G&H, “High speed detectors,” https://www.gandh.com/wp-content/uploads/2021/10/GH-DS-FO-High-Speed-Detector-EM169-7010-Rev5-1.pdf. Accessed: 2021-11-25.

47. iXblue, “DR-AN-20-MO 20 GHz Analog Driver,” https://photonics.ixblue.com/sites/default/files/2021-02/DR-AN-20-MO.pdf. Accessed: 2021-11-25.

48. S. Phang, T. M. Benson, H. Susanto, S. C. Creagh, G. Gradoni, P. D. Sewell, and A. Vukovic, “Theory and numerical modelling of parity-time symmetric structures in photonics: Introduction and grating structures in one dimension,” in Recent Trends in Computational Photonics, (SpringerNature, 2017), pp. 161–205.

49. S. Phang, A. Vukovic, G. Gradoni, P. D. Sewell, T. M. Benson, and S. C. Creagh, “Theory and numerical modelling of parity-time symmetric structures in photonics: Boundary integral equation for coupled microresonator structures,” in Recent Trends in Computational Photonics, (SpringerNature, 2017), pp. 207–233.

50. C. W. Cryer, “Numerical methods for functional differential equations,” in Delay and Functional Differential Equations and their Applications, (Elsevier, 1972), pp. 17–101.

51. A. Bellen and M. Zennaro, Numerical Methods for Delay Differential Equations (Oxford university press, 2013).

52. R. Reed and R. J. MarksII, Neural Smithing: Supervised learning in Feedforward Artificial Neural Networks (MIT Press, 1999).

53. V. Janyani, A. Vukovic, J. D. Paul, P. Sewell, and T. M. Benson, “Time domain simulation in photonics: A comparison of nonlinear dispersive polarisation models,” Opt. Quantum Electron. 37(1-3), 3–24 (2005). [CrossRef]  

54. L. Weicker, T. Erneux, O. D’Huys, J. Danckaert, M. Jacquot, Y. Chembo, and L. Larger, “Slow–fast dynamics of a time-delayed electro-optic oscillator,” Philos. Trans. R. Soc., A 371(1999), 20120459 (2013). [CrossRef]  

55. L. Weicker, T. Erneux, O. d’Huys, J. Danckaert, M. Jacquot, Y. Chembo, and L. Larger, “Strongly asymmetric square waves in a time-delayed system,” Phys. Rev. E 86(5), 055201 (2012). [CrossRef]