Add full cut optimization as introduced in pyirf 0.13

[LukasBeiske]

This changes the PointSourceSensitivityOptimizer to use the full cut optimization introduced in pyirf 0.13. The previous EventDisplay-like optimization can now be used via the PointSourceSensitivityGhOptimizer.

I did a quick comparison of the three optimizers using prod6 files (multiplicity >= 2 for gh opt and percentile cuts, as the HillasReconstructer is used):

I am not sure, why the EventDisplay-like optimization results in a better sensitivity at high energies.

Fixes #2771

[maxnoe]

How fine did you make the scanning of the cuts? In principle, the full cut opt should always be at least as good as the one that is restricted,.if it is allowed to find the same cuts.

[LukasBeiske]

How fine did you make the scanning of the cuts? In principle, the full cut opt should always be at least as good as the one that is restricted,.if it is allowed to find the same cuts.

I kept everything at the default values, as we have them in here right now. So, if I'm not mistaken, the only difference would be that the EventDisplay-like optimization has a theta cut with 68% efficiency, while the full optimization only tries 60% and 70%.
But I doubt that this makes such a difference. I'll test that and take a closer look again tomorrow.

And, I think, there is no check for a minimum number of events per bin in the full cut optimization, while the 68% theta cut for the EventDisplay-like optimization has a minimum of 10 events per bin. However, this only seems to play a role for the two lowest energy bins, if at all.

[LukasBeiske]

There was still an error with the application of the cuts in the irf tool. It should be correct now and this improved the sensitivity situation a lot, but there are still some bins where the EventDisplay-like optimization outperforms the full optimization, even though the first uses a theta cut with 70% efficiency for this plot, which also gets tested in the full optimization.

[ctao-dpps-sonarqube]

Analysis Details

1 Issue

• 0 Bugs
• 0 Vulnerabilities
• 1 Code Smell

Coverage and Duplications

• 89.10% Coverage (94.30% Estimated after merge)
• 0.00% Duplicated Code (0.70% Estimated after merge)

Project ID: cta-observatory_ctapipe_AY52EYhuvuGcMFidNyUs

View in SonarQube

[maxnoe]

If I remember correctly, these percentiles cannot be compared directly, as the eventdisplay-like optimization computes the percentile on the events surviving the initial cut, whereas the full optimization computes it on all events.

[LukasBeiske]

If I remember correctly, these percentiles cannot be compared directly, as the eventdisplay-like optimization computes the percentile on the events surviving the initial cut, whereas the full optimization computes it on all events.

Ah, good point, I forgot about that. I'll run the full optimization again with a finer gridding. I guess, these underperformance for some bins will disappear then.

[LukasBeiske]

Running both optimizations with finer grids:

PointSourceSensitivityOptimizer:
  gh_cut_efficiency_step=0.02
  theta_cut_efficiency_step=0.02
  
PointSourceSensitivityGhOptimizer:
  gh_cut_efficiency_step=0.02

gets the performance of the full optimization closer, but some bins are still worse then with EventDisplay-like optimization.

[kosack]

gh_cut_efficiency_step=0.02
theta_cut_efficiency_step=0.02

If a smaller step is important to get a good sensitivity, please also update the ctapipe-quickstart configurations to provide the best values for users, and also the default values in the tool itself.

[LukasBeiske]

gh_cut_efficiency_step=0.02
theta_cut_efficiency_step=0.02

If a smaller step is important to get a good sensitivity, please also update the ctapipe-quickstart configurations to provide the best values for users, and also the default values in the tool itself.

I just noticed this (similar for the full optimization):

The config is exactly the same (default values) besides the stepsize. This doesn't make sense. I'm re-running everything now and, if this behavior is still there, I'll convert this PR to draft until I figure out whats going on here.

[maxnoe]

There was still an error with the application of the cuts in the irf tool.

Is this an error introduced here? Or does it affect main? If it affects also main, could you open a PR just with the bugfix please?

[LukasBeiske]

There was still an error with the application of the cuts in the irf tool.

Is this an error introduced here? Or does it affect main? If it affects also main, could you open a PR just with the bugfix please?

This does not affect main, it was related to the application of the multiplicity cut introduced here.

[Hckjs]

gh_cut_efficiency_step=0.02
theta_cut_efficiency_step=0.02

If a smaller step is important to get a good sensitivity, please also update the ctapipe-quickstart configurations to provide the best values for users, and also the default values in the tool itself.

I just noticed this (similar for the full optimization):

The config is exactly the same (default values) besides the stepsize. This doesn't make sense. I'm re-running everything now and, if this behavior is still there, I'll convert this PR to draft until I figure out whats going on here.

Did you use the same dataset here to calculate the sensitivities as for the cut optimization?

[Hckjs]

