## 1.

## Introduction

Mueller matrix polarimetry has received considerable recent attention for the characterization of a wide range of optical media for numerous practical applications in various branches of science and technology.^{1}^{–}^{3} The Mueller matrix (a $4\times 4$ matrix) represents the transfer function of an optical system in its interactions with polarized light and contains complete information about the medium polarization properties.^{1}^{–}^{3} Mueller matrix polarimetry thus offers the possibility of quantifying the polarimetry characteristics (which contain wealth of morphological, structural and functional information) of the sample. However, complexities arise when more than one polarization effects occur within the sample (the three basic medium polarization properties being diattenuation, retardance and depolarization^{2}) as these contribute in a complex interrelated way to the Mueller matrix elements, masking potentially interesting polarization metrics and hindering their unique interpretation.^{4} Quantitative polarimetry is further compromised due to the presence of strong multiple scattering effects in optically thick (turbid) media such as biological tissues.^{4}^{,}^{5} Multiple scattering within a tissue not only causes extensive depolarization, but also alters the polarization state of the residual polarization-preserving signal in a complex fashion.^{4}^{,}^{6} Each of the intrinsic polarization properties, if separately extracted from the ‘lumped’ system Mueller matrix, can potentially serve as a useful biological metric. For example, the anisotropic organized nature of many tissues stemming from their fibrous structure leads to linear retardance (or linear birefringence) effect. Changes in this anisotropy resulting from disease progression or treatment response alter the linear retardance properties, making this a potentially sensitive probe of tissue status.^{7}^{,}^{8} Similarly, measurement and quantification of circular retardance (or optical rotation) may offer an attractive approach for noninvasive monitoring of tissue glucose levels.^{9}

Various kinds of Mueller matrix decompositions have therefore been proposed in the recent past for inverse polarimetry analysis.^{10}^{–}^{13} Among these, the product decomposition proposed by Lu and Chipman has been the most widely explored method.^{10} This method consists of decomposing a given Mueller matrix into sequential product of the three ‘basis’ matrices describing the three basic medium polarization properties, depolarization, retardance (linear and circular) and diattenuation (linear and circular). Since the product decompositions^{10}^{–}^{13} model the resultant polarization processes as a preferred sequential product of the ‘basis’ matrices describing the constituent polarization effects, they are potentially ambiguous with respect to the order of the basis matrices, due to the noncommuting nature of matrix multiplication. Moreover, in complex systems like tissues, no unique order (or sequence) can be assigned to the polarimetry effects; rather, these are exhibited simultaneously. The decomposition-derived polarization parameters are thus expected to be influenced by the choice of the ordering of the basis matrices. Although our recent studies have shown that under certain conditions, the ambiguity in the decomposition-derived polarization parameters can be minimal,^{14} a more generally applicable inverse analysis model is certainly required for polarized light assessment of media exhibiting simultaneous polarization effects (such as tissues). In order to address this issue, a more general kind of decomposition based on differential matrix formalism of the Mueller calculus, has recently been developed.^{15}^{,}^{16} This formalism, historically stemming from Jones $N$-matrix formalism,^{17} was originally proposed by R. M. A Azzam for nondepolarizing optical media;^{18} it has been extended and generalized to incorporate depolarization effects.^{15} Within the differential matrix formalism, all elementary polarization and depolarization properties of the medium are contained in a single differential matrix representing them simultaneously. Thus, this approach is expected to be more suitable for Mueller matrix polarimetry analysis of tissues, exhibiting simultaneous polarization effects. Implementation of these approaches for polarized light assessment of complex systems like tissues would, however, require comprehensive validation. Specifically, it is necessary to perform a comparative evaluation of the polarization parameters derived via the differential matrix decomposition and the product decompositions (specifically, the widely used Lu-Chipman polar decomposition), and to establish a link between them. We have therefore investigated this issue by performing a detailed comparative evaluation of the polar decomposition and the differential matrix decomposition of Mueller matrices from complex tissue-like turbid media exhibiting simultaneous scattering and polarization effects. We have paid particular attention on linear retardance ($\delta $), optical rotation ($\mathrm{\Psi}$) and depolarization ($\mathrm{\Delta}$) parameters, because these are the most prominent tissue polarimetry effects, which also hold promise as useful biological metrics for applications involving tissue characterization/diagnosis.^{4}

In particular, we have studied the effect of simultaneous or sequential exhibition of linear retardance and optical rotation polarization events on the derived values for the $\delta $ and $\mathrm{\Psi}$ parameters using either the retarder polar decomposition (which assumes preferred sequential occurrence of these effects) or the differential matrix decomposition. In order to establish self-consistency, we have also derived analytical relationships between the polarization parameters ($\delta $, $\mathrm{\Psi}$) extracted from both the retarder polar decomposition and the differential matrix decomposition for either simultaneous or sequential occurrence of the linear retardance and optical rotation effects. The self-consistency of both decompositions is validated on experimental Mueller matrices recorded from tissue-simulating phantoms (whose polarization properties are controlled, known *a-priori*, and exhibited simultaneously) of increasing biological complexity. Additional theoretical validation tests were performed on Monte Carlo (MC)-generated Mueller matrices from analogous turbid media exhibiting simultaneous depolarization ($\mathrm{\Delta}$), linear retardance ($\delta $) and optical rotation ($\mathrm{\Psi}$) effects. Finally, after successful evaluation, the potential advantages of the differential matrix decomposition over the previously used polar decomposition was explored for monitoring of myocardial tissue regeneration following stem cell therapy (through quantification of linear retardance obtained from the measured Mueller matrices). The results of these investigations are reported here.

The paper is organized as follows. In Sec. 2, we briefly summarize the mathematical methodology for polar decomposition and differential matrix decomposition of Mueller matrix for extraction of the constituent polarimetry characteristics. The derived conversion relations for the polarization parameters extracted from the two different decomposition formalisms are presented in this section. The analytical extension for ‘self-consistent’ extraction of linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ from both the decomposition formalisms for either simultaneous or sequential occurrence of the polarization events is also presented in this section. Section 3 briefly describes our experimental turbid polarimetry systems, optical phantoms and forward polarization sensitive Monte Carlo (MC) simulation model. The results of the comparative evaluation of polar decomposition and differential matrix decomposition and validation of the self-consistency of our extended formalisms on experimental Mueller matrices from phantoms and MC-generated Mueller matrices are presented in Sec. 4. Selected experimental results on the use of differential matrix decomposition from myocardial tissue samples are then presented, interpreted and discussed. The paper concludes with a discussion on the potential advantages of the differential matrix decomposition in quantitative polarimetric tissue characterization/diagnosis.

## 2.

## Theory

## 2.1.

### Polar Decomposition of Mueller Matrix

Polar decomposition method consists of decomposing a given Mueller matrix $M$ into the sequential product of three ‘basis’ matrices,

