3D holographic flow cytometry measurements of microalgae: strategies for angle recovery in complex rotation patterns

Abstract

Marine ecosystems are in the spotlight, because environmental changes are threatening biodiversity and ecological functions. In this context, microalgae play key ecological roles both in planktonic and benthic ecosystems. Consequently, they are considered indispensable targets for global monitoring programs. However, due to their high spatial and temporal variability and to difficulties of species identification (still relying on microscopy observations), the assessment of roles played by these components of marine ecosystems is demanding. In addition, technologies for a 3D assessment of their complex morphology are scarcely available. Here, we present a comprehensive workflow for retrieving 3D information on microalgae with diverse geometries through holographic microscopy operating in flow-cytometry mode onboard a lab on a chip device. Depending on the rotation patterns of samples, a tailored approach is used to retrieve their rolling angles. We demonstrate the feasibility of measuring 3D data of various microalgae, contingent on the intrinsic optical properties of cells. Specifically, we show that for quasi-transparent and low-scattering microorganisms, the retrieved angles permit quantitative 3D tomographic refractive index (RI) mapping to be achieved, providing full characterization of the alga in terms of its inner structure and outer shape. Moreover, even in the most challenging scenarios, where microalgae exhibit high light absorption or strong scattering, quantitative 3D shape reconstructions of diatoms and dinoflagellates can be at least achieved. Finally, we compare our direct 3D measurements with 2D inferences of 3D properties, obtained using a commercially available microscopy system. The ability to non-invasively obtain 3D information on microalgae marks a fundamental advancement in the field, unlocking a wealth of novel biological insights for characterizing aquatic ecosystems.

---

Introduction

Concerns about environmental changes trigger major challenges, in the modern era, because human activities significantly impact all ecosystems on Earth, making their preservation a top priority. In this context, particular attention is given to the functioning of marine ecosystems, which depends on their structure and biodiversity.¹ Among the biological components that are considered of primary importance, microalgae play a crucial ecological role, contributing to nearly half of the planet's primary productivity, forming the foundation of food webs, and affecting the biogeochemical cycles of key organic compounds.² Human-induced processes greatly affect the health and the community organization of primary producers, making them valuable bioindicators for monitoring the status of aquatic ecosystems.³ Due to their prominent ecological role, microalgae are considered indispensable targets for global monitoring programs⁴ and have traditionally been worldwide included in routine assessments by environmental agencies, as descriptors of the ecological status of aquatic ecosystems.⁵ However, the high spatial and temporal variability of microalgae and the difficulties intrinsic into species identification activities synergistically contribute to rendering reiterated assessments of these ecological components particularly demanding.⁶ Species identification is based mainly on inverted light microscopy,⁷ electron microscopy and DNA sequence analysis.^8,9 These techniques are typically time-consuming, costly, and limiting in efficiency.¹⁰ In addition, the classical manual microscopy is often necessarily integrated with scanning or transmission electron microscopy in order to obtain the diagnostic features of the cell structure and ultrastructure. The increasing need for rapid, cost-effective tools to characterize microalgae communities in a time-efficient way (collection, identification and quantification) is gradually going towards high throughput imaging applied in situ or on collected samples.¹¹ Although the application of machine learning algorithms gives encouraging results in classifying the heterogeneous microalgal community, there is still a lack of technology for the three-dimensional (3D) investigation of the complex morphology of these organisms. As a consequence, some morphological descriptors of microalgae, such as biovolume – a critical parameter for biomass calculation – remain rough estimations. This is a serious gap considering the crucial roles morphology plays in the survival and ecological interactions of microalgae. For instance, morphology influences predator avoidance, nutrient uptake, light absorption, buoyancy and motility.¹² Furthermore, anthropogenic activities can have a deep impact on the morphology of microalgae, potentially making morphology a good bioindicator for pollution.

