Accurate evaluation of size and refractive index for spherical objects in quantitative phase imaging

Measuring the average refractive index (RI) of spherical objects, such as suspended cells, in quantitative phase imaging (QPI) requires a decoupling of RI and size from the QPI data. This has been commonly achieved by determining the object's radius with geometrical approaches, neglecting light-scattering. Here, we present a novel QPI fitting algorithm that reliably uncouples the RI using Mie theory and a semi-analytical, corrected Rytov approach. We assess the range of validity of this algorithm in silico and experimentally investigate various objects (oil and protein droplets, microgel beads, cells) and noise conditions. In addition, we provide important practical cues for future studies in cell biology.


Introduction
Quantitative phase imaging (QPI) is a collective term for interferometric techniques that quantify the phase retardation of otherwise transparent objects. The establishment of QPI as a swift tool in single-cell analysis has given rise to a broad spectrum of applications in cell biology. For instance, QPI has been used to quantify cell dry mass [1,2], cell dynamics [3,4,5,6], bacterial infection [7], parasitic infection [8], or cellular differentiation [9]. An important quantity that can be measured with QPI is the cellular refractive index (RI). Besides characterizing cells based on RI and the associated local protein concentration [1,2,10], knowing the RI is also essential for related applications such as the optical stretcher to quantify optical forces [11,12,13], or for Brillouin microscopy to compute cell elasticity [14,15].
The accurate determination of the cellular RI with QPI is difficult, because the optical path difference (OPD) measured in QPI needs to be separated into integral RI and cell thickness. Thus, the integral RI can only be computed if the cell thickness is measured, which has been achieved using scanning probe microscopy [16,17], confocal microscopy [18], deliberate variation of the OPD [3,19], or spatial confinement [20,21]. However, RI computation with such OPD approaches is based purely on geometric considerations, neglecting light-scattering. A more general approach to this problem is optical diffraction tomography (ODT) [22,23], which yields a complete 3D map of the intracellular RI. However, ODT depends on the acquisition of multiple phase images, which is accompanied by an elaborate experimental effort. Moreover, state-of-the-art ODT techniques approximate light propagation with the Rytov approximation [24,25], which can cause an underestimation of the intracellular RI for large RI gradients [26,27].
However, if a cell is spherical and sufficiently homogeneous in RI, then its average RI can be inferred from a single phase image. Several studies have exploited cell sphericity by computing the cell radius from the cell area visible in the phase image using OPD approaches [28,29,30,31]. However, a rigorous treatment of this problem by exactly modeling light propagation with Mie theory has not been presented so far due to the large computational effort required.
Here, we present a QPI phase fitting algorithm that reliably retrieves RI and size of spherical objects. We address the computational challenge with efficient implementations of the Mie-and Rytov-scattered fields. Additionally, we derive a correction factor for the Rytov approximation, result-ing in accuracies similar to the results of Mie theory. Using in silico simulations, we compare our scattering approaches to OPD approaches for homogeneous spheres with RI values from 1.334 to 1.440. In the experimental part, we demonstrate 2D QPI phase fitting for several test targets and for various noise conditions. Our approach allows for a faster and better interpretation of QPI data for spherical objects and gives valuable insights into single-cell light-scattering.

