1. Introduction

The sorting of biological cells is of key importance in several biomedical applications, like diagnostics, therapeutics and cell biology. However, an accurate classification and separation of different cell types is usually expensive, time consuming and often requires alterations of the samples due to the use of labels, e.g. fluorescent tags, that may hinder subsequent analyses [1]. For these reasons, the development of label-free, high-speed, automated and integrated cell sorting solutions is of particular interest. Among several options, the employment of digital holographic microscopy in microfluidic flow cytometry is a promising candidate. In this technique, the classification is carried out through the analysis of the interference pattern (hologram) projected by the cell when illuminated by monochromatic light. The hologram is acquired by an image sensor and contains information on the 3D refractive index structure of the cells [2]. The large amount of information contained in a cell hologram enables nontrivial analysis and classifications. On the other hand, the computational cost of elaborating such a complex source of information by reconstructing the image from the hologram is a major hindrance to an increase in the cell sorter throughput, e.g. by parallelization of the process.

An important reduction in the required computing power was achieved by bypassing the reconstruction of the cell image and directly processing the acquired hologram with a machine learning algorithm that carries out the classification task [2, 3]. However, a further improvement of the classification simplicity and performance is desired in order to fully exploit the potential of this implementation.

Pillar scatterers for an extreme learning machine (ELM) implementation

In this work, a passive, linear, integrated photonic stage is proposed as an interface between the hologram projection and the image sensor, in order to simplify the classification process in the relatively slow electric domain. In particular, the forward scattered light coming from a cell and the corresponding background are made to propagate through a collection of silica pillar scatterers of elliptic cross section embedded in silicon nitride (Fig. 1). For computational time reasons, this process has been investigated via 2D finite-difference time-domain (FDTD, [4]) simulations as a proof-of-concept, approximating the 3D case of a cell flowing in a microfluidic channel interfaced with a photonic chip. The far-field intensity of the light exiting the scatterers cluster is collected by an array of virtual pixels that approximately simulate a line scan image sensor. The pixel outputs are fed into a linear classifier that can be, for example, implemented in the electric domain.

 

Fig. 1 Schematic of the classification process. From right to left: a monochromatic plane wave impinges on a microfluidic channel containing a cell in water $(n_{H_{2}O}~1.34)$, which has a low refractive index contrast (n_cytoplasm = 1.37, n_nucleus = 1.39); the forward scattered light passes through a collection of silica scatterers $(n_{Si{O}_2}~1.461)$ embedded in silicon nitride$(n_{S{i}_3{N}_4}~2.027)$ and organized in layers; the radiation intensity is then collected by a far-field monitor, which is divided into bins (pixels); each pixel value is fed into a trained logistic regression, which classifies the cell as a “normal” cell (small nucleus) or as a “cancer” cell (big nucleus). The logistic regression consists of a weighted sum of the pixel values. The weights are trained so that the sum exceeds a threshold value only if a certain input class is recognized.

Download Full Size | PDF

