Sparse Identification for Nonlinear Optical Communication Systems: SINO Method

TL;DR

The paper presents SINO, a low-complexity machine learning approach that identifies minimal degrees of freedom needed to mitigate nonlinear impairments in optical fiber communication systems, improving data recovery.

Contribution

Introduction of SINO, a sparse identification method that simplifies nonlinear impairment mitigation in optical systems by determining optimal degrees of freedom.

Findings

Effective in standard fiber communication links

Successful application to spatial-division-multiplexing systems

Reduces complexity of nonlinear impairment mitigation

Abstract

We introduce low complexity machine learning based approach for mitigating nonlinear impairments in optical fiber communications systems. The immense intricacy of the problem calls for the development of "smart" methodology, simplifying the analysis without losing the key features that are important for recovery of transmitted data. The proposed sparse identification method for optical systems (SINO) allows to determine the minimal (optimal) number of degrees of freedom required for adaptive mitigation of detrimental nonlinear effects. We demonstrate successful application of the SINO method both for standard fiber communication links and for few-mode spatial-division-multiplexing systems.

Tables

Bandwidth	32 GBaud
RRC Roll off	0.01
Distance	10$\times$100 km
Noise Figure	4.5 dB
Attenuation coefficient	0.2 dB/km
CD	17 ps/nm/km
Nonlinearity	1.4 1/W/km

Equations

$$\mathbf{Y}=\Theta (\mathbf{X})\Xi$$

$$\frac{\partial U_p}{\partial z}=-\frac{\alpha }{2}{U}_p-{\beta }_1^p\frac{\partial U_p}{\partial t}-i\frac{{\beta }_2^p}{2}\frac{{\partial }^2{U}_p}{{\partial }^{2}t}+i{\kappa }_{pp}{\gamma }_{pp}|U_p{|}^2{U}_p+i{\kappa }_{pq}\sum _{q\ne p}^D{\gamma }_{pq}|U_q{|}^2{U}_p+{\eta }_p(t,z),$$

$$Y_{kpi}^{\prime }=\frac{\rho }{L}{\eta }_{kpi}+\epsilon V[Y{]}_{kpi},$$

$$V[Y{]}_{kpi}=\sum _{j=1,2}{\Psi }_s(z)\sum _{m,n=-M}^M({\epsilon }_{pp}{Y}_{k+n,pj}{Y}_{k+m+n,pj}^{\ast }{\tilde{C}}_{mn}^{pp}+\sum _{q\ne p}^D{\epsilon }_{pq}{Y}_{k+n,qj}{Y}_{k+m+n,qj}^{\ast }{\tilde{C}}_{mn}^{pq})Y_{k+m,pi}$$

$$C_{mn}^{pq}=i\int dz{\Psi }_s(z)\int \int \int d\omega d{\omega }_{1}d{\omega }_2{e}^{i{\omega }_1({\beta }_1^p-{\beta }_1^q)z+i({\omega }_1^2+2{\omega }_1\omega )({\beta }_2^p-{\beta }_2^q)z-i{\omega }_1{\omega }_2{\beta }_2^{q}z-i{\omega }_{1}mT-i{\omega }_{2}nT}\times$$

$$f^{\ast }(\omega )f({\omega }_1+\omega )f({\omega }_2+\omega )f^{\ast }({\omega }_1+{\omega }_2+\omega )$$

$$Y_{kpi}^0=X_{kpi}+\rho {\zeta }_{kpi},{\zeta }_{kpi}=\int dz{\eta }_{kpi}(z)$$

$$Y_{kpi}^{1}=\sum_{m,n=-\infty}^{\infty}X_{k+n,pj}\Big{(}\varepsilon_{pp}C_{mn}^{pp}X_{k+m,pi}X_{k+n+m,pi}^{*}+\sum^{D}_{q\neq p}\varepsilon_{pq}C_{mn}^{pq}X_{k+m,qi}X_{k+n+m,qi}^{*}\Big{)}$$

