Performance of Maximum Likelihood estimation of Mueller matrices taking into account physical realizability and Gaussian or Poisson noise statistics

Haofeng Hu; Razvigor Ossikovski; François Goudail

doi:10.1364/OE.21.005117

1. Introduction

Mueller matrix polarimetry, which allows one to determine the full polarimetric response of a medium, has many applications in machine vision, remote sensing, biomedical imaging and industrial control [1–7]. Generally speaking, a Mueller matrix polarimeter is composed of a polarization state generator, which generates several independent states of polarization to illuminate the medium, and a polarization-state analyzer, which analyzes the state of polarization modified by the medium.

In practice, Mueller matrix measurements are always perturbed by noise, which causes the retrieved Mueller matrix to deviate from the true value, and may even make it non physically realizable [8,9]. One thus has to find an estimator that is as close as possible from the observed Mueller matrix, while being physically realizable. Up to now, extensive works have been conducted, and different estmators taking into account physical realizability in reasonable ways have been proposed [10–14].

However, little attention has been paid to characterizing the accuracy of these estimators, and to taking into account the statistics of the noise. In the present paper, we address this issue by using the tools of estimation theory. We propose a Maximum Likelihood constrained(MLC) Mueller matrix estimation method that takes into account both physical realizability and noise statistics. We illustrate this method on two types of noise sources frequently encountered in optical imaging systems: additive Gaussian noise and Poisson shot noise. In both cases, we demonstrate reduction of estimation error by enforcing the physical realizability constraint, and superiority of the MLC solutions compared to empirically constrained ones.

The paper is organized as follows. In Section 2, we introduce the problem of Mueller matrix estimation and enforcement of physical realizability. We describe the two most conventional estimation methods: unconstrained and empirically constrained estimators. In order to find an estimator with better performance in terms of estimation error, we derive in section 3 the general theory of MLC estimation and apply it to the cases of additive Gaussian noise and Poisson shot noise. Using Monte-Carlo simulations, we lead a quantitative study of the gain obtained with the MLC estimator in terms of estimation error and depolarization index (DI) of the estimate, in the presence of additive Gaussian noise (Section 4) and of Poisson shot noise (Section 5).

2. Estimation under physical realizability constraints

Measurement of the Mueller matrix involves a series of N intensity measurements made by illuminating the sample with different polarization states and analyzing them with different analyzers. The k^th measurement of this series, k ∈ [1, N], writes:

i_{k} = I_{0} t_{k}^{T} M s_{k}, k \in [1, N]

where I₀ is the intensity of the light source, which is assumed unpolarized, s_k is the 4-dimensional vector corresponding to first column of the Mueller matrix of the optical component used for generating the polarized light in illumination (polarization state generator), M is the Mueller matrix of the measured sample, and t_k is the 4-dimensional vector representing the first line of the Mueller matrix of the optical component used for analyzing the output light (polarization state analyzer). Let us denote m the 16-dimensional vector obtained by reading the Mueller matrix M in the lexicographic order. In the following, we will consider equivalently the matrix M or the vector m to represent the polarimetric properties of a material, since these two mathematical objects contain the same information. Eq. (1) can also be written as:

i_{k} = I_{0} {[t_{k} \otimes s_{k}]}^{T} m = q_{k}^{T} m, k \in [1, N]

where ⊗ denotes the Kronecker product and

q_{k} = I_{0} [t_{k} \otimes s_{k}], k \in [1, N]

Let us denote i the N-dimensional vector obtained by stacking the values of i_k, k ∈ [1, N], and Q the N × 16 matrix whose rows lines are the q_k, k ∈ {1, N}. With this notation, Eq. (2) writes:

i = Q m .

One must have N ≥ 16 to be able to reconstruct the Mueller matrix m from the measurement vector i. If N = 16, the Mueller matrix is simply obtained by inverting the matrix Q:

{\hat{m}}^{u} = Q^{- 1} i .

The estimator m̂^u will be called in the following the unconstrained (UC) estimator of m. In practice, the measurement vector i does not strictly verify Eq. (4) since it is perturbed by noise, which can be of different origins, such as for example detector noise, Poisson shot noise, or speckle noise. In this case, the estimator m̂^u is not exactly equal to the true value of the observed Mueller matrix, that we will denote m⁰. It may even correspond to a non physically realizable Mueller matrix [8, 9]. Indeed, a Mueller matrix M is physically realizable if and only if its coherency matrix, defined as

C [m] = \sum_{i, j = 0}^{3} M_{i j} σ_{i} \otimes σ_{j},