with $\Leftarrow $ symbol used to signify the decomposition process. Here, the matrix ${M}_{\mathrm{\Delta}}$ describes the depolarizing effects of the medium, ${M}_{R}$ accounts for the effects of linear and circular retardance (or optical rotation), and ${M}_{D}$ includes the effects of linear and circular diattenuation. This product decomposition with multiplication order of Eq. (1) was proposed by Lu and Chipman and is accordingly known as Lu-Chipman decomposition.^{10}Once calculated, these constituent basis matrices are further analyzed to derive individual polarization medium properties. Specifically, diattenuation ($D$), depolarization coefficient ($\mathrm{\Delta}$), linear retardance ($\delta $), and circular retardance/optical rotation ($\mathrm{\Psi}$) ($\mathrm{circular}\text{\hspace{0.17em}}\mathrm{retardance}=2\text{\hspace{0.17em}}\times \text{\hspace{0.17em}}\mathrm{optical}\text{\hspace{0.17em}}\mathrm{rotation}$), can be determined from the decomposed basis matrices as shown below.

The magnitude of diattenuation ($D$) can be determined from the diattenuation matrix ${M}_{D}$ as^{10}

## Eq. (2)

$$D=\frac{1}{{M}_{D}(1,1)}\sqrt{{M}_{D}{(1,2)}^{2}+{M}_{D}{(1,3)}^{2}+{M}_{D}{(1,4)}^{2}}\mathrm{.}$$^{10}

## Eq. (3)

$$\mathrm{\Delta}=1-\frac{|\mathrm{tr}({M}_{\mathrm{\Delta}})-1|}{3},\phantom{\rule[-0.0ex]{2em}{0.0ex}}0\le \mathrm{\Delta}\le 1.$$Assuming sequential effects of linear retardance and optical rotation the retarder matrix ${M}_{R}$ can be further decomposed into what can be termed “retarder polar decomposition,”

where ${M}_{C}$ and ${M}_{L}$ (${M}_{L}$) are another two retarder matrices exhibiting only circular (i.e., optical rotation) and linear retardance, respectively.^{6}

^{,}

^{19}By inserting the standard Mueller matrix forms for linear and circular retarders

^{2}

^{,}

^{19}into [Eq. (4)], the magnitudes of linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ can be calculated directly from the elements of the retarder matrix ${M}_{R}$ as

^{6}

^{,}

^{19}

## Eq. (5)

$$\delta ={\mathrm{cos}}^{-1}\{\sqrt{[{\mathbf{M}}_{\mathbf{R}}(2,2)+{\mathbf{M}}_{\mathbf{R}}{(3,3)]}^{2}+{[{\mathbf{M}}_{\mathbf{R}}(3,2)-{\mathbf{M}}_{\mathbf{R}}(2,3)]}^{2}}-1\}$$## Eq. (6)

$$\mathrm{\Psi}={\mathrm{tan}}^{-1}\left[\frac{{\mathbf{M}}_{\mathbf{R}}(3,2)-{\mathbf{M}}_{\mathbf{R}}(2,3)}{{\mathbf{M}}_{\mathbf{R}}(2,2)+{\mathbf{M}}_{\mathbf{R}}(3,3)}\right].$$Note that the magnitudes $\delta $ and $\mathrm{\Psi}$ of linear retardance and optical rotation do not depend on their order of appearance in the sequence although the matrix product is generally noncommutative. The total retardance $R$ (a parameter that represents the combined effect of linear and circular retardance) can also be determined from ${M}_{R}$^{10}^{,}^{19}

An algorithm for decomposing Mueller matrices with multiplication order of the basis matrices reverse to that of Eq. (1) has also been developed (reverse decomposition^{11}), and it has been shown that all other possible decompositions (a total number of six depending upon the order of multiplication) can be obtained from these two decompositions by using similarity transformations.^{11} Note, although, that the processes for decomposition using the order of Eq. (1) or its reverse order are qualitatively similar; there do exist subtle differences in the construction of the basis matrices and the derived polarization parameters.^{11} Never-the-less all families of product decompositions model the polarized light–matter interaction as a sequential product (of unique order or sequence) of the basis matrices describing the constituent polarization effects, and are thus potentially ambiguous with respect to the multiplication order of these matrices.

## 2.2.

### Differential Matrix Decomposition

Within the differential matrix formalism, the elementary properties of the medium are contained in the differential matrix $m$ which is related to the Mueller matrix $M$ and its spatial derivative along the propagation of light as

The physical picture behind Eq. (8) assumes that the medium is laterally homogeneous (at least within the spot dimension of the probing light) and that the polarization and depolarization effects take place simultaneously.^{18}The general form of the differential matrix $m$ for depolarizing anisotropic media in terms of basic properties of the medium is given by Ref. 15

## Eq. (9)

$$m=\left(\begin{array}{cccc}\alpha & \beta & \gamma & \partial \\ {\beta}^{\prime}& {\alpha}_{1}& \mu & -\nu \\ {\gamma}^{\prime}& -{\mu}^{\prime}& {\alpha}_{2}& \eta \\ {\partial}^{\prime}& {\nu}^{\prime}& -{\eta}^{\prime}& {\alpha}_{3}\end{array}\right).$$Here, $\beta $ and $\gamma $ is the linear dichroism along $xy$ (horizontal- vertical) and $\pm 45\xb0$ laboratory axes, respectively; $\partial $ is the circular dichroism; $\eta $ and $\nu $ are linear birefringence along $xy$ and $\pm 45\xb0$ axes, respectively; $\mu $ is the circular birefringence; $\alpha $ and ${\alpha}_{1}$, ${\alpha}_{2}$, ${\alpha}_{3}$, respectively are the isotropic absorption and the anisotropic absorptions (or depolarizations) along $xy$ $\pm 45\xb0$ and circular axes, respectively. Note that in absence of depolarization, $\alpha =\phantom{\rule{0ex}{0ex}}{\alpha}_{1}={\alpha}_{2}={\alpha}_{3}(=0)$. Moreover, in such case, $({\beta}^{\prime},\text{\hspace{0.17em}}{\gamma}^{\prime},\text{\hspace{0.17em}}{\partial}^{\prime},\text{\hspace{0.17em}}{\eta}^{\prime},\text{\hspace{0.17em}}{\nu}^{\prime},\text{\hspace{0.17em}}{\mu}^{\prime})=\phantom{\rule{0ex}{0ex}}(\beta ,\gamma ,\partial ,\eta ,\nu ,\mu )$. In the general case of a depolarizing medium, the six elementary polarization properties are in fact, the mean values of the corresponding pairs $(\beta ,{\beta}^{\prime})$, $(\gamma ,{\gamma}^{\prime})$, $(\partial ,{\partial}^{\prime})$, $(\eta ,{\eta}^{\prime})$, $(\nu ,{v}^{\prime})$ and $(\mu ,{\mu}^{\prime})$.^{15}

