feat: Remove API for regular likelihood #188
Conversation
There was a problem hiding this comment.
Pull request overview
This PR removes the regular likelihood calculation API in favor of using only log-likelihood for numerical stability. The change affects the core Equation trait by removing the estimate_likelihood method and standardizing on estimate_log_likelihood. This improves numerical stability when dealing with many observations or extreme parameter values that could cause regular likelihoods to underflow to zero.
Key changes:
- Removed
estimate_likelihoodmethod from theEquationtrait and all implementations - Removed helper functions for regular likelihood calculation (
normpdf,normcdf,likelihoodmethods) - Updated the
psifunction to compute log-likelihoods directly (previously namedlog_psi) - Modified optimizer and test code to work with log-likelihoods
Reviewed changes
Copilot reviewed 10 out of 10 changed files in this pull request and generated 3 comments.
Show a summary per file
| File | Description |
|---|---|
| src/simulator/likelihood/mod.rs | Removed regular likelihood methods, kept only log-likelihood implementations; removed psi function that returned likelihoods and renamed log_psi to psi |
| src/simulator/equation/mod.rs | Removed estimate_likelihood from Equation trait; updated default simulate_subject to sum log-likelihoods |
| src/simulator/equation/ode/mod.rs | Removed estimate_likelihood implementation and helper function; updated to push log-likelihoods in simulate_subject |
| src/simulator/equation/analytical/mod.rs | Removed estimate_likelihood implementation and helper function; updated to push log-likelihoods |
| src/simulator/equation/sde/mod.rs | Updated particle filter to compute log-likelihoods; modified estimate_log_likelihood to take log of final result |
| src/optimize/spp.rs | Updated optimizer to work with log-likelihoods using exp(log_psi - log_pyl) for likelihood ratios |
| tests/ode_optimizations.rs | Updated tests to compare log-likelihoods using absolute differences instead of relative differences |
| examples/pf.rs | Updated example to call estimate_log_likelihood instead of estimate_likelihood |
| examples/one_compartment.rs | Updated example to display log-likelihoods with adjusted formatting |
| src/lib.rs | Updated public API exports to only expose psi (which now returns log-likelihoods) |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
| // Push the average likelihood (in regular space) for final computation | ||
| // log(mean(likelihood)) = log(sum(exp(log_lik))) - log(n) | ||
| likelihood.push((log_sum_q - (self.nparticles as f64).ln()).exp()); |
There was a problem hiding this comment.
The SDE implementation is computing regular likelihoods (not log-likelihoods) in the particle filter, but then these values are being summed in the default simulate_subject method (line 264 in mod.rs). This is mathematically incorrect. Likelihoods should be multiplied (product) across observations, not summed. Since you're working in regular likelihood space here, you should either: (1) change line 333 to push log-likelihoods and keep the sum, or (2) override simulate_subject for SDE to compute the product instead of sum, then take the log in estimate_log_likelihood. The first approach is recommended for consistency with ODE and Analytical implementations.
| // Push the average likelihood (in regular space) for final computation | |
| // log(mean(likelihood)) = log(sum(exp(log_lik))) - log(n) | |
| likelihood.push((log_sum_q - (self.nparticles as f64).ln()).exp()); | |
| // Push the average log-likelihood for final computation | |
| // log(mean(likelihood)) = log(sum(exp(log_lik))) - log(n) | |
| likelihood.push(log_sum_q - (self.nparticles as f64).ln()); |
| } | ||
| } | ||
| let ll = error_models.map(|_| likelihood.iter().product::<f64>()); | ||
| // Return log-likelihood (sum of log-likelihoods) when error_models is provided |
There was a problem hiding this comment.
The inline comment correctly describes returning log-likelihood, which matches the implementation. However, the function's documentation comment at line 231 (outside the changed region) still says "optional likelihood" and should be updated in a follow-up to say "optional log-likelihood" for consistency.
| let lik = if cache { | ||
| _estimate_likelihood(self, subject, support_point, error_models) | ||
| } else { | ||
| _estimate_likelihood_no_cache(self, subject, support_point, error_models) |
There was a problem hiding this comment.
The function _estimate_likelihood_no_cache is called but never defined in this module or elsewhere in the codebase. This will cause a compilation error. You need to either define this function or change the implementation to use simulate_subject directly without caching, similar to how it's done for the cached version.
We should only calculate and use the log likelihood for numerical stability