The photonic stage containing the scatterers is intended to exploit the nonlinearity of the transfer function that relates the phase shift accumulated by the light through the cell to the corresponding interference pattern measured by an image sensor. Basically, the system can be seen as an hardware implementation of an extreme learning machine (ELM), i.e. a feedforward neural network whose hidden nodes are randomly generated and kept fixed [5] (see subsection Nonlinear phase sensitivity in Results and discussion for a further argumentation). In the past decade, this machine learning paradigm has been the object of an increasing interest in various research field due to its remarkable efficiency, simplicity, and generalization performance [6]. Generally, the main advantages of ELMs with respect to other machine learning techniques are that only a linear readout (in this case a linear classifier) needs to be trained and that it is easily implemented in hardware. In this case, the pillar scatterer stage determines the ELM hidden node structure by projecting onto the far-field intensity a very intricate nonlinear mapping (based on sinusoidal functions of the phase information). This parallel processing is carried out nearly instantaneously with respect to both the cell movement and the operating speed of an electronic computer, providing an important advantage over other machine learning solutions in the electric domain. It should be stressed that the phase-to-intensity nonlinearity is already expressed by the interference pattern projected by the cell alone, without scatterers. However, the complexity of such a nonlinear mapping can be enhanced and controlled by the use of scatterers in order to increase the performance of a subsequent linear classification. The exploitation of light interference in order to fabricate a passive integrated reservoir computing implementation (which is based on similar principles of ELMs but is applied to time-dependent signals [7]) using linear optical media was demonstrated in [8]. However, the time-dependent input information was encoded in the intensity of a laser signal and therefore the transfer function to the output detector was quadratic (amplitude to intensity). Here, the input information is encoded in the phase of a laser signal, such that the readout transfer function can be sinusoidal with respect to that input (see subsection Nonlinear phase sensitivity in Results and discussion for a further argumentation). The sine, for example, can be profitably employed as activation function in feedforward neural networks under suitable conditions [9].

2. Simulation details

Concerning the optical simulations, the FDTD method was chosen because of its reliability in dealing with complex dielectric structures (provided that a fine enough mesh size is chosen) and because it allows for an intuitive understanding of the computed physical process. In particular, visualizing the light field as a function of time provided a qualitative idea of how much the light signal coming from the cell was mixed by the dielectric scatterers.

Both a cell in a microfluidic channel and the proposed silica scatterers were modeled in the same 2D FDTD simulation (Fig. 1) employing Lumerical’s FDTD Solutions software. A monochromatic plane wave (vacuum wavelength λ = 532 nm as in [2]) of constant intensity impinges transversely on a microfluidic channel (15 µm wide) filled with water. The channel interfaces with a region containing layers of elliptic scatterers (0.5 µm wide and 1 µm long) made of SiO₂ embedded in Si₃N₄. At the end of the scatterer region, a vertical far-field monitor (with an angular resolution of 55.6 points/deg) covers the total length of the simulated space. The simulation region is 28 µm long along the vertical direction and from 20 µm to 30 µm along the horizontal direction, depending on the number of scatterer layers. The FDTD mesh size is ∼ λ/29.

It should be stressed that the scope of this work is to provide a proof-of-concept of a new approach to photonics machine learning, that can be generalized to many other implementations besides the discussed examples. Therefore, the dimensionality of the simulation, the size of the structures and the cell model are a consequence of a trade-off between closeness to reality, saving of computational time and the search of a sufficiently complex (but not overly so) task. In any case, all the simulated objects, aside from the scatterers, were designed independently of the classification results.

Cell model

The cell model is composed of a cytoplasm region (n_cytoplasm = 1.37) surrounding a nucleus region (n_nucleus = 1.39). An example is shown in Fig. 1. In order to generate a different cell shape for each simulation, a 2D randomized cell model was employed, based on distorted ellipses. Considering the ellipse equation in polar coordinates (ρ is the distance from the origin and θ the angle with respect to the horizontal axis)

In this work, two different classification tasks are considered. The first is based on average nucleus size and aims to distinguish between “normal” cells (small nucleus) and “cancer” cells (bigger nucleus). The names in quotation marks were chosen because of the common tendency of cancer cells to show evident irregularities in nucleus size [11]. The second task is based on nucleus shape and aims to distinguish between “lymphocytes” (big quasi-spherical nucleus) and “neutrophils” (nucleus divided in 3 lobes). The names in quotation marks refer to two among the most common white blood cells that are present in human blood. These two cell models are, physically and biologically speaking, only rough representations of their real counterparts when flowing in a liquid medium [12]. The scope of this work, though, is to demonstrate the differences between two classification approaches applied on models that are sufficiently representative of real cell cases from a machine learning perspective.

“Normal” and “cancer” cells