$$\Theta (\mathbf{X})=\left[\begin{array}{ccccccc}...& ...& ......& ...& ......& ...& ...\\ ...& x_k& x_k|x_k{|}^2& x_k|x_{k+1}{|}^2& ...& x_{k+m}{x}_{k+n}{x}_{k+m+n}^{\ast }& ...\\ ...& ...& ......& ...& ......& ...& ...\end{array}\right]$$

$${\mathbf{Y}}^{\mathbf{r}}={\Theta }^r{\Xi }^r-{\Theta }^i{\Xi }^i,{\mathbf{Y}}^{\mathbf{i}}={\Theta }^i{\Xi }^r+{\Theta }^r{\Xi }^i$$

$$\min_{\beta_{0},\beta}\Big{(}\frac{1}{2m}\sum_{i=1}^{m}(y_{i}-\beta_{0}-(\Theta^{r}\Xi^{r}-\Theta^{i}\Xi^{i})_{i}^{T}\beta)^{2}+\lambda P_{\alpha}(\beta)\Big{)}$$

$$P_{\alpha }(\beta )=\frac{1-\alpha }{2}||\beta |{|}_2^2+\alpha ||\beta |{|}_1$$

$$SNR=\frac{S}{N+aS^{3}+bS^{2}N},\ \$$

$$a_{PSE}=2\sum _{q=1}^D\sum _{m,n\ne 0}{\epsilon }_{pq}^2(2-{\delta }_{pq})|C_{m,n}^{pq}{|}^2$$

$$a_{PSE}=\frac{8}{27}(1+2(D-1))\frac{N_s{\gamma }^2{L}_{eff}\ln ({\pi }^2{\beta }_2{L}_{eff}{B}^2)}{\pi {\beta }_2{B}^2},$$

$$a_{SINO}=2\sum _{q=1}^D\sum _{m,n\ne 0}{\epsilon }_{pq}^2(2-{\delta }_{pq})|C_{m,n}^{pq}-\widetilde{C_{m,n}^{pq}}{|}^2$$

$$b=\frac{4}{9}(1+2(D-1))\frac{N_s{\gamma }^2{L}_{eff}\ln ({\pi }^2{\beta }_2{L}_{eff}{B}^2)}{\pi {\beta }_2{B}^2}$$

Full text

††preprint: APS/123-QED

I Introduction

There is an enormous pressure on the fiber-optic communication industry to deal with the exponentially increasing capacity demands from data traffic DJR ; CLEO driven by the existing and constantly emerging internet and broadband services as well as the fast growing machine-to-machine traffic. Current technological solutions suggest a combination of advanced modulation formats and space division multiplexing to achieve substantial enhancement of the spectral efficiency Fontaine ; Rademacher ; Randel ; PW , while fiber remains a nonlinear medium leading to signal being strongly degraded by nonlinear impairments. In future systems, potentially using different types of fibers, nonlinear effects still will be one of the fundamental limiting factors. Therefore, compensation of the nonlinear signal distortions is a critical challenge for development of the next generation of communication systems operating at higher transmission rates.

Up to date, non-linear impairment mitigation has been considered mostly for single mode fiber links with a number of techniques being proposed both in the optical and electronic domain OPC ; regen ; Radic1 ; Ip ; DSP_Seb . Most of the research in the field of electronic compensation has been focused on back-propagation algorithms that emulate optical transmission in the digital domain through an inverse fiber link (with reversed order segments and opposite sign parameters), realized either by means of a Split-Step Fourier method Ip ; DSP_Seb or Volterra Series Transfer Functions V1 ; V2 ; V3 . Despite various simplifications that have been proposed, both approaches are considered to be highly complex because they require multiple computational steps along the link and at least two samples per symbol. This has heavily discouraged any effort for future commercial deployment even for legacy SMF fiber systems.

