2.1.

Monte Carlo Algorithm

The detailed description of the polarized MC algorithm and its implementation can be found in prior publications.^13,14 Here, we recall the main steps of the statistical modeling approach while omitting some details for brevity. A point light source, placed at a fixed position in space, emits a given number (${10}^7$ to ${10}^8$) of monoenergetic photons with the preset states of polarization. Those photons impinge a top surface of the sample at a given incident angle. Every photon travels a certain distance within the sample before being scattered on a sphere or cylinder. For each collision event, this distance is determined statistically using a mean free path parameter calculated from the scattering cross-sections of scatterers and their number density. The sizes, refractive indices, and number densities of sphere and cylinder scatterers as well as refractive index of isotropic host medium or ordinary and extraordinary refractive indices and spatial orientation of the optical axis of linear birefringent host medium are the input parameters of the optical model and can be adjusted to mimic the conditions of a real biological sample.¹⁴ The photon changes its polarization state and direction of propagation after each scattering event. The angles of deflection and rotation of polarization plane are calculated using the rejection method.¹⁸ A transfer matrix for scattering is determined by Mie theory for spherical scatterers or scattering matrix theory for infinite-long cylinder.¹⁹ The host medium may also be absorbing (not implemented in this study). The random walk of a photon continues within a scattering medium until it is either absorbed within the sample or moved outside the sample volume, where it can be lost or hit a detector. The first version of the software was developed on C language and run on a CPU platform.¹³ The latest version of the software is accelerated using a GPU platform.²⁰

2.2.

Optical Models of Skin Dermal Layer

To choose an appropriate optical model of a dermal layer of skin tissue cuts, one needs to account for both fibroblasts and well-aligned collagen fibers, which form the dermal equivalent of a skin model.^21,22 While light scattering on cells and fibers produces depolarization, the optical anisotropy of the dermal layer results in retardance due to form birefringence.¹ Thus, we used the monodisperse spherical scatterers in the optical model of the dermal layer to reproduce an isotropic scattering on cells. Infinitely long cylindrical scatterers were added to the optical model to simulate the effect of form birefringence due to the presence of aligned collagen fibers in dermis. The refractive indices of spherical and cylindrical scatterers ($n_s$, $n_c$) and isotropic medium ($n_m$) were set to 1.45 and 1.33, respectively.

We also explored the validity of replacing a form birefringence by an intrinsic birefringence of a uniaxial linear anisotropic host medium with in-plane optical axis, ordinary index $n_m^o=1.33$, and extraordinary index $n_m^e=1.33+\mathrm{\Delta }n$, varying parameter $\mathrm{\Delta }n$ from ${10}^{-5}$ to ${10}^{-3}$. We took the values of refractive indices as for bulk fresh tissue,¹ while noting that those values may be somewhat different for the studied fixed unstained tissue cuts. This is most probably not so important for our consideration because the refractive index of a scatterer and its size are highly correlated parameters in the Mie electromagnetic scattering problem. With our choice of refractive index values for both scatterers and host medium, the optical contrast $n_s/n_m$ (or $n_c/n_m$) is more than 1. Keeping constant the value of optical contrast, we varied the size of scatterers to reproduce the general trends in polarization and depolarization parameters in our simulations. We believe this is a reasonable assumption for performing the parametric numerical studies to reproduce the experimental trends.