If i understand correctly, the gh cut optimization:

• first calculates (initial) theta cuts based on inital gh cuts
• optimizes gh cuts based on this inital theta cut
• calculates optimize theta cut based on optimized gh cuts

When the minimization of relative sensitivity is calculated on the initial theta cuts, it doesn't mean that its also minimized on "optimized" theta cuts based on optimized gh cuts, right? So the discrepancy should be valid by definition...

The full cut optimization should not have that problem since its doing a full grid search:

[maxnoe]

@Hckjs This is correct yes. One could probably get around this by optimizing again at least once or until it converges. But I think we should just use the global optimization scheme instead.

[maxnoe]

I just noticed this (similar for the full optimization):

I didn' realize that plot was for the old optimization, as @Hckjs points out it doesnt hold true for the gh opt that finer steps should alsway result in lower sensitivity.

For the full optimization, it only holds if

• the coarse steps are part of the finer steps (otherwise the best step could be a step in the coarser sample not contained in the finer sample)
• the senstivity is computed on the same dataset used for the cut optimization, otherwise the statistical differences could dominate the difference in sensitivity, not actually better cuts. (I.e. "overtraining", the reason why we should use a separate dataset in the first place).

[ctao-sonarqube]

Quality Gate passed

Issues
0 New issues
0 Fixed issues
0 Accepted issues

Measures
0 Security Hotspots
89.2% Coverage on New Code
0.0% Duplication on New Code

See analysis details on SonarQube

[LukasBeiske]

I just noticed this (similar for the full optimization):

I didn' realize that plot was for the old optimization, as @Hckjs points out it doesnt hold true for the gh opt that finer steps should alsway result in lower sensitivity.

For the full optimization, it only holds if

* the coarse steps are part of the finer steps (otherwise the best step could be a step in the coarser sample not contained in the finer sample)

* the senstivity is computed on the same dataset used for the cut optimization, otherwise the statistical differences could dominate the difference in sensitivity, not actually better cuts. (I.e. "overtraining", the reason why we should use a separate dataset in the first place).

Thanks @Hckjs, I missed that. However, doing it again for the full optimization using the same dataset for cut optimization and sensitivity calculation, the same problem is visible:

I'll look into this again. Maybe there is something I'm still missing in the grid search within pyirf.

[kosack]

I'll look into this again. Maybe there is something I'm still missing in the grid search within pyirf.

Looking at only the final sensitivity makes it a bit hard to debug since it has so many factors that effect it. Those fluctuations could be due to low stats if one of the cuts is too tight. Might be good to compare the cut efficiencies, background rates, PSF, and effective areas separately

[LukasBeiske]

Looking at only the final sensitivity makes it a bit hard to debug since it has so many factors that effect it. Those fluctuations could be due to low stats if one of the cuts is too tight. Might be good to compare the cut efficiencies, background rates, PSF, and effective areas separately

I did some more plots (see below) and checked the code again, but I did not find anything new.
However, since Jonas did the plot above (where the finer grid search is always as good or better as the coarser one) based on files he re-processed from dl1, my current suspicion is that something changed between ctapipe 0.23.1 (with which the dl2 files I am using where processed) and now that fixes this problem.
I will re-process the same files Jonas used and check whether this is actually the case.

[ctao-sonarqube]

Quality Gate passed

Issues
0 New issues
0 Fixed issues
0 Accepted issues

Measures
0 Security Hotspots
89.7% Coverage on New Code
0.0% Duplication on New Code

See analysis details on SonarQube

[maxnoe]

Checking with @Hckjs, I think we understand what is happening.

The cut optimization is done in a single FoV range, by default from 0° to 5° degrees offset.

If you optimize the cuts in this full range but then compute sensitivity in different, smaller offset ranges, the assumption smaller step = better or equal sensitivity does not hold.

[Hckjs]

Small update after checking again with @maxnoe: We think the root of the problem is the DL2EventLoader.make_event_weights. It keeps the default event weights of 1 for events outside the fov_offset bins which are used downstream by the cut optimization, but not by the Sensitivity2dMaker. This should be fixed by #2927

[LukasBeiske]

That fixed it (See plots below)! Thanks for looking into this and sorry that it took me this long to respond.

[LukasBeiske]

Would it make sense to already refactor the irf code here to use Karls EventPreprocessor added in #2928? Or would it be better to wait for #2927 and then refactor the irf code to use both in a separate PR?

[maxnoe]

Let's do independent things in independent PRs, this is already a great new feature

[ctao-sonarqube]

Quality Gate passed

Issues
0 New issues
0 Fixed issues
0 Accepted issues

Measures
0 Security Hotspots
90.3% Coverage on New Code
0.0% Duplication on New Code

See analysis details on SonarQube