manoharan-lab / holopy

Hologram processing and light scattering in python
GNU General Public License v3.0
131 stars 50 forks source link

Calculation of noise_sd should be revised #430

Open vnmanoharan opened 7 months ago

vnmanoharan commented 7 months ago

In holopy.core.io, we estimate the pixel noise from a stack of images as follows. We first calculate the standard deviation for each pixel in each channel over the stack, then calculate the coefficient of variation (CV) for each pixel in each channel, then average the CV over all pixels for each channel. We end up with an estimate of the pixel noise in each channel. Let $I_{ijkl}$ be the stack of images, where $i$ and $j$ are the 2D coordinates of each image, $k$ is the channel, and $l$ is the image number in the stack. The pixel noise calculation we currently do is

\mathrm{noise\_sd}_k = \mu_{ij} \left(\frac{\sigma_{l}( I_{ijkl} )}{\mu_{l} (I_{ijkl})}\right),

where $\sigma_l$ is the standard deviation calculated over the stack, $\mu_l$ is the mean calculated over the stack, and $\mu_{ij}$ is the mean calculated over image coordinates.

There are other ways to calculate the pixel noise. An alternative approach is to take the pixel-wise mean of the numerator and denominator above, and then take the ratio:

\mathrm{noise\_sd}_k = \left(\frac{\mu_{ij} \left[\sigma_{l}( I_{ijkl} )\right]}{\mu_{ij} \left[\mu_{l} (I_{ijkl})\right]}\right).

We used to do the second method, until https://github.com/manoharan-lab/holopy/commit/89450e00f21cef98afd4e93605ffba13beef9679, when we changed to the current method. The current estimator, though, has some undesirable properties. In particular, the estimate of the pixel noise will be dominated by pixels with high CV, and the CV can be infinite if the mean value of the pixel is zero. So in general the first estimator should give a higher noise level than the second, and it will be more sensitive to outliers.

We should probably go back to the second method, but since the estimated pixel noise is used in the likelihood function, changing the calculation now would change inference results. Targeting this change for version 4.