Phasor approach of Mueller matrix optical scanning microscopy for biological tissue imaging.

Mueller matrix microscopy is an advanced imaging technique providing a full characterization of the optical polarization fingerprint of a sample. The Lu-Chipman (LC) decomposition, a method based on the modeling of elementary polarimetric arrangements and matrix inversions, is the gold standard to extract each polarimetric component separately. However, this models the optical system as a small number of discrete optical elements and requires a priori knowledge of the order in which these elements occur. In stratified media or when the ordering is not known, the interpretation of the LC decomposition becomes difficult. In this work, we propose a new, to our knowledge, representation dedicated to the study of biological tissues that combines Mueller matrix microscopy with a phasor approach. We demonstrate that this method provides an easier and direct interpretation of the retardance images in any birefringent material without the use of mathematical assumptions regarding the structure of the sample and yields comparable contrast to the LC decomposition. By validating this approach through numerical simulations, we demonstrate that it is able to give access to localized structural information, resulting in a simple determination of the birefringent parameters at the microscopic level. We apply our novel, to our knowledge, method to typical biological tissues that are of interest in the field of biomedical diagnosis.

Significance

Mueller matrix microscopy provides a full characterization of the optical polarization fingerprint of a sample without using any dyes. The gold standard method for extracting the elementary optical components is based on the Lu-Chipman decomposition and has proven its efficiency in quantifying the main polarimetric parameters as the dichroism, birefringence, and scattering. In this work, we demonstrate that combining Mueller matrix microscopy with a phasor approach provides a direct interpretation of the retardance images in any birefringent material without modeling. We believe that this work is a starting point for the interpretation of the localized polarimetric fingerprint inside tissues in a direct approach, leading to future applications for in situ diagnosis of pathologies and diseases that could assist histopathological evaluation.

Introduction

Mueller matrix (MM) microscopy is a label-free method based on the control of generated and analyzed polarization states, providing a full characterization of the optical polarization fingerprint of a sample (1). At least 16 intensities are needed to completely determine the 4 × 4 matrix of the sample, in which each independent element m can be associated with a specific physical phenomenon, such as linear or circular dichroism, birefringence, and scattering (2, 3, 4, 5). MM microscopy has been successfully used in numerous research fields such as materials science (6), biomedical diagnostics (7, 8, 9), and microscopic imaging (10, 11, 12, 13), showcasing its versatility and easy implementation on multimodal systems (14). To highlight a particular physical process, decomposing the MM is a mandatory step for tracking the amplitude of a specific polarimetric effect. Among all the decomposition methods, the gold standard Lu-Chipman (LC) polar decomposition is based on a priori knowledge of the sample and models the propagation of the light through a successive arrangement of elementary optical features composed of a dichroic, a birefringent element, and a depolarizer. This method, based on the Stokes-Mueller formalism, uses multiple matrix inversions that drastically increase the estimated systematic errors, so the physical meaning of the MM could be misinterpreted. Additionally, aside from this decomposition, numerous works propose manual modification of the order of the elementary optical components to match the most realistic conditions (15,16). This step requires a priori knowledge of the sample and can be arduous, as finding the proper decomposition corresponding to the physical polarimetric fingerprint in the case of complex and mixed biological media can be particularly challenging. Based on these issues, modifying the order of the optical features that model the sample properties results in wrong polarimetric values, disconnected from the real composition of the sample.

The phasor plot approach has been extensively applied to graphically analyze spectra in a variety of advanced microscopy techniques. It was originally introduced in the framework of fluorescence lifetime imaging (FLIM) analysis, as an alternative to multiexponential decay fitting (17). In the phasor FLIM approach, the fluorescence decay at each pixel is analyzed in the frequency domain via a Fourier transform algorithm that generates, for each pixel (x, y) and for a given frequency, a complex number (g,s) (17). The values of (g, s) associated with all the pixels of an image are usually represented in a two-dimensional histogram, called the phasor plot. This graphical representation of FLIM data has gained increasing popularity among FLIM users, as it allows identifying multiple molecular species (18) and quantifying Förster resonance energy transfer interactions (19) without any prior fluorescence decay model assumptions (20,21). More recently, the phasor plot approach has been extended to the analysis of other types of microscopy data, including spectral images (22), image correlation spectroscopy data (23,24), and super-resolution images (25,26). In the context of nonlinear, label-free microscopy, the phasor plot has been used to analyze polarization-resolved second harmonic generation data and map the features of collagen architecture in tissues (27). Moreover, a recent theoretical work has proven the phasor space based on the analysis of circular polarized light scattering to be an intuitive tool, sensitive to the chiral organization of biopolymers at the submicroscopic level (28).

