Image-Based flow cytometry and angle-resolved light scattering to define the sickling process

(1)

Image-Based Flow Cytometry and Angle-Resolved

Light Scattering to De

ﬁne the Sickling Process

Polat Goktas,

1,2

* Ilya O. Sukharevsky,

3

Sandra Larkin,

4

Frans A. Kuypers,

4

Ozlem Yalcin,

2

Ayhan Altintas

1

Abstract

Red blood cells (RBCs) from sickle cell patients exposed to a low oxygen tension reveal highly heterogeneous cell morphologies due to the polymerization of sickle hemoglobin (HbS). We show that angle-resolved light scattering approach with the use of image-basedflow cytometry provides reliable quantitative data to define the change in mor-phology of large populations of RBCs from sickle cell patients when the cells are exposed for different times to low oxygen. We characterize the RBC morphological profile by means of a set of morphological and physical parameters, which includes cell shape, size, and orientation. These parameters define the cell as discocyte, sickle, elon-gated, as well as irregularly or abnormal RBC shaped cells, including echinocytes, holly-leaf, and granular structures. In contrast to microscopy, quick assessment of large numbers of cells provides statistically relevant information of the dynamic process of RBC sickling in time. The use of this approach facilitates the understanding of the pro-cesses that define the propensity of sickle blood samples to change their shape, and the ensuing vaso-occlusive events in the circulation of the patients. Moreover, it assists in the evaluation of treatments that include the use of anti-sickling agents, gene therapy-based hemoglobin modifications, as well as other approaches to improve the quality of life of sickle cell patients. © 2019 International Society for Advancement of Cytometry

Key terms

computationalﬂow cytometry; red blood cells; sickle cell disease; HbS polymerization; light scattering; precision and personalized patient-oriented medicine

S

ICKLECell Anemia, an inherited red blood cell (RBC) disorder, was deﬁned as the

first molecular disease (1). The simple point mutation in beta globin generates a hemoglobin molecule that polymerizes under low oxygen tension changing the dis-coid morphology of the RBC to the shapes that gave the disorder its name, an altered the cell rigidity, altered cell surface, adhesion, a short life span, and intravas-cular hemolysis and inflammation (2–8). Hence, the process of morphology change (sickling) of RBC from sickle cell patients (SS-RBCs) is the fundamental pathophysi-ological hallmark of the vasculopathy, and organ damage that defines sickle cell dis-ease (SCD). Perhaps more importantly, the rate by which sickling takes place is significant. A decrease in the rate of polymerization-induced sickling is thought to have a profound effect on the clinical outcome of patients. Rapid deoxygenation combined with morphology assessment of large numbers of cells to provide high-throughput data assessments of the heterogeneous cell population is therefore important. Microscopy does not provide the information to describe large numbers of highly heterogeneously shaped RBC observed in a typical blood film (9). The blood smearing technique might also lead to alter morphological RBC structures especially for irregularly or abnormally shaped cells, and induces error-prone proce-dure for the classification of the “true” morphology profile of RBCs. Hence, fixing of

1

Department of Electrical and Electronics Engineering, Bilkent University, Ankara 06800, Turkey

2

Department of Physiology, School of Medicine, Koc University, Istanbul 34450, Turkey

3

Department of Electrical and Computer Engineering, Chair of High-Frequency Engineering, Technical University of Munich, Munich 80333, Germany

4

Children’s Hospital Oakland Research Institute, Oakland, California, 94609 Received 3 January 2019; Revised 15 March 2019; Accepted 21 March 2019 Grant sponsor: Alexander von Humboldt Foundation; Grant sponsor: IEEE Anten-nas and Propagation Society (AP-S) Doctoral Research Grant

Additional Supporting Information may be found in the online version of this article.

*Correspondence to: Polat Goktas, Department of Electrical and Electronics Engineering, Bilkent University, Ankara 06800, Turkey. Email: pgoktas@ee. bilkent.edu.tr

Published online 13 April 2019 in Wiley Online Library (wileyonlinelibrary.com) DOI: 10.1002/cyto.a.23756

(2)

the morphology of the cell population is meaningful. In con-trast to light microscopy, conventionalflow cytometers (FCs) can detect cells of interest that are fluorescently labeled, and provide automatic cell sorting by capturing the scattered light and fluorescence intensity values (10). However, the optical signature of conventional FCs is based on the light scattered from individual particles at fixed observation angles only (e.g., forward scattering [FSC-A] and side scattering [SSC-A]), and yields limited information about morphology and biophysi-cal properties of cells. Cell segmentation methods such as Cel-lProfiler (11), Fiji (12), CellSegm (13), and deep convolutional neural networks (dCNNs) (14), as well as clustering algorithms based on k-means, t-SNE, and SPADE techniques (15) have been used to analyze cellular structures in the data obtained by light microscopy orflow cytometry. However, the performance of these algorithms to define morphology depends on proper initiation of labeled parameters by hand, and extensive algo-rithm training, of a large number of bright- and dark-field images. This has not been shown to be practical for the assess-ment of large populations of RBC with a highly heterogeneous morphology that is different from sample to sample. Imaging flow cytometry (IFC) combines the single cell image-based screening of microscopy with the high-throughput speed of conventional FCs (16–20). IFC can provide information regarding cellular structures, such as cell size, shape, orienta-tion, signal intensity, locaorienta-tion, and texture for each cell with high statistical validity through bright- and dark-field chan-nels, which require no staining procedures, and fluorescence channels. In a recent publication, van Beers et al. (21) revealed the first use of IFC to measure normal and abnormal RBCs obtained from bright-field image data under hypoxic condi-tion (2% oxygen) at one time sickling point after the onset of deoxygenation.

