Chaos is an irregular, stochastic-like linearity that occurs in a deterministic system; in addition, the time series with chaotic characteristics generated by a chaotic model are chaotic time series, which contain rich dynamical information [1]. Since the 1990s, chaotic science has been interpenetrating with other sciences, and has been widely applied in many fields such as mathematics, physics, life science, earth science, economics, and astronomy [2–7]. Therefore, chaotic time series are a bridge from chaos theory to reality, and the prediction of chaotic time series has very broad application prospects in short-term load prediction of electric power systems, stock market quotation prediction, sunspot prediction, short-term traffic flow prediction, weather prediction, early warning of geomagnetic storms [8,9], etc.

In recent years, the generation of chaotic signals has attracted great attention, especially optical chaos generation has been widely studied in several fields, such as physical random bit generation (RBG), secure optical communication, and optical logic [10–17]. In optical chaotic communication systems, chaotic signals are used as optical carriers for hiding information to enhance the security of the physical layer. For example, in 2020, Cai et al. investigated the effect of a fiber Bragg grating on the chaotic output of the vertical-cavity surface-emitting laser (VCSEL) coupled networks through experimental methods and generated high-speed physical random numbers using XOR and other post-processing methods [14]. In 2020, Jiang et al. proposed and demonstrated a chaos generation scheme based on external feedback semiconductor lasers, which used a parallel coupled ring resonator as a reflector, supporting both bandwidth enhancement and suppression of excellent time delay characteristics, which was closed to zero [15]. In 2022, they further proposed a new scheme for generating optical chaos, using the output of an external-cavity semiconductor laser as the driving signal of a phase modulator to modulate the output of a continuous wave laser, and conducted experimental verification. This scheme supported the simultaneous generation of two low-correlation chaotic signals with a broadband spectrum and time delay suppression characteristics [16]. In 2023, Li et al. proposed and numerically demonstrated the generation of broadband and high-dimensional chaotic signals based on optically pumped spin polarized VCSEL. The simulation results demonstrated for the first time that spin VCSEL with simple auxiliary configurations allows for the generation of chaos with desired characteristics, including effective bandwidth up to 30 GHz and above, time delay suppression not exceeding 0.2, flatness of 0.75 and above, and high complexity over a wide parameter range [17].

On the other hand, the conventional chaotic event prediction models are usually classified into two main categories: one is the dynamics method based on nonlinear mathematical models, which simplifies and numerically solves the already known physical models for the purpose of prediction. In addition, the other is the phase space reconstruction method based on the actual observational data, which assumed that the complete observational data of the modeled physical system would be obtained, the event sequence is established accurately through the fitting of the data, fixed and invariant model, and the predictive model no longer modifies the model parameters during the prediction process. In recent years, a variety of chaotic prediction models have been proposed, such as moving average model [18], autoregressive model [19], autoregressive moving average model [20], feed-forward neural network [21], recurrent neural network, and support vector machine [22]. As mature models widely used in linear time series forecasted, the moving average model, autoregressive model, and autoregressive moving average model have the advantages of high prediction accuracy and speedy prediction but have low prediction accuracy for nonlinear time series and limitations in chaotic time series forecasting. Neural networks with strong nonlinear computing ability and high prediction accuracy are powerful tools for nonlinear time series prediction, e.g., Shen et al. [23] improved the extreme learning machine and forecasted chaotic time series, which verified the strong robustness of the proposed model. However, due to the limitation of neural networks to the number of available training samples, there are problems of overfitting and falling into local optimal solutions, whereas the computation time increases significantly with the increase of network complexity, thus capping the applicability of neural networks. To illustrate, Guo [24] used a genetic algorithm to optimize a support vector machine for time series forecasting and achieved a higher forecasting accuracy than before. At the same time, Yan et al. [25] used the support vector machine model optimized by a genetic algorithm for short-term wind speed prediction and obtained the predicted value of short-term wind speed, which is of great practical significance. Nevertheless, compared with other models, it has the disadvantages of slow computation speed, too much dependence on the selection of parameters, and no unified method for determining the model parameters.