On the other hand, perturbation analysis of the Manakov equations has led to the development of efficient equalization methods that can mitigate accumulated intra-channel impairments in a single computational step and one sample-per-symbol Tao ; Karlsson . Central to this approach has been the identification of the perturbation coefficients that describe the interaction of each symbol with its preceding and succeeding symbols in the transmission channelGao ; Zhuge ; Peng ; quant . For static connections and specific pulse shapes, such as Gaussian or sinc, the perturbation coefficients can be derived analytically and stored in a look-up table Tao ; ECOC_Essiambre . Since the total number of terms can be excessively high for dispersion un-managed links, where the channel memory is long, exploiting common symmetries and quantizing the coefficients are two of the techniques that have been recently proposed for complexity reduction. Furthermore, to achieve operation in reconfigurable network environments, an adaptive version of the method that uses training sequences and decision-directed least-mean squares algorithms has been introduced in McGill . This enables to identify the perturbation coefficients before establishing any new connection and without prior knowledge of the corresponding transmission link parameters. The latter progress places the perturbation method in a broader context and signifies its practical importance in future optical networks.

With this paper, we expand the application of perturbation-based nonlinear compensation in few mode fiber transmission systems by introducing a novel channel model, which captures nonlinear interplay between different modes in weak and strong coupling regime: Sparse Identification for Nonlinear Optical communication systems: SINO method. Contrary to the aforementioned approaches that deal only with the intra-channel nonlinearities, here, the proposed SINO method takes into account also the nonlinear interaction between the co-propagating modes by introducing additional inter-channel perturbation terms. The associated complexity scaling was addressed by adapting sparse identification method Kutz , which makes use of the Lasso algorithm LASSO1 ; LASSO2 , thus, enabling computation of perturbation coefficients with inherent principal component analysis. The latter minimizes a mean squared error (MSE) estimation based on the training sequence and removing redundant predictors to improve model accuracy. The method does not require knowledge of the transmission line and is applicable for multi-dimensional multichannel systems.

II Sparse identification for nonlinear optical systems

II.1 Problem formulation

To demonstrate the application of the SINO method we consider the transmission of a $D$ spatial modes and two polarizations along a few mode fiber link of $N_{s}$ amplified spans, see Fig. 1. At the receiver, a demultiplexer performs an initial separation of the spatial super-channel by projecting its tributaries onto a fixed mode basis. After coherent detection, matched filtering and down-sampling, the signals are fed into a linear MIMO equalizer, which compensates dispersion and remaining mode mixing effects. The transmission nonlinearities are treated separately by a non-linear MIMO processor (see Fig. 2). Unlike the linear equalizer, the processor creates a nonlinear matrix $\Theta(\mathbf{X})$ (Fig. 3), which we call library, by forming mixing products of the input symbols. Depending on the channel memory the library will encompass different combinations - nonlinear mixing - of input symbols, thus, the library will have more elements that the original modulation format. See how library expands by different combinations of 16-QAM symbols in Fig. 3: without memory in linear scenario the library has 16 symbols, taking into account nonlinearity (self-phase modulation) results in additional 16 symbols (formed from the original 16-QAM by a simple nonlinear transformation - $x_{k}|x_{k}|^{2}$, while adding memory expands the library further - for $M=1$ there are 48 new symbols formed as $x_{k}|x_{k+1}|^{2}$, further increase of memory expands the library even more. Subsequently, it combines together the resulting terms through a coupling matrix $\Xi$, which contains information about the transmission line parameters and the pulse-shape of the signal. As symbols interact with different weights, the impact of some of these terms can be negligible, hence the matrix $\Xi$ is sparse and we can employ machine learning tools in identifying it. We may write:

[TABLE]

where $\mathbf{Y}$ is a vector of length $m$ representing the output signal. The nonlinear system of Eq. 1 can be viewed as a representation of the output $\mathbf{Y}$ via a set of coordinates $\Xi$ on a complex functional space $\Theta(\mathbf{X})$. The sparsity of $\Xi$ allows us to simplify significantly the computational complexity of the nonlinear filtering process by removing redundant coefficients. As we will see below, this is an inherent characteristic of the lasso method, which has been adopted here for identifying the optimum $\Xi$.

On the other hand, the exact form of the library matrix $\Theta(\mathbf{X})$ depends on the type of non-linearity the equalizer has to address. Given that in our case the fiber channel has memory we should expect each row of the matrix to involve not only elements, such as, $x(t_{1})...x^{2}(t_{1})...$, but also the non-instantaneous interactions between the symbols in different time slots, such us $|x(t_{1})|^{2}x(t_{2})$, therefore $\Theta(\mathbf{X})$ will be a dynamic matrix. Being able to identify its exact form is expected to improve significantly the accuracy and the convergence of the equalization process. This was achieved through a perturbation analysis on the Manakov equations that govern signal propagation along the FMF link and the derivation of a novel discrete-time multivariate channel model, which is described in detail below. An important result of this derivation was the fact that we could analytically define the memory effects characterizing nonlinear interference within each channel and among the co-propagating modes and incorporate this dynamic behavior into $\Theta(\mathbf{X})$.

II.2 SDM Channel Model

For deriving the discrete-time channel model, and consequently the library matrix $\Theta(\mathbf{X})$, we assume that the propagation of an optical signal via an few mode optical fiber is governed by the Manakov equation for the $D$-modes SDM_eq ; Mecozzi_SDM ; Xiao :

[TABLE]

here deterministic distortions are described by fiber losses $\alpha$, second-order dispersion $\beta_{2}$, and nonlinearity coefficient $\gamma$, whereas noise is zero-mean AWGN with variance $\langle\eta(z,t),\eta^{*}(z^{\prime},t^{\prime})\rangle=\frac{N_{D}}{L}\delta(z-z^{\prime})\delta(t-t^{\prime})$ where $N_{D}$ and $L$ are notations for noise spectral density and transmission length correspondingly. The formula captures both a) weak and b) strong coupling regimes, when a) $\kappa_{pp}=8/9,\kappa_{pq}=4/3$ and b) $\kappa_{pq}=\kappa_{pq}=8D/3/(2D+1)$.

