# Double differential distribution using set of data

Consider a set of data having the form $\{\{x_{1}, y_{1}\}, …,\{x_{n},y_{n}\}\}$. Each point from the set corresponds to 1 event with coordinate $x, y$. I'm looking for a procedure in Mathematica that allows me to automatically create bins within given intervals in $x$ and $y$ and then to compute double differential events distribution $d^{2}N/dxdy$. Or, equivalently, a procedure that returns a table with three columns in the form $$\{\{\bar{x}_{1}, \bar{y}_{1},N_{\text{events}|_{\bar{x}_{1},\bar{y}_{1}}}\},\{\bar{x}_{2}, \bar{y}_{2},N_{\text{events}}|_{\bar{x}_{2},\bar{y}_{2}}\},...\}$$ Here $\bar{x}_{i}, \bar{y}_{i}$ denotes $i$th bin.

I expect that somehow this can be done by using Histogram function, but I am not sure about this. Could you please tell me how to do this?

• Can you please explain what you mean by "double differential events distribution"? Do just what the number of events in each bin? – Henrik Schumacher Sep 8 '18 at 15:19

I don't know what double differential events are but if a histogram is acceptable which implies some sort of random sample from an underlying smooth distribution, how about a nonparametric density estimate?

xy = RandomVariate[BinormalDistribution[{0, 0}, {1, 3}, 0.6], 1000];
skd = SmoothKernelDistribution[xy];
Show[ListPlot[xy, PlotRange -> All],
ContourPlot[PDF[skd, {x, y}], {x, -4, 4}, {y, -8, 8}, ContourShading -> None]]


Plot3D[PDF[skd, {x, y}], {x, -4, 4}, {y, -8, 8},
BoxRatios -> {1, 2, 1}, PlotRange -> All]


Maybe something like this?

xbins = {-1, 1, 0.1};
ybins = {-1, 1, 0.2};
pts = RandomReal[{-1, 1}, {1000, 2}];
tab = Join[
Tuples[{Most[Range @@ xbins], Most[Range @@ ybins]}],
Partition[Flatten[BinCounts[pts, xbins, ybins], 1], 1],
2
]


HistogramList will also get you the frequency counts (or the associated PDF value if desired) but there's a bit of work to get the midpoints of the bins (which I assume is what you mean by $(\bar{x}_1, \bar{y}_i)$):

* Generate some data *)
pts = RandomReal[{-1, 1}, {1000, 2}];

(* Get histogram counts *)
xbins = {-1, 1, 0.1};
ybins = {-1, 1, 0.2};
h = HistogramList[pts, {xbins, ybins}]

(* Midpoints of bins *)
nx = Length[h[[1, 1]]]
ny = Length[h[[1, 2]]]
xmid = (h[[1, 1, 2 ;; nx]] + h[[1, 1, 1 ;; nx - 1]])/2
ymid = (h[[1, 2, 2 ;; ny]] + h[[1, 2, 1 ;; ny - 1]])/2

(* {xmidpoint, ymidpoint, count} *)
counts = Flatten[Table[{xmid[[i]], ymid[[j]], h[[2, i, j]]}, {i, nx - 1}, {j, ny - 1}], 1]