The above-mentioned reported deep neural networks as well as support vector machines for chaotic neural networks achieve the prediction of chaotic sequences through software programs, which consume a large amount of computer hardware in the prediction of more complex event sequences. Therefore, it is necessary to further investigate new models with simple structure, good generalization performance, high accuracy, and efficiency. Reservoir computing (RC) is a special kind of recurrent neural network. The input weights and internal weights of the reservoir are randomly generated and fixed. Only the output weights need to be obtained through the training algorithm, reducing a significant amount of training and computational costs [26]. Time delayed RC uses a nonlinear node to realize the exchange of space and events through the method of event multiplexing, effectively simplifying the hardware structure [27–29]. For prediction of time series data, time delayed RC based on nonlinear semiconductor lasers has the following advantages: fast, efficient, and parallel computing capability. Numerous research results have shown that this method is effective in predicting chaotic tracks in circuits or nonlinear optical systems [30–36], and it can predict the synchronization of chaotic tracks reasonably well. However, the works of using RC systems to predict chaotic time series mainly focus on the numerical simulation stage, and the chaotic time series used are generally standard datasets. Therefore, predicting the time series with different complexities generated from actual optical chaotic systems in the RC experimental system is worth further research.

In this paper, a multi-delay mutual coupling VCSEL ring network structure is used to obtain laser chaotic time series of different complexities by adjusting the parameters such as injection current, coupling strength, and detuning frequency. Furthermore, the prediction of short-term laser chaotic time series is realized by the photonic RC system based on VCSEL with feedback loops in experiment and simulation. The main contents of the remaining sections of this paper are summarized as follows. In Section 2, the experiment and simulation schemes for chaotic time series generation and prediction are presented. The results for chaotic time series generation and prediction are revealed in Section 3. The conclusion is presented in Section 4.

The RC, as a brain-inspired computational paradigm, includes input, reservoir, and output layers. As Fig. 1(b) shows, the input weight and reservoir connection weight of reservoir are generated randomly, which reduces the computing complexity. Here, in this work, we explore methods capable of achieving short-term prediction for chaotic time series. The entire system can be divided into two parts: the multi-delay mutual coupling VCSEL ring network used for chaotic time series generation shown in Fig. 1(a) and the photonic RC system revealed in Fig. 1(c).

The chaotic time series are generated by a chaotic system based on VCSEL with a mutual coupling network as shown in Fig. 1(a). The optical outputs of VCSEL1 and VCSEL3 are injected into VCSEL2, while the chaotic optical output from VCSEL2 is injected into VCSEL1 and VCSEL3, forming the simplest structure of star network topology. Subsequently, the optical outputs of VCSEL1 and VCSEL3 are injected into the feedback loop through a 50:50 fiber coupler (FC) from the 2-port of an optical circulator (OC1) to the 3-port. After being combined with the erbium-doped optical fiber amplifier (EDFA) and variable optical attenuator (VOA) for power adjustment, it is injected into the 1-port of OC2 to the 3-port, while the output of VCSEL2 is injected into the feedback loop from the 2-port of OC2. In addition, the feedback information is added into the system from the 1-port of OC1. Here, polarization controllers (PC1, PC2, and PC3) are employed to adjust the polarization state of light. Higher-stability and low-noise diode controllers (LACs) could be used to provide stable temperature and pump current for VCSELs.

Besides, in this work, the spin-flip model is also employed to model the mutual coupling network based on VCSEL for simulation, which can be described as the following rate equations [37]: $$\begin{gathered} dEmxdt=k(1+iα)(NmEmx−Emx+inmEmy)−(γa+iγp)Emx \\ +∑n=1,n≠m3krEnx(t−τnm)e−i(ωnτnm−Δωnmt), \end{gathered}$$(1)$$\begin{gathered} dEmydt=k(1+iα)(NmEmy−Emy+inmEmx)−(γa+iγp)Emy \\ +∑n=1,n≠m3krEny(t−τnm)e−i(ωnτnm−Δωnmt), \end{gathered}$$(2)dNmdt=γN[μ−Nm(1+|Emx|2+|Emy|2)+inm(EmxEmy*−EmyEmx*)],(3)dnmdt=−γsnm−γn[n(|Emx|2+|Emy|2)+iNm(EmyEmx*−EmxEmy*)].(4)