In this study, we investigated the changes in time as cells undergo deoxygenation. We aimed to show that a quick and well-defined deoxygenation, morphological fixation, together with a computational Nystrom-type meshless discretization morphology framework, can be used to follow the rate of sickling with an imaging FC. This approach improves the automated, high-throughput classification performance of the heterogeneous morphological RBC shapes as they develop in time in SCD blood. The combination of Thickness Minimum, Length, and Area features as a Sickle Index shape feature obtained from a user-defined custom mask in the IFC data enhances the potential of the optimal discrimination between normal, intermediate, and sickle cells of SS-RBCs. Impor-tantly, the main merit of the paper lies in showing for thefirst time that the angle-resolved light scattering analysis based on boundary shape structures is correlated with the measured SSC-A pattern realized by IFC to provide the refractive index distribution for each SS-RBCs. Moreover, our proposed shape feature criteria with the angle-resolved light scattering approach when the cells are exposed for different times to low oxygen would provide better identification of “true” nor-mal, intermediate, and sickle cell region boundaries in IFC data. With this approach, we are able to notably reduce the data analysis procedure to identify normal and sickle cells

from a cell population in a given experiment through IFC with angle-resolved light scattering method. An automated cell morphological/physical classiﬁcation method provides an effective and accurate tool in monitoring of SCD status in dif-ferent clinical conditions, resulting in a rapid evaluation for treatment scenarios in SCD.

M

ATERIALS AND

M

ETHODS

Sample Preparation and Image Data Collection Whole blood (50μl) was added at time zero to 1 ml of buffer (HEPES-buffered saline [HBS] with 2% bovine serum albumin [BSA], pH 7.4) equilibrated in a Tonometer (model 1237, Instrumentation laboratory, Bedford, MA) at a set oxy-gen tension of 1% oxyoxy-gen (PO2 = 8 mmHg). A regular spin

(10 s on, 10 s off) cycle was applied to rapidly equilibrate the RBC with the oxygen tension applied. At different time points, 10μl fractions were automatically collected and fixed by addi-tion to HBS with 2% paraformaldehyde and 0.5% glutaralde-hyde. Samples were stored and subsequently measured using the ImageStream®XMark IIflow cytometer (Amnis Corpora-tion, Seattle, WA, part of Millipore Sigma). This study was approved by the UCSF Benioff Children’s Hospital Oakland’s IRB. Data representing at least 10,000 intact single RBC were collected based on bright-field (BF) cell images in focus (gated on 30–70 gradient root mean square [RMS]) for events 50–200 Area and 0.3–1.0 Aspect Ratio in Amnis INSPIRE software. IFC Data Analysis in IDEAS®

IFC data exploration and analysis software (IDEAS® ver-sion 6.2.188.0, Amnis Corporation, part of Millipore Sigma) was used for the determination of parent populations, such as normal, intermediate, and sickle cell region boundaries in the deoxygenated SS-RBCs. Ten thousand events of SS-RBCs deoxygenated for 0, 2, 5, 10, and 20 min. were analyzed with the use of a cell classifier set at a minimum area of 30 μm2to exclude 1μm nonfluorescent polystyrene beads (speed beads). First, focused cells from all events were selected from the BF channel by the Gradient RMS feature as a criterion of BF Gradient RMS feature above 30 arbitrary units. Focused cells were further masked using a user-defined custom mask, which is a combination of the Range, System, and Object function masks, to reduce noise and artifacts arising from objects that change feature properties resulting from more than one cell in a frame. This user-defined custom mask was visually determined to betterfit the actual contour of the BF region as compared to the default M01 mask (Fig. 2), and this approach was used for all subsequent feature calculations. For the determination of a parent population, a bivariate plot of BF Area versus BF Aspect Ratio was used and a region was created to select a single cell population in order to exclude small particles, doublets, and aggregates. This parent population was then used for subsequent assessment of mor-phological SS-RBC subpopulations and determination of dis-tinguishing feature sets. Details about masks and criterion of features for the classification of deoxygenated SS-RBCs used throughout this article are given in“Results” section.

(3)

Angle-Resolved Light Scattering Technique, Muller Boundary Integral Equation

The analysis of a light scattering pattern by an individual RBC has been reported by using various electromagnetic numerical methods (10, 22–29) and integral equation methods. In this article, a boundary integral equation method (Muller Boundary Integral Equation [MBIE]), which enables fast and efﬁcient procedures for spherical and nonspherical biological particles with a Nystrom-type meshless discretization, was used. When the boundary is smooth and has a regular 2π-periodic parametrization, the Nystrom-type discretization (30–33) using trigonometric approximation of integrand functions is con-venient, as it leads to exponential convergence of numerical solution. The time factor e−jwt is assumed and suppressed throughout the article. The numerical code developed in the work is available upon request from the corresponding author.

The shapes of SS-RBCs are modeled as a homogeneous dielectric particle with boundaries by means of the Gielis for-mula (34): rð Þ =θ 1 acos mθ 4 β+ 1 bsin mθ 4 γ !−1₌ α ; a, b > 0: ð1Þ

where r is the radial coordinate;θ is the polar angle; a > 0 and b > 0 are the size coefficients of the shape structure; m is the radial symmetry number; andα, β, and γ are shape coeffi-cients. Although parameters in the formula have a substantial meaning on the modeling of the shapes, they are not simple numbers to provide more insight into our understanding of each individual coefficient (34) combined to create the shape. However, broadly, the radial symmetry (or the number of vertices) is adjusted by the parameter m, inside the sine and cosine functions. By changing the value of m, the number of sides of a shape can be adjusted. The parameterα determines sharpness or flatness of corners and convexity of sides in the cell shape. The parameters β and γ determine whether the structure is inscribed or circumscribed in the shape (e.g., β, γ < 2.0 for sub-polygon or β, γ > 2.0 for super-polygon structures) (35). Thus, the modeling of the cell struc-ture is depending on the combination of the a, b, m,α, β, and γ coefficients in the formula. The MBIEs for unknown equiv-alent currents at the boundary of scatterer (36) are formulated as follows: ð L K11ðr ! , r!0ÞU r !0 ds0− ð L K12ðr ! , r!0ÞV r !0 ds0= Uinc r!: ð2Þ pi+ pe 2pe V r ! + ð L K21 r!, r!0 U r !0 ds0 −ð L K22 r ! , r!0 V r !0 ds0= Vinc r!, r!2 L: ð3Þ

where r!= x, yð Þ and r!0= xð 0, y0Þ are the integration and obser-vation points, respectively. U corresponds to theﬁeld compo-nent Ezin the inner domain Di. Moreover, V r

!₀

is the limit value of the normal derivative of the totalﬁeld on the closed