Given the expansion of the signal over pulses (i.e. the field component of polarization state $i$ $U_{pi}(t,0)=\sum_{k=-\infty}^{\infty}x_{kpi}f(t-kT))$, after matched filter the continuous-time signal $U_{p}(t,L)=\{U_{p1},U_{p2}\}$ undergoes dispersion compensation and is sampled at $t=kT$: $Y_{kpi}(z)=P^{-1/2}\int dt\mathrm{D}[U_{pi}(t,z)]f(t-kT)$ results in a discrete-time channel model:

[TABLE]

with nonlinear potential

[TABLE]

with AWGN noise term having correlation $\langle\eta_{kpi}(z),\eta_{k^{\prime}p^{\prime}i^{\prime}}^{*}(z^{\prime})\rangle=\delta(z-z^{\prime})\delta_{kk^{\prime}}\delta_{pp^{\prime}}\delta_{ii^{\prime}}$ and where $\Psi_{s}(z)=e^{-\alpha mod(z,L_{s})}$ is the signal power profile and $L_{s}$ is the span length. The coupling coefficients define the memory effects in the channel, it depends on pulseshape and system parameters:

[TABLE]

[TABLE]

Thus, we employ multiple scale analysis with small parameters: $\varepsilon_{pp}=\kappa_{pp}\gamma_{pp}Pz_{d},\ \ \ \varepsilon_{pq}=\kappa_{pq}\gamma_{pq}Pz_{d},\ \ \ \rho=\sqrt{N_{D}B/P}$ with $B$ denoting signal bandwidth. Assuming solution in the form: $Y=\sum_{n=0}^{\infty}\varepsilon^{n}Y_{n}$, in the main order we have linear channel with AWGN noise $\zeta_{k}$ having statistics $\langle\zeta_{k},\zeta^{*}_{m}\rangle=\delta_{km}$:

[TABLE]

Then in the first order over nonlinearity we have the nonlinear distortion:

[TABLE]

Equation 5 reveals the memory property of the channel, which is that the interference between symbols decays exponentially with the symbol distance. Having exact knowledge of this property enables the development of a more efficient equalization process.

The impact of simplifications on coupling matrix was studied in quant , where $C_{mn}$ was approximated by a finite number of quantized levels and the degradation was estimated in terms of Q2-factor. To calculate the latter one needs to calculate numerically a $(2M+1)^{2}$-number of 4-dimensional integrals (for a single channel case or $D(2M+1)^{2}$ for D-mode case). For a specific case of sinc-pulses single-channel WDM the problem has been reduced to 1-dimensional integrals ECOC_Essiambre and has experimentally demonstrated the benefits of the model application for nonlinearity mitigation of intra-channel interference in multichannel-WDM.