The n and m (m,n=1,2,3) in Eqs. (1)–(4) represent VCSEL1, VCSEL2, and VCSEL3, respectively. Emx and Emy represent slowly varying complex electric field amplitudes of two polarization modes called X-PC and Y-PC, respectively. Nm stands for the total carrier reversal between the conduction and valence bands, while nm means the carrier reversal with opposite spins. τnm is the time delay between VCSELn and VCSELm. The descriptions and typical values of other parameters are shown in Table 1 [38].Table 1.

Parameter Values for Three VCSELs in the Numerical Simulation

In addition, the Lyapunov exponent [39,40], Kolmogorov-Sinai entropy [41], sample entropy, and permutation entropy (PE) [42–46] are several popular methods used to quantify the chaotic complexity of laser systems. Compared to other methods, PE is widely employed to measure the complexity of chaotic time series in the laser chaotic systems due to its flexible and fast computation [44–46]. PE, proposed by Bandt and Pompe in 2002 based on information theory, is applied to quantify the complexity of the time series [38].

Generally, a set of measured values that forms the time series X={X(i),i=1,2,…,n} could be obtained by observing the nonlinear systems. A D-dimensional vector Xj=[X(j),X(j+τ),…,X(j+(D−1)τ)] and j=1,2,…,T−(D−1)τ, where the parameters τ and D are the embedding delay and the embedding dimension, respectively, could be constructed by adopting phase space reconstruction. Subsequently, each vector Xj is arranged in an ascending order, and the subscript order is denoted by a combination π=(r1,r2,…,rD). D! possible permutations should be considered. The probability distribution of the combination is recorded as p=p(π) revealed in Eq. (5); thus the PE could be defined as Eq. (6). H(p)=0 stands for all vectors Xj having the same array, making the series completely predictable, while H(p)=1 means each permutation appears with different probabilities: p(π)=#{j|j≤T−(D−1)τ;Xjhas typeπ}T−(D−1)τ,(5)H(p)=−∑p(π)logp(π)log(D!).(6)

However, for multi-delay systems, a single-point PE value is not sufficient to reflect the complexity of the system. Therefore, Guo et al. proposed to use the average of multi-point PE values to measure the complexity of multi-delay systems, averaging the PE of embedded delays from one to the maximum delay [46]. The MPEs are defined in Eq. (7), where τmax is the max embedded delay [46]: HM(p)=∑τ=1τmaxH(pτ)τmax.(7)

The photonic RC based on VCSEL with feedback loops can be employed to predict short-term laser chaotic time series, as illustrated in Fig. 1(c). The experimental setup can be summarized in the following steps. At first, a commercial 1550 nm VCSEL (Seoul Viosys) without isolation is regarded as the nonlinear node, which is employed to achieve nonlinear mapping from low to high dimensions and driven by a low-noise laser diode controller (ILX-Lightwave, LAC-3724C). Additionally, the external optical injection generated by the tunable laser (TL) with isolation is sent into a Mach–Zehnder modulator (MZM). The signal is modulated with the injection signal S(t), generated by an arbitrary waveform generator (AWG, 70002A) and amplified by an RF amplifier. Here, the S(t) could be obtained by multiplying the input signal u(t) by the mask signal m(t). The modulated optical signal is divided into two parts using a 20:80 FC2. The 20% port is used to send signal into a photodetector (PD, HP 11982A) for detection and the oscilloscope (OSC, Agilent DSOV334A) for collection. Signal from the 80% port passes through the EDFA and VOA, which could adjust the power of the injection signal. The modulated injection signal and feedback signal from PC6 are combined by 50:50 FC3 and then injected into VCSEL through OC3. Furthermore, the reservoir consists of VCSEL as a nonlinear node and a feedback loop with a delay line (DL), PC6, FC4, and OC3 for achieving time division multiplexing (TDM). Lastly, 10% of the output of the VCSEL is fed back to the VCSEL by OC3 after passing through the delay line (DL) and PC6, which is employed to adjust the polarization of the feedback signal. The remaining 90% of the output of VCSEL is sent to PD2 for photoelectric conversion and the oscilloscope (OSC) for collecting information from the reservoir.