If the medium properties (both polarizing and depolarizing) are uniformly distributed along the propagation direction $Z$, i.e., the differential matrix $m$ is $Z$-independent, Eq. (8) can be integrated to yield

with $L$ being the matrix logarithm of $M$ and $l$ being the optical pathlength, i.e., the equivalent length light would travel in the medium if optical properties and depolarization occurred truly simultaneously and were perfectly uniformly distributed longitudinally. If these conditions are met in practice, then $l$ can be identified either with the sample thickness in case of forward scattering experiment or with the penetration depth in backscattering geometry. (It should be noted for completeness that Mueller matrix logarithm decomposition, Eq. (10), is formally equivalent^{20}to the matrix roots decomposition that has been introduced by Chipman,

^{21}thus initiating the differential (simultaneously occurring properties) decomposition approach, in contrast to the product, (sequentially occurring properties) one.

The elementary properties of the medium can be determined by constructing the Lorentz antisymmetric, ${L}_{m}$, and symmetric, ${L}_{u}$, components of $L(L={L}_{m}+{L}_{u})$ as follows^{15}

## Eq. (11)

$${L}_{m}=\frac{1}{2}(L-G{L}^{T}G)\phantom{\rule[-0.0ex]{1em}{0.0ex}}\text{and}\phantom{\rule[-0.0ex]{1em}{0.0ex}}{L}_{u}=\frac{1}{2}(L+G{L}^{T}G),$$For a depolarizing medium, the off-diagonal elements of ${L}_{m}$ represent the mean values of the six elementary properties [as per Eq. (9)] accumulated over pathlength $l$, while the off-diagonal elements of ${L}_{u}$ express their respective uncertainties. Further, the main diagonal elements of the matrix ${L}_{u}$ (${\alpha}_{1}$, ${\alpha}_{2}$ and ${\alpha}_{3}$) represent the depolarization coefficients (after the subtraction of the isotropic absorption $\alpha $ from the diagonal) along the $xy$, $\pm 45\xb0$ and circular axes.^{15} The accumulated polarization parameters (intrinsic properties scaled by pathlength $l$) can be determined from the elements of the ${L}_{m}$ and ${L}_{u}$ matrices as

## Eq. (12)

$$\mathrm{Net}\text{\hspace{0.17em}}\mathrm{dichroism}\text{\hspace{0.17em}}d=\sqrt{{[{L}_{m}(1,2)]}^{2}+{[{L}_{m}(1,3)]}^{2}+{[{L}_{m}(1,4)]}^{2}}$$## Eq. (13)

$$\mathrm{Net}\text{\hspace{0.17em}}\mathrm{depolarization}\text{\hspace{0.17em}}{\mathrm{\Delta}}_{m}=\frac{1}{3}|\alpha 1+\alpha 2+\alpha 3|$$## Eq. (14)

$$\mathrm{Linear}\text{\hspace{0.17em}}\mathrm{retardance}\text{\hspace{0.17em}}{\delta}^{\mathrm{log}\text{-}M}=\sqrt{{[{L}_{m}(2,4)]}^{2}+{[{L}_{m}(3,4)]}^{2}}$$## Eq. (15)

$$\mathrm{Optical}\text{\hspace{0.17em}}\mathrm{rotation}\text{\hspace{0.17em}}{\mathrm{\Psi}}^{\mathrm{log}\text{-}M}=\frac{1}{2}{L}_{m}(2,3)$$## Eq. (16)

$${\theta}^{\mathrm{log}\text{-}M}=\frac{1}{2}{\mathrm{tan}}^{-1}\frac{\{{L}_{m}(2,4)\}}{\{{L}_{m}(3,4)\}}.$$In this convention, the total retardance can be determined as

## Eq. (17)

$${R}^{\mathrm{log}\text{-}M}=\sqrt{{({\delta}^{\mathrm{log}\text{-}M})}^{2}+{4({\mathrm{\Psi}}^{\mathrm{log}\text{-}M})}^{2}}.$$## 2.3.

### Conversion Relations for the Polarization Parameters

An important issue while relating the medium polarimetry characteristics defined in the polar decomposition formalism [Eqs. (2)–(7)] and those in the differential matrix formalism (the accumulated polarization parameters defined in Eqs. (12)–(17), is the conversion from one set of parameters to the other. Note that the Mueller matrix corresponding to each accumulated polarization parameter (elements of the differential matrix $m$ scaled by pathlength $l$) can be obtained from the solution of Eq. (8) and the eigenvalues of $m$.^{18}^{,}^{22} From these derived Mueller matrices (for each of the polarization effects) and the corresponding known forms of Mueller matrices in the polar decomposition formalism,^{2}^{,}^{4} one can derive the following relations.

## Eq. (18)

$${D}^{\mathrm{log}\text{-}M}=\mathrm{tanh}(d)=\mathrm{tanh}\phantom{\rule{0ex}{0ex}}\left\{\sqrt{{[{L}_{m}(1,2)]}^{2}+{[{L}_{m}(1,3)]}^{2}+{[{L}_{m}(1,4)]}^{2}}\right\}\leftrightarrow {D}^{\text{polar}},$$## Eq. (19)

$${\mathrm{\Delta}}^{\mathrm{log}\text{-}M}=1-\frac{1}{3}({e}^{-{\alpha}_{1}}+{e}^{-{\alpha}_{2}}+{e}^{-{\alpha}_{3}})\leftrightarrow {\mathrm{\Delta}}^{\text{polar}},$$## 2.4.

### Analytical Extension for “Self-Consistent” Extraction of Linear Retardance $\delta $ and Optical Rotation $\mathrm{\Psi}$

Having discussed the normalization criteria for diattenuation and depolarization effects, we now turn our attention to the important issue of self-consistent determination of linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ using differential matrix decomposition and the retarder polar decomposition [(Eq. 4)], for either simultaneous or sequential effects. The form of a Mueller matrix exhibiting simultaneous linear retardance and circular retardance can be derived, from the differential matrix (${L}_{R}$) representing a retardance effect alone, as ${M}_{R}=\mathrm{exp}({L}_{R})$ yielding^{18}^{,}^{22}

## Eq. (20)

$${M}_{R}=\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& \mathrm{cos}\text{\hspace{0.17em}}R+{\eta}^{2}\rho & -\eta \nu \rho +\mu \sigma & -\eta \mu \rho -\nu \sigma \\ 0& -\eta \nu \rho -\mu \sigma & \mathrm{cos}\text{\hspace{0.17em}}R+{\nu}^{2}\rho & \nu \mu \rho -\eta \sigma \\ 0& -\eta \mu \rho +\nu \sigma & \nu \mu \rho +\eta \sigma & \mathrm{cos}\text{\hspace{0.17em}}R+{\mu}^{2}\rho \end{array}\right),$$^{18}

## Eq. (21)

$${L}_{R}=\left(\begin{array}{cccc}0& 0& 0& 0\\ 0& 0& \mu & -\nu \\ 0& -\mu & 0& \eta \\ 0& \nu & -\eta & 0\end{array}\right).$$^{2}

^{,}

^{18}However, as is obvious from Eqs. (7), (17), and (20), the value of total retardance ($R$) would remain the same. This discrepancy in determination of $\delta $ and $\mathrm{\Psi}$ parameters can be rectified, if these are directly estimated from the elements of the retardance vector introduced below.

The elements of the retardance vector $\overrightarrow{R}=R{({r}_{1},{r}_{2},{r}_{3})}^{T}$ can be derived from the $3\times 3$ submatrix ${m}_{R}$ of the $4\times 4$ retardance matrix ${M}_{R}$ obtained via the polar decomposition (first column) as^{10}

## Eq. (22)

$${r}_{i}=\frac{1}{2\text{\hspace{0.17em}}\mathrm{sin}\text{\hspace{0.17em}}R}\sum _{}^{}{\epsilon}_{ijk}{({m}_{R})}_{jk},$$Using Eqs (20) and (22), the elements of the retardance vector can be determined as

This shows the equivalence between the two approaches when the retardance effects (linear retardance and optical rotation) occur simultaneously. The values for net linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ can thus be determined with the help of Eqs. (14) and (15) as

rather than from Eqs. (5) and (6) based on the retarder polar decomposition.The validity of the determination of $\delta $ and $\mathrm{\Psi}$ parameters using Eqs. (24) and (25) in combination with polar decomposition has been tested on experimental Mueller matrices from optical phantoms exhibiting simultaneous linear retardance, optical rotation and depolarization effects and on analogous MC-generated Mueller matrices. These results are presented subsequently.

Conversely, when the medium exhibits sequential linear and circular retardance effects the use of differential matrix formalism—Eqs. (14) and (15)–would yield erroneous estimation of $\delta $ and $\mathrm{\Psi}$ parameters for the same reason described above. In fact, the differential matrix formalism can also be ‘corrected’ to yield the ‘true’ values for $\delta $ and $\mathrm{\Psi}$, even for sequential occurrence of these effects. Here, the correction process is based on determination of these parameters from the eigenvalues and the eigenvectors of the differential matrix for the retardance alone, ${L}_{R}$. [${L}_{R}$ is derived from the Lorentz antisymmetric component ${L}_{m}$ of the matrix logarithm $L$ of the experimental Mueller matrix; see Eqs. (10) and (11)] The eigenvalues of the matrix ${L}_{R}$ are (0, 0, $iR$, $-iR$), whereas one of the eigenvectors corresponding to the doubly degenerate zero eigenvalue turns out to be the vector ${(1,{r}_{1},{r}_{2},{r}_{3})}^{T}$ i.e., its last three components are proportional to the retardance vector $\overrightarrow{R}$. Using the known forms of linear retarder and optical rotation matrices^{2}^{,}^{4} and assuming sequential effects (of either order), the analytical formulations of the components of $\overrightarrow{R}$ can be derived. Following a series of algebraic manipulations, sequentially occurring linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ can be related to the components of the retardance vector (i.e., to the eigenvector of ${L}_{R}$) as

## Eq. (26)

$$\delta (\text{sequen}\text{tial})=2\text{\hspace{0.17em}}{\mathrm{sin}}^{-1}\sqrt{({r}_{1}^{2}+{r}_{2}^{2})\frac{1\pm \sqrt{1-{\mathrm{sin}}^{2}R}}{2}}$$## Eq. (27)

$$\mathrm{\Psi}(\text{sequen}\text{tial})=\frac{1}{2}\text{\hspace{0.17em}}{\mathrm{sin}}^{-1}\frac{2{r}_{3}\mathrm{sin}\text{\hspace{0.17em}}R}{1+\mathrm{cos}\text{\hspace{0.17em}}\delta},$$The analytical relationships provided above allows one to self consistently extract the values of linear retardance $\delta $ and optical rotation $\mathrm{\Psi}$ from both decomposition formalisms for either simultaneous or sequential occurrence of these polarization events. While we have focused our attention on the validation of the two formalisms for media exhibiting simultaneous polarization effects (because this is the situation typically encountered in tissues^{4}), the validity of the determination of $\delta $ and $\mathrm{\Psi}$ parameters using Eqs. (26) and (27) with differential matrix decomposition for the sequential occurrence of linear and circular retardance effects, is also warranted. It is also worth mentioning here that, in the differential matrix formalism, the total retardance (or birefringence) is represented in a basis using two linear ($xy$, horizontal–vertical, birefringence and $\pm 45\text{\hspace{0.17em}}\mathrm{deg}$ birefringence) and one circular axis^{18}^{,}^{22} whereas, in the polar decomposition formalism, retardance is generally represented by linear retardance $\delta $, its orientation angle $\theta $ and circular retardance ($2\mathrm{\Psi}$).^{10} By establishing a link between these two formalisms, the parameters of the differential matrix formalism are thus transformed to the latter representation basis, as is apparent from Eqs. (12)–(17).

## 3.

## Experimental Methods

## 3.1.

### Experimental Polarimetry Systems

The experimental Mueller matrix polarimetry system employed in our study consists of both point measurement and imaging systems. The point polarimetry system employs polarization modulation with synchronous detection and thus provides a high signal-to-noise ratio (SNR) polarimetric signal. The imaging system employs $dc$ measurements involving sequential measurements with different combinations of source polarizers and detection analyzers, albeit with lower SNR than the synchronous $ac$ detection schemes. The imaging system, although suffering from a lower SNR, allows for spatial mapping of the polarization parameters of a sample.

The details of our point measurement polarimetry system have been reported previously.^{19}^{,}^{23} Briefly, our turbid-polarimetry system employs polarization modulation using a photoelastic modulator (Hinds Instruments IS-90) and synchronous lock-in detection (Stanford Research Systems, SR830), allowing sensitive low-noise measurements of the Stokes vector of light interacting with a turbid sample in various geometries.^{4}^{,}^{19}^{,}^{23} In our experimental embodiment, we employ polarization modulation on the sample-emerging light by placing the polarization modulator (along with a removable quarter-wave plate and a linear analyzer) between the sample and the detector.^{4}^{,}^{19} By cycling the polarization of the incident beam (He-Ne laser, 15 mW, $\lambda =632.8\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{nm}$) using a linear polarizer and quarter wave-plate combination, and measuring the output Stokes vectors ($I$, $Q$, $U$, $V$), the Mueller matrix of the sample can be constructed.^{19} The input polarization was sequentially cycled between four states (linear polarization at 0, 45, and 90 deg, and right circular polarization) and the output Stokes vector for each respective input states were measured. The Mueller matrix was constructed from the measured set of Stokes vectors as^{19}

## Eq. (28)

$$\mathbf{M}(i,j)=\left[\begin{array}{cccc}\frac{1}{2}({I}_{H}+{I}_{V})& \frac{1}{2}({I}_{H}-{I}_{V})& {I}_{p}-\mathbf{M}(1,1)& {I}_{R}-\mathbf{M}(1,1)\\ \frac{1}{2}({Q}_{H}+{Q}_{V})& \frac{1}{2}({Q}_{H}-{Q}_{V})& {Q}_{p}-\mathbf{M}(2,1)& {Q}_{R}-\mathbf{M}(2,1)\\ \frac{1}{2}({U}_{H}+{U}_{V})& \frac{1}{2}({U}_{H}-{U}_{V})& {U}_{p}-\mathbf{M}(3,1)& {U}_{R}-\mathbf{M}(3,1)\\ \frac{1}{2}({V}_{H}+{V}_{V})& \frac{1}{2}({V}_{H}-{V}_{V})& {V}_{p}-\mathbf{M}(4,1)& {V}_{R}-\mathbf{M}(4,1)\end{array}\right].$$Here, the four input states are denoted with the subscripts $H$ (0 deg), $P$ (45 deg), $V$ (90 deg), and $R$ (right circular). The indices $i$, $j=1$, 2, 3, 4 denote rows and columns, respectively.

In the imaging polarimetry system, light from a 635-nm diode laser (ThorLabs) was used as the excitation source. The polarization state of the incident light is controlled using a removable quarter-wave plate and a linear polarizer (forming a polarization state generator, PSG). As the beam exits the sample, a removable quarter-wave plate and linear analyzer (polarization state analyzer, PSA) select a single polarization state of the outgoing light, which is finally detected with a CCD camera (CoolSNAP K4, Photometrics, $2048\text{\hspace{0.17em}}\times \text{\hspace{0.17em}}2048\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{pixels}$, $7.4\text{\hspace{0.17em}}\times \text{\hspace{0.17em}}7.4\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mu \text{m}}^{2}\text{\hspace{0.17em}}\phantom{\rule{0ex}{0ex}}\mathrm{pixel}$ size). Four input states (horizontal 0 deg, vertical 90 deg, $+45$, right-circularly polarized) and six output states (horizontal, vertical, $\pm 40$ deg, right- and left circularly polarized light) were recorded for a total of 24 combinations of measurements. The output Stokes parameters for each input state were measured from the detected intensities as^{24}

