Repository navigation
Derive log CDF and log CCDF of transformed distributions from std_log_lcdf #3364
Description
Activity
Most of this issue is implemented in #3406. There,
normal_lcdf,lognormal_lcdf/lognormal_lccdfandexp_mod_normal_lcdf/exp_mod_normal_lccdfcompute the standard normal log CDF and its gradient through one shared internal implementation,internal::std_normal_lcdf_value_gradinprim/fun/std_normal_lcdf_impl.hpp. #3406 also answers the OpenCL question: it adds OpenCL device functions for the same computation, and the OpenCLnormal,lognormalandexp_mod_normalfunctions use them.normal_lccdfandstd_normal_lccdfalready go through thelcdffunctions ondevelopsince #3363.The skew normal item does not fit the same pattern. The skew normal CDF is
Phi(z) - 2 T(z, alpha), with Owen's T function, anddevelopcomputes it witherfcandowens_t. In the tails the precision is lost in that subtraction, so a call tonormal_lcdfwould not make it stable. It needs its own representation of the skew normal tail.I think this issue could be closed when #3406 is merged, and then open a separate issue for stable
skew_normal_lcdfandskew_normal_lccdfif there is interest.
Issue #1985 highlighted how changes improving numerical accuracy do not always propagate their benefits to closely related distributions. In the case of standard transformations, the Stan community could capture some DRY benefits from connecting as many derived functions as possible to a reference function that keeps up with the latest numerical improvements (or at least tries to).
While I don't have the overview to make a super-issue of this – there are surely some places where it does make sense to duplicate even lengthy and complicated numerical code – I think there is some low-hanging fruit for the normal distribution. Along the lines of what I tried to do in pull request #3363, I would propose specifically:
normal_lcdfto rely onstd_normal_lcdf(with some extra partials lines to avoid autodiffing, as discussed here);lognormal_lcdfandlognormal_lccdfto rely onnormal_lcdf(likewise with some extra partials lines);exp_mod_normal_lcdfandexp_mod_normal_lccdf; andskew_normal_lcdfandskew_normal_lccdf.The
wiener5_lpdffunction already uses an approach somewhat like this.On the OpenCL side, there is also
std_normal_lcdf_scaled_impl, but I don't understand enough about OpenCL to know whether it would make sense to tie that tostd_normal_lcdfalso.If there is agreement that this kind of structure would be useful, I am happy to take a stab at it.