Modeling light-scattering by a sphere
The OPD approach resembles a highly simplified model for light-scattering, treating the propagation of light as a line integral through the sphere's constant RI. To estimate the average RI of a cell from a single quantitative phase image with the OPD approach, the cell radius must be determined. The radius can be determined with a series of 1D phase profile fits [28], or by fitting a circle to a contour along the cell edge (found using e.g. the Canny edge detection algorithm) [32]. The former approach is equivalent to a direct fit of the projected RI of a sphere to the measured phase image (referred to as "OPD projection approach" in the remainder of the manuscript), while the latter approach only uses the contour data to determine the radius (referred to as "OPD edge-detection approach"). Note that the OPD edge-detection approach is somewhat erratic, because the contour found with an edge detection algorithm depends on image resolution. Both OPD approaches do not take into account diffraction and thus only describe optically thin objects correctly.
The Rytov approximation gives a better estimation of light propagation [33]. It is state-of-the-art in ODT applications, because it provides a simple linear model for diffraction and because it is sufficiently accurate for many cell types [26]. In general, to compute the Rytov approximation, a full 3D model of the specimen is required and its Fourier transform must be computed according to the Fourier diffraction theorem [22]. In case of a sphere, however, we reduced this problem to two dimensions by using the analytical solution of the Fourier transform of a sphere (see appendix A.1 for a detailed description). Our semianalytical approach enables a fast and accurate computation of the Rytov approximation for homogeneous spheres.
Mie theory yields an exact solution for the scattered field. Here, we used the software BHFIELD 1 [34] and refocused the obtained fields to the plane at the sphere center using the Python package nrefocus 2 . Generating the complete 1 versioned October 5th 2012, available at https://seafile.zfn.
uni-bremen.de/f/2d6fc70841/?dl=1 (accessed August 2017) 2 version 0.1.5, available at https://pypi.python.org/pypi/ nrefocus scattered field of a sphere using Mie theory is computationally expensive. Hence, we approximated the full 2D field by radially averaging two 1D fields, perpendicular and parallel to the polarization axis. As discussed in appendix A.2, the error made by this approximation is negligible for the purpose of the present study.
A comparison of these scattering models with Mie theory as a reference is shown in figure 1 for an exemplary homogeneous sphere. While the error made by the Rytov approximation is small, the error of the OPD projection approach scores high values, mostly due to the discontinuity in the gradient of the simulated phase image. To allow a comparison to the OPD edge-detection approach as well, the lower right quadrants of figures 1a and 1b show the OPD projection and resulting error based on RI and radius as obtained with the OPD edge-detection approach applied to the Mie data. The OPD edge-detection approach exhibits a large error which can be explained by an underestimation of the radius due to the fact that the contour found by the edge detection algorithm does not coincide with the lateral perimeter of the sphere. This example illustrates that the determination of the RI and the radius of spherical objects should be addressed with models that take into account diffraction (Mie theory, Rytov approximation), rather than plain OPD approaches (OPD projection, OPD edgedetection). Figure 1: Employed scattering models. a) Simulated quantitative phase images of a sphere with a radius of 2.5 µm and a refractive index (RI) of 1.36 embedded in a medium with an RI of 1.333 at a wavelength of 550 nm. Each quadrant shows a phase image for one of the scattering models described in the text. The phase image of the optical path difference (OPD) edge-detection approach (lower right quadrant) was modeled with the OPD projection approach using the parameters obtained from the detected edge (edge shown in white in the upper left quadrant). b) Difference between the scattering models and the Mie model (upper left quadrant in (a)) in % of the maximum OPD.

Fitting scattering models in-silico
The determination of the RI and the radius of a homogeneous sphere from a single quantitative phase image using one of the scattering models introduced above requires a phase image fitting algorithm. Previous work on this topic employed a set of 1D OPD fits to the 2D phase image [28]. For Mie theory or the Rytov approximation, this 1D approach is inefficient, either because fitting requires an enormous number of field computations (Mie) or because each iteration of the 1D fit requires the computation of a 2D field (Rytov, our implementation). To take advantage of the performance-enhancing procedures described in the previous section (Mie and Rytov), an actual 2D phase image fit is, in fact, necessary.
Here, we propose an image fitting algorithm that we specifically designed for homogeneous spheres in QPI. Our algorithm iteratively fits RI, radius, and lateral position of a sphere as well as phase offset to a quantitative phase image. In contrast to commonly used fitting techniques that call the modeling function for every parameter set during fitting, our algorithm interpolates the modeled phase image in a predefined interval and thus requires less calls to the modeling function. This approach, which is described in detail in appendix A.3, reliably retrieves RI, radius, lateral position, and phase offset for homogeneous, spherical phase objects 3 .