The parameters for “normal” and “cancer” cell models are chosen as follows (subscript c stays for “cytoplasm” and n for “nucleus”): a_c = b_c = 5 × (1 + 0.1ε) µm, A_c = 0.1 × (1 + 0.9ε), ω_c = 3 × (1 + 0.9ε)rad⁻¹, B_c = 0.2 for cytoplasm and

“Lymphocytes” and “neutrophils”

The main difference between the “lymphocyte” and the “neutrophil” models is that the first has a big quasi-spherical nucleus while the latter has 3 quasi-spherical nuclei whose average total area is equal to the average area of the “lymphocyte” nucleus. The corresponding parameters are chosen as follows (subscript c stays for “cytoplasm” and n for “nucleus”): a_c = b_c = 5.5 × (1 + 0.1ε) µm, A_c = 0.03 × (1 + 0.67ε), ω_c = 3 × (1 + 0.9ε)rad⁻¹, B_c = 0.2 for cytoplasm and

In order to provide the reader with an intuitive idea of the shapes and the randomness of the employed cell models, some examples are represented in Fig. 2.

 

Fig. 2 Examples of cells automatically generated by the employed randomized models. a) Comparison between generated examples of “normal” cell and “cancer” cell. b) Comparison between generated examples of “lymphocyte” and “neutrophil”.

Download Full Size | PDF

Scatterer configuration

The scatterer configuration has a large number of degrees of freedom and its complete exploration, even using smart methods such as evolutionary algorithms, would be computationally quite expensive. In fact, for each tested configuration, hundreds or thousands of simulations have to be performed in order to provide the classifier with a sufficient number of training and test samples. Therefore, only a few general parameters that control the complexity of the collected interference pattern were explored, looking for a maximum in the classifier performances. The scatterers are placed in vertical layers (e.g. Fig. 1) with an average vertical distance of 1µm between their centers. The center of each scatterer is randomly displaced with respect to their unperturbed center in the layer, both along the vertical and the horizontal directions. All the random displacements belonging to the same architecture are sampled from the same uniform probability distribution. The considered parameters for the structure optimization are:

 

Fig. 3 Far-field intensity profiles of the light scattered by a cell: a) without the presence of scatterers, the interference pattern is relatively simple and smooth, most of the intensity is confined between −6° and 6°; b) with 1 layer of scatterers, the far-field intensity is distributed around periodically placed peaks, most of the field stays between −40° and 40°; c) with 4 layers of scatterers, the far-field intensity is distributed in a complex pattern mostly between −60° and 60°.

Download Full Size | PDF

3. Machine learning aspects

The far-field intensity profile was divided into N_pix bins (or pixels) and the integration over each bin was fed into a logistic regression. The Scikit-learn Python library [10] was employed, using the “liblinear” solver. For each tested scatterer configuration a number N_samp of simulations was performed randomly varying the cell shape, as described in the previous section. In particular, the classification results reported in the next section are obtained from sample sets of N_samp = 3200 samples each. In half of the N_samp simulations a “normal” cell was considered, while in the remaining half a “cancer” cell (with bigger nucleus) was used. 75% of these two sets was employed in the training of the logistic regression, while the rest was used as test set.

Gaussian white noise was added a post)eriori on the interference patterns before they were divided into bins. Different values for the noise standard deviation were chosen: 1%, 5%, 10%, 20%, 30%, 50% and 100% of the average intensity over the sample set.

A study on the dependence of the classification test error on the regularization strength (L1 and L2) and on N_samp was carried out in the two cases with and without the use of scatterers, on a set of 2000 samples each. For this investigation, a 4 layer scatterer configuration was considered, with A_r = 150nm and D = 1.846µm (this particular configuration was one of the best performing in a preliminary exploration based on a lower (800) number of samples). The study pointed out that regularization had no significant positive effect on the performances of the two classification tasks. Furthermore, it showed that the learning curve (test error v.s. number of samples) converged around N_samp = 800 when the scatterers were used. Thus, for the classifications presented in the next section no regularization was considered, but the investigation on the effect of different noise levels partially covers the topic at least with regard to the L2 regularization. The classification dependence on the number of pixels and training epochs was instead kept under direct control by performing sweeps for each tested configuration.