Here, we combine MM microscopy with the phasor approach to facilitate the interpretation of MM imaging of biological samples. In our method, named MM-phasor, the full MM is first converted into four Stokes vector images, and then a Fourier transform algorithm is applied to each of these images. In this way, the coupling of both techniques allows an alternative interpretation of the multiparametric polarimetric contrast at a molecular level, independent of the sample position, giving similar values to the LC decomposition. Moreover, we demonstrate numerically that MM-phasor is able to retrieve the retardance value in the case of randomly layered, oriented birefringent materials. This new, to our knowledge, approach is a model-free method for accessing some of the polarimetric information, whereas LC fails because of the requirement for a priori ordering of the successive optical features of the sample for such MM decomposition. For this purpose, we discuss the sensitivity of such an approach for discriminating between structures of interest that exhibit similar contrast. As a proof of concept, we apply this approach to imaging of biological specimens (collagen and myosin fibers) showing mostly birefringence. Such materials exhibit negligible dichroism and/or depolarization at the microscopic level because of their high transparency and weak scattering. For this reason, in this work we briefly consider the influence of these two effects on the phasor map. Finally, to understand the influence of the polarized light on sample orientation, we propose decoupling the effect of both the input and output polarized light from the sample optical orientation.

Materials and methods

MM optical scanning microscope

The experimental setup was a modified inverted commercial laser scanning microscope Olympus FV1000 (Tokyo, Japan), reported in our previous work (13). The excitation source is a laser diode at 808 nm. Polarization information is encoded by four polarization states (0, 90, and 45° and right circular polarization) from a motorized polarization states generator (PSG) placed before the microscope body, composed of a linear polarizer, a half-wave plate, and a quarter-wave plate. To image the samples in a reflection configuration, we placed a silver plane reflective mirror at the focal plane of a 10×/0.4 NA objective. The passive polarization states analyzer (PSA) analyzed the transformation of the polarization after the interaction of the polarized light with the sample. It consisted of a thin glass slice, a Fresnel rhomb, two Wollaston prisms, and four photodiodes to collect the signal pixel by pixel for four polarization states. The four input polarization states generation and the four photodiodes were synchronized with the scanning unit via a custom LabVIEW routine via a data board acquisition. Finally, a MATLAB (The MathWorks, Natick, MA) program was used to extract the full MM image of the sample from the 16 acquired polarization-resolved 256 × 256 images.

As we have shown in our previous work (13), the calibration procedure includes taking into account the remaining optical orientation misalignment of the PSG and PSA, as well as the polarimetric signature of the optical devices in the microscope body, composed mainly of lenses and reflective mirrors. It is based on the eigenvalue calibration method to calibrate the PSG and PSA using reference elements such as air, linear polarizer, and quarter-wave plate (29). The accuracy of acquiring the full MM is determined by the condition number of the polarimeter (PSG and PSA without microscope); in this study, this corresponds to condition number = 2.6, presenting uncertainty less than 5% on the Mueller coefficients. Next, we image the MM of the whole system to take into account the polarimetric fingerprint of the microscope pixel by pixel over the entire field of view. Thereby, MM images of the sample are obtained by directly inverting the microscope fingerprint by applyingwhere [M(x, y)], [M(x, y)], and [M(x, y)] are the MM images of the sample, of the total system microscope + sample, and of the microscope only, respectively. This equation neglects the influence of optical components that follow the sample in the reflection path, as discussed in our previous work (13). Furthermore, we have demonstrated that the double pass of light through the sample in the epigeometry imaging causes a π-dephasing effect at the reflection on the mirror surface. Thus, the retardance value R from the double-pass optical path is defined as 2π, so that the exact value R is recovered by applying R = 2π − R. Considering that the samples are nondepolarizing, the increase of the scattering through the sample after the double pass is considered in the standard deviation (SD).

MM image acquisition and analysis

The physical effects are extracted from an implemented LC decomposition algorithm (1). It is based on the representation of the sample by the three main optical interactions successively crossed by the excitation light: dichroism, birefringence, and scattering. These interactions correspond to the product of the three matrices associated with the polarimetric parameters:where [MΔ], [M], and [M] are the sample matrices describing the total depolarization, diattenuation (linear and circular), and birefringence (linear and circular), respectively. After normalization by the total transmitted light, corresponding to the m00 element, the physical parameters denoted in this work are retrieved byandwhere P, R, and D correspond to the depolarization index (30), retardance, and diattenuation; m are the MM elements; tr(M) is the trace of the retardance matrix M; and (r1, r2) are the elements of the retardance vector (31). In this way, each of these three physical parameters varies between 0 ≤ P ≤ 1, 0° ≤ R ≤ 180°, and 0 ≤ D ≤ 1 and is used for coding the post-treated images pixel by pixel. It is worth noting that the amplitudes of these three parameters are directly related to the sample properties. The values of orientation corresponding to a weak signal of diattenuation and retardance are not physical and are coded in black in the α images.

The MM and the associated LC decomposition data are presented in Table 1, averaged over the region of interest (ROI) 1 and ROI 2 for air and a transparent adhesive film area, respectively, from Fig. 6.

Table 1

Experimental normalized MM of homogeneous air and transparent adhesive film samples and associated LC polarimetric parameters averaged over the ROIs