Error of the scattering models
To assess the accuracy of our 2D fitting approach, we performed a series of Mie simulations and retrieved RI and radius using the light-scattering models described above. Figure 2 shows the relative errors in RI and radius made using the OPD edge-detection approach and the 2D phase fitting approach (OPD projection, Rytov approximation). The relative errors for RI (n) and radius (r) were computed using ∆n = n measured − n exact n exact , The black cross and the red square label the positions of the exemplary spheres used in figure 1 and 3. The white asterisk labels a region where the phase image resolution cannot resolve steep phase gradients (>2π phase jump between two pixels) [31]. As a result, all approaches failed to determine the correct RI and radius in this region, which we also observed with the otherwise error-free 2D Mie fit shown in figure 8 in the appendix. As discussed above, the OPD edge-detection approach underestimates the radius of Figure 2: Quantification of modeling errors. Left column: error for refractive index (RI) and right column: error for radius when using (a,b) the optical path difference (OPD) edge-detection approach or the proposed 2D fitting algorithm with either (c,d) the OPD projection approach or (e,f ) the Rytov approximation. The ground truth data were generated from 9600 Mie simulations with an RI of the medium of 1.333 (water), and a grid size of 128×128 px. The lateral simulation size was 4r for r < 5λ (vacuum wavelength λ, sphere radius r) to capture diffraction effects and 3r for r > 10λ to enable Mie simulations for larger spheres, with a linear transition in-between. Initial values for the 2D fit were obtained using the OPD edgedetection approach. The cross and the square depict the spheres shown in figures 1 and 3. The asterisk marks a region where 2D fitting fails due to large phase gradients (see text). the sphere and thus overestimates the RI, which is clearly visible in figures 2a (green) and 2b (brown). The 2D fit with the OPD projection (figs. 2c and 2d) exhibited a noisy error signal, which can be attributed to the discontinuity in the gradient of the modeled phase. In comparison, the error of the 2D Rytov fit shown in figures 2e and 2f is smooth and has lower values than the OPD projection in both RI and radius for RI values below 1.39. The data suggest that the comparatively faster Rytov approximation can be preferred over Mie theory, if an error in the radius below 2% and an error in the RI below 0.1% is acceptable and if the imaged object has a radius above three wavelengths (3λ) and an RI below 1.36. . c) For different initial parameters of RI and radius (dots at the border) and for several noise levels (color code), the phase images were fitted with our new fitting algorithm using the Rytov approximation. For each case, the fitting process, which took on average seven iterations, is shown as a dotted line. The red square indicates the exact value. The green triangle indicates the value obtained with the OPD edge-detection approach. d) Radial residual profiles for the 5 % noise case (profile indicated in (b)) for selected iterations.

Convergence with noise
In practice, quantitative phase images are subject to phase noise which can be described as a background pat-tern that is varying over the range of several pixels. Here, we modeled phase noise using 2D Perlin noise as implemented in the Python package noise 4 . Figure 3 illustrates the convergence of the proposed fitting algorithm for various starting parameters and for realistic noise conditions ranging from 0% to 5% standard deviation measured relative to the maximum OPD. The data show that the algorithm converges to the same value after an average of seven iterations, independent of noise or initial conditions. Unexpectedly, the RI and radius obtained with the OPD edge-detection approach, labeled with a green triangle, is closer to the correct value. This is a coincidental result that can be explained by the fact that the OPD edge-detection approach is resolution-dependent, i.e. the resolution chosen lead to better results for the OPD edge-detection approach (see fig. 2a,b). Note that the initial conditions shown in figure 3c are rather extreme. In practice, we determined the initial guess for radius, RI, and position of the sphere with the OPD edge-detection approach. Thus, experimental phase noise does not affect the convergence of the proposed 2D phase image fitting algorithm.

