Model PDF and PDF slices

[CristianAntochi]

I was thinking that we do not provide any easy to use tool to retrieve the PDF of a model in Flamedisx.

As is, it seems really difficult (or is it even possible?) to retrieve the PDF of a model. Say one fits an ER model with some data and wants to analyse the model in detail, look at the extreme tails etc, and to do this the best way would be to retrieve the best fit model pdf in some coordinates, say, s1-s2. With Flamedisx at the moment one can get rather easily the differential rate of the model at the given data points.
Does anyone have any idea how to maybe get a pdf evaluation of the model in the flamedisx observables (s1,s2,r,z,t)?
On the other hand multidimensional PDFs are weird objects, and not easy to evaluate. Could maybe a tool that could evaluate and interpolate the differential rate in some dimensions while summing(integrating) on the others, be a solution?

[JelleAalbers]

Hi Cristian, if you're asking for the full PDF in all observables, you could create events at the points you're interested in, ask flamedisx for the differential_rate, then divide by the mu estimate to normalize to a proper density. We could certainly make a utility function to simplify that.

But for most applications (e.g. plotting something), I imagine you'll want the PDF integrated along some dimensions instead -- e.g. the (S1, S2) or perhaps the (cS1, cS2, r) PDF. For that, the easiest way is to simulate many events and histogram them. Actually integrating differential_rate is also possible in theory, e.g. for an (S1, S2) PDF, you could compute something like

   mu^-1 * Integrate dx dy dz d[x_observed] d[y_observed] ... differential_rate(S1, S2, x, y, ...) * delta^(3)(field_distortion(x,y,z) - [x_obs, y_obs, z_obs]) ...

with possible additional delta functions for the relations between observables. But this is going to be super slow. For high-dimensional integrals, Monte Carlo methods are often good, e.g. you could solve the integral above by feeding the integrandto a (Hamiltonian) MCMC. This will then do a lot of work just to get something that samples from the (integrated/marginalized/howeveryoucallit) PDF, while we already have a way of sampling directly: simulate! Simulate is ~1000x faster (~us/event) than even a single differential_rate call (~ms/event).