In this work, for simulation, the reservoir is implemented by employing the VCSEL with optical feedback as a nonlinear node, which has two modes called X-PC and Y-PC. Furthermore, the output of X-PC is fed back into X-PC while the output of Y-PC is fed back into Y-PC, indicating parallel-polarized optical feedback (PPOF) [38,39].

Additionally, the photonic RC could be simulated by using a spin-flip model with feedback to model the VCSEL with a feedback loop; the rate equations are described in Eqs. (8)–(12) [36]: $$\begin{gathered} dErxdt=kr(1+iαr)(NrErx−Erx+inrEry)−(γra+iγrp)Erx \\ +kdErx(t−τ)e−i(ωxτ)+kinjε(t), \end{gathered}$$(8)$$\begin{gathered} dErydt=kr(1+iαr)(NrEry−Ery+inrErx)−(γra+iγrp)Ery \\ +kdEry(t−τ)e−i(ωyτ), \end{gathered}$$(9)dNrdt=γrα[μr−Nr(1+|Erx|2+|Ery|2)+inr(ErxEry*−EryErx*)],(10)dnrdt=−γrsnr−γrn[nr(|Erx|2+|Ery|2)+iNr(EryErx*−ErxEry*)],(11)ε(t)=|ε0|2{1+ei[u(t)×mask(t)]}ei2πΔft,(12)where Erx and Ery represent slowly varying complex electric field amplitudes of X-PC and Y-PC modes in VCSEL used in the photonic RC system, respectively. Additionally, Nr stands for the total carrier reversal between the conduction and valence bands, and nr describes the carrier reversal with opposite spins. μr represents the normalized bias current of the VCSELr. In Eqs. (8) and (9), the feedback terms can be found in the third term. kd represents feedback strength, and τ denotes feedback time delay. The injected term is described in the last term in Eq. (8). kinj stands for the injected strength. ε(t) is described as Eq. (8) [39]. In addition, ε(t) stands for the output of MZM, which achieves the input weight in the optical domain. Therefore, u(t) and mask(t) represent the input signal and mask signal used in the system, respectively. The system of equations is resolved numerically using a second order Runge-Kutta with 2 ps per step. The simulation parameters are shown in Table 2 [47]. Table 2.

Parameter Values for the Photonic RC System in Our Numerical Simulation^a

It is noteworthy that the VCSEL with self-feedback light is utilized as the reservoir. The number of virtual nodes is N=220. The interval between the adjacent virtual nodes is θ=1ns. The sampled period of input signal (equal to the feedback time delay) is τ=n×θ (the information processing rate is R=1/τ). In the post-processing phase, the virtual node states extracted from the photonic RC system are placed into matrix X, which can be multiplied by the output weight to obtain the prediction or recognition results.

The minimum normalized mean square error (NMSE) can be used to evaluate the performance of the system: NMSE=1L∑j=1L[y¯(j)−y(j)]2σ2,(13)where y¯(j) is the target value, y(j) is the predicted value, L is the total number of experimental data, and σ is the standard deviation of the target value. For the chaotic time series prediction task, when NMSE=1, it represents that the photonic RC system is completely unable to predict the next output of the chaotic sequence; when NMSE=0, it represents that the photonic RC system can accurately predict the next output of the chaotic sequence. Generally speaking, when NMSE≤0.1, it can be considered that the system can complete chaotic time series effectively.