cross-section contour L of the scattered object from the inner side of it, ds0 is the elementary arc length, and the constants are pi,e= 1/μi,e in the E-polarization case. The primaryﬁeld

and its normal derivative are deﬁned as Uinc and Vinc, respectively.

The kernel functions of the MBIE are deﬁned as follows: K11ðr ! , r!0Þ = ∂Giðr ! , r!0Þ=∂n0−∂Geðr ! , r!0Þ=∂n0: ð4Þ K12ðr ! , r!0Þ = Giðr ! , r!0Þ− pð i=peÞGeðr ! , r!0Þ: ð5Þ K21ðr ! , r!0Þ = ∂2Giðr ! , r!0Þ=∂n∂n0−∂2Geðr ! , r!0Þ=∂n∂n0: ð6Þ K22ðr!, r!0Þ = ∂Giðr!, r!0Þ=∂n− pð i=peÞ∂Geðr!, r!0Þ=∂n: ð7Þ

where Gð Þi, e = Gð Þi, eðr!, r!0Þ = i=4ð ÞH0ð Þ1kð Þi, eρ! are the Green

functions of the corresponding media; ke= ko ffiffiffiffiffiffiffiffiεeμe

p _and

ki= ko ffiffiffiffiffiffiffiεiμi

p

, where ko is the free-space wavenumber, and

ρ!

= r!−r!0 with ρ = ρ ! is the distance between the points r

!

and r!0, and H₀ð Þ1ð Þ is the Hankel function of the ﬁrst kind: and zero-order.

For the case of piece-smooth contour (37,38), we approxi-mated unknown equivalent currents in Eqs. (4)–(7) with a step-constant function. Following Smotrova et al. (32), we apply the following quadrature formula:

ð2π 0 ln 4sin2t−et 2 f t,et det≈ X 2N−1 l = 0 PN_l ð Þf t,tt ð lÞ: ð8Þ

where f has a logarithmic singularity, and

PN l ð Þ = −t 2π N XN−1 m = 1 cosfm tð −tlÞg m − π N2cos N tf ð −tlÞg; tl=πl=N, l = 0,1,…,2N −1: ð9Þ ð2π 0 f t,et det≈ π=Nð ÞX 2N−1 l = 0 f t, tð lÞ: ð10Þ

We obtain the following matrix equation block after such a discretization: I + K ð Þ U V = U inc 2pe pi+ peV inc ! : ð11Þ

where the vectors of unknowns are U = {U(tl)}l = 0, 2N− 1and

V = {V(tl)}l = 0, 2N− 1, and I is a unit 4N × 4N matrix, and

the matrix K has a block structure:

K = eK11ðtk, tlÞ eK12ðtk, tlÞ 2pe pi+ pe eK21ðtk, tlÞ 2pe pi+ pe eK22ðtk, tlÞ 0 @ 1 A: ð12Þ where k, l = 0, 1,…, (2N−1). In Eqs. (4)–(11), the quantities

ρ!

:n!

, ρ!:n!0, and n!0:n! are the scalar products of the corresponding vectors. Assuming that the variables t and et in

(4)

parametric form correspond to r!and r!0, the distance between two points on the boundary contour isρ = x t ð Þ−x et 2+

y tð Þ−y et

2

Þ1=2, and the outer normal vector in the bound-ary is n!= 1 I tð Þ dy dt− dx dt

, where I(t) is the Jacobian matrix as I tð Þ = ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi dx=dt ð Þ2_{+ dy=dt}_ð _Þ2 q : Each eKijis a 2N× 2N matrix: eKij= PkNð ÞKtl ij1ðtk, tlÞ + π NK 2 ijðtk, tlÞ h i I tð Þ:l ð13Þ

Functions Kij1, 2take the following values:

K₁₁1 = K₂₂1 = 0: ð14Þ K121 = 1 4π JoðkiρÞ− pi pe JoðkeρÞ , r!6¼ r!0 0, o:w: 8 > < > : ð15Þ K₂₁1 = 1 4π k 2 iJ2ðkiρÞ−k2eJ2ðkeρÞ ρ! :n! ρ! :n!₀ ρ2 −_4π1½kiJ1ðkiρÞ−keJ1ðkeρÞ n!:n!0 ρ , r ! 6¼ r!₀ 1 8π k 2 e−k 2 i , o:w: 8 > > > > > > > < > > > > > > > : ð16Þ K112 = i 4 kiH 1 ð Þ 1 ðkiρÞ−keH1ð Þ1ðkeρÞ h i , r!6¼ r!0 0, o:w: 8 < : ð17Þ K₁₂2 = i 4 pi pe Hð Þ1 o ðkeρÞ−Hð Þo1ðkiρÞ + K1 12ln 4sin2 tk−tl 2 n o , r!6¼ r!0 − 1 4π ln 4 k2 i −pi pe ln4 k2 e , o:w: 8 > > > > > > > < > > > > > > > : ð18Þ K221¼ −i 4 k2_iHð Þ21ðkiρÞ−k2eH 1 ð Þ 2 ðkeρÞ h _{i ρ}! :n! ρ! :n!₀ ρ2 − kiHð Þ11ðkiρÞ−keHð Þ11ðkeρÞ h _{i n}! :n!₀ ρ −K1 21ln 4 sin2 tk−tl 2 n o , r! 6¼ r!0 k2_i−k2_e 8π ðiπ −2C þ 1Þþ 1 8π k2iln 4 k2_iI2_{ð Þ}_t −k 2 eln 4 k2_eI2_{ð Þ}_t " # , o:w: 8 > > > > > > > > > > > > > > > > > > < > > > > > > > > > > > > > > > > > > : ð19Þ K₂₂2 = i 4 kiH 1 ð Þ 1 ðkiρÞ− pi pe keH1ð Þ1ðkeρÞ ρ! :n! ρ , r ! 6¼ r!₀ 0, o:w: 8 > < > : ð20Þ

where C is the Euler-Mascheroni constant, and J0, J1, J2and

Hoð Þ1, H1ð Þ1, H 1 ð Þ

2 are the complex-argument Bessel and Hankel

functions, respectively.

For the far-ﬁeld computation, we use the asymptotic expansion of the Hankel function for large arguments (39), leading to the following expression:

Ufarðθ,φÞ =