## Eq. (29)

$$I={I}_{{0}^{\xb0}}+{I}_{{90}^{\xb0}}\phantom{\rule{0ex}{0ex}}Q={I}_{{0}^{\xb0}}-{I}_{{90}^{\xb0}}\phantom{\rule{0ex}{0ex}}U={I}_{+45\xb0}-{I}_{-45\xb0}\phantom{\rule{0ex}{0ex}}V={I}_{R}-{I}_{L},$$Solid polyacrylamide phantoms were used as tissue simulating phantoms in our study. These solid optical phantoms were developed using polyacrylamide as a base medium, with sucrose-induced optical activity, polystyrene microspheres-induced scattering and mechanical stretching to cause linear birefringence.^{19} These phantom systems mimic the complexity of biological tissues, in that they exhibit simultaneous effects of linear birefringence, optical activity and depolarization due to multiple scattering. Transmission Mueller matrices from these phantoms were recorded using the point measurement polarimetry system.^{19}^{,}^{23} To investigate the utility of linear retardance measurements and quantification (using differential matrix decomposition) for monitoring regenerative treatments, samples from a rat model of myocardial infarction and regeneration were used.^{23}^{,}^{24} All animal studies were carried out under institutional approval from the University Health Network, Toronto, Canada. Myocardial infarctions were induced in Lewis rats through coronary artery ligation. A group of animals was used as a control and mesenchymal stem cells were administered to a treatment groups by intra-myocardial injection into the site of infarction, two weeks post ligation. Both control and treatment groups were sacrificed four weeks post infarction (two weeks post stem-cell injections for the two treatment groups), at which time the hearts were removed, fixed in 10% formalin, and sectioned to one-mm-thick axial sections. Mueller matrix images were recorded in the transmission geometry, from these 1-mm-thick *ex vivo* tissue samples.^{24}

