Problem getting DECODE to work with an sCMOS camera (no EM gain)

[cleterrier]

Hi, I'm trying to train DECODE for acquired sequences from a microscope with sCMOS cameras rather than EMCCD.
I took care of changing the camera parameters in SMAP for the Z calibration and processing a sample of frames to obtain a Training parameter yaml file, which is the following (not that as I defined that there was no EM gain on the sCMOS camera, the em_gain value is 0.0):

Camera:
  baseline: 100.0
  e_per_adu: 0.48
  em_gain: 0.0
  px_size: [97.0, 97.0]
  read_sigma: 74.4
  spur_noise: 0.002
InOut: {calibration_file: 'D:/data SMLM/20210303 Calib 3D Karo/210209_2C_2Cameras_100nm/zStackCamera1_3dcal.mat',
  experiment_out: 'D:/data SMLM/20210303 Calib 3D Karo/210209_2C_2Cameras_100nm'}
Simulation:
  bg_uniform: [24.927, 38.2988]
  emitter_av: 25.0
  emitter_extent:
  - [-0.5, 39.5]
  - [-0.5, 39.5]
  - [-62.0, 507.0]
  intensity_mu_sig: [2284.5032, 456.90063]
  lifetime_avg: 1.0438

When using this yaml file as input in the training Jupyter notebook, first I get an error at the step that makes:

tar_em, sim_frames, bg_frames = simulator.sample()
sim_frames = sim_frames.cpu()

frame_path = 'D:/data SMLM/20210303 Calib 3D Karo/C1_MT-AF647_Clath-CF568_dual_astig_25ms_47um_27k_frames/ROI1_5K_crop40.tif'  # change if you load your own data
data_frames = decode.utils.frames_io.load_tif(frame_path).cpu()

print(f'Data shapes, simulation: {sim_frames.shape}, real data: {data_frames.shape}')
print(f'Average value, simulation: {sim_frames.mean().round()}, real data: {data_frames.mean().round()}')

The error being:

---------------------------------------------------------------------------
ZeroDivisionError                         Traceback (most recent call last)
<ipython-input-10-4d2b30743787> in <module>()
----> 1 tar_em, sim_frames, bg_frames = simulator.sample()
      2 sim_frames = sim_frames.cpu()
      3 
      4 frame_path = 'D:/data SMLM/20210303 Calib 3D Karo/C1_MT-AF647_Clath-CF568_dual_astig_25ms_47um_27k_frames/ROI1_5K_crop40.tif'  # change if you load your own data
      5 data_frames = decode.utils.frames_io.load_tif(frame_path).cpu()

C:\Users\chris\.conda\envs\decode_env\lib\site-packages\decode\simulation\simulator.py in sample(self)
     54 
     55         emitter = self.em_sampler()
---> 56         frames, bg = self.forward(emitter)
     57         return emitter, frames, bg
     58 

C:\Users\chris\.conda\envs\decode_env\lib\site-packages\decode\simulation\simulator.py in forward(self, em, ix_low, ix_high)
     93 
     94         if self.noise is not None:
---> 95             frames = self.noise.forward(frames)
     96 
     97         return frames, bg_frames

C:\Users\chris\.conda\envs\decode_env\lib\site-packages\decode\simulation\camera.py in forward(self, x, device)
     94         """Gamma for EM-Gain (EM-CCD cameras, not sCMOS)"""
     95         if self._em_gain is not None:
---> 96             camera = self.gain.forward(camera)
     97 
     98         """Gaussian for read-noise. Takes camera and adds zero centred gaussian noise."""

C:\Users\chris\.conda\envs\decode_env\lib\site-packages\decode\simulation\noise_distributions.py in forward(self, x)
     56 
     57     def forward(self, x):
---> 58         return torch.distributions.gamma.Gamma(x, 1 / self.scale).sample()
     59 
     60 

ZeroDivisionError: float division by zeroge value, simulation: {sim_frames.mean().round()}, real data: {data_frames.mean().round()}')

The "division by zero" suggested that DECODE had a problem with em_gain=0, so I replaced it by em_gain = 1 in the yaml file from SMAP:

Camera:
  baseline: 100.0
  e_per_adu: 0.48
  em_gain: 1.0
  px_size: [97.0, 97.0]
  read_sigma: 74.4
  spur_noise: 0.002