In our experiments, a dermal layer of all skin model cuts demonstrated higher circular depolarization compared with the linear one ($|{\alpha }_{44}|>|{\alpha }_{22}|$, $|{\alpha }_{44}|>|{\alpha }_{33}|$). This suggests the dominance of the Rayleigh scattering regime over the Mie scattering regime²³ and justifies the use of subwavelength spherical and cylindrical scatterers in an optical model of the dermal layer. The wavelength of probing light was fixed at $0.533\mu \mathrm{m}$, so we tested spherical and cylindrical scatterers whose diameter ranged from 0.01 to $0.5\mu \mathrm{m}$. Their concentrations $c_{\mathrm{s}}$ and $c_{\mathrm{c}}$ were described by the scattering coefficients ${\mu }_s$ and ${\mu }_c$ [${\mu }_{s,c}=1/(c_{s,c}{\sigma }_{s,c})$, where ${\sigma }_{s,c}$ is the scattering cross-sections of sphere and cylinder, respectively), which were varied from 5 to $5000{\mathrm{cm}}^{-1}$. The parameter $\mathrm{\Delta }n$ for uniaxial birefringent host medium was adjusted to fit the experimental results¹² for total linear retardance parameter $R_T$ (see Sec. 3.1). The optical axis of linear birefringent host medium was always oriented parallel to the sample surface, reflecting the arrangement of collagen fibers in a dermal layer of histological cuts. The GPU acceleration allowed us to carry out the simulations in a wide range of parameters to find the best-fit values.

Histological cuts of skin tissue models of varying thicknesses (nominal values 3 to $30\mu \mathrm{m}$) were mounted on 1-mm-thick microscopy glass slides in our experiments.¹² During the calibration of the Mueller polarimetric microscope with the eigenvalue calibration method,²⁴ the measurements of air (one of the reference samples) were performed through a bare microscopy glass slide. Hence, the contribution of glass was excluded from Mueller matrices of all measured histological cuts. To model our experimental setup, we performed MC simulations in transmission configuration for the range of histological cut thicknesses defined by profilometer measurements¹¹ without adding a 1-mm-thick glass layer to our optical model. A spatially uniform light beam was normally incident onto the flat front surface of a sample. No back surface roughness was taken into account in our optical model. The simulated images of forward scattering Mueller matrix elements were spatially averaged over a central circle of $600\mu \mathrm{m}$ in diameter to reproduce the experimental conditions, and the resulting Mueller matrices were decomposed using the LMMD method.⁹

3.1.

Rotation Invariants of Logarithmic Decomposition

The set of polarization and depolarization parameters obtained from the LMMD includes the values of linear $(LB,L{B}^{\prime })$ and circular ($CB$) retardance, linear $(LD,L{D}^{\prime })$ and circular ($CD$) dichroism, and linear (${\alpha }_{22}$, ${\alpha }_{33}$) and circular $({\alpha }_{44})$ depolarization coefficients.⁹ The parameters $LB$, $CD$, and ${\alpha }_{22}$ are defined with respect to the framework of 0 deg to 90 deg (linear polarization); the parameters $L{B}^{\prime }$, $L{D}^{\prime }$, and ${\alpha }_{33}$ are defined with respect to the framework of $\pm 45\mathrm{deg}$ (linear polarization).

Neither optical activity nor circular dichroism was detected in the polarimetric measurement data for skin model histological cuts ($CB=0$, $CD=0$). The presence of well-aligned collagen fibers in a dermal layer of skin model cuts indicated the direction of the optical axis of uniaxial linear birefringent medium. In our experiments,^11,12 the orientation of histological cuts in a sample holder plane was performed manually, thus, producing an uncertainty in azimuth of the optical axis from sample to sample. Therefore, the measured values of polarization parameters $LB$, $L{B}^{\prime }$ for skin model cuts of different thicknesses were affected not only by the different optical paths of the detected photons but also by the in-plane rotation of the samples. The experimental results for histological cuts of skin models have also demonstrated the existence of anisotropy of linear depolarization (${\alpha }_{22}\ne {\alpha }_{33}$).^11,12 Both parameters ${\alpha }_{22}$ and ${\alpha }_{33}$ are not invariant under the in-plane rotation of an anisotropic sample. Therefore, we also derived a rotation-invariant parameter for the linear retardance and linear depolarization anisotropy parameters to eliminate the impact of sample orientation with respect to the laboratory coordinate frame.

The logarithm of Mueller matrix $\mathbf{M}$ is calculated from the eigenvalue decomposition as

3.2.

Choice of an Appropriate Optical Model of Dermis

Several optical models were tested to reproduce the optical effects observed in a dermal zone of samples. Some models (e.g., the SC model) were ruled out for the reasons discussed in Sec. 3.2.1. Finally, the SCB model was adopted in these studies, and the results are presented in Sec. 3.2.2.