The presented results are obtained through a validation process in which the simulated data samples are randomly shuffled before they are split into training and test sets. From the results generated by repeating this procedure 20 times a mean value and a confidence interval (chosen as ± twice the standard deviation) were calculated and plotted. Note that a different noise vector is added to the intensity profiles after each shuffling.

4. Results and discussion

In this section, the results concerning the classification improvements due to the employment of silica scatterers are discussed in 4 parts:

Green light

Let us consider a green laser source (λ = 532nm) and let us compare the classification error on the test samples when no scatterers are present and when, instead, 4 scatterer layers are employed (considering the random displacement amplitude A_r = 150nm and the layer distance D = 1.846µm). In the first case, when no scatterers are used, the angle range for which the far-field intensity is not negligible is estimated to be between −6° and 6° (Fig. 3(a)) and this is the range to which the number of pixels refers. In the latter case, when the 4 scatterer layers are present, the chosen angle range is between −60° and 60° (Fig. 3(c)). Let us stress that the error rate expected value and the confidence intervals drawn in all the following graphs are extracted from the results obtained from 20 random permutations of the simulated samples and correspond to an optimal choice of the number of training epochs. The resulting error rate values for different numbers of pixels and for different noise levels (Fig. 4) show that the use of scatterer layers allows for a significant error rate reduction (up to ∼ 50%), provided that a sufficient number of pixels and a low enough noise level are considered. The increased sensitivity of classification performance towards added noise level when scatterers are used is ascribed to the fact that the scatterers presence unfolds the cell diffraction pattern into a higher number of components (Fig. 3(c)) that may be important for classification. Thus, it is probable that some of these components have low intensity with respect to the average pattern intensity and are therefore easily overcome by high relative levels of noise.

 

Fig. 4 Comparison between the test error rates of “normal” and “cancer” cell classification, corresponding to the absence (in red) and the presence (in blue) of scatterers (4 layers, A_r = 150nm, D = 1.846µm). A green laser source (λ = 532nm) is employed. a) Test error rate as a function of the number of employed pixels, with 5% added white noise. The darker and the lighter versions of the two line colors respectively represent the mean value and the confidence interval (of 2 standard deviations) over the 20 sample sets generated for validation. b) Test error rate (averaged on the values obtained considering N_pix = 250, 260, …, 300) as a function of the added noise percentage. In order to avoid error bar overlap, some of the blue points are slightly shifted to the right. Both the plots show that the scatterers’ presence allows for an error rate reduction up to ∼ 50%, provided that a sufficient number of pixels and a low enough noise level are considered.

Download Full Size | PDF

The chosen scatterer configuration (N_layers = 4, A_r = 150nm and D = 1.846µm) was the best performing among the tested combinations (800 simulated samples each) of N_layers = 1, 2, 3, 4, 5, A_r = 50nm, 100nm, 150nm, 200nm, 250nm and D = 1.846µm, 2.85µm, 3.40µm. The two last values of D respectively give a maximum and a minimum in the far-field transmission, considering four layers of scatterers without any random displacement. Performances similar to the ones presented in this paper were obtained for a broad number of cases. Therefore, the main conclusion of this (non-exhaustive) geometrical exploration is that the employment of scatterer layers of different kinds can decrease the classification error rate significantly, up to 50%.

Nonlinear phase sensitivity

If there seems to be a limit to the improvements obtainable via the use of scatterers when a green coherent source is used, this limit is easily overcome by decreasing the source wavelength. An explanation of this effect can be provided by the following simplified argumentation.