Digital holography (DH) has emerged as a promising alternative to conventional microscopy. Among the advantages offered, digital holography (DH) is a quantitative, label-free technique with the ability to numerically refocus images after acquisition and to capture the entire complex field modulated by the sample. This makes it very flexible and powerful, as it can image samples in correct focus post-holographic acquisition.¹³ The reconstruction process of digital holograms furnishes at the same time amplitude as well as the so-called quantitative phase maps (QMPs) that represent the quantitative integral measurements of the optical path delay introduced by the sample, when a coherent laser wavefront is transmitted through it. One possibility offered by DH is also the ability to reconstruct the 3D RI profile via tomography. Among others, In-Flow Holographic Tomography (FHT) is a breakthrough that exploits the sample rotation induced in microfluidic channels, so that the illumination is kept fixed while the object rolls and gets probed from multiple directions.¹⁴ However, one challenging issue is the accurate retrieval of rotation angles to be associated with each QPM. Recently, it has been demonstrated that the retrieval of rotation angles can be obtained for spherical or quasi-spherical samples (such as suspended cells) under the assumption of uniform rotation along a single axis¹⁵ (see SI paragraph 2). Applying this approach to microalgae is challenging because of the difficulties associated with their non-spherical morphology, interclass heterogeneity and much more complex rotation patterns experienced by non-spherical objects. The flow speed is generally adjusted to achieve rotation optimal for tomographic reconstruction, meaning that a sufficient portion of the prepared sample reaches the microfluidic channel and most of the specimens rotate in the FOV in a manner useful for tomography (e.g., the rotation was slow enough to acquire at least 30–50 QPMs for each cell). A too high flow rate would lead the samples to simply slide in the channel without rotating, while a too low flow would result in sample loss in the syringe tubes. However, within the optimal range (10–20 nL s⁻¹), the flow rate does not significantly influence rotation; instead, the rotational behaviour is primarily dictated by the sample's shape and its entry point into the channel. Furthermore, among the extremely vast diversity of algae, some exhibit low light scattering and absorbance, while others exhibit non-negligible scattering and light absorbance. Thus, the former can be characterized by FHT, while the latter cannot, as the QPMs cannot be extracted. Nonetheless, even in the cases of algae having strong scattering and/or non-negligible light absorbance, the holographic technology combined with flow cytometry permits retrieving 3D information that is meaningful for ecosystem studies. In fact, in the latter cases, valuable information for classification, morphometric evaluations and biovolume estimation can be recovered by using algorithms such as the Shape From Silhouette (SFS).^16,17 Essentially, DH permits the amplitude and phase of the diffracted wavefront to be reconstructed and we used refocused amplitude maps (AMs) for the SFS recovery. Regardless of the technique exploited to recover the 3D shape (FHT or SFS), the in-flow 3D cytometry requires an automated and reliable approach to recover the rolling angle sequence for each microalga.

Here, we show a novel workflow designed for tracking the angular sequences of non-spherical objects in continuous flow onboard a lab on a chip (LoC) device. The block-scheme is sketched in Fig. 1. We show that a different approach needs to be developed to reconstruct its three-dimensional structure, depending on the specific shape of the microalga, its optical properties and the type of rotation it performs. We demonstrate, as well, that a single holographic optical instrument can provide quantitative 3D measurements of the tested sample, regardless of the microalgae's intrinsic optical properties. Moreover, we discuss the discrepancy between the proposed approach and estimations obtained with an Imaging FlowCytoBot (IFCB). Such commercial imaging flow cytometers are typically used to image and characterize microalgae in 2D and thus have inherent limitations in retrieving 3D morphological features from raw 2D measurements. It is important to note that this work introduces a novel strategy for measuring angle rotations that allows for retrieving FHT as well as SFS. However, a full analysis of the 3D distribution of a species from tomographic data is beyond the scope of this study. Here, we restrict our analysis to the 3D outer shape of the microalgae, focusing on morphometric evaluations and biovolume estimation. Such parameters are common data for techniques like FHT and SFS. Future studies will explore the analysis of inner structures from the 3D RI tomograms of transparent and weakly scattering microalgae, made viable after having established a robust angle retrieval workflow.

Fig. 1 Workflow for the proposed angle retrieval pipeline for non-spherical samples. Depending on the specific rotation pattern followed by the sample during its rotation and by its optical properties, different strategies are used for the recovery of its 3D shape and internal structures.

Methods

Rotation of non-spherical samples

To develop a correct angle retrieval method for non-spherical samples, we previously analysed their rotation pattern in the microfluidic setup in use, described in detail in the SI, paragraph 5. We observed the rotation of all samples in a laminar flow, which is complex and varies among algal species in relation to their shape and entering orientation. Accurate mathematical descriptions have been proposed in ref. 18–20 for a spheroid, which is a rotational solid. According to the pivotal study,¹⁸ the rotation has components along all three Cartesian axes x, y, z, which complicates the 3D recovery. Therefore, it is necessary to quantify and possibly compensate the spurious rotations and then proceed with reconstruction. The possibility to compensate the spurious rotations depends on whether the sample is rotating around its major or minor semiaxis, i.e., which profile it offers with respect to the acquisition plane. Since this differentiation is of major importance in the proposed angle retrieval method, throughout this paper we will refer to the rotation happening around the major or second minor semiaxes as “rotation pattern 1” and “rotation pattern 2”, respectively, as sketched in Fig. 2. The reason why only these two patterns can be assumed lies in geometrical considerations made over the samples, explained in Fig. 2b and c. Decomposing the rotation along its three Cartesian components, the component of our interest for tomography is the one happening in the direction perpendicular to the flow, as the flowing samples align to the rotation axis with one of their semiaxes. Spheroidal samples are characterized by two semiaxes of the same length, therefore, they can only rotate along two directions, since two of them coincide; accordingly, they are completely described by rotation patterns 1 and 2. Ellipsoidal samples are characterized by three different semiaxes and, in principle, could exhibit three distinct rotational patterns. However, we experimentally noticed that the rotation around the smallest semiaxis is an unstable condition. In fact, it is only kept for few frames, before any asymmetry or perturbation breaks these conditions and leads the sample to assume one between rotation pattern 1 and 2. Therefore, in our discussion we will only consider these two patterns.