3.2.1.

Ruled out optical models

Experimental studies of tissue phantoms²⁵ and tissue models¹² with transmission Mueller polarimetry had confirmed that retardance values calculated using LMMD depend linearly on thickness, while the depolarization parameters ${\alpha }_{22}$, ${\alpha }_{33}$, and ${\alpha }_{44}$, calculated using LMMD show quadratic dependence on layer thickness. To reproduce the experimental values of polarization and depolarization parameters of a dermal layer of histological cuts measured with the Mueller microscope in transmission configuration,¹² we tested the SC, SB, and SCB optical models.

It was already demonstrated that at normal incidence an isotropic medium with spherical scatterers does not exhibit any retardance effect.^26,27 A phase shift in the detected signal can be induced by scattering of polarized light by cylindrical scatterers as well as by light passing through a birefringent medium. Our modeling results demonstrated that the SC model with isotropic host medium fell short of reproducing the experimental values of retardance for low values of scattering coefficient ${\mu }_c$. The volume density (or concentration) of cylindrical scatterers had to increase significantly to fit the experimental trends in retardance, but with increase of parameter ${\mu }_c$ this model produced very high values of dichroism and depolarization which by far exceed the corresponding experimental values. Hence, we concluded that a uniaxial linear birefringent host medium is a necessary component of our optical model. It will increase the simulated retardance values without pushing up the dichroism and depolarization parameters of a simulated medium. Therefore, we ruled out the SC model for further consideration.

We also examined the SB optical model of the skin dermal layer. The values of $\mathrm{\Delta }n$ parameter, radius of spherical scatterer $R_s$, and scattering coefficient ${\mu }_{\mathrm{s}}$ were varied to find the best-fit to the experimental data. The SB optical model fits well the experimental values of retardance¹¹ with optimal values of $\mathrm{\Delta }n=0.009$ and $R_s=0.05\mu \mathrm{m}$ [see Fig. 1(a)]. It is worth mentioning that the optimal value of $\mathrm{\Delta }n$ for the fixed tissue cuts was found to be about 2 orders of magnitude larger compared with the values reported for the fresh biological tissue.¹ The simulation results with the SB optical model confirmed the linear dependence of total linear retardance $R_T$ (calculated from LMMD) on layer thickness. The simulated values of depolarization parameters ${\alpha }_L$, ${\alpha }_{44}$ have demonstrated a quadratic dependence on thickness within the experimental range of errors [see Figs. 1(b) and 1(c)]. No anisotropy of linear depolarization was observed with the SB model, as simulated values of ${\alpha }_{LA}=0$ for all layer thickness [Fig. 1(d)].

Fig. 1

Results of MC simulations with the SB optical model for different layer thicknesses: (a) total linear retardance $R_T$ (radians) and dimensionless depolarization parameters (b) ${\alpha }_L$, (c) ${\alpha }_{44}$, and (d) ${\alpha }_{LA}$. Simulated data are shown by open symbols corresponding to different scattering coefficients (i.e., different concentrations of spherical scatterers). Open boxes with error bars represent the experimental data. Solid lines show the results of (a), (d) linear and (b), (c) quadratic fit of the simulated data.

The experimental values of depolarization coefficients ${\alpha }_{22}$, ${\alpha }_{33}$, and ${\alpha }_{44}$ for a dermal layer of histological cuts obey the relation $|{\alpha }_{22}|<|{\alpha }_{33}|<|{\alpha }_{44}|$.¹² However, the simulations with the SB model could not reproduce the effect of anisotropy of linear depolarization ($|{\alpha }_{LA}|\ne 0$) observed experimentally [see Fig. 1(d)]. Moreover, no linear dichroism can be simulated with the SB model, while the nonzero values of the linear dichroism were measured in our experiments.¹² Therefore, we conclude that the optical SB model of the dermal layer of skin model histological cuts has to be modified to reproduce experimental trends.

3.2.2.

SCB model