Systematic correction for the Rytov approximation
The errors in RI and radius made by the Rytov approximation shown in figures 2e and 2f become large for objects with RI values above 1.36. Thus, 2D fitting with the Rytov approximation would yield poor accuracies for many cell types with RI values reaching up to 1.40 and above. In terms of efficiency, this would be a drawback, because high accuracies could only be achieved by falling back to the computationally more expensive Mie theory. To address this issue, we derived a systematic correction for the Rytov approximation, making it possible to access large RI values with high accuracy. We found that our implementation of the Rytov approximation (see appendix A.1) exhibits a systematic error that is independent of the image resolution used. The major fraction of this systematic error is dependent only on the RI for object radii of three wavelengths and above. This allowed us to derive a systematic correction for the Rytov approximation (Rytov-SC), considerably reducing the error made. We derived the correction formulas for RI (n Ryt-SC ) and radius (r Ryt-SC ) by fitting a polynomial function to the Rytov error, yielding where n med is the RI of the medium and n Ryt is the RI obtained using our 2D fitting algorithm and our implemen-4 version 1.2.1, available at https://pypi.python.org/pypi/noise/ tation of the Rytov approximation. Note that we chose the variable x such that the correction is independent of the RI of the medium n med (see figure 9 in the appendix). As we used the error maps shown in figure 2 to derive equations 3 and 4, the systematic correction is valid for spheres with RI n sph and radius r sph of at least 3λ ≤ r sph ≤ 20λ (7) with the vacuum wavelength λ of the light used. This systematic correction extends the applicability of the Rytov approximation in 2D fitting to cells with high RI values (up to 1.44 and above), yielding high theoretical accuracies for RI (<0.1%) and radius (<1%).

Results
To compare the five approaches for the retrieval of RI and radius (OPD edge-detection, OPD projection fit, Rytov fit, Rytov-SC fit, Mie fit), we applied them to quantitative phase data of lipid droplets, microgel beads, and cells.

2D Mie fits to doplets, beads, and cells
Mie simulations in combination with our 2D fitting algorithm served as a benchmark for the other four approaches. Figure 4 shows a representative set of phase images, the corresponding Mie fits, and the fit residuals, which are discussed in the following. Liquid droplets are ideal test samples for the investigation of scattering by spheres, because they are homogeneous and assume a spherical shape due to the surface tension at the droplet-medium interface. Figure 4a shows a quantitative phase image of a silicone oil droplet embedded in phosphate buffered saline (PBS). The oil droplet was produced by vortexing a two phase solution made of 1 mL of silicon oil (Sigma-Aldrich, 10cSt) and 10 mL of de-ionized water containing 2% w/v poly(ethylene glycol) monooleate (Sigma-Aldrich). The image was recorded with quadriwave lateral shearing interferometry (QLSI) [35] using a commercial QPI camera (SID4Bio, Phasics S.A.) attached to an inverted microscope (AxioObserver Z1, Zeiss) with a 40× objective (NA 0.65, 421060-9900, Zeiss). The illumination wavelength was confined to an average of 647 nm using a bandpass filter (F37-647, 647/57, AHF analysentechnik). We separately measured the RI of the silicone oil using an Abbe refractometer (2WAJ, Arcarda), yielding a value of 1.402 that matches the value of 1.405 from the 2D Mie fit shown in figure 4b. The resulting relative error of 0.04% is small which is reflected by residuals below 5% of the maximum OPD shown in figure 4c. Figure 4d shows the quantitative phase image of a protein droplet that was recorded with the Experimental quantitative phase images (1st column), corresponding 2D fits with Mie theory (2nd column), and resulting fit residuals (3rd column) for liquid droplets (a-f), microgel beads (g-m) and HL60 cells (n-t). The scale bar is 5 µm. same setup as above, except for the bandpass filter. Here we used an average imaging wavelength of 550 nm for the RI analysis. The droplet consisted of the RNA/DNA-binding protein FUS and was embedded in a protein buffer that is described in detail in reference [36]. FUS has been linked to neurodegenerative diseases such as amyotrophic lateral sclerosis (ALS) (see e.g. [37,38]) which is correlated to a liquid to solid phase transition of the protein with time [36]. Here, we obtained an average RI of 1.435 for a liquid FUS protein droplet (fig. 4e) using a protein buffer RI of 1.3465, measured with the Abbe refractometer named above. The fit residuals (fig 4f) are largely below 5%, indicating good agreement of theory and experiment. Note that the liquid droplets shown had comparatively high RI values and were imaged at low resolution, conditions that both were handled well by the proposed 2D fitting algorithm.
Microgel beads offer an additional, convenient possibility to test scattering by a sphere. Here, we used polyacrylamide (PAA) beads that were produced in a flow focusing microfluidic system described elsewhere [39,40]. We imaged two types of PAA beads with different density (named A and B for simplicity) using two different imaging setups. The microgel bead of type A was imaged with the QLSI setup used for the FUS droplet above with the addition of a telescope (f 1 =−15 mm, LD2020-A and f 2 =75 mm, LBF254-075-A, Thorlabs, Germany) to increase the magnification, and thus the sampling of the recorded phase image, by a factor of five. The microgel bead of type B was imaged with digital holographic microscopy (DHM). A detailed description of the DHM setup used can be found in reference [32]. A comparison of the resulting quantitative phase images in figures 4g and 4k shows that the DHM data have a higher background noise than the QLSI data (3.4% vs. 0.7%). The convergence of the 2D fitting algorithm is not affected by such noise magnitudes (see also fig. 3) as can be seen in the fit residuals that are either below 5% ( fig. 4j) or reproduce the background phase noise magnitudes ( fig. 4m).
HL60 cells in suspension have a spherical shape. Even though the internal RI distribution of a cell is generally nonuniform, a 2D Mie fit can be used to estimate its average RI. Figures 4n and 4r show quantitative phase images of two HL60 cells, recorded with the DHM and QLSI setups described above; To match the wavelength of the laser used for DHM imaging (HeNe 633 nm, HNL050L-EC, Thorlabs), the illumination light of the QLSI setup was confined by a bandpass filter (BP 640/30, 488050-8001, Zeiss) in addition to filter F37-647 named above. Compared to liquid droplets or beads, the fit residuals shown in figures 4q and 4f exhibit heterogeneous structures that represent more and less dense regions within the cell. The similar RI values fitted to the representative cells, 1.3677 (QLSI) and 1.3673 (DHM), indicate consistency across the two different imaging setups and confirm resilience and accuracy of the 2D Mie fit with regard to phase noise, which is further discussed in the next section.

