## 3 Numerical Setup and Results
### Parameters Optimization  
In the parameter optimization part, the signal-to-noise (SNR) ratio was set to 16 dB, the training length was \(3\times 10^{3}\) PAM-4 signals, and the testing length was \(5\times 10^{4}\) PAM-4 signals. Here, we mainly show the effects of the virtual node interval (\(\theta\)), mask size (\(\gamma\)), feedback coefficient (\(\beta\)), and the ridge regression value (_reg_) in the RC training process. As mentioned in Section 2.3, the larger \(\theta\) contributes to the reduced hardware requirements for implementing the RC scheme. However, if the selected \(\theta\) is too large, it may exceed the characteristic time scale of the nonlinear node. In that case, the influence of neighboring nodes on the state of the virtual node will be weakened, and the performance of RC will be drastically reduced. As shown in Figure 5, when \(\theta\) was set between 10 and 80 ps, the BER was relatively low, but when \(\theta\) was set between 50 and 80 ps, the error was more stable, so \(\theta\) = 80 ps was superior here. As shown in Figure 6, when the mask size \(\gamma\) was set to a value of approximately 1, BER was relatively small. As \(\gamma\) increases from 0, the BER decreases first and then increases when \(\gamma\) is larger than 1. This phenomenon is straightforward. When \(\gamma\) is close to 0, the input multiplied by \(\gamma\) is equivalent to almost zero in the system. When \(\gamma\) is too large, the signal will deviate from the nonlinear region of the transfer function. Both situations will deteriorate the system's transmission performance. Here, we chose \(\gamma\) = 0.70. Figure 5: Schematic diagram of the input signal pre-processing process. The trend in the relationship between the normalized feedback coefficient \(\beta\) and BER was similar to that of \(\gamma\). The trend is shown in Figure 7. It was necessary to make the value of the function as close as possible to the nonlinear region of the transfer function, so \(\beta\) also needed a value region that was not too big or too small. When \(\beta<1\), the internal dynamics of the reservoir were stable. When \(\beta>1\), it may become divergent. Therefore, for \(\beta>1\), the performance decreased as \(\beta\) increased. Finally, the optimal value of \(\beta\) was determined to be 0.11 by dividing into smaller intervals as shown in the inset of Figure 7.  
In the training process of the output connection weight matrix (\(W_{out}\)), the optimization of the regularization parameter (_reg_) was also involved. The output layer collects \(x(n)\) at all times, and the node state of each cycle was a column to form matrix \(B\). This matrix \(B\) was multiplied by the output connection weight matrix, \(W_{out}\), to obtain the output. The output  
Figure 6: Relationship between \(\gamma\) and BER. Figure 7: Relationship between \(\beta\) and BER. matrix \(T\) of the training process was the expected output \(d(n)\). \(T\) is a known quantity so that \(W_{out}\) can be solved as the following formula:  
\[W_{out}=T*B^{-1} \tag{6}\]  
Considering that matrix \(B\) may be a singular matrix and enhance the model's generalization ability, the ridge regression algorithm was used to satisfy these two conditions. Equation (6) was modified to the following formula [20]:  
\[W_{out}=TB^{T}(BB^{T}+reg*E)^{-1} \tag{7}\]  
where \(E\) is the identity matrix of the same size as \(B^{*}B^{T}\). The value of \(reg\) directly affects \(W_{out}\), affecting the performance of the RC. As shown in Figure 8, when \(reg\) was very close to 0, the curve fluctuated sharply, and the error was large. This phenomenon occurred because the regularization parameter was too small, causing the model to transform into an over-fitting state. At the same time, because the training length was small, regularization was needed to improve the network's generalization ability. When the \(reg\) value was relatively large, especially ranging from 1E-6, the BER also increased.