Fig. 2 Illustration of the geometry of our microfluidic system and schematization of the nomenclature used for the rotation patterns discussed in this work. (a) Velocity profile in the microfluidic channel. (b) and (c) 3D schematization of spheroidal and ellipsoidal samples. (d) Schematization of the terminology introduced in this work. In rotation pattern 1, the sample rotates around its major semiaxis, while in rotation pattern 2, the sample rotates around the second minor direction, resulting in the peculiar sliding–tumbling alternation.

Another consideration is due for approximate ellipsoidal samples, characterized by three different semiaxis values. In this case, the angular velocity changes whether the sample faces the flow with its major or minor cross-section. A frequent condition is that the sample exhibits a sharp rotation when its major cross-section is perpendicular to the flow. In contrast, when its major cross-section is aligned to the flow, the sample experiences a sliding along the microfluidic channel before completing its rotation, because the torque is insufficient to induce a rotation (tumbling–sliding alternation; Fig. 2).

In fact, we observed that the sliding of the sample stops, resulting in a quick 180° flip, when a sufficient cross-section is exposed to the flow. This occurs due to the induced torque, which depends on the sample's minor cross-section and the velocity profile within the microfluidic channel. Furthermore, we frequently observed that the samples did not undergo a complete 360° rotation in the holographic microscope's imaging region, as no perturbation was present to modify the sliding phase. This phenomenon tends to occur more frequently for samples following rotation pattern 2, as the two cross sections have hugely different areas. Lastly, the intrinsic optical properties of microalgae warrant consideration. Indeed, as already said above, certain species characterized by elevated absorption or scattering cannot be treated as pure phase, weakly scattering samples, and thus are more appropriately described by their amplitude images rather than their quantitative phase-contrast maps. However, DH gives access to the entire complex field of the object. For each microalgal species, the most reliable component can be selected and extracted (i.e., QPMs or AMs), by monitoring an opportune metric used for the automatic refocusing, the Tamura coefficient (TC). As explained in ref. 21, the TC exhibits an absolute maximum in the refocusing range in the case of amplitude samples, or an absolute minimum in the case of phase samples, allowing us to discern between them. As the component to be used cannot be known a priori, our proposed angle retrieval pipeline bases its estimation solely on geometrical considerations made over the binarized projection of the QPMs/AMs, without utilizing their internal values. In the following, we name such areas binary maps (BMs), for proceeding with the most appropriate method for 3D reconstruction (i.e., RI tomography or SFS).

Angle retrieval in rotation pattern 1