Comparison of the scattering models
Each of the methods presented to determine the average RI of spherical objects has benefits and drawbacks. Mie theory yields the best theoretically possible result, but it is computationally expensive. The Rytov approximation is less expensive and, with a systematic correction (equations 3 and 4) it can achieve noteworthy accuracy. The OPD projection approach is faster than the Rytov approximation, but its accuracy is limited because it does not take into account diffraction. The same applies to the OPD edge-detection approach which is even faster but exhibits a resolution-dependent overestimation of the RI 5 . To extend the theoretical insights obtained with the error maps shown in figure 2, a comparison of the scattering models based on four representative experimental data sets is shown in figure 5. Figure 5a shows the RI values determined for a single FUS protein droplet from a time series recorded with the setup introduced in the previous section. Note that, due to the low resolution of the phase images ( fig. 4d), the RI values computed using the OPD edge-detection approach are severely overestimated and spread-out. Figure 5b shows RI values for microgel beads of type A imaged with the corresponding QLSI setup described in the previous section ( fig. 4g). Here, the difference between the scattering models is not as large as for the protein droplet, because the resolution is higher and because the average RI of the beads is much closer to the surrounding medium (PBS, n PBS =1.335). Figures 5c and 5d show the RI values fitted to two different HL60 cell populations recorded with QLSI and DHM (see figures 4n and 4r). Note that each method produces consistent RI values across both imaging modalities (cell populations). A comparison of the cell radii, shown in figure 10 in the appendix, indicates slight differences between the two populations, which is supposedly caused by biological variation.
In all examples shown, the systematically corrected Rytov approximation presents the best trade-off between accuracy and computational time. The high accuracy becomes most evident for the case of the FUS protein droplet which had a large RI and thus caused a prominent systematic error in the Rytov approximation. In case of the microgel beads, with RI values of about 1.34, the OPD projection approach 5 In fact, a re-assessment of the average RI of a set of microgel beads presented in a previous publication shows that a 2D fit with the Rytov approximation (n hydro2D-Ryt =1.354), when compared to the OPD edge-detection approach used in that publication (n hydro2D-edge =1.356), is a better match for the microgel RI of an artificial cell phantom obtained using ODT (n=1.354) [27]. and the non-corrected Rytov approximation still yield accurate results. However, at RI values that are observed in HL60 cells, the systematic correction of the Rytov approximation becomes important, producing RI values that are considerably closer to Mie theory than all other approaches. The total fitting time of the Rytov approximation is more than thirteen times shorter than the total fitting time of Mie simulations and less than three times longer than the fitting time of the OPD approach (e.g. for figure 5c: 17 minutes OPD projection, 48 minutes Rytov, 11 hours Mie) on a single CPU (Intel Core i7-4600U @2.10GHz, microcode version 0x21). Therefore, to accurately determine size and RI for spherical objects in QPI, we suggest the systematically corrected Rytov approximation in combination with a 2D fitting algorithm as presented here.