Experimental MM	LC
[1.0000.003±0.009−0.009±0.0110.012±0.022−0.005±0.0240.999±0.022−0.002±0.025−0.008±0.0490.011±0.0190.007±0.0181.002±0.025−0.011±0.0420.007±0.014−0.009±0.0100.004±0.0121.014±0.031]	D = 0.027 ± 0.014
	R = 1.9 ± 0.9°
	Pd = 1.007 ± 0.01
[1.0000.012±0.009−0.005±0.0110.071±0.029−0.012±0.0291.023±0.032−0.222±0.026−0.327±0.0630.007±0.0540.376±0.0240.678±0.0261.022±0.078−0.023±0.0160.055±0.019−0.342±0.0180.562±0.036]	D = 0.072 ± 0.026
	R = 55.5 ± 1.9°
	αR = −85.7 ± 0.9°
	Pd = 1.051 ± 0.029

MM values taken from Fig. 6, ROIs 1 and 2.

Figure 6

Phasor plots of ROIs from a transparent adhesive film sample surrounded by air. (a) Phasor plots for the air ROI 1 (red), (b) phasor plots for the transparent film ROI 2 (blue), and (c) phasor plots for both air and transparent film tape ROI 3 (green) for each of the four input polarization states. (Top panel) Images coded in total collected light intensity m00. (Left panel) Raw images coded in H and RCP intensity.

As shown in Table 1, the calculated values match with the expected values for air and transparent film. Indeed, the expected values for a nondepolarizing medium (P = 1) such as air are D = 0, R = 0°, and P = 1, whereas for the transparent film they are D = 0, 0° ≤ R ≤ 180°, and P = 1. It is worth noting that the MM of transparent film presents higher MM coefficient values related to the dichroism induced by the interface of the layer. The wider SD for the transparent film, especially for α, comes from the inhomogeneity of the surface for such a sample.

The phasor plot analysis described in Result and discussion: MM-phasor analysis method section was implemented in MATLAB (The MathWorks). To remove residual biological components and clean the phasor plots from background contributions, a ROI was defined selecting only the fingerprint of the fiber.

Sample preparation

Muscular tissue was harvested from the leg of a raw rabbit purchased in a supermarket and was cleaned from other connective tissue. The muscular tissue was rinsed several times with phosphate-buffered saline solution (pH 7.4) and fixed with 4% paraformaldehyde in phosphate-buffered saline solution overnight. The day after, the tissue was dehydrated by using a series of ethanol solutions of increasing alcoholic concentration. Next, the sample was embedded in paraffin before sectioning using a microtome (Leica, Wetzlar, Germany). The thickness of the muscle sections was of 5 μm. The sections were placed and fixed on a microscope slide and imaged.

Results and discussion

MM-phasor analysis method

In our previous work, we developed a full MM optical scanning microscope based on the sequential generation of four polarization states 0, 90, and 45° and right circular polarization states (RCPs) through the PSG (13). The polarized radiation can be described by a 1 × 4 Stokes vector  = (S0, S1, S2, S3), where the coefficients S0, S1, S2, and S3 are obtained from the combination of polarized light intensities I0° (horizontal), I90° (vertical), I45°, I135°, IRCP (right circular), and I (left circular), giving the following set of equations (32):and

The transformation of each input polarization state after interaction with the sample is recorded pixel by pixel by the PSA allowing the determination of the 16-element MM of the sample. The link between the input and output polarization states, and respectively, is established using the MM of the sample [M(x, y)]:

In the case in which an optical block is added between the PSG and PSA, the polarimetric signature of this feature must be taken into account. More particularly, the output Stokes vector can be directly retrieved as a function of input experimental polarization states [W] = byandwhere , , , are the experimental input Stokes vectors linear horizontal, linear vertical, linear 45°, and RCP polarization states.

First, to convert the MM into a phasor plot, we reconstruct the linear polarization-resolved intensity from the output Stokes vectors of Eqs. 12, 13, 14, and 15. Indeed, the linear output polarization states intensities I (θ = 0, 45, 90, or 135°) can be extracted directly from each previous Stokes vector defined by Eqs. 7, 8, 9, and 10 byand

Analyzing the polarization-resolved signal I is similar to rotating the polarization analyzer after the sample, with the aim of following the anisotropic polarization emission of molecules under illumination, as widely realized for studying nonlinear anisotropic molecular orientation (33). Because we follow the linear polarization dependence, the circular polarization intensity is omitted but can be extracted from the combination of the S0(x, y) and S3(x, y) Stokes coefficients.

The phasor P(x, y) associated with the polarization-resolved image I(x, y) has coordinates g(x, y) and s(x, y) given byandwhere Θ = 180° and 1/Θ is the frequency used for the Fourier transform. Thus, every single pixel can be identified as a vector in phasor space with two phasor coordinates (g, s) mapped to (x, y) coordinates. The step-by-step protocol for converting the MM signal into a phasor, for each pixel, is summarized in Fig. 1. In general, from the phasor approach, it is possible to analyze multiple frequencies in the signal. In this work, we have limited our analysis to a single frequency. The phasor is a vector of modulus M and angle φ given byand