For chaotic time series generation, the setup shown in Fig. 1(b) was employed. On the one hand, for experiment, the optical outputs of VCSEL1 and VCSEL3 injected into VCSEL2 are considered. Here, chaotic time series are collected from the output of VCSEL2. In addition, the chaotic time series with different MPEs are obtained by adjusting the detuning frequency between VCSEL2 and VCSEL1, which could be achieved by operating VCSELs at different temperatures. The VCSELs provide a longer central wavelength under higher operating temperatures. The frequency detunings between VCSEL1 and VCSEL2, and between VCSEL2 and VCSEL3 are Δf21=fvcsel2−fvcsel1 and Δf23=fvcsel2−fvcsel3, respectively. The three coupling delays are set as 2 ns between VCSEL1 and VCSEL2, 2.2 ns between VCSEL2 and VCSEL3, and 2 ns between VCSEL1 and VCSEL3. Different lengths of fiber extensions can be utilized to achieve different delays.

During the experiment, the pump currents of VCSEL1, VCSEL2, and VCSEL3 are maintained at 4.53 mA, 5.97 mA, and 5.05 mA with the temperatures fixed at 15.92°C, 26.16°C, and 24.09°C when there is no frequency detuning among VCSEL1, VCSEL2, and VCSEL3. The detuning frequency is adjusted by varying the operating temperature of VCSEL1. The time-domain signal, power spectrum, PE, and ACF are revealed in Fig. 2. From Figs. 2(a)–2(d), chaotic signals with MPE=0.9678, 0.9789, 0.9972, and 0.9977 are acquired under the VCSEL1 temperature of 17.00°C, 16.72°C, 16.52°C, and 16.30°C, which correspond to center wavelengths of 1557.064 nm, 1557.030 nm, 1557.010 nm, and 1556.986 nm, individually. Thus, the detuning frequencies |Δf21|=|fvcsel2−fvcsel1| are 15 GHz, 10 GHz, 7.5 GHz, and 5 GHz, respectively.

On the other hand, three chaotic time series with similar complexity are generated through simulation, maintaining the frequency detuning Δf21=Δf21=Δf=0GHz, normalized injection current μ=2.7, and coupling strength kr=20ns−1. The time-domain signal, power spectrum, PE, and autocorrelation function (ACF) curves of VCSEL1, VCSEL2, and VCSEL3 output signal as shown in Fig. 3, reflect the obvious time delay characteristics in the system. Note the valley in the PE curve and periodic pulses in the ACF curve. After calculation, the MPEs for the chaotic time series collected from VCSEL1, VCSEL2, and VCSEL3 are 0.9933, 0.9913, and 0.9911, respectively, which meets the requirements of chaotic time series. By adjusting the kr, chaotic time series MPEs are 0.9839, 0.9839, 0.9760, 0.9762, 0.9760, and 0.9753. Here, the permutation entropy is calculated using D=6, T=10,000, and τmax=6.1ns.

On the one hand, in simulation, chaotic time series with different and similar complexities obtained in Fig. 3 by simulation are considered to explore the difficulty of the chaotic time series with different complexities. The NMSE as a function of injection strengths is shown in Fig. 4(a) for different complexity signals and in Fig. 4(b) for similar signals. The NMSE for signals of MPE=0.9760 (0.9839, 0.9839, 0.9913, or 0.9933) decreases when kinj≤10ns−1 (kinj≤25ns−1, kinj≤20ns−1, or kinj≤25ns−1) and fluctuates at 0.0047 (0.0094, 0.02942, or 0.0361) in Fig. 4(a). Subsequently, in Fig. 4(b) the NMSE for signals of MPE=0.9753 (0.9760 or 0.9762) increases when kinj≤15ns−1 and fluctuates at 0.0047 (0.0054 or 0.0058). Comparing the results presented in Figs. 4(a) and 4(b), it could be observed that for the more complex chaotic time series, the more difficult it is to predict, as reflected in the figure that the NMSE of the system’s prediction result is higher for the chaotic time series with larger MPE, the NMSE is insignificantly different.