As mentioned in the previous section, the preliminary step to recover the angular projections is to identify the spurious rotation contributions. Referring to the reference system presented in Fig. 3a, the rotation component of our interest for the 3D reconstruction is the one happening around the x axis. The additional components around y and z are spurious. We define the BM planes as the family of planes perpendicular to the z axis. The motion around the y axis of an angle α happens in the xz plane, perpendicular to the BM planes, as sketched in Fig. 3b. This movement reflects in the BM planes as a variation of the length of the projection of the major semiaxis, a, resulting in a change in the focusing distance. Fig. 3c plots the value of a for the diatom Navicula sp. The maximum value of a (a_max) corresponds to the condition for which the oscillation angle, α, is equal to 0, i.e., the sample lies in the BM planes and a corresponds to the real major semiaxis length of the sample. The minimum values of a (a_min) correspond to the maximum of α and represent the points in which the sample has the biggest tilt with respect to the BM planes. The value of α can be estimated recalling that a_min = a_max cos α. Fig. 3d represents the evolution of α for the diatom Navicula sp., and the maximum value of α is 18.7°. The maximum relative variation of a measured as is equal to 0.053 μm. Since this oscillation occurs in the xy plane, it is not possible to compensate it, thus leading to a certain approximation in the 3D reconstruction. The second spurious rotation component occurs around the z axis. This component determines an oscillation of the BMs with respect to the horizontal direction, a phenomenon that we herein name as apparent kayaking, sketched in Fig. 3e. It is possible to estimate the apparent kayaking angle λ evaluating the angle between a and the x direction. In Fig. 3g, we plot the evolution of λ over time for Navicula sp., while in Fig. 3f, we show the effects of the apparent kayaking over the retrieved BM. The angular sequence is estimated starting from the BM once the compensation of the spurious rotation is performed. As mentioned previously, the compensation of the spurious rotation around the y axis cannot be performed in our experimental setup. Thus, the quantification step is performed to assess this error component, comparing the variation of length of the major semiaxis with respect to the resolution of the setup. Differently, the apparent kayaking angle can be compensated in the BMs evaluating λ for each available frame and back rotating the BM towards the x axis to align it. The rotation is performed clockwise for positive λ values, and anticlockwise for negative λ. Once the two spurious rotation contributions have been accounted for, the angular projections can be estimated starting from the apparent kayaking compensated BMs. To do so, we propose an approach inspired by ref. 22 and 23 presented in Fig. 4 for a sample of Navicula sp. Fig. 4a, top, shows two BMs of Navicula sp. during its rotation after the correction of the apparent kayaking. Starting from the BMs, the dimensions of the sample in the three Cartesian directions are measured. A spheroid of the corresponding dimensions is then synthesized, which models the specimen, sketched in Fig. 4a (bottom); following, a complete and exhaustive set of numerical phase maps BMs_β is generated by projecting the spheroid along the angular sequence [small beta, Greek, macron] = [β₁ … β_N], obtained choosing an angular step of 0.01 degrees and sketched in Fig. 4a, centre. The synthetic step is strongly inferior to the expected real one (linked to the framerate of acquisition). From all the BMs_β, the value of the minor semiaxis D is measured, from which the angular frequency of the rotation will be calculated. Plotting D for all the simulated projections yields a harmonic curve, depicted in blue in Fig. 4b. The corresponding curve measured from the real data is superimposed in red as a comparison. Starting from the curve representing the variation of D, its sinusoidal-fitted corresponding curve y(t) = A sin 2πf̂t is obtained, depicted also in Fig. 4b with a dashed yellow line. In this equation, A is the amplitude of the fitted sinusoid, while f̂ is its oscillation frequency. It must be noted that the oscillation frequency f̂ is double compared to the rotation frequency of the spheroid, as D reaches its maxima twice during a full rotation of the spheroid, respectively, for 180° (half rotation) and 360° (full rotation). In formulas, naming f and T the frequency and the period of the rotation of the spheroid, respectively, we obtain . By observing the amplitude Fourier spectrum of this fitted curve after its mean value, m, is subtracted, two peaks appear, located respectively at ±f̂, as reported in Fig. 4c. Accordingly, the angular pulsation for the spheroid rotation can be expressed as rad s⁻¹, which corresponds to 180f̂ deg s⁻¹. Finally, the sought-after angular step is obtained as , where F_s is the camera framerate. Once the angular step is known, the entire angular sequence can be retrieved by assigning to the first BM the angle 0°. The sinusoidal interpolation of the variation of D along the BMs corresponds to a frequential filtering of the angular speed. Indeed, due to the slight asymmetry of real samples which are not ideal spheroids, the angular speed might not be homogeneous during the period, resulting in an altered sinusoid and therefore spurious frequential component. These contributions are filtered by our method. However, as experimentally demonstrated, since the samples discussed in this case exhibit small asymmetries, this approximation does not alter the final reconstruction.

Fig. 3 Illustration of the spurious rotation components in rotation pattern 1. (a) Visualization of the rotation movement in the 3D space. (b) Sketch of the y-spurious component rotation in the 3D space, defining α as the angle between the major semiaxis and the x axis. (c) Measurement of the length of the major axis for a sample of Navicula sp., where the change due to the y-spurious rotation component can be appreciated. (d) Evaluation of the angle α for a sample of Navicula sp. (e) Definition of the apparent kayaking angle λ in the BM planes. (f) Identification of λ over a real BM. (g) Evaluation of λ for a sample of Navicula sp.

Fig. 4 Proposed angle retrieval method for non-spherical samples for rotation pattern 1. (a) From the BMs, the dimensions of the sample are measured and used to retrieve the spheroid representing it. From the numerical rotation and integration of the spheroid, the synthetic BMs_β are synthetized. (b) Value of the minor axis length D along the real BMs (red line), the BMs_β (blue line) and the interpolated sinusoid (yellow line). (c) Spectral component of the interpolated sinusoid, from which the angular step is computed.

Angle retrieval in rotation pattern 2