Let us neglect for a moment the light deflection due to the cell refractive index structure and let us consider only the phase shift of the light along all the possible fixed paths (here labeled with n) through the cell to one pixel on the screen. At the top of Fig. 5 a schematic representing 3 of these paths is shown. Let us state that the light along a path n has unitary initial amplitude and null initial phase (the reasoning is independent from the initial conditions) and that it accumulates a phase shift θ_n through the cell. Moreover, let us say that its amplitude is reduced by a factor A_n and its phase is increased by ϕ_n along the path to the pixel excluding the path inside the cell. Thus, the complex amplitude of the radiation impinging on the pixel is ${\mathrm{\Sigma }}_n{A}_n{e}^{i({\theta }_n+{\varphi }_n)}$ and the acquired intensity I is proportional to:

 

Fig. 5 Two equivalent schematics of the proposed classifying system. At the top, a physical schema shows an example of amplitude and phase evolution along 3 optical paths that end up impinging on the same pixel of the image sensor. The acquired light intensity is then weighted and summed by a linear classifier. At the bottom, a diagram (under the form of a neural network architecture) represents the corresponding mathematical operations on the light phase accumulated through the cell refractive index structure (see Eq. (7)). For simplicity’s sake, the light deflection due to the cell presence is neglected and thus the amplitudes A_n and the factors A_nm, B_nm and C are considered as constants with respect to the inputs θ_n.

Download Full Size | PDF

Let us now consider, for instance, the difference Δθ between phase shifts corresponding to a path through the nucleus and a neighboring path that instead does not intersect the nucleus. Let us call this phase shift difference Δθ_n in the case of a “normal” cell (smaller nucleus) and Δθ_c in the case of a “cancer” cell (bigger nucleus). We can intuitively say that if the readout linear classifier is able, for example, to detect the difference ΔI between the intensity contributions produced by Δθ_n and Δθ_c respectively among the other intensity contributions, the system can be successfully trained in carrying out the classification task. From Eq. (7) follows that an estimate of this critical intensity difference is given by:

UV light

Considering a system comprising 4 scatterer layers (A_r = 150nm and D = 2.85µm) with source wavelengths λ = 532nm, 400nm, 300nm and 200nm, the overall change in the acquired diffraction patterns was calculated for small modifications (of 130 nm) of the nucleus size (Fig. 6), keeping the rest of the simulation parameters fixed. The change between two interference patterns has been calculated by summing the absolute values of the elements of their point-wise difference vector. The dispersion of the employed materials was accounted for and the same absolute contrast was kept between the refractive indexes in the cell model and the water. As expected, the impact of such nucleus modifications on the acquired interference pattern increases significantly as the laser wavelength is decreased. On the other hand, the same happens to the transmission losses due to the scatterers presence (from ∼ 29% at λ = 532nm to ∼ 54% at λ = 200nm).

 

Fig. 6 Change in the acquired diffraction pattern due to small increases (130 nm) of the nucleus size as a function of the starting nucleus size. The change between two interference patterns has been calculated by summing the absolute values of the elements of their point-wise difference vector. It can be noted that the smaller the employed wavelength is, the larger the pattern modifications become, implying an easier classification task.

Download Full Size | PDF

In order to consider a plausible UV laser source, the classification performance was investigated using λ = 337.1nm. In particular, a comparison between the case without scatterers and the case with 4 scatterer layers (A_r = 150nm and D = 2.85µm) is discussed (Fig. 7). In comparison with the results obtained using a green laser source (Fig. 4), the UV laser implementation substantially enhances the beneficial effects of the scatterers employment. In fact, not only the achieved error rate reduction (up to an order of magnitude) is remarkably larger, but also the robustness to noise is significantly improved. This confirms the predictions arisen from the previous discussion on the nonlinear phase sensitivity of the cell diffraction pattern: as expected, even a slight nonlinearity increase of the hidden nodes processing determines a considerable enhancement of the classification performance.

 