## 3.2.

### Forward Polarization Sensitive Monte Carlo Simulation Model

We have used the polarization sensitive Monte Carlo (PSMC) model developed and validated previously by our group,^{25} to generate Mueller matrices from the simulated turbid media exhibiting simultaneous linear birefringence, optical activity and depolarization effects (analogous to the optical phantoms). An important aspect of our PSMC model is that it has been extended^{25} to simulate the simultaneous effects of linear birefringence and optical activity in the presence of multiple scattering, through the use of Jones’ $N$-matrix formalism.^{17} This model thus mimics the experimental optical phantoms in that the polarization effects are exhibited in a simultaneous fashion over the path length of the simulated medium. The specifics of the model can be found in references.^{25}^{,}^{26}

Simulations were run for a set of input optical parameters of the scattering medium (scattering coefficient ${\mu}_{s}$, absorption coefficient ${\mu}_{a}$) exhibiting simultaneous linear and circular birefringence effects. In the simulations, circular and linear birefringence was modeled through the optical activity $\chi $ in degrees per centimeter, and through the anisotropy in refractive indices ($\mathrm{\Delta}n$), respectively. Here, $\mathrm{\Delta}n(={n}_{e}-{n}_{0})$ is the difference in refractive index along the extraordinary axis (${n}_{e}$) and the ordinary axis (${n}_{o}$). For simplicity, it was assumed that the direction of the extraordinary axis and the value for $\mathrm{\Delta}n$ is constant throughout the scattering medium. In each simulation, ${n}_{e}$ and ${n}_{o}$ were given as input parameters and a specific direction of the extraordinary axis was chosen. The Mueller matrices were generated for a slab of scattering medium having varying optical properties (${\mu}_{s}$, $\mathrm{\Delta}n$ and $\chi $), for light exiting the medium through any user-selected direction.^{26} The absorption was taken to be small (${\mu}_{a}=0.001\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{cm}}^{-1}$) for all simulations. The photon collection geometry was chosen to have a detection area of $1\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{mm}}^{2}$ and an acceptance angle of 20 deg.

## 4.

## Results and Discussion

Figure 1 shows the results of the decomposition analysis on the $4\times 4$ Mueller matrix experimentally recorded in the forward (transmission) detection geometry from a birefringent ($\mathrm{extension}=\phantom{\rule{0ex}{0ex}}4\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$ for strain applied along the vertical direction), chiral (optical activity $\chi =1.96\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg}/\mathrm{cm}$, corresponding to 1 M concentration of sucrose), nonscattering phantom (${\mu}_{s}=0\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{cm}}^{-1}$).^{19} The Mueller matrix logarithms ${L}_{m}$ and ${L}_{u}$ [derived via Eq. (10)] are shown and the estimated polarization parameters using both the differential matrix and polar decomposition approaches are also listed in the accompanying Table 1.