Rotation pattern 2 represents the most challenging condition in our experimental setup due to the presence of the tumbling–sliding sequence. As we did for the previous class, the first step is to correct the spurious rotation contributions around the y and z axes. The rotation around the z axis, which determines the apparent kayaking, can be evaluated in the sliding phase, still defining it as the angle between the horizontal direction and the major semiaxis. However, due to the tendency of the samples to align with the flow direction in these rotation cases, this angle tends to be larger than 45°. Therefore, the correct compensation should be made by rotating the sample towards the y direction, to keep the rotation happening around the second minor semiaxis, of a quantity equal to 90 − λ. This definition applies only to the BMs in the sliding phase, as when the sample is tumbling, the major semiaxis of the recorded BMs does not coincide with the specimen's one. Indeed, the definition stated above for the apparent kayaking angle does not correspond to λ for these BMs. The correct compensation should be made in the xz plane, which is perpendicular to the BM planes and therefore not possible in our experimental setup. However, since the tumbling phase is generally fast, this correction can be neglected.

For the angular sequence estimation, a serious concern is also related to the automatic identification of a periodicity in the rotation. The approach we propose for the automatic identification of the rotation frames is summarized in Fig. 5. Fig. 5a represents the application of the method in ref. 15 to a typical spherical sample. The clear alternation between maxima and minima in the correlation coefficient helps to assign the 360° frame and to identify the BMs assigned to 90°, 180° and 270° for all the periods recorded. The correlation coefficient leverages the phase values to find correspondences, as no consideration can be made on the geometry. The angle retrieval pipeline is completely automatic, provided that the full rotation BM can be identified. However, the application of this procedure to a non-spherical sample (in this example, we used Cocconeis sp. in rotation pattern 2) would result in an ambiguity, as reported in Fig. 5b. According to this plot, multiple BMs maximise the similarity metrics with respect to the first orientation, all relative to the sliding phase, thus making impossible to apply the automated procedure typically adopted for spherical samples. According to the method we propose, the starting point is the evaluation of the area of the BMs, as plotted in Fig. 5c. From this graph, the points of minimum area are clearly recognizable, corresponding to the sample being normal to the flow and, therefore, rotating. The regions of maximum area, instead, correspond to the sliding phase. According to these considerations, we assign the 0° angle to the first point minimizing the area of the specimen, as this point is always clearly identifiable, instead of assigning it to the first BM available. All the BMs before the first and successive to the last minimum value of the area correspond to an ambiguous assignation and are therefore discarded. By examining the other points of minimum area, the rotation halves can be identified, i.e., the projections corresponding to a flipping of the sample (180°) or a full rotation (360°).

Fig. 5 Proposed angle retrieval method for non-spherical samples for rotation pattern 1. Application of the method in ref. 15 (a) to a typical spherical object and (b) to a non-spherical sample (specifically a monocyte). (c) Proposed method for the elimination of the redundancies in the sliding phase and the identification of the rotation quarters. (d) Resulting BM sequence for the evaluation of the projection angles.

In order, to identify the BMs corresponding to rotation quarters (90°, 270°, etc.), the ambiguity corresponding to the maximum area projections should be solved. The ensemble of BMs corresponding to the maximum area is identified by fixing a threshold over the area plot equal to the 90% of the maximum area value. Among all the candidates for the rotation quarter projections, only one BM is kept, corresponding to the middle point between two consecutive rotation halves, marked in black in Fig. 5c. All the other candidate BMs are eliminated as they are redundant. The resulting exploited projections are indicated in red in Fig. 5d, where the rotation quarters are now evident. It is worth underlying that this procedure is completely automatic, just like the one for spherical samples in ref. 15. This is remarkable in terms of repeatability and throughput of the tomographic reconstruction of flowing samples. To finally construct the angular sequence, the same method presented in section 2.2 is used for each quarter of rotation. This method is particularly robust because it considers different cases. Indeed, in the case of symmetrical samples, the same angular speed is obtained for each quarter, while for asymmetrical samples (for example, replicating samples consisting of two asymmetric halves), the rotational speed detected for each quarter might be different. As a final consideration, the proposed angle tracking method does not rely on the phase values, as this would result in a method not applicable to samples for which, due to the severe scattering and or absorption, the phase estimation is uncertain.

Results and discussion

Exploiting the angle retrieval approaches presented above, we reconstructed the 3D shape of several microalgae, specifically the diatoms Cocconeis sp., Skeletonema pseudocostatum, Skeletonema marinoi, and Navicula sp. and the dinoflagellates Scrippsiella acuminata, Heterocapsa sp. and Prorocentrum sp., all belonging to the rotation pattern 1 class, and the diatoms Cocconeis sp., Thalassiosira eccentrica and Skeletonema marinoi belonging to the rotation pattern 2 class. Details about the sample's preparation can be found in the SI, paragraph 4. For most diatoms, we could perform tomography due to a sufficient quality of the acquired QPMs, therefore the entire 3D RI profile of the algae was reconstructed. For all the dinoflagellates and Thalassiosira eccentrica, due to excessive absorption and excessive scattering by the frustule, respectively, we computed SFS. The angle retrieval approach presented, indeed, is suitable both for tomography and any other 3D reconstruction method like SFS. The differentiation between the two rotation patterns described in section 2.1 was performed automatically evaluating the area of the sample in the sequence of QPMs. The samples for which a sharp decrease in area was followed by several frames with almost no area variation were assigned to rotation pattern 2, while the samples for which a relatively small area variation frame to frame was calculated were assigned to rotation pattern 1. Additional details regarding the effects of the application of the angle retrieval pipeline for spherical samples for rotation pattern 2 are given in the SI, paragraph 3, along with the effects of a missing compensation of the apparent kayaking.