In the next section we show how to obtain the matrix numerically using sparse identification and generalization of the model for higher nonlinearity.

II.3 Stage1: Building library

As we know the type of nonlinearity that takes place we can construct the library (see Fig. 4a) as follows

[TABLE]

Also, as all vectors are complex we will separate real and imaginary parts (see Fig. 4b):

[TABLE]

Finally, we can apply this method directly to output signal $A(t)$, while it is very important to receive discrete-time representation by using appropriate filter, e.g. matched filter. The latter is particularly fit for identification of discrete-time model within information-theoretic treatment.

For real symmetric pulseshapes $f(\omega)=-f(\omega)=f^{*}(\omega)$ and $\beta_{1,2}^{p}=\beta^{q}_{1,2}$ we can reduce complexity by exploiting matrix symmetries: $C_{mn}=C_{nm}=C_{-m,-n}=-C^{*}_{-m\neq 0,n}$, which reflects in the appropriate reduced-dimensionality of the library $\Theta(\mathbf{X})$.

II.4 Stage2: LASSO application

As the problem is formulated, next we can apply any of the machine learning algorithms for a sparse solution of an overdetermined system. So that LASSO inherently performs principal component analysis, determining the minimum number of non-zero elements in $\Xi$, consequently, the minimum number of required operations. Once we have the function basis we can use various methods to calculate the sparse matrix $\Xi$. In particular, here we applied least absolute shrinkage and selection operator (LASSO) LASSO1 ; LASSO2 using mean-squared error as a cost function.

As data and library are prepared we are looking to solve the problem to find a minimum:

[TABLE]

where the penalty term:

[TABLE]

In our simulations we applied the lasso fit with ten-fold cross validation (see Fig. 4c). The resulted matrix $\Xi$ comprises coefficients of the $C_{mn}$ matrix in Eq. 4 with simultaneous principle component analysis, which determines the minimum complexity operation. Once, the system is fully identified -the matrix $\Xi$ is determined, one can use it to compensate nonlinear impairments (see Fig. 4d).

III Results and discussion

The simulations parameters are summarized in Table 1. We investigated the transmission of a single wavelength spatial super-channel, along a 10-span FMF link of total length 10x100-km=1000-km. We considered $D$ spatial modes, each one of them including two polarization states, with $D$ taking values of 1, 3, and 6. On each of the 2D subchannels we launched a 4096-symbol, 16-QAM modulated stream of root-raised cosine pules (0.01 roll-off) at a symbol rate of 32GBaud and sampling rate of 16 samples-per-symbol (SpS). A FM-EDFA of 4.5 dB noise figure was considered after each span for compensating the propagation losses. Our focus was on investigating the equalization performance on the Kerr nonlinear effects only. Therefore, for the signal transmission we considered for every mode a Manakov-type of the propagation equation, see 3, where the nonlinear effects are already averaged over all polarization states due to the randomly changing birefringence. Also, for simplicity and without loss of generality, the modes belonged in the same group and they were strongly coupled, which implies a single nonlinear coefficient both for intra- and inter-modal effects (i.e. 1.4 1/W/km), and common group velocity and chromatic dispersion parameters. Finally, the performance was compared in terms of Q2-factor which as calculated from the error vector magnitude (EVM), according to $Q^{2}=1/EVM^{2}$.

In Fig. 5, we compare the equalization performance of our proposed technique (dotted lines) with the case of equalizing linear impairments only (dashed lines), as well as, with the case of having ideal digital back propagation (DBP) (solid lines). Ideal DBP endures the full compensation of deterministic non-linearity and it is limited by nonlinear signal-noise interactions (here we ignored modal dispersion effects). Thus, the numerical estimations can be captured in a straightforward way by a simple analytic expression:

[TABLE]

the coefficient governing signal-signal interactions for ideal DBP scenario $a_{DBP}=0$, while for phase shift equalization

[TABLE]