InOut: {calibration_file: 'D:/data SMLM/20210303 Calib 3D Karo/210209_2C_2Cameras_100nm/zStackCamera1_3dcal.mat',
  experiment_out: 'D:/data SMLM/20210303 Calib 3D Karo/210209_2C_2Cameras_100nm'}
Simulation:
  bg_uniform: [24.927, 38.2988]
  emitter_av: 25.0
  emitter_extent:
  - [-0.5, 39.5]
  - [-0.5, 39.5]
  - [-62.0, 507.0]
  intensity_mu_sig: [2284.5032, 456.90063]
  lifetime_avg: 1.0438

Now it runs up to the comparison between the simulated frames and the experimental frames. However, the simulated frame show only noise and no visible blinking event using the parameters from the SMAP file, with an intensity span from 0-300 on the simulation vs 0-80 on the experimental frame:

So there seem to be a problem here. I don't understand why I don't see any blinking event on the simulated frame with this yaml parameter file. If you want to reproduce using the yaml file above, I also attach the calibration file used:
http://www.neurocytolab.org/up/DECODE_sCMOS_Zcalib.mat

Thank you for your help!

PS.
A side issue:
The text just above the comparison step in the notebook says:

Looking at the mean brightness, we see that there is a large mismatch. The reason is that we simulate the photons that are emitted by the fluorophores, while the tiff file shows the photons that are recorded by the camera (after amplification). To enable direct comparison we convert the real data frames into photon numbers.
Be aware that camera.forward() and camera.forward() are not each others inverse. forward() performs noise sampling, while backward() simply rescales by the em_gain

Indeed, the mean intensity of the simulation is 44.0, while the average intensity of the experimental frames as reported in the notebook is 157.0 before the camera.backward rescaling and 27.0 after the rescaling (which is consistent with [(157-baseline)/(e_per_adu•em_gain)], with the values from the yaml file above (baseline 100, e_per_adu 0.48, em_gain 1). As an aside I think "rescaled by the em_gain" in the notebook test is imprecise as it's rescaled by (e_per_adu•em_gain) if I understand correctly what happens.

[Haydnspass]

Hey @cleterrier,

for no EM-Gain, it should be set to None not 0 or 1, but you have to explicitly set it to nothing, you must not leave it out.
Therefore this part should look like:

Camera:
  baseline: 100.0
  e_per_adu: 0.48
  em_gain: 
  px_size: [97.0, 97.0]
  read_sigma: 74.4
  spur_noise: 0.002

I'll have a look at your calibration file. If you want to debug you can overwrite your camera in the notebook with

# ...
camera = decode.simulation.camera.PerfectCamera(device='cpu')

which will only result in shot noise.

[cleterrier]

I added the perfect camera in the training notebook, just before:
# finally we derive some parameters automatically for easy use param = decode.utils.param_io.autoset_scaling(param)

and the result is the same pattern of pure shot noise. So I guess the Z calibration is the culprit. This is a spline calibration made by SMAP using a larger window (31x31 pixels) because there is a lot of astigmatism on this setup, so PF are extended. Here's the SMAP Matlab Figure corresponding to the calibration:
http://www.neurocytolab.org/DECODE_sCMOS_CalibFig.fig

and the corresponding Z stack tif (34 steps, 50 nm Z steps, 97 nm/pixel, 0.48 e-/ADU):
http://www.neurocytolab.org/DECODE_sCMOS_ZcalibStack.tif

[Haydnspass]

To me the PSF looks fine.
I sampled a few emitters and plotted a single frame. See the example below.
Good catch with the notebook. camera.forward sample from the noise distribution, camera.backward simply converts camera units back to photon units

[cleterrier]

Thanks @Haydnspass! Seems that the problem is not in the calibration then?

Can you try to execute the Training notebook with

• this calibration file (same as above)
• this yaml file from SMAP (with em_gain modified to none)
• these sample frames?

I'd be interested to see if you get something that looks like pure noise also... You can see the experimental data doesn't have a high SNR, maybe the SMAP parameters inferred make it even worse, and it just looks like noise?

Because with the above files, to me right now the comparison looks like this:

[Haydnspass]

Just tried, it is just the low photon count. As soon as I used a factor of 10 for intensity_mu_sig the emitters became visible to the eye (so no bug).
I checked the normalization of the PSF throughout the z range which looks fine as well. (see below, first normalization check, second one ROI without noise). Notice that the max. intensity is far less than the summed intensity so I guess intensity_mu is just off by a factor of 10-100.

Where did you get the intensity and bg values from? From the SMAP plugin or by guess?