Tomography for rotation pattern 1 samples

We retrieved the sequences of rolling angles of various microalgae exhibiting weak scattering, and thus we were able to reconstruct successfully the respective tomograms. Such diatoms experienced rotation pattern 1 according to our nomenclature. Specifically, we report the 3D RI tomograms of Navicula sp., a chain of two Skeletonema pseudocostatum samples, a chain of two Skeletonema marinoi cells, a single cell of Cocconeis sp. and a cluster of two cells of Cocconeis sp. samples rotating jointly. To compute the tomograms, the high order total variation (HOTV) of order = 2 was used as described in ref. 24. As a result, we present five examples of 3D RI tomograms of five microalgae species exhibiting rotation pattern 1 (Fig. 6). The central xz slice of the tomograms is displayed, together with the isolevel visualization obtained with three thresholds set as the 10%, 20% and 50% of the maximum value of RI, respectively. The peculiar internal structure of each diatom is retrieved (Fig. 6). For Navicula sp., the two distinctively shaped chloroplasts are visible along the borders of the diatom, together with the low-RI central nucleus. Silicate structures such as the raphe cannot be visualized in the transmission RI measurement. The peculiar morphology of Cocconeis sp. is reconstructed, characterized by a thin elliptical shape with a high internal RI region associated with the chloroplast. The tomograms in Fig. 6 have been calculated through the correct angle sequence estimation, which compensate the compensation of the apparent kayaking, and in fact the 3D tomograms have the expected morphology, while in the case of missing the correction for kayaking, the morphology appears to be clearly aberrated (see Fig. S3a–d). In the SI paragraph 3, we discuss the severe effects of the lack of compensation of the apparent kayaking. Fig. S6 and Video S1 summarize the rotation behavior of some of the microalgae shown in this section, while Fig. S7 shows the BMs associated with these samples.

Fig. 6 Typical tomographic reconstruction of diatoms following rotation pattern 1. Up, isolevel visualization; down, central slice of the tomogram. (a) Navicula sp., (b) Skeletonema pseudocostatum, (c), Skeletonema marinoi, (d) Cocconeis sp. cluster, and (e) Cocconeis sp. (Video S1).

Tomography for rotation pattern 2 samples

Fig. 7 illustrates the tomograms computed for diatoms belonging to the rotation pattern 2 class. In this case, the diatoms reconstructed are again three cells of Cocconeis sp. and one chain of Skeletonema marinoi. However, due to a different orientation with which they entered the microfluidic channel, the rotation pattern exhibited is different from the case described above. This means that, unless precise engineering of the fluid in the microfluidic channel is made, the orientation of the samples entering the channel is crucial to determine their rotation behavior, and therefore it is not possible to use the same angle retrieval approach for each sample without a previous assessment. The detrimental effects of a non-correctly retrieved angular sequence are highlighted in the SI in Fig. S3e–h. Due to their thin disk-shaped structure, the cells of Cocconeis sp. typically rotate with rotation pattern 2, exhibiting the characteristic sliding–tumbling behavior. Their motion is often turbulent and non-stable among the consecutive 360° rotations in the field of view (FOV). The chains of Skeletonema marinoi tend to show the tumbling–sliding motion pattern too in most of the cases, as we assessed experimentally, often performing one single 180° flip in the available FOV. However, differently from Cocconeis sp., they tend to tumble in a less abrupt way, while remaining in the sliding phase for a bigger number of frames. Due to the above considerations, regarding the possibility to correct the spurious rotation components, in this case, only tomograms of samples not presenting the apparent kayaking rotation were computed. Generally, the approximation of the tomograms obtained is worse compared to the rotation pattern 1 case. However, using the proposed correct angle retrieval approach, it is still possible to compute tomograms that retrieve the actual sample's shape and where the internal structures and the mapped expected RI range are in agreement with expectations. We have experimentally isolated a dead Cocconeis sp. cell and reconstruct its tomogram, in which the cytoplasm extrusion was visible (Fig. 7a). Only the external membrane was left, and the chloroplasts collapsed together forming a round structure. These results further highlight the capability of the method to discern in diatoms signs of stress occurring, e.g., as a consequence of non-ideal environmental conditions. Fig. S6 and Video S2 summarize the rotation behavior of some of the samples shown in this section, while Fig. S7 shows the BMs associated with these samples.