ð

L

e−ikeðxcosθ + ysinθÞ _ik

e n ! :r! U r ! +pi pe ∂U r ! ∂n ( ) dl: ð21Þ whereθ is a scattering angle and φ is an angle of incidence.

An angle-resolved light scattering pattern is calculated as follows: Ithðθ,φÞ = Ufarðθ,φÞ j j2 4ke : ð22Þ Since the area of a side scattering pulse (SSC-A) in FC is detected at approximately 90 to the laser beam at large angles (13–26), a side scattering pattern is calculated as:

SSC−A φð Þ = 10log X θ2 90_½ o₊_φ−13o_;₉₀o₊_{φ + 13}o Inormthðθ,φÞ 8 < : 9 = ;, ð23Þ where Inormthðθ,φÞ = Ithðθ,φÞ min θ Ithðθ,φÞ , φ 2 0½ o,360oÞ: ð24Þ

R

ESULTS Optical Modeling of SS-RBCs

From the image data collected, six representative cell morphological states were chosen from a SS-RBC population obtained from IFC bright- and dark-ﬁeld data, reproducing them with Eq. 1 (see Fig. 1A). RBC shape models ranging from biconcave to elongated, and to typical sickle cell shapes are modeled as a homogenous dielectric particle with losses, whose shape is determined by choosing the corresponding a, b, m,α, β, and γ coefﬁcients as follows: Discocyte (a = 0.99, b = 0.99, m = 2,α = 1, β = 8, and γ = 8), Granular (a = 1.2, b = 2.4, m = 5,α = 3, β = 24, and γ = 4), Holly-leaf (a = 1.2, b = 2.4, m = 7, α = 1.3, β = 1.3, and γ = 2.1), Echinocyte (a = 2.1, b = 2.1, m = 7, α = 3, β = 6, and γ = 6), Elongated (a = 1.07, b = 1.07, m = 2,α = 24, β = 128, and γ = 128), and Sickle (a = 1, b = 1, m = 2, α = 36, β = 128, and γ = 512) (see Fig. 1B). We performed numerical electromagnetic simu-lations based on the light scattering technique in order to investigate the relationship between cell morphology (e.g., cell shape, size, and orientation) and refractive index on the scattering pattern. Numerical simulations were carried out in the wavelength of 632.8 nm. The refractive index spectra of SS-RBCs were used in the range of ni= 1.6× 10−3 and

(5)

1.36≤ nr≤ 1.41 whereas the refractive index of the host

medium, ne was considered to be 1.335. We assume that

dielectric function, εi,e(λ) has a nonzero positive imaginary

part (pffiffiffiffiεi= nr+ ini). Cell rotations were taken in the range of

0and 360with 4interval range for each SS-RBC structure. In Figure 1C, we calculated the angular far-ﬁeld distributions of

the scattered E-polarized light for three sample cell orientations: φ = 0

(face-on incidence), φ = 45(oblique incidence), and φ = 90

(rim-on incidence). The simulation results indicate that the scattering cross-section values have a remarkable change depending on cell types and orientation in the range of side scattering (SSC-A) pattern.