where the σ_i are the Pauli matrices, is nonnegative [11]. There is a one-to-one relation between the matrix M (or, which is equivalent, the vector m) and the matrix C[m]. Thus, in order to check if the estimated Mueller matrix m̂^u is physical, one can compute its associated coherency matrix with Eq. (6), and diagonalize it so as to obtain:

C [{\hat{m}}^{u}] = V diag (λ) V^{T},

where λ = {λ₁, λ₂, λ₃, λ₄} denotes the eigenvalues and V is a 4 × 4 unitary matrix. If some of these eigenvalues are negatives, the estimate m̂^u is a nonphysical Mueller matrix. One can easily define an estimator that is physically realizable in the following way, proposed by Cloude and Pottier [11]. One defines a new vector λ′ as follows:

λ_{k}^{'} = {\begin{matrix} λ_{k} & if & λ_{k} \geq 0 \\ 0 & otherwise \end{matrix}

from which a new coherency matrix can be constructed [11]:

C [{\hat{m}}^{e}] = V diag (λ^{'}) V^{T} .

The Mueller matrix m̂^e corresponding to this coherency matrix will be called in the following the empirically constrained (EC) estimator of m. It is easily shown that this estimator is the physically realizable matrix which is closest to the UC estimator in terms of Euclidean distance for a given realization of the noise, that is [10, 11],

{\hat{m}}^{e} = arg min_{m \in ℛ} {{‖ m - {\hat{m}}^{u} ‖}^{2}},

where ℛ denotes the set of physically realizable Mueller matrices. This estimator is the standard way to obtain a physically realizable Mueller matrix from noisy measurements. However, as far as estimation precision is concerned, the optimal estimator would rather be the one that minimizes the distance with the true value m⁰ of the Mueller matrix, averaged over the realizations of the noise, that is:

{\hat{m}}^{opt} = arg min_{m \in ℛ} {〈 {‖ m - m^{0} ‖}^{2} 〉},

where 〈.〉 denotes ensemble averaging over the noise realizations. There is no systematic way to determine a closed-form expression for this estimator, that is likely to depend on the statistical distribution of the noise affecting the measurement. However, in the next section, we will propose an estimator based on the Maximum Likelihood principle, that can perform better than the UC or the EC estimator in terms of the mean-square error criterion defined in Eq. (11).

3. Maximum Likelihood constrained estimation

In the presence of noise, the measurement vector i in Eq. (4) is a random vector whose statistical properties depend on the true Mueller matrix m⁰ and of the data acquisition process. Let us denote its joint probability density function (PDF) $P_{i}^{m^{0}} (i)$ . The Maximum Likelihood estimator with physicality constraints, that will be called in the following Maximum Likelihood Constrained (MLC) estimator, is given by [15]

{\hat{m}}^{mlc} = arg max_{m \in ℛ} {log [P_{i}^{m} (i)]} .

In order to solve efficiently this constrained optimization problem, one needs a way to describe the elements of ℛ. It is presented in the following section.

3.1. Parametrization of the Mueller matrix

A simple way to describe the elements of ℛ is introduced in Ref. [12, 13]. It consists in representing the coherency matrix by its Choleksy decomposition, C[m] = AA^†, with

A = [\begin{matrix} a_{1} & 0 & 0 & 0 \\ a_{5} + i a_{6} & a_{2} & 0 & 0 \\ a_{11} + i a_{12} & a_{7} + i a_{8} & a_{3} & 0 \\ a_{15} + i a_{16} & a_{13} + i a_{14} & a_{9} + i a_{10} & a_{4} \end{matrix}],

where the elements of the parameter vector a = (a₁,..., a₁₆) are real-valued. Knowing the vector a is equivalent to knowing the coherency matrix C[m], which is equivalent to knowing the Mueller matrix m. Consequently, by inspecting all the possible values of a in ℝ¹⁶, one inspects the whole set ℛ of physically realizable Mueller matrices. The MLC estimator defined in Eq. (12) can thus be written as:

{\hat{m}}^{mlc} = arg max_{a \in ℝ^{16}} {log [P_{i}^{m (a)} (i)]},

where m(a) explicitly denotes the Mueller matrix m associated to the vector a. Standard numerical optimization algorithms can be used to solve this unconstrained optimization problem.

The general Cholesky decomposition with matrix A in Eq. (13) ensures that the Mueller matrix is physically realizable. In practice, one may have some prior knowledge about the measured Mueller matrix, that make it possible to strengthen the constraints by modifying the parametrization of matrix A. Let us assume for example that the observed matrix is a pure Mueller-Jones matrix. This means that its coherency matrix has only one nonzero eigenvalue, and can thus be written as C[m] = BB^† with:

B = [\begin{matrix} b_{1} & 0 & 0 & 0 \\ b_{2} + i b_{3} & 0 & 0 & 0 \\ b_{4} + i b_{5} & 0 & 0 & 0 \\ b_{6} + i b_{7} & 0 & 0 & 0 \end{matrix}],

which contains only 7 variables. This method is equivalent to constructing the Mueller matrix based on the Jones matrix [16].

Furthermore, in addition to the constraint enforced by the Hermitian positive semidefinite coherency matrix, which make sure the degree of polarization (DOP) of the Stokes vector transformed by the Mueller matrix is no greater than one, we also enforce the constraint that the intensity of the input Stokes vector cannot be amplified by the Mueller matrix, which can be expressed as

M_{00} + \sqrt{M_{01}^{2} + M_{02}^{2} + {M_{03}}^{2}} \leq 1.

This constraint is enforced by subtracting a great value from the likelihood as a penalization term when Eq. (16) is violated.

3.2. Maximum Likelihood constrained estimators in the presence of additive Gaussian and Poisson noise sources

In this paper, we will evaluate the MLC method on two examples of noise sources frequently encountered in optical imaging systems: white additive Gaussian noise, which can represent detector noise, and Poisson shot noise, that represents the fundamental precision limit of well-designed systems. With these two types of noise sources, each measurement is statistically independent, and the likelihood can be written:

P_{i}^{m^{0}} (i) = \prod_{k = 1}^{N} P_{i_{k}}^{m^{0}} (i_{k}) .

Moreover, to simplify discussions, we will consider that only N = 16 measurements are made. In this case, the UC estimator is given by Eq. (5) whatever the noise source.

In the presence of white additive Gaussian noise, the expression of the loglikelihood is [17]:

ℓ_{G} = - \frac{1}{2 σ^{2}} {‖ i - Q m ‖}^{2} - ln (\sqrt{2 π} σ),

where σ² is the variance of the noise. The MLC estimator in the presence of additive white Gaussian noise, denoted GMLC estimator in the following, is thus:

{\hat{m}}^{gmlc} = arg min_{m \in ℛ} {{‖ i - Q m ‖}^{2}} .

It is interesting to notice that since the UC estimator is m̂^u = Q⁻¹i, this equation can also be written:

{\hat{m}}^{gmlc} = arg min_{m \in ℛ} {{‖ Q ({\hat{m}}^{u} - m) ‖}^{2}} .

Comparing this expression with Eq. (10), we see that it differs from the EC estimator only by a weighting of squared error with the matrix Q.

If the measurements are perturbed with Poisson shot noise, the measurement vector i is a Poisson random vector, with mean equal to Qm and probability law equal to

P (i_{k}) = e^{{[- Q m]}_{k}} \frac{{[Q m]}_{k}^{i_{k}}}{i_{k}!}, k \in [1, N] .

The loglikelihood in the presence of Poisson shot noise is given by [17]

ℓ_{P} (m) = \sum_{k = 1}^{16} [- {[Q m]}_{k} + i_{k} ln {[Q m]}_{k} + ln (i_{k}!)] .

The MLC estimator in the presence of Poisson shot noise, denoted PMLC estimator in the following, is thus given by

{\hat{m}}^{pmlc} = arg min_{m \in ℛ} {\sum_{k = 1}^{16} {[Q m]}_{k} - i_{k} ln ({[Q m]}_{k})} .

By comparing this equation with Eq. (19), we observe that the MLC estimators for the two considered types of noise sources are different. It is thus important to take into account the real statistical distribution of the noise in the design of the estimator.

3.3. Virtual experiment

In some practical cases, one has only access to the measured Mueller matrix m and not the intensities i. If the measurement matrix Q used to obtain the measurements is unknown, the previously described MLC methods cannot be directly applied. In this case, a “virtual experiment” has to be conducted [16], for which we need to set a matrix Q′, and consequently calculate the “virtual intensities” i′. The estimator in the “virtual experiment” for additive white Gaussian noise is thus given by

{\hat{m}}_{G}^{virtual} = arg min_{m \in ℛ} {{‖ i^{'} - Q^{'} m ‖}^{2}},

and that for Poisson shot noise is given by

{\hat{m}}_{P}^{virtual} = arg min_{m \in ℛ} {\sum_{k = 1}^{16} {[Q^{'} m]}_{k} - i_{k}^{'} ln ({[Q^{'} m]}_{k})} .

If the matrix Q′ used in the “virtual experiment” is the same as the matrix Q used in the real experiment, the estimators m̂^virtual will be identical to the MLC estimators. However, if Q′ is different from Q, the estimators m̂^virtual are different from the MLC ones, since in this case, the likelihoods in Eqs. 24 and 25 do not correctly represent the statistics of the noise. This may lead to a larger estimation error.

4. Estimation in the presence of additive white Gaussian noise

To model the acquisition process (see Eq. (1)), we define the measurement matrices S = [s₁, s₂, s₃, s₄] and T = [t₁, t₂, t₃, t₄] that contain the illumination and analysis vectors that are combined to obtain the 16 intensity measurements. In order to minimize the variance of the UC estimator in the presence of white Gaussian noise, it has been shown that the columns of matrices S and T have to form a regular tetrahedron on the Poincaré sphere [18–20]. The same result is obtained in the presence of Poisson noise, but in order to equalize the variances of each coefficient of the Mueller matrix, the following particular tetrahedron has to be chosen [21, 22]:

S = T = \frac{1}{2} [\begin{matrix} 1 & 1 & 1 & 1 \\ 1 / \sqrt{3} & - 1 / \sqrt{3} & - 1 / \sqrt{3} & 1 / \sqrt{3} \\ 1 / \sqrt{3} & - 1 / \sqrt{3} & 1 / \sqrt{3} & - 1 / \sqrt{3} \\ 1 / \sqrt{3} & 1 / \sqrt{3} & - 1 / \sqrt{3} & - 1 / \sqrt{3} \end{matrix}] .

We will use these measurement matrices in the following of the paper. We will also assume that the true Mueller matrix is

M^{0} = [\begin{matrix} 0.5 & 0.14 & 0 & 0 \\ 0.14 & 0.5 & 0 & 0 \\ 0 & 0 & 0.48 & 0 \\ 0 & 0 & 0 & 0.48 \end{matrix}],

It corresponds to a pure partial linear polarizer. The eigenvalues of its coherency matrix is [0, 0, 0, 1], of which three eigenvalues are equal to 0. Therefore, in the presence of noisy measurements, the constraint that the eigenvalues should not be below 0 is significant.

4.1. Comparison of estimation accuracies

We assume the measurements are disturbed by additive white Gaussian noise. Obviously, when the UC estimate is physically realizable, the UC, EC and GMLC estimates are identical. On the other hand, when the UC estimator is not physically realizable, enforcement of the physical realizability condition comes into play and the UC, EC and GMLC estimators take different values. We have represented in Table 1 the different estimates obtained for one realization of the noise. The estimates in the matrix forming in Table 1 are obtained by arraying the vectorized estimates m̂ described in Section 4 in 4 × 4 Mueller matrices. It is seen that the coherency matrix of the UC estimate contains negative eigenvalues, which indicates that it is not physically realizable. The EC estimate only puts to zero the negative eigenvalues of the UC estimate, while keeping the positive eigenvalues, whereas the GMLC estimate not only puts the negative eigenvalues to zero, but also modifies the positive eigenvalues of the UC estimate. If we further check the eigenvectors of the coherency matrix (not shown in Table 1), we can find that the GMLC estimate also modifies the eigenvectors of the UC estimate, whereas the eigenvectors of the non-zero eigenvalues of the EC estimate are identical to those of the UC estimate.

Table 1. UC, EC and GMLC estimates of one realization as well as the corresponding eigenvalues of the coherency matrix. The measurements are perturbed with additive Gaussian noise.

View Table | View all tables in this article

Using the results in Table 1, we have calculated in Table 2 the Frobenius distance between the different estimators. In the first row of Table 2, it is seen that EC estimate is closer to the UC estimate than the GMLC. It was expected since the EC estimate indeed minimizes this distance [10, 11]. On the other hand, it is seen in the second row of Table 2 that the GMLC is closer to the true Mueller matrix M⁰ than the EC estimate, which indicates a better estimation accuracy.

Table 2. The norms between the estimates and the UC estimate or the true value M⁰.

View Table | View all tables in this article

These results are an illustration of the behavior of the different estimators on a single noise realization. Let us now characterize their average performances. We will characterize this performance with the total root mean square error (RMSE), defined as,

RMSE = \sqrt{\sum_{i, j = 0}^{3} 〈 {({\hat{M}}_{i j} - {M^{0}}_{i j})}^{2} 〉},

We have represented in Fig. 1 the total RMSEs of the UC, EC and GMLC estimators as a function of signal to noise ratio (SNR) defined as:

{SNR}_{G} = 10 {log}_{10} [\frac{I_{0}^{2}}{σ^{2}}]

This SNR is expressed in the unit of dB. The RMSE has been estimated with Monte-Carlo simulations on 2000 noise realizations.

Fig. 1 Total RMSE as a function of SNR for the UC, EC and GMLC estimators in the presence of the additive Gaussian noise. The RMSE is estimated from 2000 noise realizations and the error bars correspond to the standard deviation of this estimated RMSE.

Download Full Size | PDF

To obtain the results presented in the present paper, we have compared a global optimization algorithm, the Shuffled Complex Evolution (SCE-UA), and a local optimization method, the Nelder-Mead Simplex algorithm, with a starting point equal to the empirically constrained estimator. We have checked that optimisation of the ML constrained likelihood with these two methods lead to the same RMSE. This is probably due to the fact that our starting point for the local optimization is a good guess, close to the global optimum. Since SCE-UA Method is time consuming, we have rather used local optimization in the simulation presented in this paper.

It is seen in Fig. 1 that the RMSEs of EC and GMLC estimators are lower than that of the UC estimator. This gain can be attributed to the fact of taking into account the a priori knowledge that the estimate must be physically realizable. Furthermore, it is seen that the GMLC estimator has a lower RMSE, and thus a better accuracy, than the EC estimator. It has to be noted that in Fig. 1, we have considered values of SNR_G between 20 and 45 dB since it is in this noise range that polarimeters are most often operated. However, we have checked that even when the SNR_G is as low as 1 dB, the RMSE of the GMLC estimator remains lower than that of the EC estimator.

In addition, we also tried to employ various measurement matrices (including different regular tetrahedron-based measurement matrices and different arbitrary measurement matrices), and we found that the total RMSE of the GMLC estimator is always lower than those of the EC and UC estimators. Therefore, it can be concluded that the quality of estimation can be improved by taking into account the statistics of the noise and the measurement matrix through the GMLC estimator.

Let us now assume that we know that the true Mueller matrix is a Mueller-Jones Matrix. We can represent its coherency matrix as in Eq. (15). The corresponding estimator, which is called Mueller-Jones constrained (MJC) estimator in this paper, can be obtained by employing the ML method. The MJC estimator is compared with the GMLC in Fig. 2, and it is seen that the MJC estimator further decreases the RMSE, because it integrates more a priori knowledge on the observed matrix. Furthermore, it also decreases the computation time since less variables have to be estimated.

Fig. 2 Total RMSE as a function of SNR for the GMLC and MJC estimators in the presence of the additive Gaussian noise. The RMSE is estimated from 2000 noise realizations

Download Full Size | PDF

4.2. Comparison of depolarization indices

In this section, we will compare the different estimators in terms of their depolarization index (DI) [24], defined as:

DI = \frac{\sqrt{\sum_{i, j = 0}^{3} M_{i j}^{2} - M_{00}^{2}}}{\sqrt{3} M_{00}},

This parameter characterizes the average behavior of the DOP of the Stokes vector transformed by the Mueller matrix. We have plotted in Fig. 3 the mean value of DI of EC and GMLC estimators as a function of SNR. Since the true Mueller matrix is a pure partial polarizer (Eq. (27)), its DI is equal to 1. It can be seen in Fig. 3 that the DI of GMLC estimator is statistically closer to the true value of 1 than the DI of EC estimator. This indicates that the EC estimator is too “cautious”, that is, it overly decreases the DI in order to avoid that the DOP of the Stokes vector transformed by the Mueller matrix be greater than 1. However, this decrease makes the DI of EC estimator farther away from the true value. This result is coherent with the superiority of the GMLC in terms of RMSE.

Fig. 3 Mean of DI of EC and GMLC estimators as a function of SNR. The average DI is estimated from 2000 noise realizations and the error bars correspond to the standard deviation of this estimated mean.

Download Full Size | PDF

In order to verify the universality of the superiority of GMLC estimator, we conducted Monte-Carlo simulations for various Mueller matrices with different values of DI and at various values of SNR, and it is found that the total RMSE of the GMLC estimator is always the lowest. In addition, we also tried to employ various measurement matrices (including different regular tetrahedron-based measurement matrices and different arbitrary measurement matrices), and we found that the DI of the GMLC estimator is always greater than that of the EC estimator and is closer to the true value of DI under any condition.

5. Estimation in the presence of Poisson shot noise

Let us now assume that the measurements are perturbed with Poisson shot noise. For the purpose of illustration, we employ measurement matrices S and T given by Eq. (26) and the tested Mueller matrix given by Eq. (27). According to our simulation results for one realization of the noise shown in Table 3, we observed that the PMLC estimator (defined in Eq. (23)) not only sets the negative eigenvalues of coherency matrix to 0, but also modifies the positive eigenvalues and the eigenvectors of the coherency matrix as the GMLC estimator does in the presence of additive Gaussian noise. In addition, if we look at the Frobenius distance between the different estimates shown in Table 4, we can find that PMLC estimate is closer to the true Mueller matrix than the EC estimate.

Table 3. UC, EC and PMLC estimates of one realization as well as the corresponding eigenvalues of the coherency matrix. The measurements are perturbed with Poisson shot noise.

View Table | View all tables in this article

Table 4. The norm between the EC or PMLC estimate and the UC estimate or the true value M⁰.

View Table | View all tables in this article

In order to compare the average estimation accuracies of the different estimators, we have plotted in Fig. 4 the total RMSEs of the UC, EC and PMLC estimators as a function of SNR. In the presence of Poisson shot noise, the SNR is defined as

{SNR}_{P} = 10 {log}_{10} [I_{0}]

It can be seen in Fig. 4 that the total RMSE of the EC estimator is lower than that of the UC estimator, which is attributed to taking into account the physical realizability constraint. It is also noticeable that the total RMSE of the PMLC estimator is lower than that of the EC estimator, which is attributed to additionally taking into account the statistics of the noise.

Fig. 4 Total RMSE as a function of SNR for the UC, EC and PMLC estimators in the presence of Poisson shot noise. The number of realizations is 2000.

Download Full Size | PDF

Besides, it is seen in Fig. 5 that the mean value of DI of the PMLC estimate is greater than that of the EC estimate, and closer to the true value of DI than that of the EC estimator. We have thus a property similar to that observed in the presence of additive Gaussian noise.

Fig. 5 Mean of DI as a function of SNR for the EC and PMLC estimators in the presence of Poisson shot noise. The number of realizations is 2000.

Download Full Size | PDF

In the presence of Poisson shot noise, we also conducted simulations for estimating various Mueller matrices with different values of DI at various values of SNR and with various measurement matrices. We found that the total RMSE of the PMLC estimator is always lower than those of the UC and EC estimators. In addition, it was found that the DI of PMLC estimator is always statistically greater than that of the EC estimator. It can thus be concluded that in the presence of Poisson noise, the MLC technique also provides better accuracy than classical unconstrained or constrained methods.

6. Conclusion

We have investigated the estimation performance of a Maximum Likelihood-based estimator of the Mueller matrix that enforces physical realizability through the use of the nonnegative Hermitian coherency matrix. Using numerical simulations, we showed that Maximum Likelihood Constrained estimator has a better accuracy than unconstrained and empirically constrained estimators, either in the presence of additive white Gaussian noise or Poisson shot noise. It thus presents a statistically justified and efficient alternative to the classical constrained estimator proposed in [10, 11], and it can be applied in principle to any kind of noise source as soon as the likelihood is known. Interesting perspectives to this work thus include application of the MLC technique to other types of noise sources. Besides, Mueller polarimetry experiment for comparing the performances of different estimators for real measured matrices is also a promising perspective.

Acknowledgments

The authors wish to thank Guillaume Anna, Enric Garcia-Caurel and Alex Kostinski for fruitful discussions.

References and links

1. J. E. Solomon, “Polarization imaging,” Appl. Opt. 20, 1537–1544 (1981). [CrossRef] [PubMed]

2. J. S. Tyo, M. P. Rowe, E. N. Pugh, and N. Engheta, “Target detection in optical scattering media by polarization-difference imaging,” Appl. Opt. 35, 1855–1870 (1996). [CrossRef] [PubMed]

3. S. Breugnot and P. Clémenceau, “Modeling and performances of a polarization active imager at λ = 806 nm,” Opt. Eng. 39, 2681–2688 (2000). [CrossRef]

4. S. L. Jacques, J. C. Ramella-Roman, and K. Lee, “Imaging skin pathology with polarized light,” J. Biomed. Opt. 7, 329–340 (2002). [CrossRef] [PubMed]

5. Y. Y. Schechner, S. G. Narasimhan, and S. K. Nayar, “Polarization-based vision through haze,” Appl. Opt. 42, 511–525 (2003). [CrossRef] [PubMed]

6. J. M. Bueno, J. Hunter, C. Cookson, M. Kisilak, and M. Campbell, “Improved scanning laser fundus imaging using polarimetry,” J. Opt. Soc. Am. A 24, 1337–1348 (2007). [CrossRef]

7. A. Pierangelo, B. Abdelali, M.-R. Antonelli, T. Novikova, P. Validire, B. Gayet, and A. De Martino, “Ex-vivo characterization of human colon cancer by mueller polarimetric imaging,” Opt. Express 19, 1582–1593 (2011). [CrossRef] [PubMed]

8. J. Cariou, B. L. Jeune, J. Lotrian, and Y. Guern, “Polarization effects of seawater and underwater targets,” Appl. Opt. 29, 1689–1695 (1990). [CrossRef] [PubMed]

9. B. J. Howell, “Measurement of the polarization effects of an instrument using partially polarized light,” Appl. Opt. 18, 809–812 (1979). [CrossRef] [PubMed]

10. D. G. M. Anderson and R. Barakat, “Necessary and sufficient conditions for a mueller matrix to be derivable from a jones matrix,” J. Opt. Soc. Am. A 11, 2305–2319 (1994). [CrossRef]

11. R. S. Cloude and E. Pottier, “Concept of polarization entropy in optical scattering,” Opt. Eng. 34, 1599–1610 (1995). [CrossRef]

12. A. Aiello, G. Puentes, D. Voigt, and J. P. Woerdman, “Maximum-likelihood estimation of mueller matrices,” Opt. Lett. 31, 817–819 (2006). [CrossRef] [PubMed]

13. J. Zallat, C. Heinrich, and M. Petremand, “A bayesian approach for polarimetric data reduction: the mueller imaging case,” Opt. Express 16, 7119–7133 (2008). [CrossRef] [PubMed]

14. J. E. Ahmad and Y. Takakura, “Estimation of physically realizable mueller matrices from experimentsusing global constrained optimization,” Opt. Express 16, 14274–14287 (2008). [CrossRef] [PubMed]

15. S. M. Kay, Fundamentals of statistical signal processing - Volume I : Estimation Theory (Prentice-Hall, Englewood Cliffs, 1993).

16. R. Ossikovski, “Retrieval of a nondepolarizing estimate from an experimental mueller matrix through virtual experiment,” Opt. Lett. 37, 578–580 (2012). [CrossRef] [PubMed]

17. H. H. Barrett and K. J. Myers, Foundations of Image Science (Wiley, New York, 2004).

18. A. Ambirajan and D. C. Look, “Optimum angles for a polarimeter: part ii,” Opt. Eng. 34, 1656–1658 (1995). [CrossRef]

19. D. S. Sabatke, M. R. Descour, E. L. Dereniak, W. C. Sweatt, S. A. Kemme, and G. S. Phipps, “Optimization of retardance for a complete stokes polarimeter,” Opt. Lett. 25, 802–804 (2000). [CrossRef]

20. J. S. Tyo, “Noise equalization in stokes parameter images obtained by use of variable-retardance polarimeters,” Opt. Lett. 25, 1198–1200 (2000). [CrossRef]

21. F. Goudail, “Noise minimization and equalization for stokes polarimeters in the presence of signal-dependent poisson shot noise,” Opt. Lett. 34, 647–649 (2009). [CrossRef] [PubMed]

22. G. Anna and F. Goudail, “Optimal mueller matrix estimation in the presence of poisson shot noise,” Opt. Express 20, 21331–21340 (2012). [CrossRef] [PubMed]

23. Q. Y. Duan, V. K. Gupta, and S. Sorooshian, “A shuffled complex evolution approach for effective and efficient global minimization,” J. Optim. Theory Appl. 76, 501–521 (1993). [CrossRef]

24. J. J. Gil and E. Bernabeu, “Depolarization and polarization indices of an optical system,” Optica Acta: International Journal of Optics 33, 185–189 (1986). [CrossRef]

	Estimates	Eigenvalues
UC	$[\begin{matrix} 0.5009 & 0.1368 & 0.0022 & 0.0025 \\ 0.1433 & 0.5103 & - 0.0082 & - 0.0009 \\ - 0.0025 & - 0.0053 & 0.4956 & 0.0073 \\ 0.0049 & - 0.0049 & 0.0027 & 0.4824 \end{matrix}]$	$[\begin{matrix} - 0.0148 \\ - 0.0065 \\ 0.0089 \\ 1.0143 \end{matrix}]$
EC	$[\begin{matrix} 0.5116 & 0.1386 & - 0.0007 & 0.0024 \\ 0.1399 & 0.5052 & - 0.0047 & 0.0009 \\ - 0.0005 & - 0.0025 & 0.4886 & 0.0044 \\ 0.0044 & 0.0030 & - 0.0002 & 0.4840 \end{matrix}]$	$[\begin{matrix} 0 \\ 0 \\ 0.0089 \\ 1.0143 \end{matrix}]$
GMLC	$[\begin{matrix} 0.5034 & 0.1390 & - 0.0005 & 0.0030 \\ 0.1398 & 0.5017 & - 0.0005 & 0.0011 \\ - 0.0006 & - 0.0006 & 0.4830 & 0.0018 \\ 0.0038 & - 0.0007 & - 0.0009 & 0.4819 \end{matrix}]$	$[\begin{matrix} 0 \\ 0 \\ 0.0020 \\ 1.0047 \end{matrix}]$

	EC	GMLC
‖M̠ − M̂^u‖	0.0146	0.0194
‖M̂ − M⁰‖	0.0130	0.0064

	Estimates	Eigenvalues
UC	$[\begin{matrix} 0.5084 & 0.1286 & - 0.0051 & 0.0010 \\ 0.1321 & 0.5040 & - 0.0064 & 0.0050 \\ - 0.0001 & - 0.0021 & 0.4912 & - 0.0534 \\ 0.0062 & 0.0225 & - 0.0111 & 0.4622 \end{matrix}]$	$[\begin{matrix} - 0.0354 \\ 0.0125 \\ 0.0386 \\ 1.0011 \end{matrix}]$
EC	$[\begin{matrix} 0.5261 & 0.1287 & - 0.0059 & 0.0003 \\ 0.1316 & 0.4878 & - 0.0021 & 0.0107 \\ - 0.0003 & 0.0013 & 0.4835 & - 0.0379 \\ 0.0052 & - 0.0163 & 0.0042 & 0.4685 \end{matrix}]$	$[\begin{matrix} 0 \\ 0.0125 \\ 0.0386 \\ 1.0011 \end{matrix}]$
PMLC	$[\begin{matrix} 0.5084 & 0.1342 & - 0.0032 & 0.0023 \\ 0.1371 & 0.4860 & - 0.0002 & 0.0150 \\ 0.0015 & 0.0006 & 0.4818 & - 0.0304 \\ 0.0077 & - 0.0144 & 0.0119 & 0.4697 \end{matrix}]$	$[\begin{matrix} 0 \\ 0.0007 \\ 0.0235 \\ 0.9926 \end{matrix}]$

	EC	PMLC
‖M̂ − M̂^u‖	0.0188	0.0268
‖M̂ − M⁰‖	0.0415	0.0363

	Estimates	Eigenvalues
UC	$[\begin{matrix} 0.5009 & 0.1368 & 0.0022 & 0.0025 \\ 0.1433 & 0.5103 & - 0.0082 & - 0.0009 \\ - 0.0025 & - 0.0053 & 0.4956 & 0.0073 \\ 0.0049 & - 0.0049 & 0.0027 & 0.4824 \end{matrix}]$	$[\begin{matrix} - 0.0148 \\ - 0.0065 \\ 0.0089 \\ 1.0143 \end{matrix}]$
EC	$[\begin{matrix} 0.5116 & 0.1386 & - 0.0007 & 0.0024 \\ 0.1399 & 0.5052 & - 0.0047 & 0.0009 \\ - 0.0005 & - 0.0025 & 0.4886 & 0.0044 \\ 0.0044 & 0.0030 & - 0.0002 & 0.4840 \end{matrix}]$	$[\begin{matrix} 0 \\ 0 \\ 0.0089 \\ 1.0143 \end{matrix}]$
GMLC	$[\begin{matrix} 0.5034 & 0.1390 & - 0.0005 & 0.0030 \\ 0.1398 & 0.5017 & - 0.0005 & 0.0011 \\ - 0.0006 & - 0.0006 & 0.4830 & 0.0018 \\ 0.0038 & - 0.0007 & - 0.0009 & 0.4819 \end{matrix}]$	$[\begin{matrix} 0 \\ 0 \\ 0.0020 \\ 1.0047 \end{matrix}]$

Performance of Maximum Likelihood estimation of Mueller matrices taking into account physical realizability and Gaussian or Poisson noise statistics

Abstract

1. Introduction

2. Estimation under physical realizability constraints

3. Maximum Likelihood constrained estimation

3.1. Parametrization of the Mueller matrix

3.2. Maximum Likelihood constrained estimators in the presence of additive Gaussian and Poisson noise sources

3.3. Virtual experiment

4. Estimation in the presence of additive white Gaussian noise

4.1. Comparison of estimation accuracies

4.2. Comparison of depolarization indices

5. Estimation in the presence of Poisson shot noise

6. Conclusion

Acknowledgments

References and links

Cited By

Figures (5)

Tables (4)

Equations (31)

Optics Express