Fig. 7 Tomographic reconstruction of diatoms following rotation pattern 2. Right, isolevel visualization; left, central slice of the tomogram. (a), (c) and (d) Cocconeis sp. (b) Skeletonema marinoi (Video S2).

Shape From Silhouette for rotation pattern 1 and 2 samples

We performed the SFS reconstruction of four additional species exhibiting rotation pattern 1 or 2 and that were not suitable for FHT due to their optical properties. Specifically, we reconstructed Scrippsiella acuminata, Heterocapsa sp., Prorocentrum sp. and Thalassiosira eccentrica, shown in Fig. 8a. From all the available 3D reconstructions of microalgae shown in this work, we extracted morphological parameters of interest (Table 1). The parameters were measured directly from the computed SFS (Fig. 8a), or from the outer shells of the tomographic reconstructions, when available (Fig. 8b). The intrinsic 3D nature of our reconstructions allows access to the entire profile of the samples, accounting also for their anisotropy (Fig. 8). Accordingly, the measurements of the parameters of interest and in particular of the biovolume do not rely on approximations of the 3D shape of the sample, differently from the state-of-the-art imaging tool for microalgae, namely the Imaging FlowCytobot (IFCB)²⁵ (see SI paragraph 1). Video S3 summarizes the rotation behavior of the samples shown in this section.

Fig. 8 (a) SFS reconstructions of scattering or absorbing samples. (b) Outer shells of tomographic reconstructions for a set of samples discussed in the previous sections (Video S3).

Table 1 Parameters measured from the computed 3D reconstructions. When available, the parameters were measured over the outer shells of the tomographic reconstruction; otherwise, they were extracted from the SFS

Sample	Volume [μm^3]	Eccentricity	Equivalent diameter [μm]
Cocconeis sp.	97.61	0.66	6.71
Navicula sp.	443.91	0.97	13.41
Skeletonema marinoi	135.74	0.92	7.64
Skeletonema pseudocostatum	90.70	0.91	9.80
Scrippsiella acuminata	3079.02	0.61	18.42
Heterocapsa sp.	354.16	0.69	8.66
Prorocentrum sp.	10–577.40	0.38	33.14
Thalassiosira eccentrica	14–058.46	0.51	31.76

---

Conclusions and discussion

This work presented a workflow for the 3D reconstruction of non-spherical samples via in-flow holography onboard a LoC, embedding a novel angle retrieval approach. The workflow considers the optical properties of the samples, reconstructing the 3D shape, starting from amplitude or phase-contrast measurements. The angle retrieval approach correctly considers the movement of the samples in the microfluidic environment and allows the projection sequence to be effectively recovered under different rotation conditions, thus allowing the 3D reconstruction that instead is missed by commercial systems working in 2D mode. As a proof of concept, this method was employed for the 3D reconstruction of microalgae cells. Microalgae are marine microorganisms of prime interest due to their critical ecological roles, but they are also a challenging group for in-flow 3D imaging because of their wide range of morphologies, which induce complex rotational patterns. Although flowing cells of microalgae have been imaged in 3D before,^14,26 we are the first to show that 3D quantitative measurements can be obtained for a wide variety of microalgae presenting different shapes and exhibiting various uncontrolled rotation patterns (summarized in Fig. S6 and Videos S1–S3). Achieving 3D imaging of microalgae is an important step forward in their monitoring. An improved visualization of the 3D complexity of the cells can help the taxonomist to better identify, for example, the species belonging to the same genus. Better characterization of the 3D structure of the cells will also provide more accurate biovolume measurements and, consequently, biomass calculations.^27,28 For now, the highest degree of accuracy (taxonomically and geometrically) is achieved using light microscopy. However, if light microscopy is used, it is very difficult to measure the third linear dimension of microalgae cells. Consequently, biovolume is often calculated from estimations of the third dimension based on the width of the cell and by using species-specific conversion factors with theoretical shape assignations. Manually measuring microalgae cells using a light microscope is time consuming and in certain cases, the biovolume calculation remains a rough estimation. For these reasons, a technology that can provide fast, high throughput measurements of the 3D morphology of microalgae is a significant breakthrough. The creation of extensive databases of high-resolution 3D rotatable images of species could serve as a standard reference source for machine learning-based automatic taxonomic identifications, biovolume calculations and biomass estimations. Here, we demonstrated that a single holographic system can be used to retrieve the 3D morphology of microalgae independent of their shape, size, and intrinsic optical properties. The proposed workflow extends the variety of objects that can be measured in 3D in a LoC. It could be used in the future as an extensive tool to characterize in 3D the structure of microalgae by extracting more accurate morphological parameters such as more precise biovolume measurements. This represents a significant step forward in comparison with the actually existing methods for microalgae measurements including the state-of-the-art 2D imaging systems.