Summary
The present study provides an efficient method to accurately compute the RI and the radius for spherical objects from a single quantitative phase image. We investigated OPD approaches (OPD edge-detection, OPD projection), which were used in previous publications, and introduced efficient implementations of diffraction-based models (Mie theory, systematically corrected Rytov approximation) in combination with a 2D phase fitting algorithm, yielding improved accuracy and stability. Figure 6 illustrates the validity ranges of the models used (shaded regions) and provides orientation by indicating the exemplary samples shown in figure 4. We identified the systematically corrected Rytov approximation (Rytov-SC) as a feasible model for everyday-use in basic research, offering high accuracy at low computational costs over a broad range of sample sizes and RI values. Our 2D phase fitting approach can be applied to most biological cells (taking cell sphericity as given), is resilient to phase noise (see fig. 3), and operates well even at low image resolution (see fig. 4a-f). To provide orientation to the reader, the RI values (normalized to water by multiplication of n sph with n water /n med with n water =1.333) and radii (converted to wavelengths) of the samples shown in figure 4, are displayed as individual labels.

Outlook
The presented combination of the 2D fitting algorithm and the systematically corrected Rytov approximation for a sphere is efficient, but still offers room for improvement. First, graphical processing units (GPUs) could significantly speed-up the Fourier transform and the image interpolation steps of our current implementation. Second, the choice for the sampling of the radius with 42 points to compute the Rytov approximation, as discussed in appendix A.1, is conservative and could be reduced to about 20 points. In this case, however, the correction formulas (eqns. 3 and 4) need to be updated. With these modifications, the fitting algorithm could in principle achieve real-time performance for live QPI analysis.
The fit residuals shown in figure 4 exhibit small systematic phase errors, visible as blue and red circles around the object perimeter. These can be attributed to the fact that the point spread function (PSF) of the imaging setups used is not taken into account in the simulation. At the cost of more computation time, the simulation images could in principle be convolved with the PSF which would yield lower residuals. However, we do not expect higher fitting accuracies from this approach, because the residuals shown are already low and because the number of pixels in the described region is small compared the total number of pixels that resemble the imaged object.
To resolve large intracellular compartments or to analyze cells that have an approximately ellipsoidal shape, our algorithm could be extended to support a superposition of spheres [41] and ellipsoids. This would allow to resolve subcellular compartments such as the nucleus and the nucleolus from a single phase image or enable a frame-by-frame analysis of cell volume and RI of elongated cells in an optical stretcher experiment. Such extensions would yield an enhanced picture of the imaged cell at the cost of additional fitting parameters.
An accurate measurement of the RI of spherical objects is an important prerequisite for emerging topics in adjacent fields of research. For instance, our approach could allow valuable insights into the aging process of FUS [36] by tracking and characterizing the fit residuals during the liquid to solid phase transition of individual protein droplets. This would allow to determine whether the aging process starts at the center of the drop, at its surface, or homogeneously throughout. Furthermore, our algorithms could be used to assess the quantitative imaging quality of ODT and QPI techniques, both of which frequently employ microgel or polymer beads as reference samples (e.g. [35,23,32,31,27]). Finally, complementing the optical analysis of microgel beads as presented here with a mechanical description using, for instance, atomic force microscopy, will allow to establish a well-characterized microgel bead toolbox [40] which would be an invaluable reference for the optomechanical analysis of cells using techniques such as Brillouin microscopy or the optical stretcher. Thus, the optical characterization of spheres is fundamentally important to address topical questions in cell biology and biophysics.