Figure 1

Schematic diagram of the protocol for converting the MM image into a phasor plot. (a) The MM image conversion into Stokes-Mueller images, then giving linear polarization images for input polarization states H, V, 45°, and RCP. (b) For a given input polarization state, a stack of images of the modulated linear polarized light signal is built pixel by pixel (x, y), which can be plotted in the phasor space at the coordinates (g, s).

The parameters M and φ are usually called modulation and phase and, as we shall see, are important in the interpretation of the output polarimetric signal.

The phasor coordinates (g, s) can be calculated back from M and φ fromand

In the field of FLIM, the phasor approach defined M to be between 0 and 1 and φ between 0 and π (17). Here, from the Stokes-Mueller formalism, the orientation of the optical fast axis of the sample is defined between −π and π, modifying the definition range of M ∈ [0; 0.5] and φ ∈ [−π; π] (1,34). Finally, the (g, s) coordinates are calculated pixel by pixel from the polarization-resolved images using Eqs. 20 and 21. In this work, the two-dimensional histogram of the values g(x, y) and s(x, y) associated with each pixel are then represented in the phasor plot, in which the phasor spot is coded from red to blue, for the highest to the lowest bin. A pixel-by-pixel phasor representation could also be used, but then imaging contrast information is lost, and it becomes arduous distinguishing the location of each component of the image.

Numerical simulation on birefringent materials

To study the phasor plot changes as a function of the polarization input, we simulate the behavior in phasor space of well-known reference samples that are easy to assemble experimentally, for instance, air and birefringent materials.

Air sample

First, we simulate the phasor signature from the air sample, i.e., without placing a specimen in the point spread function (PSF). To have comparable results with experimental data, we add the statistical SD measured experimentally and reported in Table 1 in Materials and methods. Thus, the corresponding air MM is given bywhere I4 is the 4 × 4 identity matrix and Δ(x, y) is the additional error matrix corresponding to the experimental SD of the microscope. By applying Eqs. 20 and 21, the phasor plots for the air sample for each input polarization states are presented in Fig. 2.

Figure 2

Phasor plot of a simulated air sample for each of the four polarization state inputs: H (linear horizontal), V (linear vertical), 45° (linear 45°), and RCP (right circular polarization).

The simulated air sample exhibits the properties of the input excitation light because of the absence of any interaction with polarimetric objects. The H and V polarization state inputs reveal a phasor in the positive (φ = 0°) and negative (φ = π) x axis, respectively. Considering that the electric field oscillates in the perpendicular orientation between these two polarization states, this means that the phasor defined the polarization orthogonality at π. This is confirmed by the 45° input, which exhibits a phasor at φ = π/2 from the H polarization states input in the positive y axis. This behavior depends on the period used for the Fourier transform (Θ = 180°). Indeed, the observed polarization angle ψ (i.e., the angle at which the polarimetric signal has a peak) is the angle of the output polarized light propagation plane with the (x, y) laboratory axis, equal to ψ = 0°, 45°, and 90° for the H, 45°, and V output polarization states, respectively. It is related to the phase φ in the phasor space by the relationship

Thus, the phase of the MM-phasor has a range −π ≤ φ ≤ π, corresponding to the range −π/2 ≤ ψ ≤ π/2 for the linear polarization angle. Here, the modulus in the absence of any anisotropic materials is equal to M = 0.5 for all inputs. In addition, the RCP input shows a phasor centered at the origin, indicating that the polarized light interacts homogeneously and has no preferential orientation, as predicted.

Birefringent sample

In this work, we study only the MM-phasor of a simulated birefringent material. Indeed, birefringence is the main polarimetric parameter, presenting the highest imaging contrast in biological tissues. We simulated the phasor plot of a birefringent material of retardance R and optical fast axis orientation α, using the following MM:withwhere Δn = (n − n) is the birefringence, e is the thickness of the retarder, and λ is the excitation wavelength.

The phasor plots are presented in Fig. 3 for different R-values as a function of α for the linear and circular input polarization states. To improve the visibility of each spot on the phasor plot, we kept the same signal/noise ratio (SNR) as for air, and the orientation angle corresponds to 0, 10, 20, 30, 40, and 45°. To simplify the phasor representation, we show only the H and RCP phasor plots. We discuss in the Supporting materials and methods how to merge these into a single phasor plot by averaging the three linear polarization inputs (H, V, and 45°).

Figure 3

Phasor plots of a simulated retarder for linear and circular polarization input states as a function of the optical fast axis orientation α (between 0 and 45° in steps of 10°) for a retardance of (a) π, (b) 3π/4, (c) π/2, and (d) π/4. (e) Schematic diagram of the phasor plot pattern for a retarder as function of the R and α parameters for linear and circular polarization input states.