Figure 1. Representative IFC images and simulated optical signature of SS-RBCs upon deoxygenation. (A) Representative images of deoxygenated SS-RBCs acquired by bright-field (BF/Ch-01) and dark-field (SSC/Ch-06) channels in IFC, obtained using the IDEAS® software. BF objects were masked using an investigator-de_{fined custom mask (System 80, Object Tight, Range-Area 350–1,200, and} Range-Aspect Ratio 0.3–1). Scale bar = 7 μm. (B) Two-dimensional cross-sections of SS-RBC shape models obtained from Eq. (1). (C) The simulated light scattering pro_{file of individual SS-RBCs with the use of MBIE method at n}r= 1.38 and E-polarization conditions. [Color

(6)

Finding an Optimal Shape Quantiﬁcation Factor for Deoxygenated SS-RBCs in IFC

In order to determine the optimal settings for normal, intermediate, and sickle cell region boundaries in deoxygen-ated SS-RBCs, a user-deﬁned custom mask including the Range, System, and Object function masks was used to increase the possibility of analyzing one cell per frame in BF channel whereas a default M01 mask contains a small particle in the masked area, resulting in distorted feature values (Fig. 2A). The Range mask was created using the System 80, Range Area 350–1,200 pixels and Range Aspect Ratio 0.3–1 masks to reach a tighterﬁt to the actual contour of the BF region, and thus providing an accurate measurement of cellular information for the individual cell. This range is slightly different as com-pared to the RangeSystem mask at van Beers and coworkers (18) which indicates that Area varies between 350 and 5,000 pixels, and Aspect Ratio is in the range of 0–1. The combina-tion of Shape Ratio (that is, minimum thickness of the masked object divided by the length) and Area features as a Sickle Index shape feature (= Shape Ratio×1,000/Area) enhances the selection of the“true” sickle cell region bound-ary from other delineated SS-RBC phenotypes (Fig. 2B). The well-known highly heterogeneous morphologies of SS-RBC under low oxygen can be distinguished as sickle, elongated, intermediate, and normal cells in the range of 0–20 Sickle

Index whereas the intermediate cells would be located in the sickle cell region boundary at Shape Ratio versus Circularity (that is, the degree of the deviation of mask from a circle) fea-ture (Fig. 2C). The quantitation of percent morphological RBC subpopulations were analyzed with the use of“Shape Ratio ver-sus Circularity” feature obtained from the default M01 mask and the user-deﬁned custom mask as Shape RatioDefault (M01) Maskand Shape RatioCustom Masklegends, and Shape Ratio feature

with a cut-off range (21) as a Shape RatioCut-off level, Custom Mask

legend and “Sickle Index versus Circularity” feature obtained from the user-deﬁned custom mask as a Sickle IndexCustom Mask

legend (Fig. 3A). Candidate sickle cell region boundaries were selected at Shape Ratio < 0.6 and Circularity < 10 on the default M01 and user-defined Range mask, and at Shape Ratio feature with a cut-off of below 0.5 on the combination of default and custom BF mask (21) and at Sickle Index < 7 and Circularity < 7 on the custom BF mask as our proposed shape feature criteria for deoxygenated SS-RBCs over a period of 0, 2, 5, 10, and 20 min. The“Ratio Sickled to Normal” is calculated as the percentage of sickled cells divided by the percentage of normal cells on the defined mask. We observed that the Sickle Index fea-ture enables better identification of sickle cells especially for SS-RBCs after 5 min. of deoxygenation. Moreover, we performed the“Ratio Abnormal to Normal” analysis, which is calculated as the percentage of abnormal cells (thus, sickle + intermediate

Figure 2. Masking and gating strategy to identify normal, intermediate, and sickle cell region boundaries in deoxygenated SS-RBCs. (A) Representative BF images of SS-RBCs overlaid with the respective boundary masks (cyan color). The left column shows the BF image with no mask. The second column shows the default M01 mask, which contains a small particle in the masked area. The third column reveals the user-de_{fined custom mask, Range(System(Object(M01, Ch01, Tight), Ch01, 81), 350–1,200, 0.3–1). This custom BF mask} reaches a tighterfit to the actual contour of the BF region, resulting in an accurate measurement of cellular information for the individual cell. (B) Normal, intermediate, and sickle cell region boundaries are identi_{fied with the use of “Shape Ratio versus Circularity” and “Sickle} Index versus Circularity” features obtained from the user-defined custom BF mask. (C) The representative sample of Shape Ratio and Sickle Index feature values are shown in the upper left side of each BF image. [Color_{figure can be viewed at wileyonlinelibrary.com]}

(7)

cells) divided by the percentage of normal cells obtained from the deﬁned mask, in Figure S1. It is noted that at time zero, when the cells are not exposed to low oxygen, around 10% and 25% of the cells have sickle and abnormal (sickle + intermediate) cell structures obtained from the criterion of our proposed fea-ture, respectively. However, the difference between sickle and abnormal cell region boundaries obtained from Shape Ratio fea-ture with a cut-off range (21) is smaller than of our proposed feature since sickle cells would be seen in the undeﬁned region

boundaries, resulting in lower estimates of the percentage of “true” sickle cells in SCD. When cells were exposed to air after 20 min., this number is regained (not shown). These abnormal cells have the typical morphology of irreversible sickled cells (ISC), generated in the circulation of the patient. For complete-ness, the quantitation of“number of sickled and abnormal cells” and “percentage of sickled and abnormal cells” that are mea-sured for each reaction time of deoxygenation is provided in Table S1.

Figure 3. Morphological/physical response of normal and sickle RBCs of SCD as a function of exposure time to low oxygen acquired by IFC. (A) The quantitation of“percent sickled cells” and “ratio sickled to normal cells” using Shape RatioDefault (M01) Mask, Shape RatioCustom Mask,

Shape RatioCut-off level, Custom Mask, and our proposed shape feature (Sickle IndexCustom Mask) criteria. (B) Morphological/physical parameters

obtained with the criterion of Sickle Index feature on the user-deﬁned Range mask as a function of reaction time of deoxygenation. The dots in each case represent the mean value, whereas the vertical error bars represent the standard deviation of the corresponding cellular information. [Colorﬁgure can be viewed at wileyonlinelibrary.com]

(8)

We examined additional features used in the BF analysis such as Cell Diameter, Area, Orientation, and SSC-A (that is, 10log(Side Scattering Intensity)) in IFC to assess the effect of exposure time during the deoxygenation procedure on the morphological/physical properties obtained with the crite-rion of our proposed shape feature (Fig. 3B). The cell diame-ter, area, and side scattering intensity of sickle cells are significantly higher than that of normal cells after 5 min. of deoxygenation. However, among the first 5 min. of deoxy-genation, there is no statistically significant difference in the cellular features. We observed that sickle cells are character-ized by limited angular width, and are growing in a certain

direction much longer than the diameter of the normal cell while normal cells appear to have originated from multiple domains with homogenous growth directions. Moreover, the comparison plots of the mean morphological/physical features obtained with the criterion of our proposed shape feature and Shape Ratio feature with a cut-off range (21) as a function of reaction time of deoxygenation, and also as a function of normal and sickle RBCs are provided in Figure S2 and S3, respectively. Speciﬁcally, mean values of cellular information for normal and sickle RBCs obtained from our proposed method are higher than of van Beers et al.’s (21) approach. For completeness, we also report the quantitation of the

Figure 4. Angle-resolved light scattering approach in IFC allows for label-free, noninvasive classi_{ﬁcation of morphological RBC} subpopulation region boundaries in the deoxygenated SS-RBCs. (A) The scatter histogram plots of measured SSC-A patterns as a function of cell diameter and orientation for SS-RBCs at t = 0 and 20 min. of deoxygenation obtained from bright- and dark-_{ﬁeld images in} IFC. (B) The difference plots of the simulated SSC-A patterns between (Left) Discocyte and Elongated, and (Right) Discocyte and Sickle shaped RBC models. (C) The comparison between simulated and measured SSC-A patterns for the case of 1.36_{≤ n}r≤ 1.41 and average

radius, d = 8μm in normal and sickle cell region boundaries of SS-RBCs. [Color ﬁgure can be viewed at wileyonlinelibrary.com]

(9)

measured cellular information by means of mean standard deviations (SD) range obtained with the criterion of features used in this study in Table S2.

Optimization of Normal and Sickle Cell Region Boundaries in High-Throughput IFC with Angle-Resolved Light Scattering Approach

We further considered the angular information on the light scattering pattern to improve the accuracy of classifying “true” normal, intermediate, and sickle RBCs since it is uncer-tain tofind a specific range of a shape quantification feature extracting the exact population of sickle cell and other delin-eated RBC phenotypes. Additional features used in the BF analysis were coupled with SSC intensity measurements using the default M06 mask. Figure 4A shows measured SSC-A intensity values as a function of cell diameter and orientation for SS-RBCs at t = 0 and 20 min. of deoxygenation. There is an increase in the cell diameter as the pattern of the SSC-A increased at t = 20 min. of deoxygenation. In order to pro-vide more insight into our understanding of normal and sickle cell region boundaries from measured SSC-A pattern, we calculated the angular far-field distribution obtained from MBIE, and modified the predicted light scattering dataset to consider minimal and maximal changes on the scattering pat-tern as a normalization procedure, and then converted the normalized pattern into SSC-A pattern to make the compari-son with the measured IFC data. First, we performed the light scattering analysis to investigate the effect of cell type and refractive index on the SSC-A pattern for the case of 1.36≤ nr≤ 1.41 and average radius, d = 8 μm in the

Dis-cocyte, Elongated, and Sickle shape RBC models (Fig. 4B). The minimum difference of SSC-A patterns between Discocyte and Elongated shaped RBC models was observed at around nr= 1.39 whereas that of between Discocyte and Sickle shaped

RBC models was found at around nr= 1.36 or 1.39 depending

on the cell orientation. We compared the predicted and mea-sured SSC-A patterns for the case of average radius, d = 8μm in normal and sickle cell region boundaries. Then, we applied the statistical analyses performed using GraphPad Prism 7.0a in order toﬁnd optimal refractive index values for normal or sickle cell region boundaries from SSC-A pattern. For illustra-tive purposes, the data are shown in the mean SD range in Figure 4C. For consistency, all data are treated as nonparamet-ric. Between the data differences are tested with the Kruskal-Wallis or Dunn’s multiple comparisons test. One-sample tests were performed with the Wilcoxon signed-rank test. From the multiple comparison tests, we selected the refractive index value with respect to the criteria of“nonsigniﬁcant” and “mini-mum value of mean rank difference.” The refractive index values were found as nr= 1.39 and nr= 1.37 for normal and

sickle cell region boundaries of SS-RBCs, respectively (Fig. 4C).

D

ISCUSSION

The polymerized state of hemoglobin S (HbS) has been characterized by a double nucleation mechanism at the molec-ular scale (40). The process by which the polymer distorts the

normal RBCs of SCD patients is complex and leads to diverse morphological shapes depending on the configuration of the intracellular HbS polymer affected by physiological conditions, such as temperature, pH, rate of the deoxygenation, and mean corpuscular hemoglobin concentration (41–43). The presence of molecules such as HbF or HbA as the result of hydroxyurea treatment or transfusion will further affect the morphology diversity of the RBC population when exposed to low oxygen. In addition, bloods collected from patients and not exposed to low oxygen contain a fraction of abnormally shaped RBC clas-sified as ISC. Previous studies have clasclas-sified the morphology of SS-RBCs as“discocyte” (cells with a nearly normal appearance), “intermediate cells” (irregularly or abnormal shaped cells, such as echinocyte, holly-leaf, and granular structures),“sickle”, and “elongated” types under deoxygenation states (41,42).

The analysis of a light scattering pattern by an individual RBC has been reported by using various electromagnetic numerical methods such as Mie-series (22), physical-geomet-ric-optics (23), discrete-dipole (10) approximations, anomalous diffraction (24), Lorentz-Mie (25), and discrete sources method (26) theories, ﬁnite-difference time-domain (FDTD) (27), T-matrix (28), multilevel fast multipole algorithm (29) and integral equation methods. In this article, a boundary inte-gral equation method (MBIE), which enables fast and efﬁcient procedures for spherical and nonspherical biological particles with a Nystrom-type meshless discretization, was used to solve biological cells in health and disease conditions. In our case, when the boundary is smooth and has a regular 2π-periodic parametrization, the Nystrom-type discretization (30–33) using trigonometric approximation of integrand functions is more convenient, as it leads to exponential convergence of numerical solution.

Although it has been previously reported that the mean values of Hb content for sickle RBCs are measured higher than that for normal RBCs (44), Byun et al. (45) reported that the refractive index value of HbS is not different from the refrac-tive index of Hb. However, Finch et al. (46) mentioned that the refractive index and absorption coefﬁcient exhibit low value when the incident light is parallel to the deoxygenated HbS molecules in theﬁber compared with normal condition. More-over, the refractive index spectra of the malaria parasite Plas-modium-falciparum RBCs are provided in the range of 1.36–1.41 (47), but the refractive index spectra for individual SS-RBCs from sickle cell patients exposed to a low oxygen ten-sion are unknown. Hence, we used the refractive index spectra of SCD in the range of 1.36–1.41 in this study in order to determine biophysical properties of individual SS-RBCs by means of cell diameter, size, and orientation on SSC-A pattern. Our results showed that the morphological/physical features of the cell are correlated with the cell orientation, shape, and diameter parameters on the predicted and measured scattering patterns.

Conventional FCs is limited to scatter measurements at ﬁxed observation angles in the forward and side scattering detectors which in turn limits the cellular information that could be gained. Additional optical detectors in FCs could be utilized to enhance cellular speciﬁcity, and to achieve a better

(10)

understanding of the scattering cross-section values to pro-vide more quantitative results for the biological questions.

High-throughput IFC allows for simultaneous phenotypic and function studies in the customized creation features and masks, which provides a cell-by-cell prediction of cellular infor-mation and highly sensitive selection of parent populations from all events. A curved object such as a sickle shaped cell is more accurately obtained using Thickness, Length, and Area fea-tures instead of Height and Width feafea-tures, which can be used to separate rectangular shaped objects. It is often not known which morphological shape quantification feature is invaluable to distinguish heterogeneous shapes of SS-RBCs in IFC. Since small particles in the same frame cannot be removed from the default M01 mask, normal and intermediate cells can be classi-fied in the range of sickle cell region boundary, resulting in dis-torted sickle cell percentage of SS-RBCs. Hence, a default M01 mask cannot be optimal for determining accurate cellular fea-tures in the deoxygenated SS-RBCs. Moreover, intermediate cells would be located in the criterion of Shape Ratio < 0.6 and Circularity < 10 as a sickle cell region boundary on the user-defined custom mask, thus leading to overestimating the per-centage of sickle cells. In a recent publication, van Beers et al. (21) reported that Shape Ratio feature with a cut-off range improves the accuracy of the classified cell structures as abnor-mal, norabnor-mal, and undefined region boundaries in SS-RBCs. However, sickle cells can be seen in the undefined region boundaries, resulting in lower estimates of the percentage of sickle cells. Thus, the combination of Thickness Minimum, Length, and Area features as a Sickle Index shape feature created with a user-defined custom mask enables better automated quantification of SS-RBCs partially following deoxygenation and hemoglobin polymerization, resulting in improvement of image analysis in mask generation and feature selection in IFC data.

In this study, we show that a novel label-free, computa-tional sickle cell classification framework based on the bound-ary integral equation method with the use of image-basedflow cytometry data improves the accuracy of classification perfor-mance of“true” normal, intermediate, and sickle RBCs in deox-ygenated SS-RBCs. By using a rapid deoxygenation allows the measure of sickling in time. Our results show the rate of sick-ling when cells are exposed to 1% oxygen (PO2= 8 mmHg),

but such measurements can be expanded to any oxygen tension chosen. Addition offluorescent probes would further enhance the analysis of the population, an approach readily available as IFC allows measurement of fluorescence in addition to the parameters indicated in this study. As an example, labeling with a fluorophore for reticulocytes allows the distinction between the rates of sickling of the youngest cells in the population as compared to the rest. Defining the percentage of sickled cells in time upon deoxygenation is a rapid way to define the propen-sity sickle RBC in the circulation of the patient to sickle under specific conditions. The morphological/biophysical information obtained from our proposed shape feature criteria with angle-resolved light scattering approach in IFC data offers more insights into our understanding of SCD pathophysiology as well as the effect of different treatment options. Small blood samples collected from patients before and after treatment with, for

example, transfusion, anti-sickling agents, or hemoglobin gene therapy, provide the ability monitoring of the effects of treat-ment on sickling of the RBC population in relation to clinical outcome. Hence, our method could offer a unique route to be applied for a tool in the clinic to facilitate robust and objective characterization, and follow-up care of SCD status, and to pro-vide a rapid action for the conditions that would lead to vaso-occlusion in the capillaries of SCD patients causing to a reduced lifespan, organ damage or painful crisis.

A

CKNOWLEDGMENTS

P.G. was supported by the IEEE Antennas and Propaga-tion Society (AP-S) Doctoral Research Grant and Turkish Sci-entiﬁc Technical Research Council (TUBITAK-BIDEB 2211/ A) fellowship. I.O.S. was supported by the Alexander von Humboldt Foundation Fellowship. We thank the Flow Cytometry Core (Gladstone Institutes) for experimental and technical support, and Dr. Vakur Behcet Erturk (Bilkent Uni-versity) for helpful discussions regarding this study.

C

ONFLICT OF

I

NTEREST

The authors declare no competingﬁnancial interests.

L

ITERATURE

C

ITED

1. Pauling L, Itano HA, Singer SJ, Wells IC. Sickle cell anemia, a molecular disease. Science 1949;110:543–548.

2. Kato GJ, Piel FB, Reid CD, Gaston MH, Ohene-Frempong K, Krishnamurti L, Smith WR, Panepinto JA, Weatherall DJ, Costa FF, et al. Sickle cell disease. Nat Rev Dis Primers 2018;4:18010.

3. Kaul DK, Xue H. Rate of deoxygenation and rheologic behavior of blood in sickle cell anemia. Blood 1991;77:1353–1361.

4. Ballas SK, Mohandas N. Sickle red cell microrheology and sickle blood rheology. Microcirculation 2004;11:209–225.

5. Barabino GA, Platt MO, Kaul DK. Sickle cell biomechanics. Annu Rev Biomed Eng 2010;12:345–367.

6. Kuypers FA. Hemoglobin S polymerization and red cell membrane changes. Hematol Oncol Clin North Am 2014;28:155–179.

7. Li X, Du E, Lei H, Tang Y-H, Dao M, Suresh S, Karniadakis GE. Patient-speciﬁc blood rheology in sickle-cell anaemia. Interface Focus 2016;6:20150065.

8. Papageorgiou DP, Abidi SZ, Chang H-Y, Li X, Kato GJ, Karniadakis GE, Suresh S, Dao M. Simultaneous polymerization and adhesion under hypoxia in sickle cell dis-ease. Proc Natl Acad Sci U S A 2018;115:9473–9478.

9. Eaton WA, Hofrichter J. Sickle cell hemoglobin polymerization. In: Anﬁnsen CB, Edsall JT, Richards FM, Eisenberg DS, editors. Advances in Protein Chemistry. Volume 40. New York: Academic Press Inc.; 1990; p. 63–279.

10. Brown M, Wittwer C. Flow cytometry: Principles and clinical applications in hema-tology. Clin Chem 2000;46:1221–1229.

11. Carpenter AE, Jones TR, Lamprecht MR, Clarke C, Kang IH, Friman O, Guertin DA, Chang JH, Lindquist RA, Moffat J, et al. CellProﬁler: Image analysis software for identifying and quantifying cell phenotypes. Genome Biol 2006;7:R100. 12. Schindelin J, Arganda-Carreras I, Frise E, Kaynig V, Longair M, Pietzsch T,

Preibisch S, Rueden C, Saalfeld S, Schmid B, et al. Fiji - An open source platform for biological image analysis. Nat Methods 2012;9:676–682.

13. Hodneland E, Kögel T, Frei D, Gerdes H-H, Lundervold A. CellSegm - A MATLAB toolbox for high-throughput 3D cell segmentation. Source Code Biol Med 2013;8:16. 14. Xu M, Papageorgiou DP, Abidi SZ, Dao M, Zhao H, Karniadakis GE. A deep con-volutional neural network for classiﬁcation of red blood cells in sickle cell anemia. PLoS Comput Biol 2017;13:e1005746.

15. Aghaeepour N, Finak G, Hoos H, Mosmann TR, Brinkman R, Gottardo R, Scheuermann RH. Critical assessment of automatedﬂow cytometry data analysis techniques. Nat Methods 2013;10:228–238.

16. Basiji DA, Ortyn WE, Liang L, Venkatachalam V, Morrissey P. Cellular image anal-ysis and imaging byﬂow cytometry. Clin Lab Med 2007;27:653–670.

17. Samsel L, McCoy JP. Imagingﬂow cytometry for the study of erythroid cell biology and pathology. J Immunol Methods 2015;423:52–59.

18. Fertrin KY, Samsel L, van Beers EJ, Mendelsohn L, Kato GJ, McCoy JP. Sickle cell imagingﬂow cytometry assay (SIFCA). Methods Mol Biol 2016;1389:279–292. 19. Blasi T, Hennig H, Summers HD, Theis FJ, Cerveira J, Patterson JO, Davies D,

Filby A, Carpenter AE, Rees P. Label-free cell cycle analysis for high-throughput imagingﬂow cytometry. Nat Commun 2016;7:10256.

(11)

20. Eulenberg P, Köhler N, Blasi T, Filby A, Carpenter AE, Rees P, Theis FJ, Wolf FA. Reconstructing cell cycle and disease progression using deep learning. Nat Commun 2017;8:463.

21. van Beers EJ, Samsel L, Mendelsohn L, Saiyed R, Fertrin KY, Brantner CA, Daniels MP, Nichols J, McCoy JP, Kato GJ. Imagingﬂow cytometry for automated detection of hypoxia-induced erythrocyte shape change in sickle cell disease. Am J Hematol 2014;89:598–603.

22. Pedersen GD, McCormick NJ, Reynolds LO. Transport calculations for light scatter-ing in blood. Biophys J 1976;16:199–207.

23. Mazeron P, Muller S. Light scattering by ellipsoids in a physical optics approxima-tion. Appl Opt 1996;35:3726–3735.

24. Streekstra GJ, Hoekstra AG, Heethaar RM. Anomalous diffraction by arbitrarily oriented ellipsoids: Applications in ektacytometry. Appl Opt 1994;33:7288–7296. 25. Ulicný J. Lorenz-Mie light scattering in cellular biology. Gen Physiol Biophys 1992;

11:133–151.

26. Eremina E, Eremin Y, Wriedt T. Analysis of light scattering by erythrocyte based on discrete sources method. Opt Commun 2005;244:15–23.

27. Lu JQ, Yang P, Hu X-H. Simulations of light scattering from a biconcave red blood cell using theﬁnite-difference time-domain method. J Biomed Opt 2005;10: 024022.

28. Bi L, Yang P. Modeling of light scattering by biconcave and deformed red blood cells with the invariant imbedding T-matrix method. J Biomed Opt 2013;18: 055001.

29. Ergül Ö, Arslan-Ergül A, Gürel L. Computational study of scattering from healthy and diseased red blood cells. J Biomed Opt 2010;15:045004.

30. Colton D, Kress R. Inverse Acoustic and Electromagnetic Scattering Theory. New York: Springer Science & Business Media, 2012.

31. Wang L, Cox JA, Friedman A. Modal analysis of homogeneous optical waveguides by the boundary integral formulation and the Nyström method. J Opt Soc Am A 1998;15:92–100.

32. Smotrova EI, Tsvirkun V, Gozhyk I, Lafargue C, Ulysse C, Lebental M, Nosich AI. Spectra, thresholds, and modalﬁelds of a kite-shaped microcavity laser. J Opt Soc Am B 2013;30:1732–1742.

33. Tsalamengas JL. Exponentially converging Nyström methods applied to the integral-integrodifferential equations of oblique scattering/hybrid wave propagation in presence of composite dielectric cylinders of arbitrary cross section. IEEE Trans Antennas Propag 2007;55:3239–3250.

34. Gielis J. A generic geometric transformation that uniﬁes a wide range of natural and abstract shapes. Am J Bot 2003;90:333–338.

35. Gielis J, Natalini P, Ricci PE. A note about generalized forms of the Gielis formula. In: Gielis J, Ricci PE, Tavkhelidze I, editors. Modeling in Mathematics. Paris, France: Atlantis Press, 2017; p. 107–116.

36. Müller C. Foundations of the Mathematical Theory of Electromagnetic Waves. Berlin, Germany: Springer-Verlag, 1969.

37. Sukharevsky IO, Shapoval OV, Altintas A, Nosich AI. Validity and limitations of the median-line integral equation technique in the scattering by material strips of sub-wavelength thickness. IEEE Trans Antennas Propag 2014;62:3623–3631. 38. Sukharevsky IO, Nosich AI, Altintas A. Manipulation of backscattering from a

dielectric cylinder of triangular cross-section using the interplay of GO-like ray effects and resonances. IEEE Trans Antennas Propag 2015;63:2162–2169. 39. Watson GN. A Treatise on the Theory of Bessel Functions. Cambridge, England:

Cambridge University Press, 1995;p. 814.

40. Ferrone FA, Hofrichter J, Eaton WA. Kinetics of sickle hemoglobin polymerization. I. Studies using temperature-jump and laser photolysis techniques. J Mol Biol 1985; 183:591–610.

41. Kaul DK, Fabry ME, Windisch P, Baez S, Nagel RL. Erythrocytes in sickle cell ane-mia are heterogeneous in their rheological and hemodynamic characteristics. J Clin Investig 1983;72:22–31.

42. Asakura T, Mattiello JA, Obata K, Asakura K, Reilly MP, Tomassini N, Schwartz E, Ohene-Frempong K. Partially oxygenated sickled cells: Sickle-shaped red cells found in circulating blood of patients with sickle cell disease. Proc Natl Acad Sci U S A 1994;91:12589–12593.

43. Christoph GW, Hofrichter J, Eaton WA. Understanding the shape of sickled red cells. Biophys J 2005;88:1371–1376.

44. Jones CR, Jhonson CS Jr, Penniston JT. Photon correlation spectroscopy of hemo-globin: Diffusion of oxy-HbA and oxy-HbS. Biopolymers 1978;17:1581–1593. 45. Byun H, Hillman TR, Higgins JM, Diez-Silva M, Peng Z, Dao M, Dasari RR,

Suresh S, Park Y. Optical measurement of biomechanical properties of individual erythrocytes from a sickle cell patient. Acta Biomater 2012;8:4130–4138. 46. Finch JT, Perutz MF, Bertles JF, Döbler J. Structure of sickled erythrocytes and of

sickle-cell hemoglobinﬁbers. Proc Natl Acad Sci U S A 1973;70:718–722. 47. Park Y, Diez-Silva M, Popescu G, Lykotraﬁtis G, Choi W, Feld MS, Suresh S.

Refractive index maps and membrane dynamics of human red blood cells parasit-ized by Plasmodium falciparum. Proc Natl Acad Sci U S A 2008;105:13730–13735.