We then added cylindrical scatterers to the optical model of dermis to simulate the effects of anisotropy of linear depolarization and linear dichroism. First, we used the same set of parameter values as for the SB optical model, but added a group of cylindrical scatterers of radius $R_c=0.05\mu \mathrm{m}=R_{\mathrm{s}}$, aligned along the $X$ axis. The scattering coefficient for spherical scatterers ${\mu }_s$ was fixed at $1500{\mathrm{cm}}^{-1}$ and scattering coefficient for cylindrical scatterers ${\mu }_c$ was varied from 500 to $1500{\mathrm{cm}}^{-1}$.

The results obtained with the SCB optical model after applying LMMD to the simulated Mueller matrices of layers of varying thicknesses are shown in Fig. 2. The presence of cylindrical scatterers had very limited influence on the values of linear retardance, but it had a significant impact on the values of linear dichroism and anisotropic depolarization effect. The SCB model yields the values of total linear retardance that also match well the experimental data [see Fig. 2(a)] The linear increase of total linear retardance $R_T$ and linear dichroism values on thickness are shown in Figs. 2(b) and 2(a), respectively. The nonzero intercept of linear regression curve for the experimental linear dichroism values with the $Y$ axis [Fig. 2(b)] was explained by scattering of transmitted light on the rough surface of tissue.¹² It was shown that for an anisotropic media a surface scattering does not affect the retardance values but has an impact on linear dichroism values.²⁸ The quadratic dependence of ${\alpha }_L$, ${\alpha }_{44}$, and ${\alpha }_{LA}$ on thickness is also confirmed within the experimental range of errors. Moreover, an anisotropic depolarization effect ($|{\alpha }_{LA}|\ne 0$) is well reproduced with the SCB optical model [see Figs. 2(c)–2(e)].

Fig. 2

Results of MC simulations using the SCB optical model. Dependence on layer thickness of (a) total linear dichroism $D_T$ (dashed line is a linear regression curve for the experimental data), (b) total linear retardance $R_T$, (c)–(e) depolarization parameters ${\alpha }_L$, ${\alpha }_{44}$, and ${\alpha }_{LA}$, respectively. Simulated data are shown by open symbols corresponding to different scattering coefficients ${\mu }_c$ (i.e., different concentrations of cylindrical scatterers), ${\mu }_s=1500{\mathrm{cm}}^{-1}$, $R_s=R_c=0.05\mu \mathrm{m}$. Open boxes with error bars represent the experimental data. Solid lines show the results of (a), (b) linear and (c)–(e) quadratic fit of the simulated data.

The impact of the radius of a cylindrical scatterer $R_c$ on anisotropic linear depolarization parameter ${\alpha }_{LA}$ was also studied (see Fig. 3). Our simulations show that cylinders with a smaller radius produce stronger anisotropy in linear depolarization. Hence, the presence of anisotropic linear depolarization is an indication of the scattering on small-size fibroid scatterers in the studied medium (so-called form birefringence). The values of parameter ${\alpha }_{LA}$ can be used for the estimation of the characteristic size of nonspherical scatterers.

Fig. 3

Dependence of parameter ${\alpha }_{LA}$ on layer thickness for different radii of cylindrical scatterers. The parameters of the SBC optical model are: ${\mu }_s={\mu }_c=1500{\mathrm{cm}}^{-1}$, $R_s=0.05\mu \mathrm{m}$. Open symbols correspond to the different radii of cylindrical scatterers: $R_c=0.05$, 0.5, and $2\mu \mathrm{m}$, respectively. The concentration of cylindrical scatterers $c_c$ was adjusted to keep a constant value of the scattering coefficient ${\mu }_c$. Open boxes with error bars represent the experimental data. Solid lines show the results of a parabolic fit of the simulated data.

To summarize, an appropriate optical model for a dermal layer of unstained fixed histological cuts of skin model tissue should include the subwavelength spherical scatterers and well-aligned cylindrical scatterers both distributed in a uniaxial linear birefringent medium. This model can qualitatively reproduce the thickness dependence of polarization and depolarization properties obtained from LMMD of the experimental Mueller matrices of a dermal layer of skin model histological cuts.¹²