First, we consider the H phasor plot. For α = 0°, the phasor is identical to that obtained without a sample in Fig. 2. The reason for this is that there is no polarization change when the input polarized light is aligned with the fast axis of the birefringent object, in accordance with the polarization formalism (1). For R = π, the phasor rotates toward π as the value of α increases from 0 to 45° (Fig. 3 a). For R < π, from R = 3π/4 to R = π/4 (Fig. 3, b–d), the phasor follows a similar rotation in function of α, but the center of this rotation changes as a function of R (Fig. 3, b–d, top). This center translates until the limit case R = 0, in which the phasor becomes identical to the air phasor. Thus, in the H phasor plot, there is no straightforward relation between the values of R and α and the values of M and φ. If we consider the RCP phasor plot, for R = π, the RCP phasor is identical to the air phasor, whatever the value of α. For R < π, the phasor rotates from the angle φ = π/2 (α = 0°) to the angle φ = π (α = 45°) for any value of R (Fig. 3, b–d, bottom). Thus, for R < π, the retarder orientation α is linked directly to the phase angle φ as (see Supporting materials and methods)where the phasor angle φ is defined between −π/2 and 3π/2, we get values of α between −π/2 and π/2. The modulation observed at the RCP phasor is independent of the optical fast axis orientation of the retarder and depends only on R (Fig. 3, b–d, bottom). Thus, for R < π, the retardance R is linked directly to the phasor modulation M bywhere the modulation M varies between 0 and 0.5, corresponding to R-values defined between 0 and π/2. For any given value of modulation M, another valid solution is

In other words, in the RCP phasor we get the same modulation M for the retardance values R and π − R (see, for instance, Fig. 3, b and d). Notably, we can discriminate between the cases R < π/2 and R > π/2 by additionally reading the H phasor. For R < π/2, the phasor lies within the small circle of center (0.25, 0) and radius 0.25 (area inside the red line in Fig. 3 e, top left), whereas for R > π/2, the phasor lies outside the small circle and within the circle of center (0, 0) and radius 0.5 (area inside the blue line in Fig. 3 e, top right).

In summary, as presented in Fig. 3 e, we numerically demonstrated the capability of the MM-phasor to extract R and α from the RCP and H phasor plots for any type of birefringent specimen. For 0 < R < π, these polarimetric parameters are obtained by simply extracting M and φ in both H and RCP phasors based on Eqs. 30, 31, and 32). For the specific case R = π, only the H phasor is needed for obtaining α, as α = φ/4 (Fig. 3 e, top), which speeds up the analysis. However, the limitation of this approach comes when determining the values from a weak birefringent material, resulting in M-values too small compared to the SNR. The reason is that the phasor spot evolves too close around the center and is assimilated to the case in which the sample is only air. We calculated the modulation value as a function of R for determining the sensitivity of MM-phasor in Fig. 4 at the same orientation α = 0°. This is similar to finding the minimal retardance in which the phasor spot is included in the experimental SD reported for the air sample. As explained in the previous section, the phasor spot lies on the g axis with s = 0. The simulated data have been extracted from the RCP phasor because R is extracted from this space as demonstrated in Eq. 32. Hence, we have presented the modulation value for a retardance from π/2 to π/40, resulting in M-values between 0.5 (maximal M) and 0.03 as previously discussed above.

Figure 4

Evolution of the modulation value as a function of the inverse retardance as a multiple of π and the corresponding raw simulated data. The data are calculated in RCP phasor space, with retardance between π/2 (= 90°) and π/40 (= 4.5°), with the same orientation α = 0°. a.u., arbitrary unit.

In the Materials and methods, the measured SD leads to an uncertainty on the air M-value equal to 0.05. This means that below this modulation, the phasor spot provided by a weak birefringent specimen can be assimilated to the air, therefore comprising the resolution limit of the MM-phasor approach. From Fig. 4, it can be seen that M = 0.05 is reached for R = π/31 ∼5.8°. Thus, the MM-phasor method is not able to discriminate structures of interest that exhibit a contrast below this limit, comparable to the SNR reported in Table 1 for the LC decomposition.

Layered birefringent fibers

To characterize the performance of the MM-phasor approach in the case of a histopathological sample, we simulate a structure composed of a successive arrangement of elementary randomly oriented birefringent features. The model shown in Fig. 5 a is basically related to the situation one encounters in the presence of a pathology that has induced a disorganization of the connective tissue. The modeling of such a medium consists of simply multiplying the retarder MMs reported in Eq. 28 with the same R-value but at different orientations. Here, we simulate an image containing one, two, three, or four sets of simulated fibers, calculate (g, s) pixel by pixel from this image and the associated phasor plot, and then extract the retardance from this plot. In our model, shown in Fig. 5 a, the medium under illumination corresponds to a successive arrangement of birefringent components with the same thickness e and oriented at 0, 10, 20, and 30°. To make a comparison between the LC and the phasor approach in realistic conditions, we consider a type I collagen fiber with a thickness of e = 50 μm and Δn = 0.0017 (35). Thus, R is calculated from Eq. 29 and provides the ground truth of the simulation, equal to R1 = 37°. We calculate the MMs for a single zero-degree layer and then add the differently oriented layers one by one. The determination of R between the two approaches is compared to the ground truth as a function of the number of fibers in the illumination volume and is presented in Fig. 5 b. Noting that the value of R providing in this section corresponds to a single point in the birefringent fiber with a specific thickness and without dimension in the xy plane (perpendicular to the PSF volume).