## Table 1

The values for the polarization parameters extracted using Eqs. (14)–(19) (noted as Differential matrix, 2nd column). The values for the parameters extracted via retarder polar decomposition and its associated algebra [Eqs. (2)–(7)], which assumes sequential occurrence of linear retardance and optical rotation (retarder polar decomposition sequential, third column). The R, δ and Ψ-parameters extracted using Eqs. (23)–(25), for simultaneous δ and Ψ effects (retarder polar decomposition simultaneous [Eq. (4)]. The uncertainties of the parameters were calculated from the matrix lu. The uncertainty of the parameter Ψlog-M was very small (∼10−5) and hence not given.

Parameters | Differential matrix | Retarder polar decomposition (Sequential) | Retarder polar decomposition (Simultaneous) |
---|---|---|---|

$\mathrm{\Delta}$ | 0.005 | 0.005 | 0.005 |

$\delta $ (radian) | $1.345\pm 0.003$ | 1.345 | 1.345 |

$\mathrm{\Psi}$ (radian) | 0.026 | 0.031 | 0.026 |

R (radian) | 1.346 | 1.346 | 1.346 |

$D$ | $0.023\pm 0.004$ | 0.024 | 0.024 |

Since the phantom was not supposed to exhibit any intrinsic diattenuation and depolarization effects, the derived values for diattenuation ($D$) and the depolarization coefficient ($\mathrm{\Delta}$) are also accordingly quite small. Interestingly, the derived value for the optical rotation $\mathrm{\Psi}$ is higher for the retarder polar decomposition when estimated using Eq. (6) (assuming sequential occurrence of linear retardance and optical rotation) from the decomposed retardance matrix ${M}_{R}$. As noted previously, this arises because these two effects are exhibited simultaneously in the phantoms and thus the resulting retarder Mueller matrix is not of the same form as assumed for sequentially occurring effects, thus accounting for this discrepancy. Since, in the phantom, the linear retardance effect is much stronger than the circular retardance one, the discrepancy is observed to be more pronounced for the derived optical rotation $\mathrm{\Psi}$, whereas the $\delta $-parameter remains relatively unaffected. Never-the-less, this overestimation is observed to be corrected when $\delta $ and $\mathrm{\Psi}$ parameters are directly estimated from the elements of the retardance vector {shown in the [Eq. 4)]} via Eqs. (24) and (25), thus establishing self-consistency between the two formalisms. One can also notice that for this nonscattering phantom, the differential matrix-derived uncertainties in the parameters ($\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$, $\mathrm{\Delta}{\mathrm{\Psi}}^{\mathrm{log}\text{-}M}$ obtained from the matrix logarithm ${L}_{u}$) are quite small and physically insignificant.

In Fig. 2, we show the experimentally recorded Mueller matrix (in the forward scattering geometry) and the corresponding ${L}_{m}$ and ${L}_{u}$ matrices for a birefringent ($\mathrm{extension}=4\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{mm}$), chiral ($\mathrm{concentration}\text{\hspace{0.17em}}\mathrm{of}\text{\hspace{0.17em}}\mathrm{sucrose}=1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{M}$, $\chi =1.96\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg}.{\mathrm{cm}}^{-1}$), turbid (${\mu}_{s}=30\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{cm}}^{-1}$, anisotropy parameter $g=0.95$) phantom.^{19} The determined values for the polarization parameters are listed in Table 2. In absence of any intrinsic diattenuation property, once again the values for $D$ estimated from both the decomposition formalisms are quite low and comparable. The estimated values for the depolarization coefficient $\mathrm{\Delta}$ from the two formalisms are self-consistent and are also in excellent agreement with the controlled input. This confirms the validity of the conversion relation [Eq. (19)] between the depolarization coefficients of the two formalisms (note that the un-normalized value of depolarization coefficient, ${\mathrm{\Delta}}_{m}=0.24$ of Eq. (13) deviates from the corresponding estimate from polar decomposition). The accumulation of optical rotation due to an increase in optical pathlength engendered by multiple scattering is apparent.^{19} In fact, the accumulated value of $\mathrm{\Psi}$ matches well the expected one, $\mathrm{\Psi}={\mathrm{\Psi}}_{0}\times \langle L\rangle $. The significantly lower value of the linear retardance $\delta $ as compared to that expected for accumulated pathlength has been shown previously to arise due to the reduction of the apparent $\delta $ as a consequence of curved zig-zag propagation paths for a group of photon population a component of which is directed along the linear birefringence axis.^{19} Once again, the estimation of $\mathrm{\Psi}$-parameter using the retarder polar decomposition which assumes sequential occurrence of linear and circular retardance leads to an overestimation; self-consistency is established by using retarder polar decomposition for simultaneous occurrence of linear retardance and optical rotation effects [Eqs. (24) and (25)]. Interestingly, the differential matrix-derived uncertainty in linear retardance ($\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$) is higher for this scattering phantom as compared to the nonscattering phantom. This arises both due to the stronger depolarization present in this phantom and is also due to the random orientation effects of linear birefringence (as observed by different subpopulation of photons undergoing zig-zag paths in the medium).

## Table 2

The values for the polarization parameters extracted using Eqs. (14)–(19) (third column). The values for the parameters extracted via retarder polar decomposition using Eqs. (2)–(7) for sequential linear retardance and optical rotation effects. The R, δ and Ψ- parameters extracted using Eqs. (23)–(25), for simultaneous δ and Ψ effects (fifth column). The uncertainty of the parameter Ψlog-M was very small ∼10−5) and hence not given. The first column shows the control inputs. The control input for the net depolarization coefficient Δ was determined from the experimental Mueller matrix of the pure depolarizing phantom having no birefringence and optical activity. The expected value for linear retardance (δ) and optical rotation (Ψ) are estimated by using the corresponding values from the nonscattering phantom (δ0 and Ψ0) and accounting for increased pathlength [δ=δ0×〈L〉, Ψ=Ψ0×〈L〉, 〈L〉=1.17 cm MC-derived average photon pathlength).

Parameters | Control input | Differential matrix | Retarder polar decomposition (Sequential) | Retarder polar decomposition (Simultaneous) |
---|---|---|---|---|

$\mathrm{\Delta}$ | $\sim 0.188$ | 0.211 | 0.210 | 0.210 |

$\delta $ (radian) | $\sim 1.574$ | $1.386\pm 0.022$ | 1.384 | 1.386 |

$\mathrm{\Psi}$ (radian) | $\sim 0.030$ | 0.030 | 0.036 | 0.030 |

R (radian) | $\sim 1.575$ | 1.387 | 1.387 | 1.387 |

$D$ | $\sim 0$ | $0.030\pm 0.009$ | 0.032 | 0.032 |