Author contributions

Investigation: G. G., J. B., P. F., and L. M.; supervision: V. B.; methodology: F. B., P. F., and V. B.; resources: L. R., A. S., M. C., and V. Z.; software: F. B. and P. M.; validation: G. G., E. H., and L. R.; visualization: D. P., G. G., and E. H.; formal analysis: F. B., L. M., and P. M.; funding acquisition: V. B. and L. R.; writing – original draft: F. B. and E. H.; writing – review & editing: V. B., L. M., P. M., and P. F.

Conflicts of interest

There are no conflicts to declare.

Data availability

Supplementary information available: SI text contains the following supplementary sections: S1. IFCB characterization of microalgae. S2. Angle tracking for spherical particles rotating uniformly. S3. Quality enhancement of tomograms. S4. Isolation and culture conditions of microalgae. S5. Optical setup. S6. AM/QPM retrieval. S7. Video. S8. BM corresponding to the reconstructed samples. See DOI: https://doi.org/10.1039/D5LC00559K.

Data for this article, including 3D tomograms and SFS, are available at http://Figshare.com at https://doi.org/10.6084/m9.figshare.29256932.v1 and https://doi.org/10.6084/m9.figshare.29256977.v1.

Acknowledgements

This work was supported by Progetto PRIN 2022 PNRR P2022MA95R – microFluidic platfOrm fOr The toMographic phAse micRoscopy of phytoplanKton assemblageS (FOOTMARKS)-CUP B53D23023980001 and funded by European Union-NextGenerationEU. Anna Italiano and Emanuele Somma aided the culture and taxonomical identification of Cocconeis sp. diatoms, respectively. L. R. acknowledges support by the National Biodiversity Future Centre (NBFC) Program, Italian Ministry of University and Research, PNRR, Missione 4 Componente 2 Investimento 1.4 (Project: CN00000033). The authors thank Wiebe Kooistra and Urban Tillman for providing some of their cultures.

Notes and references

1. C. B. Field, et al. , Science, 1998, 281, 237–240.
2. F. Azam, et al. , Mar. Ecol.: Prog. Ser., 1983, 10, 257–263.
3. EU WFD, European Commission Directive 2000/60/EC of the European Parliament and of the Council of 23 October 2000 establishing a framework for Community action in the field of water policy, Off. J. Eur. Communities L 327, Brussels, 2000, pp. 1–72.
4. S. D. Batten, et al. , Front. Mar. Sci., 2019, 6, 395.
5. J. F. Tweddle, et al. , Mar. Policy, 2018, 97, 1–9.
6. L. Roselli, et al. , Mar. Environ. Res., 2022, 183, 105791.
7. H. Utermöhl, SIL Commun., 1958, 9, 1–38.
8. M. E. Sonmez, et al. , Micron, 2023, 172, 103506.
9. M. W. Fawley, et al. , J. Appl. Phycol., 2020, 32, 2699–2709.
10. I. Agnarsson, et al. , Syst. Biol., 2007, 56, 531–539.
11. J. J. Pierella Karlusich and E. Pelletier, et al. , Mol. Ecol. Resour., 2022, 23, 16–40.
12. V. Sonnet, et al. , Limnol. Oceanogr., 2022, 67, 1850–1864.
13. L. Miccio, et al., Optical Methods for Inspection, Characterization, and Imaging of Biomaterials VI, 2023, p. 18.
14. F. Merola, et al. , Light:Sci. Appl., 2017, 6, e16241.
15. D. Pirone, et al. , Appl. Opt., 2020, 60, A277.
16. A. Laurentini, IEEE Trans. Pattern Anal. Mach. Intell., 1994, 16, 150–162.
17. S. Prince, Computer Vision: models, learning and inference, Cambridge University Press, 2012.
18. G. B. Jeffery, Proc. R. Soc. London, Ser. A, 1922, 102, 161–179.
19. T. Tohme, et al. , Micromachines, 2021, 12, 277.
20. G. Lauricella, et al. , Phys. Fluids, 2022, 34, 082021.
21. J. Dohet-Eraly, et al. , Opt. Lett., 2014, 39, 1109.
22. F. Merola, et al. , Lab Chip, 2013, 13, 4512.
23. Z. Wang, et al. , Light:Sci. Appl., 2021, 10, 187.
24. J. Běhal, et al. , Cell, 2022, 11, 2591,  DOI:10.3390/cells11162591.
25. H. M. Sosik and R. J. Olson, Limnol. Oceanogr.:Methods, 2007, 5, 204–216.
26. K. Umemura, et al. , MethodsX, 2020, 7, 100889.
27. L. Roselli, et al. , J. Plankton Res., 2012, 35, 135–145.
28. L. Roselli, et al. , J. Microsc., 2015, 258, 200–211.

---