Conclusion
The methods presented here resemble an important foundation for the accurate characterization of spherical microobjects, including the optomechanical phenotyping of live cells, the study of protein droplets, or the classification of artificial beads. The presented approach is complementary to tomographic techniques, limited to spherical objects but able to deliver average values for RI and radius using only a single quantitative phase image. The derivation of an efficient model for light scattering by a sphere and the implementation of an automated 2D phase-fitting algorithm resemble a major advancement in accuracy, stability, and throughput. Hence, the presented approach is an attractive tool for many emerging techniques that rely on an exact optical characterization of spherical objects to address biophysical questions in basic and applied research.

Funding
This project has received funding from the European Union's Seventh Framework Programme for the Starting Grant "Light Touch" (grant agreement no. 282060) and from the Alexander-von-Humboldt Stiftung (Humboldt-Professorship to J.G.).

A.1. The Rytov near-field of a sphere
We obtained the field scattered by a sphere in the Rytov approximation by computing its analytical 2D Fourier transform followed by an inverse Fourier transform using the Fourier diffraction theorem. The use of an analytical solution in Fourier space has the advantage that it reduces artifacts that arise from the sharp boundaries of the spherical volume in real space and that the 3D problem is broken down to a 2D problem. According to the Fourier diffraction theorem, the Fourier transform of the 2D backgroundcorrected scattered field U B (k) is mapped onto a semispherical surface in the Fourier transform of the 3D object potential F (k) according to [22,42,43] where we used the notation of [43] (unitary angular frequency Fourier transform) with the Fourier transformed detector coordinates k D , the distance between the center of the object potential and the detector plane l D , the wave number k m = 2πn med /λ (refractive index of medium n med and vacuum wavelength λ), the factor M = 1 − k 2 x − k 2 y /k m , and the unit vector representing the direction of illumination, i.e. the rotational position of the sample relative to the detector normal, s 0 . The rotationally symmetric Fourier transform of a homogeneous sphere with radius R is given by where j 1 is the spherical Bessel function of the first kind of order one. Evaluating equation 9 at the frequencies that correspond to the spherical surface described by equation 8 and performing an inverse Fourier transform of the resulting U B,sph yields the background-corrected scattered field component in the Born approximation u B,sph . The Rytov approximation of the scattered field component u R,sph is then obtained by computing the exponential of u B,sph [42,44] with the inverse Fourier transform operator F −1 . Assuming a incident plane wave, the full complex field then computes to u sph (r) = 1 + u R,sph (r).
In practice, we computed the Rytov approximation in two steps. First, the scattered field is computed with a resolution that samples the sphere radius with approximately 42 points. Second, the amplitude and phase data are linearly interpolated to match the resolution of the images recorded in the experiment. This two-step approach ensures that the field computation is resolution-independent, allowing to derive a systematic error correction for the Rytov approximation to yield results comparable to Mie theory as presented and discussed in section 2.3.

A.2. Efficient Mie computation by averaging
To compute the field scattered by a sphere, an exact solution based on Mie theory can be used [34,45]. However, compared to linear models, e.g. the optical path difference (OPD) projection or the Rytov approximation, Mie theory is computationally expensive. Here, we propose an averaged Mie approach, where the cross-sectional fields along a line perpendicular and parallel to the polarization of the simulated light are averaged. This approach reduces the required Mie field computation from a 2D image to two 1D lines, resulting in an effectively unpolarized field. Figure 7 illustrates the error made using this approach, which is oscillating within the ±1% interval relative the the maximum optical path difference and can thus be neglected for the studies presented here.