In Fig. 3, we show the Monte Carlo simulation-generated Mueller matrix and the corresponding ${L}_{m}$ and ${L}_{u}$ matrices for a birefringent (anisotropy in refractive index $\mathrm{\Delta}n=1.36\times {10}^{-5}$, that corresponds to a value of $\delta =1.34\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{rad}$ for a path length of 1 cm), chiral (optical activity $\chi =1.96\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{deg}.{\mathrm{cm}}^{-1}$), turbid medium (${\mu}_{s}=60\text{\hspace{0.17em}}\text{\hspace{0.17em}}{\mathrm{cm}}^{-1}$, $g=0.95$, $\mathrm{thickness}=1\text{\hspace{0.17em}}\text{\hspace{0.17em}}\mathrm{cm}$), in the backward detection geometry. The results are shown for scattered light collected at a spatial position 3 mm away from the point of illumination in the backscattering geometry. The determined values for the polarization parameters are listed in Table 3. Once again, the depolarization coefficients ($\mathrm{\Delta}$) derived from both the formalisms agree reasonably well. The slightly larger value for the diattenuation $D$ (extracted from both formalisms) in the backscattering geometry as compared to forward scattering geometry is attributed to the scattering-induced diattenuation effect (which primarily arises due to the contribution of singly/weakly backscattered photons).^{26} Since this backscattering-induced diattenuation effect is not exhibited in a completely distributed fashion, the resulting diattenuation parameter extracted from the differential matrix deviates slightly from that extracted from the polar decomposition. Further, since the effect of the curved-propagation path is more pronounced in the backscattering geometry, the resulting reduction in the apparent linear retardance is also more prominent in the backscattering geometry as compared to forward scattering geometry.^{26} The reduction of net optical rotation (in contrast to forward scattering geometry) is also known to originate from the contribution of the backscattered photons (that suffer scattering in the backward hemisphere only) which changes the handedness of rotation, leading to cancelation of net optical rotation.^{26} As observed before, use of retarder polar decomposition leads to systematic underestimation of the value of $\mathrm{\Psi}$. Never-the-less, this discrepancy is resolved when these are determined using Eqs. (24) and (25), for simultaneously occurring $\delta $ and $\mathrm{\Psi}$ effects. Further, the differential matrix-derived uncertainty in linear retardance ($\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$) is observed to be higher in the backscattering geometry as compared to the forward scattering geometry. The stronger random orientation effects of linear birefringence (as observed in the photon’s reference frame) and enhanced depolarization effects in the backscattering geometry account for the larger uncertainty in linear retardance. Similarly, contribution of optical rotation of different handedness from two distinct subgroups of photons (photons that undergo a series of forward scattering events to eventually emerge in the backward direction, and the back-scattered photons) might lead to the observed slightly nonzero value for uncertainty of net optical rotation ($\mathrm{\Delta}{\mathrm{\Psi}}^{\mathrm{log}\text{-}M}$).

## Table 3

The values for the polarization parameters extracted using Eqs. (14)–(19) (third column). The values for the parameters extracted via retarder polar decomposition using Eqs. (2)–(7) for sequential linear retardance and optical rotation effects. The R, δ and Ψ -parameters extracted using Eqs. (23)–(25), for simultaneous δ and Ψ effects (5th column). The first column shows the control inputs. The control input for the net depolarization coefficient Δ was determined from the MC-generated Mueller matrix of the pure depolarizing phantom having no birefringence and optical activity. The expected value for linear retardance (δ) and optical rotation (Ψ) are estimated by using the intrinsic Δn and χ values and accounting for increased pathlength [δ=(2π/λ)Δn×〈L〉, Ψ=χ×〈L〉, 〈L〉=1.33 cm MC-derived average pathlength).

Parameters | Control input | Differential matrix | Retarder polar decomposition (Sequential) | Retarder polar decomposition (Simultaneous) |
---|---|---|---|---|

$\mathrm{\Delta}$ | $\sim 0.28$ | 0.293 | 0.290 | 0.290 |

$\delta $ (radian) | $\sim 1.782$ | $1.45\pm 0.033$ | 1.444 | 1.444 |

$\mathrm{\Psi}$ (radian) | $\sim 0.046$ | $0.036\pm 0.003$ | 0.044 | 0.036 |

R (radian) | $\sim 1.784$ | 1.452 | 1.446 | 1.446 |

$D$ | $\sim 0$ | $0.052\pm 0.007$ | 0.044 | 0.044 |

In order to further test the validity of the retarder polar decomposition for simultaneous effects (Eqs. (24) and (25) of Sec. 2.4), for the extraction of $\delta $ and $\mathrm{\Psi}$-parameters from media exhibiting these effects simultaneously, in the following we present the decomposition results of Mueller matrices from a liquid crystal-based twisted nematic spatial light modulator (TNSLM, LC-R 2500, Holoeye Photonics, Germany) in transmission geometry. TNSLM is an optical element whereby the linear retardance and optical rotation effects can be assumed to be occurring simultaneously over the thickness of the liquid crystal layer.^{27}^{,}^{28} The linear retardance effect arises due to the orientation of the liquid crystal molecules whereas the optical rotation arises owing to the molecular twist over the thickness of the liquid crystal layers.^{28} Figure 4(a) and 4(b) shows the variation of the derived $\delta $ and $\mathrm{\Psi}$-parameters respectively, as a function of varying gray level (address voltage) of the TNSLM. One can observe that with decreasing gray levels (increasing the address voltage), both $\delta $ and $\mathrm{\Psi}$-parameters approach zero. This occurs because, with increasing address voltage, the bulk liquid crystal directors tend to get aligned along the electric field direction, resulting in reduction of both $\delta $ and $\mathrm{\Psi}$-parameters.^{27}^{,}^{28} In agreement with previous reports,^{27} at higher gray levels (close to the maximum value 255), the value for $\mathrm{\Psi}$ reaches its maximum and the value for $\delta $ initially increases and then shows a decreasing trend. In general, the Mueller matrices of the TNSLM exhibited very low depolarization values (not shown), possibly arising due to orientation fluctuations of the molecules. Importantly, $\delta $ and $\psi $-parameters derived via differential matrix formalism and retarder polar decomposition with its associated algebra [Eqs. (5) and (6)] for sequential occurrence of linear retardance and optical rotation, show significant discrepancy. It is apparent that retarder polar decomposition [Eq. 4)] leads to systematic underestimation of $\delta $ and overestimation of the value for $\mathrm{\Psi}$. While the overestimation effect of the $\mathrm{\Psi}$-parameter was also noticed for the optical phantoms, the underestimation effect of $\delta $-parameter becomes more apparent here. This is because of the considerably higher magnitude of the optical rotation for the TNSLM as compared to that of the optical phantoms. As can be noticed, the retarder polar decomposition for simultaneous linear retardance and optical rotation effects [Eqs. (24) and (25)], results in complete agreement of the derived $\delta $ and $\mathrm{\Psi}$-parameters with the differential matrix formalism.