Figure 5

(a) Schematic of the layered oriented birefringent medium. Each layer has a thickness of e = 50 μm and is oriented at 10° compared to the previous one. (b) Retardance value, as a multiple of R1 = 37°, as a function of the number of fibers in the illumination volume. The blue, red, and orange plots correspond to R-values calculated from the ground truth, the LC decomposition, and the MM-phasor approach, respectively. (c) Depolarization index, P, calculated from the LC decomposition, as a function of the number of fibers in the illumination volume.

In the case of a single fiber oriented at 0°, both LC and MM-phasor approaches provide a good agreement with the expected value R1 = 37°. By adding a second layer and a third one, oriented at 10 and 20°, the expected values are 2 × R1 = 74° and 3 × R1 = 111°. Thus, we notice a strong dispersion by using the LC approach, which yields 93 and 141°, whereas the MM-phasor approach is still consistent with the expected values, 68 and 121°. By adding a fourth fiber oriented at 30°, the expected value is 4 × R1 = 148° and the MM-phasor approach is closer to the real value of 149°, whereas the LC fails by giving 185°. The reason for the nonconsistency of the LC decomposition when the complexity of the sample increases is that the model composed of an arrangement of a dichroic, birefringent element and a depolarizer is no longer suitable. As shown above, the ordering of the optical components is crucial for physically interpreting the MM of the sample, and inverting it drastically influences the value of the polarimetric parameters, as illustrated in the Supporting materials and methods. In this section, we illustrated the difference when the sample is composed by a retarder and a dichroic element and its inverted ordering. In the case of a sample composed of randomly oriented structures, the polarimetric signature becomes more complex to interpret in the absence of any spatial sectioning. This is shown by the P-value reported in Fig. 5 c, extracted through LC decomposition, that decreases when the number of fibers increases, from P = 0.97 with a single fiber to P = 0.54 with four fibers. Meanwhile, the MM-phasor approach is capable of retrieving the value of R and its orientation without being affected by the degree of P. This is only true if we assume that the sample is birefringent, and the influence of depolarization is discussed in the Supporting materials and methods.

In this work, dedicated to the specific case for imaging thin and microscopic birefringent objects, we neglect the presence of any dichroism and depolarization, as reported in Table 1. However, as demonstrated in the case of birefringent randomly oriented multilayers, the sample could exhibit a strong depolarization, making determination of the polarimetric parameters as arduous as it is the for the LC decomposition. We have performed numerical simulations related to the influence on the phasor for a pure dichroic (D < 0) and a pure depolarizing (P < 1) medium (P < 1 and R 0°) presented in the Supporting materials and methods. To sum up, in the case of a pure dichroic medium, the zero-degree orientation angle (α = 0°) for the linear polarization input phasor presents a modulation value smaller than for a retarder (M < 0.5) and becomes maximal on the s axis at α = 45°. Additionally, the translation in RCP phasor space as a function of α is inside the positive g region, as opposed to a pure retarder. Thus, these two pieces of information present an interesting approach for discriminating a pure retarder from a pure dichroic medium. In the case of a pure depolarizing medium, it has been shown that the phasor spot again exhibits a smaller modulation value than for a pure retarder, translating at any orientation to the origin of (g, s) coordinates for highly depolarized samples. This means that for a highly scattering sample (P ≪ 1) mixed with a retarder, the birefringent parameters cannot be retrieved properly and this phasor approach is limited, as is the case when LC decomposition is applied (36).

Validation on reference sample

We apply the MM-phasor approach on a homemade sample composed of transparent adhesive film surrounded by air, presented in Fig. 6. This sample is of interest because of its strong birefringence signal, its simple model, and the fact that it can be discriminated from air or other weak polarimetric objects. To analyze the phasor plots for the different parts of the sample, we segmented the image separately for the air and the transparent film area by the corresponding ROI 1 (red) and ROI 2 (blue), respectively, and combined the two objects in ROI 3 (green). The arbitrarily chosen ROIs are drawn on the associated m00 image of the sample, which corresponds to the total collected light image from the MM.