A.3. A 2D fitting algorithm for spheres in QPI
The 2D fitting of a phase image with a sphere-scattering model is burdened by many calls of the modeling function. In practice, this bottleneck becomes problematic when the scattering model is computationally expensive which is the case for e.g. Mie simulations. Here, we used a custom fitting algorithm that relies on interval-based phase image interpolation to reduce the number of calls to the modeling function. The algorithm reliably fits Mie, Rytov, and OPD projection models to a 2D phase image of a sphere with the five parameters φ bg background phase image offset, x x-coordinate of the center of the sphere, y y-coordinate of the center of the sphere, n sph refractive index (RI) of the sphere, and r sph radius of the sphere.
Given initial estimates for center (x 0 , y 0 ), RI n 0 , and radius r 0 of the sphere, the algorithm iteratively fits a sphere model in six steps.
1. Vary radius. First, compute three phase images using the sphere model at the radii, r i − b r , r i , and r i + b r , where initially i = 0 and b r = 0.05r 0 . Second, generate 47 phase images in the interval [r i − b r , r i + b r ] by interpolation. Third, compare the interpolated phase images with the measured phase image using the root mean square (RMS) error. The radius with the lowest RMS error r i+1 is used for the following steps.
2. Vary RI. As in step one, compute three phase images with the sphere model at the RI values, n i − b n , n i , and n i + b n , where initially b n = 0.1(n 0 − n med ) (RI of the medium n med ). Then, interpolate 47 phase images in the interval and select the one with the lowest RMS error n i+1 .
3. Vary center. Initially, the interval parameter for the center coordinate b c is set to one wavelength (λ) or 5% of the radius, depending on which is larger b c = max(λ, 0.05r 0 ).
Then, 13×13 phase images are computed for the inter- The phase image with the lowest RMS defines the center position (x i+1 , y i+1 ) for subsequent iterations.
4. Phase background estimation. Experimental phase images can contain a constant phase offset. We estimate the phase offset φ i by averaging the background phase values of the experimental phase image. To determine the pixel locations of the background phase image, the intersection of two pixel sets (a) and (b) is used: a) all pixels within a frame border of the phase image with a width of either 5 pixels or a fifth of the image size, depending on which is larger, and b) all pixels of the simulated phase whose absolute value is lower than 1% of the maximum simulated phase.

5.
Scale down interval parameters. The interval parameters b r , b n , and b c are individually divided by two, depending on the following conditions: b c : The change in the location of the center position Scaling down the interval parameters leads to a parameter refinement in each iteration of the algorithm.
6. Stopping criteria. The algorithm stops iterating when all of the following conditions are met: • The change in radius |r i+1 − r i |/r i+1 is smaller than 0.1%.
• The change in RI |n i+1 −n i | is smaller than 0.0005.
If the stopping conditions are not all met, then the algorithm proceeds with step one.
As discussed in section 2.2.2, the initial fitting parameters can be obtained with the OPD edge-detection approach. To demonstrate the correct convergence of the fitting algorithm, figure 8 shows the 2D Mie fits to the Mie simulations, as is done in figure 3 for OPD edge-detection, OPD projection, and Rytov approximation. Besides the fitting errors due to resolution-impaired unwrapping problems (region labeled with white asterisk), which is discussed in section 2.2.1, the Mie fitting error is negligible.

A.4. The corrected Rytov approximation
As discussed in section 2.3, the Rytov approximation, implemented as described above, exhibits a resolutionindependent systematic error compared to Mie theory. Thus, we could derive a correction factor (equations 3 and 4) for the Rytov approximation, yielding fitting accuracies close to Mie theory with comparatively little computational effort.
To ensure that the systematic correction for the Rytov approximation is independent of the RI of the medium, we fitted the Rytov approximation to a set of Mie simulations with varying RI of the medium. Figure 9 shows the fitting results of the Rytov approximation and its systematic correction for a sphere with a radius of ten wavelengths and a constant RI difference between the sphere and the medium of 0.05. Figure 9: Validity of corrected Ryotv approximation. The systematic correction for the Rytov approximation given by equations 3 and 4 is independent of the refractive index of the medium as shown for a) the refractive index of the object and b) the radius of the object. The exact data to which the Rytov approximation was fitted were computed using Mie theory for a sphere with a radius ten wavelengths and a resolution of five pixels per wavelength on a grid of 150×150 pixels.
As for the fitted RI shown in figure 5, figure 10 shows the fitted radii for all algorithms presented in our study. The fit with the systematically corrected Rytov approximation (Rytov-SC) consistently yields radii resembling those of the Mie theory fit when compared to the non-corrected Rytov approximation fit.

A.5. Cell refractive index values in literature
Cells can exhibit a wide range of RI values. Table 1 lists the average RI of several cell types with literature references. Note that the differences in the RI values between subsequent publications (e.g. pancreas carcinoma cells: 1.375 versus 1.380) can be explained by large standard deviations (data not shown here). Overall, most cell types have an RI value that resides in the interval above 1.35 and below 1.40, which we used to indicate the biological relevance of this study in figure 6.