Fig. 7 Comparison between the test error rates of “normal” and “cancer” cell classification, corresponding to the absence (in red) and the presence (in blue) of scatterers (4 layers, A_r = 150nm, D = 2.85µm). An UV laser source (λ = 337.1nm) is employed. a) Test error rate as a function of the number of employed pixels, with 5% added white noise. The darker and the lighter versions of the two line colors respectively represent the mean value and the confidence interval (of 2 standard deviations) over the 20 sample sets generated for validation. b) Test error rate (averaged on the values obtained considering N_pix = 250, 260, …, 300) as a function of the added noise percentage. Both the plots show that the scatterers presence allows for a considerable error rate reduction (up to an order of magnitude) in the whole investigated ranges of pixel number and noise level.

Download Full Size | PDF

Generalization to another task

It is worth wondering whether and how the proposed technique can generalize to other kinds of tasks. In fact, one of the key advantages of ELM implementations is that one network type can work for different applications [6]. In order to (non exhaustively) test this property, the classification of “lymphocyte” and “neutrophil” cells was attempted using exactly the same hyperparameters of the previously discussed classification with a green light source, i.e. without performing any kind of further optimization. The obtained classification improvement due to the use of scatterers (Fig. 8) is even greater than the one for “normal” and “cancer” cell classification (Fig. 4). Thus, the technique discussed in this paper was proven to generalize to two fundamentally different classification tasks, as one is based on nuclei with different area but same shape, the other on nuclei with same overall area but different shape.

 

Fig. 8 Comparison between the test error rates of “lymphocyte” and “neutrophil” cell classification, corresponding to the absence (in red) and the presence (in blue) of scatterers. The scatterers configuration and the light source are the same as in Fig. 4. a) Test error rate as a function of the number of employed pixels, with 5% added white noise. The darker and the lighter versions of the two line colors respectively represent the mean value and the confidence interval (of 2 standard deviations) over the 20 sample sets generated for validation. b) Test error rate (averaged on the values obtained considering N_pix = 250, 260, …, 300) as a function of the added noise percentage. In order to avoid error bar overlap, some of the blue points are slightly shifted to the right. Both the plots show that the scatterers’ presence allows for an error rate reduction greater than 50%, provided that a sufficient number of pixels and a low enough noise level are considered.

Download Full Size | PDF

5. Conclusions

In this work, a simple passive integrated photonic structure is proposed to improve the performance of a linear classifier employed for the optical classification of different cell kinds. The main advantage of this application is that the light processing performed by the photonic structures is much faster than traditional computations in the electric domain. Therefore, only the computational power needed to perform a linear combination of the intensities measured by a limited number of pixels is required for the cell classification. This allows, for example, for high throughput label-free cell sorting employing several parallel channels along which the diffraction patterns of the lighted cells are acquired. Another advantage is that a linear classifier is usually much easier to train with respect to more complex machine learning applications, such as multi-layer neural networks.

The proposed integrated structure consists of a collection of silica pillars embedded in silicon nitride that scatters the monochromatic light diffracted by a cell. When a green laser source is employed, the scatterers’ presence leads to a significant improvement in performance (up to 50%) of the classification of simulated cells with different nucleus size, provided that a sufficient number of pixels and a low enough noise level are considered. Furthermore, an intuitive reasoning on the nonlinear phase sensitivity of a cell diffraction pattern suggests that remarkably bigger improvements can be achieved, for example, by decreasing the source light wavelength. Indeed, when a UV laser source (λ = 337.1nm) is employed, the scatterers presence leads to a considerable error rate reduction (up to an order of magnitude) for all the considered pixel numbers and noise levels. Finally, the generalization capability of the proposed technique was tested by considering the classification of cells with equal overall nucleus size but different nucleus shape. Indeed, also in this case the use of scatterers was proven to be beneficial, increasing the classification performance of more than 50%.