The results presented above clearly demonstrate that the intrinsic polarization parameters (specifically, $\delta $ and $\mathrm{\Psi}$ and $\mathrm{\Delta}$) of a medium exhibiting simultaneous polarization effects, can be ‘self-consistently’ extracted from both the differential matrix and the polar decomposition approach corrected for simultaneous occurrence of linear retardance and optical rotation effects. Moreover, the results on the tissue phantoms indicate that the differential matrix formalism yields additional polarization metrics on the uncertainty of the polarization parameters (specifically, $\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$ of linear retardance), which can be of potential use for probing local organization (disorganization) of birefringent tissue structures. We have thus explored such possibility in the context to monitoring of myocardial tissue regeneration following stem cell therapy (through quantification of linear retardance from the measured Mueller matrices), as presented below.

The Mueller matrix images from myocardial tissue samples from Lewis rat, both with and without stem-cell treatment, were recorded using our imaging polarimetry system. In agreement with our previous results on the use of our point polarimetry system,^{23}^{,}^{24} the results showed a large decrease in the magnitude of linear retardance ($\delta $) in the infarcted region of the untreated myocardium. In contrast, in the infarcted region after stem-cell treatment, an increase of $\delta $ towards native levels was observed (images not shown here), indicating re-growth and re-organization of the myocardium. Moreover, the differential matrix-derived polarimetry images showed interesting spatial behavior, as shown in Fig. 5(a), 5(b), 5(c), and 5(d). Comparison of ($a$) linear retardance ${\delta}^{\mathrm{log}\text{-}M}$ and ($b$) depolarization ${\mathrm{\Delta}}^{\mathrm{log}\text{-}M}$– images show certain common features. Both show considerable spatial variation. Most of the regions associated with high depolarization exhibit higher retardance values also. Note that regions exhibiting higher depolarization are possibly associated with increased photon pathlength (due to increased ventricular thickness with stem-cell treatment^{24}), and thus consequently exhibit higher linear retardance values. Never-the-less, many of the features in the linear retardance image may be attributed to the intrinsic variation of birefringence also. The orientation angle of birefringence ${\theta}^{\mathrm{log}\text{-}M}$–image ($c$), derived via Eq. (15), also shows interesting spatial variation. In general, the spatial variations of ${\theta}^{\mathrm{log}\text{-}M}$ are more prominent towards the edges, whereas in the middle region these are weaker (i.e., show better uniformity). The variations through the myocardial wall are likely due to the change in orientation of the myocyte fibers through the wall. The orientation of the fibers undergoes a rotation from the outside of the ventricle to the inside, leading to a rapid change in orientation angle towards the edges.^{24} In contrast, near the center of the ventricular wall the fibers are nearly uniformly aligned. It is interesting to note that the derived uncertainty of linear retardance ${\mathrm{\Delta}\delta}^{\mathrm{log}\text{-}M}$—image ($d$) exhibits some kind of correlated behavior with the ${\theta}^{\mathrm{log}\text{-}M}$—image ($c$). Regions with larger variation of orientation angle ${\theta}^{\mathrm{log}\text{-}M}$ are observed to exhibit larger uncertainty in the derived linear retardance (${\mathrm{\Delta}\delta}^{\mathrm{log}\text{-}M}$). This is expected because rapid local changes in orientation angle of birefringence should in general lead to larger uncertainty in linear retardance derived via the differential matrix decomposition (from ${L}_{u}$ matrix). Moreover, such local changes in orientation angle of birefringence should also lead to additional depolarization effects. This follows because addition of the individual retardance matrices having random orientation of axes manifests as depolarization.^{4} This is likely to be the case observed here [Fig. 5(b)]. Thus, it appears that the differential matrix decomposition-derived parameters, ${\mathrm{\Delta}}^{\mathrm{log}\text{-}M}$, ${\delta}^{\mathrm{log}\text{-}M}$, ${\theta}^{\mathrm{log}\text{-}M}$, and $\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$ provide complementary and useful information on the morphology of the tissue: overall scattering properties via ${\mathrm{\Delta}}^{\mathrm{log}\text{-}M}$, anisotropic organization of tissue structures via ${\delta}^{\mathrm{log}\text{-}M}$, local organization (or disorganization) of birefringent tissue structures via ${\theta}^{\mathrm{log}\text{-}M}$ and $\mathrm{\Delta}{\delta}^{\mathrm{log}\text{-}M}$. These polarimetry parameters thus hold promise as useful biological metrics for quantitative polarimetric tissue characterization/diagnosis/monitoring treatment response.

## 5.

## Conclusions

To conclude, we have performed a detailed comparative evaluation of the polar decomposition and the differential matrix decomposition of Mueller matrices for extraction/quantification of the individual, intrinsic polarimetry characteristics (with special emphasis on linear retardance $\delta $, optical rotation $\mathrm{\Psi}$ and depolarization $\mathrm{\Delta}$ parameters) from complex tissue-like turbid media exhibiting simultaneous scattering and polarization effects. The results suggest that for media exhibiting simultaneous polarization effects, the use of retarder polar decomposition and its associated analysis (which assumes preferred sequential occurrence of linear retardance and optical rotation effects) results in systematic underestimation of $\delta $ and overestimation of $\mathrm{\Psi}$ parameters. Complete agreement between the derived $\delta $ and $\mathrm{\Psi}$ parameters from differential matrix formalism and retarder polar decomposition was obtained when the simultaneous occurrence of linear retardance and optical rotation was incorporated in the latter formalism. The self-consistency of both the decompositions was validated on the experimental and the MC-generated Mueller matrices for media exhibiting simultaneous polarization effects (linear retardance, optical rotation and depolarization). It should be noted for completeness that the influence of the remaining polarization parameter, diattenuation, on the polar decomposition and the differential matrix decomposition has recently been reported.^{29} Since the magnitude of diattenuation is typically very weak in most tissues,^{14} it has not been investigated in details here. We have instead paid attention on linear retardance, optical rotation and depolarization parameters, because these three are the most prominent tissue polarimetry effects, which also hold promise as useful biological metrics. Never-the-less, the initial exploration of the differential Mueller matrix formalism for monitoring of myocardial tissue regeneration following stem cell therapy, showed that this formalism yields additional polarization metrics, and combined usage of these multiple polarimetry parameters hold promise for quantitative polarimetric tissue characterization/diagnosis. In general, the extended formulation, which is capable of tackling either the simultaneous or sequential occurrences of polarization events, might turn out to be a useful approach, particularly for the characterization of biological tissue, where no unique order (or sequence) can be assigned to the polarimetry effects.

## Acknowledgments

This work was supported by Indian Institute of Science Education and Research (IISER)-Kolkata, an autonomous teaching and research institute funded by the Ministry of Human Resource Development (MHRD), Govt. of India.