while for the considered case $\beta_{2}^{p}=\beta_{2}^{q},\beta_{1}^{p}=0,\gamma_{pq}=\gamma_{qq}$ the above expression can be easily simplified as

[TABLE]

(where $L_{eff}$ is an effective length). Since numerical estimations for the coupling matrix $\widetilde{C_{m,n}^{pq}}$ might deviate from the precise estimation, the coefficient transforms to

[TABLE]

here $\delta_{pq}$ is the Kronecker symbol. While the coefficient governing signal-noise interactions is given as

[TABLE]

The aforementioned Eqs. 10-13 give a good estimate of the non-linearity impact for the different mode propagation scenarios (filled symbols in Fig. 5 by green, blue and red for 1, 3, 6- spatial modes - each including two polarizations states). For single mode propagation and at high launched powers, a deviation occurs suggesting the need of considering higher order nonlinear terms.

The proposed nonlinear equalization scheme, results in an above 3 dB improvement for single mode transmission, which decreases to 2 and 1 dB for the 3 and 6 mode cases, respectively. The decrease in performance is attributed to the fact that we had considered only the 1st order terms in the nonlinear expansion of the perturbation model. Our proposed method can be straightforwardly applied to higher orders by augmenting the library matrix $\Theta$ with higher order term of $\mathbf{X}$, however, this is an issue to be addressed in a future work.

The performance of the proposed algorithm depends strongly on the calculation accuracy of the coupling matrix $\mathbf{C}$ through the evaluation of the sparse matrix $\Xi$. The calculations can be further improved by using more advanced methods for sparse matrix calculation than LASSO. Yet, in this simple configuration, which operated with a single sample per symbol and required just a single matrix multiplication, we were able to outperform conventional compensation techniques such as DBP. In Fig. 6 we have compared, in terms of received signal quality, our approach with the symmetric DBP algorithm which was solved for different number of steps per span and 2 samples per symbol. For single mode systems our method slightly outperforms the DBP of 2 steps per span. As the number of modes increases the performance difference is reversed in favor of DBP, so that for 6 modes our method compares to a DBP of single step per span. Even in that case, the complexity of the achieved equalization is significantly lower.

Indeed, an important advantage of the proposed algorithm is that it ensures lowest complexity by evaluating the coupling matrix elements while removing the redundant terms. For example, in the simulated case the number of unique non-zero elements was 48 compared to the estimation of lowest complexity previous algorithm with adaptive filtering McGill 96 components $(M-1)\log(M-1)+3M-1$ and non-optimized 324 $(2M+1)^{2}$ Tao , where $M=\lfloor B^{2}\beta_{2}L/2\rfloor$ - estimation on channel memory ($\lfloor x\rfloor$ denotes a floor function of a variable $x$). Furthermore, the previous algorithms were developed only for compensation of inter-channel nonlinearities, while in multichannel operation the complexity will increase proportionally to the number of channels. And the proposed algorithm with the inherent principal component analysis is crucial.

Finally, the SINO method is scalable to different modulation formats and signal powers as the matrix $\Xi$ remains unchanged, while power and modulation format influence are naturally incorporated in $\Theta$. Thus, once the matrix $\Xi$ has been identified during the establishment of a connection, and as long as the memory properties of the channel do not change, there is no need for retraining the algorithm and the same $\Xi$ can be used for any other modulation format (see Fig. 7). This property makes the method extremely useful for flexible smart-grid network applications.

IV Conclusions

We have developed a low complexity machine learning based nonlinear impairment equalization scheme and demonstrated its successful performance in SDM transmission links achieving compensation of both inter- and intra- channel Kerr-based nonlinear effects. The method operates in one sample per symbol and in one computational step. It is adaptive, i.e. it does not require a knowledge of system parameters, and it is scalable to different power levels and modulation formats. Finally, although it has been developed for single wavelength spatial super-channels it can be straightforwardly expanded to multi-channel systems and to any other type of nonlinear impairment.

This work has been supported by the EPSRC project UNLOC EP/J017582/1. S. Sygletos acknowledges the support from the EU-FP7 INSPACE project under grant agreement N.619732. F. M. Ferreira is acknowledged for useful discussions