Additionally, we utilize the chaotic time series with MPE=0.9839 to investigate the application of photonic RC in predicting duration. In Fig. 5(a) [Fig. 5(b)], the NMSE of the system as a function of kinj (kd) is considered for 1-, 2-, 3-, or 4-steps-ahead. When kinj<35ns−1, the NMSE decreases rapidly and then stabilizes at 0.0108, 0.0732, 0.1963, or 0.4168 for 1- (black line), 2- (blue line), 3- (red line), or 4- (green line) steps-ahead prediction, respectively. A similar rhythmic pattern can be observed in the NMSE as a function of the feedback strength kd. When kd<1.3ns−1, the NMSE decreases rapidly and can be maintained at 0.0895, 0.2889, 0.5249, or 0.7513 for different step-ahead prediction tasks, with kinj=10ns−1. From Fig. 5, it is evident that the NMSE remains below 0.2 with the 3-steps prediction task.

Here, for the experiment of prediction results, the basic characteristics of the VCSEL used are shown in Fig. 6. When the temperature of the experiment environment is fixed at 25°C, the bias current is 1.7 mA. It is observed that the two peaks are defined as X-PC (on the left, the wavelength is 1554.41 nm) and Y-PC mode (on the right, the wavelength is 1554.65 nm). Correspondingly, the X-PC mode is the dominant mode and the Y-PC mode of the VCSEL is the suppressed mode. After the optical signal is injected into the Y-PC mode of VCSEL, the Y-PC mode becomes the dominant mode, and the X-PC is suppressed as shown in Fig. 6(b). In addition, the power-current curve is revealed in Fig. 6(c) and the threshold point of the VCSEL used is approximately 1.6 mA. When the bias current adds up to 5 mA, the optical power of the VCSEL output can reach 195.4 μW.

In the experiment, the chaotic time series illustrated in Fig. 2 are utilized to validate the impact of greater complexity in the chaotic time series corresponding to increased prediction challenges. Averaging across multiple experiments reveals a decrease in NMSE from 0.0604 to 0.0496 as the MPE of the signal varies from 0.9977 to 0.9678, consistent with the simulation results in Fig. 4.

Furthermore, the chaotic time series with MPE=0.9977 are used to investigate the influence of varying dataset sizes on system prediction performance. As depicted in Table 3, an increase in the training dataset size from 3000 to 7000 results in improved performance, indicating the utility of partial datasets for the scheme.Table 3.

Impact of Dataset Size on the System

Lastly, the impact of the attenuation coefficient on the system performance is revealed in Fig. 7 for 1- (black line), 2- (blue line), 3- (red line), and 4- (green line) steps-ahead prediction tasks. Here, the complexity of the chaotic time series used is fixed at MPE=0.9977. The NMSE fluctuates at 0.0487, 0.0719, 0.2543, and 0.4868 with the attenuation coefficient increasing, respectively. In addition, the predictive performance of the RC system deteriorates with the increase of prediction step size of the chaotic time series.

To further demonstrate the predicted results, the blue sequences in Fig. 8 represent the predicted values obtained from the RC system while the red sequences represent the target values. Correspondingly, Figs. 8(a)–8(d) reveal the time domain prediction results of 1-step-ahead, 2-steps-ahead, 3-steps-ahead, and 4-steps-ahead tasks, respectively. Here, for Figs. 8(a)–8(c), the values of NMSE are 0.0487, 0.0719, 0.253, and 0.4868, respectively. It means that the photonic RC system can accurately predict the output response of a chaotic system, which includes a multi-delay mutual coupling VCSEL ring network structure within 1 ns for the accuracy requirement of NMSE less than 0.1.

In conclusion, a photonic reservoir computing system based on VCSEL, along with a feedback loop, has been proposed for the short-term prediction of the chaotic time series with different levels of complexity. We generated chaotic time series from the VCSEL’s mutual coupling network system in both simulation and experimentation. Applying the photonic RC system enabled the prediction of chaotic time series with an MPE greater than 0.98. The minimum NMSE reached 0.0047 in simulation and 0.0487 in experimentation. Through comparing chaotic time series having different complexities, it was observed that the more complex the chaotic time series, or the larger the MPE, the greater the difficulty in prediction. Furthermore, the NMSE could remain below 0.2 (0.3) for the 3-steps-ahead prediction task in simulation (experiment). Due to the inherent advantages of photonic RC, issues such as high computational complexity and difficulty in hardware implementation, which exist in another scheme, are overcome.