The phasor for the generated linear and circular polarization states of the air area, ROI 1 (in red), is presented Fig. 6 a, in accordance with our simulation (Fig. 2). A small difference from the previous simulation comes from the fact that the modulation of the phasor spot is less than 0.5, induced by an R-value between 0 and π/2. Fig. 6 b presents the H and RCP phasors of the transparent film area, ROI 2 (in blue), and Fig. 6 c corresponds to the area with both air and tape, ROI 3 (in green), for comparison. The H phasors for the air and the film regions are almost identical, indicating that the optical fast axis of the film is close to 0°, parallel to the input polarized light. A superposition of the phasors of the two materials is presented in Fig. 6 c, corresponding to an ROI taking into account both air and the transparent film area. This is of particular interest because it demonstrates our capacity to provide a distinct segmentation of the different areas in the sample by directly identifying the phasor point distribution in phasor space.

Conjointly, we demonstrate the capability of MM-phasor for converting the phasor (g, s) mapping into the (R, α) images. The reconstructed R and α images with the histogram distribution of the pixel distribution is presented in Fig. 7, extracted from LC decomposition and from the phasor plots in Fig. 7, aand b, respectively.

Figure 7

(a) R and α images extracted by LC decomposition and their associated pixel histograms. (b) R and α images extracted by comparing the H and RCP phasor plots with the associated pixel histogram.

The mean and SD values shown in Fig. 7 a correspond to those reported in Table 1 for a homogeneous surface of the tape. The absolute retardance value R = 2π × Δn × e/λ of the transparent film is retrieved, knowing its birefringence Δn = 0.0077 from the coefficients of Cauchy’s equation (37) and its approximate thickness e = 0.5 mm at the excitation wavelength λ = 808 nm. Thus, the retardance obtained for such a material corresponds to 60.7 (taking into account the phase wrapping for such a thick birefringent sample), reasonably close to the measured retardance by LC decomposition (R = 55.5°). The mean values of R and α from the phasor plots are equal to R = 54.1 ± 8.4° and α = −83.4 ± 1.5°, which provide a satisfying agreement with the data from LC decomposition (Fig. 7 b). However, we note a stronger dispersion in the values obtained by the phasor plot. In fact, LC is unaffected by changes in transmission, but the phasor approach is affected by undesirable changes in overall transmission. More precisely, LC uses the trace of the retardance Mueller 3 × 3 submatrix (M) composed by the block diagonal terms, excluding any attenuation terms (38,39). By comparison, our approach extracts R from the modulation intensity performed by a virtual rotation of the analyzer affected by the total transmittance of the sample. This decreases drastically the SNR of such an approach, but this obstacle could be overcome by studying thinner and/or microscopic objects, which is the aim of using the phasor representation presented in the next section.

MM-phasor analysis for biological samples

Collagen fiber

As a first application of the MM-phasor, we study a collagen fiber, well-known to be a uniaxial birefringent material (36). Collagen fibers are commonly used as an indicator to diagnose pathologies by detecting tissue disorder using MM imaging (40,41). To understand the influence of the sample orientation relative to the output polarization states, Fig. 8 presents collagen fibers from rabbit tendon after manual physical rotation of the sample. Here, we present the associated LC images, as well as RCP phasor plots for each sample rotation. Because of the weak birefringence effect exhibited by biological samples, the R-value is less than π/2. Therefore, we present only the RCP phasor here.

Figure 8

(a) and (b) Retardance images and orientation α obtained by LC decomposition of a collagen fiber from rabbit tendon as a function of its spatial rotation (∈ [0°, 180°]) at the RCP polarization input. (c) Associated RCP phasor plots of the same collagen fiber.

From the LC images (Fig. 8, a and b), the R-values are similar for different orientations of the fiber, R = 56.3 ± 6.9°, because of the independence of the MM to sample orientation (13). As expected, the associated α images in Fig. 8 b display shifts of 45° with respect to the physical rotation of the optical fast axis of the sample. The phasor for each of these sample rotations is presented in Fig. 8 c. Because of the inhomogeneity and the weaker birefringence of the sample compared with the transparent film, the phasor plot in the practical case exhibits a more scattered spot lying around the center. It is worth noting that the apparent SDs between each rotation of the fiber are explained by an inhomogeneous structure of the sample, coupled with a different z position of the PSF volume between each measurement. However, the spot spreads clearly in a preferential direction resulting from the orientation of the optical fast axis. Additionally, as predicted by Eq. 32, the preservation of the modulation amplitude is explained by the conservation of the R-value for all orientations, as shown in LC images (Fig. 8 a). We confirm that a sample rotation of Δα = 45° is characterized by a rotation in the phasor space for the RCP input of Δφ = 90°, in keeping with the simulations reported in Fig. 3.

To summarize, Table 2 presents (R, α) extracted with the gold standard polarimetric images from LC decomposition and with the RCP phasor plots. These results are obtained by averaging the angle taken from the ROI all along the collagen fiber, taking into account the phase wrapping effect for a physical angle greater than 90°. To confirm the 45° rotation, the differences of α between each rotation are presented with Δα.

Table 2

Averaged mean value of (R, α) for collagen fiber as a function of spatial rotation ∈ [0°, 180°] obtained from LC decomposition and the RCP phasor

Orientation		0°	45°	90°	135°	180°
Parameters
LC	R	45.4°	58.7°	59.6°	60.9°	64.6°
	αR	−68.6°	−20.4°	16.3°	55.8°	104.4°
	ΔαR	N/A	48.2°	36.7°	39.5°	48.6°
RCP phasor	R	69.9°	51.2°	62.3°	56.1°	71.9°
	αR	−69.2°	−38.2°	9.3°	51.9°	95.6°
	ΔαR	N/A	31.0°	47.5°	42.6°	43.7°

Collagen fiber presented in Fig. 8, a and b. RCP phasor shown in Fig. 8c. N/A, not applicable.

Table 2 shows an adequate preservation of the R-value for both methods, demonstrating our capability to properly predict the phasor plot behavior with respect to the rotation of the optical fast axis. Indeed, for LC, the averaged R-value for each the sample rotation is R = 57.9 ± 7.3° and becomes R = 62.3 ± 8.8° when using the RCP phasor. For the optical fast axis α, the averaged difference between each angle, Δα = θ − θ, is Δα = 43.3 ± 6.0° for LC decomposition, whereas it is Δα = 41.4 ± 6.1° for the RCP phasor. It is worth noting that for such a sample, the SDs are identical between the two approaches, in contrast to the transparent tape. Indeed, the sample transmittance is much weaker because of its thickness of a few hundred micrometers.

Myosin fiber

As a second application of the MM-phasor, our study of myosin fibers is able to provide an interesting polarimetric signal coming from microscopic-localized change in the complex of actin and myosin (42,43). Fig. 9 presents a myosin tissue sample from rabbit tendon to understand the single influence at the molecular level of the input and output polarization states in the most general case of in situ organized tissue. The reconstructed phasor images, encoded in modulation and phase, are compared with LC decomposition images.

Figure 9

(a) Total transmitted light m00 image and the raw H and RCP images segmented in four associated ROIs of myosin fibers from a rabbit tendon. (b) The retardance (R) and its orientation (α) images from LC decomposition, compared with R and α images from the phasor analysis. (c) Phasor plots for RCP input polarization, corresponding to the four ROIs in (a). ROI 1 (black area) corresponds to air area, ROI 2 (green) to the whole muscular tissue area, ROI 3 (red) to the α = −20° area of (b), and ROI 4 (blue) to the α = 20° area of (b).

Fig. 9a presents the m00 image of the myosin fibers, which provides poor contrast based only on the absorption and scattering of the incoming light. However, after applying LC decomposition, Fig. 9 b, left panel presents better image contrast, coded in (R, α), in which the fibers can be discriminated according to the orthogonal orientation of the actin-myosin complex. Indeed, it can be observed that the fibers rotate between two orientations, α = −20° and α = 20°, encoded in cyan and yellow, respectively. This symmetrical organization between the actin and the myosin was previously reported in the literature from nonlinear polarization microscopy (44). As expected, images reconstructed from the RCP phasor (Fig. 9 b, right panel) are in correlation with the structures highlighted by LC decomposition MM images. The RCP phasor plot (Fig. 9 c) from ROI 1 (air region), is in accordance with our previous analysis for such reference materials, meaning that no modulation is presented for the RCP input. Nevertheless, for the other ROIs, the phasor spots are slightly wider, probably because of the heterogeneity of the biological components. In particular, in ROI 2, a relatively large region containing fibers of different orientation, the wider spot corresponds to a simple averaging of the overall fingerprint of the muscular tissue. Nevertheless, when segmented separately by ROI 3 and ROI 4, the two oriented fibers reveal clearly two phasors having different phase value at the same modulation (same R-value).

Conclusion

In this work, we demonstrate that the combination of a phasor approach with MM microscopy provides an easy and intuitive interpretation of polarization-resolved images from birefringent samples. The method is shown to work when using samples with low dichroic and depolarizing (scattering) content. The impact of using this representation, instead of applying the gold standard LC decomposition, is the absence of a priori knowledge of the sample and of complex mathematical modeling of the elementary polarimetric components. We implemented an algorithm that provides access to the retardance and its orientation in any birefringent material by simply reading the (g, s) coordinates of the linear and circular polarization phasor plots. As preliminary results, we have applied a post-treatment protocol from all the raw polarization-resolved images, resulting in conversion into the phasor map. As future work, the phasor map of each pixel could be accelerated toward a real-time application by directly calculating the (g, s) coordinates for each scanned pixel, which would give an almost immediate result, depending on the scanning time for an image. We show that the MM-phasor can perform fast analysis of tissues, opening possibilities for future applications for in situ diagnosis of pathologies and diseases that could assist histopathological evaluation.

The data sets generated and/or analyzed and the MATLAB algorithm in this work are available from the corresponding author on reasonable request.

Author contributions

A.L.G., L.L., and A.D. conceived research, conducted the experiments, and analyzed the results. A.L.G., L.L., A.B., and C.J.R.S. wrote the manuscript. P.B., C.J.R.S., and A.D. contributed to useful discussions. All the authors reviewed and approved the manuscript.