2a:null 2b:["$","title",null,{"children":"Skellam distribution | Mattias Villani | Observable"}] 31:[] 2e:{"id":"ff603c873e322a52","slug":"skellam-distribution","title":"Skellam distribution","thumbnail":"1718fd37cdda5bda8ecb1d3016af7b0c57a1a1b5e9d225fb575f2dc04e2afb7f","access_level":"public","create_time":"$D2023-07-21T13:57:56.729Z","update_time":"$D2024-01-02T11:20:03.654Z","view_time":null,"first_public_time":"$D2023-07-21T13:57:56.729Z","fork":{"id":"f01a7daa73a030c5","slug":"von-mises-distribution","title":"von Mises distribution","owner":{"id":"7592a7ab121f74a0","login":"mattiasvillani","name":"Mattias Villani"}},"stars":0,"starred":false,"owner":{"type":"team","id":"7592a7ab121f74a0","login":"mattiasvillani","name":"Mattias Villani","avatar_url":"https://avatars-next.observableusercontent.com/a440aaf39a3abefb86cbcc758fd7a5ab42555ae68d6825b52ec006f385a94d1d","tier":"pro_2024"},"version":1436,"publish_version":null,"role":null,"authored":false,"authors":[{"type":"individual","id":"56499be8cf43f0ef","login":"mattiasvillani","name":"Mattias Villani","avatar_url":"https://avatars-next.observableusercontent.com/a440aaf39a3abefb86cbcc758fd7a5ab42555ae68d6825b52ec006f385a94d1d","tier":null}],"members":[{"type":"individual","id":"56499be8cf43f0ef","login":"mattiasvillani","name":"Mattias Villani","avatar_url":"https://avatars-next.observableusercontent.com/a440aaf39a3abefb86cbcc758fd7a5ab42555ae68d6825b52ec006f385a94d1d","tier":null,"role":"editor"}],"comments":[],"custom_thumbnail":false,"trashed":false,"subscribed":null,"forks":1,"collection":{"id":"40c8ef5af1a51aff","slug":"distributions","title":"Statistical distributions","owner":{"id":"7592a7ab121f74a0","login":"mattiasvillani","name":"Mattias Villani"}},"type":"notebook","body":{"title":"Skellam distribution","cells":[{"id":38,"value":"# Skellam distribution","pinned":false,"hidden":false,"mode":"md"},{"id":629,"value":"

${tex`X\\sim \\mathrm{Skellam}(\\mu_1, \\mu_2)`}, for ${tex`X \\in \\{\\ldots,-2,-1,0,1,2,\\ldots\\}` }","pinned":false,"hidden":false,"mode":"md"},{"id":496,"value":"A random variable following the ${tex`\\mathrm{Skellam}(\\mu_1, \\mu_2)`} distribution is obtained as the difference between two independent Poisson variables with means ${tex`\\mu_1`} and ${tex`\\mu_2`}.
\nThe density function is\n\n${tex`f(x) = e^{-(\\mu_1+\\mu_2)} \\Big( \\frac{\\mu_1}{\\mu_2}\\Big)^{x/2} I_x(2\\sqrt{\\mu_1\\mu_2})`}\n\nwhere ${tex`I_x(\\cdot)`} is the [modified Bessel function](https://en.wikipedia.org/wiki/Bessel_function#Modified_Bessel_functions:_I%CE%B1,_K%CE%B1) of the first kind of order x.","pinned":false,"hidden":false,"mode":"md"},{"id":461,"value":"dist_quantile = tex`\\text{If } X \\sim \\mathrm{Skellam}(${mu1},${mu2}) \\text{ then } `","pinned":false,"hidden":false,"mode":"ojs"},{"id":555,"value":"moments = {\n \n return tex`\n \\begin{aligned}\n &E( X) = \\mu_1 - \\mu_2 = ${(mean).toPrecision(3)} \\\\[1em]\n &\\mathrm{Stdev}( X) = \\sqrt{ \\mu_1 + \\mu_2 } = \n ${(Math.sqrt(variance)).toPrecision(3)}\\\\[1em]\n &P( X\\leq ${quantile}) \\approx ${cdf.toPrecision(4)} \\\\\n \\end{aligned}\n`\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":118,"value":"viewof params = Inputs.form([\n Inputs.range([0, 5], {\n value: 2,\n step: 0.1,\n label: tex`\\mu_1`\n }),\n Inputs.range([0, 5], {\n value: 1,\n step: 0.1,\n label: tex`\\mu_2`\n }),\n Inputs.range([-10, 10], {\n value: -2,\n step: 1,\n label: \"quantile:\"\n })\n])","pinned":false,"hidden":false,"mode":"ojs"},{"id":305,"value":"plt_pdf = Plot.plot({\n color: { \n legend: false\n },\n x: {\n label: \"x\",\n axis: true,\n domain: d3.range(Math.floor(mean - 4*Math.sqrt(variance)), Math.floor(mean + 4*Math.sqrt(variance) ) )\n },\n y: {\n label: \"f(x)\",\n axis: true\n },\n tooltip: {\n fill: \"#FF1493\",\n stroke: \"#FF1493\",\n opacity: 1,\n },\n marks: [\n Plot.ruleY([0]),\n Plot.barY(pdfdata,{x: \"x\", y: \"pdf\", fill: \"pink\", strokeWidth: 0, opacity: 1,\n title: (d) => `P(X=${d.x}) = ${(d.pdf).toPrecision(4)}`}),\n Plot.barY(pdfdata, {filter: d => d.x <= params[2], x: \"x\", y: \"pdf\", fill: \"#FF69B4\", opacity: 1,\n title: (d) => `P(X=${d.x}) = ${(d.pdf).toPrecision(4)}`})\n ]\n })","pinned":false,"hidden":false,"mode":"ojs"},{"id":109,"value":"
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n
\n\n---\n## 🙈 Appendix 🙈\n\nHere is some code and data details that you probably do not want to read.","pinned":false,"hidden":false,"mode":"md"},{"id":221,"value":"#### Imports","pinned":false,"hidden":false,"mode":"md"},{"id":1375,"value":"import {Plot} from \"@mkfreeman/plot-tooltip\"","pinned":true,"hidden":false,"mode":"ojs"},{"id":26,"value":"jstat = require('jstat')","pinned":false,"hidden":false,"mode":"ojs"},{"id":750,"value":"mathjs = require('mathjs')","pinned":false,"hidden":false,"mode":"ojs"},{"id":1173,"value":"besselfunc = require(\"bessel@1.0.2/bessel.js\")","pinned":false,"hidden":false,"mode":"ojs"},{"id":223,"value":"#### Data","pinned":false,"hidden":false,"mode":"md"},{"id":654,"value":"mu1 = params[0]","pinned":false,"hidden":false,"mode":"ojs"},{"id":656,"value":"mu2 = params[1]","pinned":false,"hidden":false,"mode":"ojs"},{"id":662,"value":"quantile = params[2]","pinned":true,"hidden":false,"mode":"ojs"},{"id":706,"value":"mean = mu1-mu2","pinned":false,"hidden":false,"mode":"ojs"},{"id":766,"value":"variance = mu1 + mu2","pinned":false,"hidden":false,"mode":"ojs"},{"id":1199,"value":"gridsize = 0.01","pinned":false,"hidden":false,"mode":"ojs"},{"id":830,"value":"function skellampdf(x, mu1, mu2) {\n return Math.exp(-(mu1+mu2)) * (mu1/mu2)**(x/2) \n * besselfunc.besseli(2*Math.sqrt(mu1*mu2), Math.abs(x));\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":845,"value":"function skellamcdf(quantile, mu1, mu2){\n return d3.sum(pdfdata.filter(d => d.x <=quantile).map(d => d.pdf))\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":17,"value":"pdfdata = {\n const x = d3.range(Math.floor(mean - 10*Math.sqrt(variance)), Math.floor(mean + 10*Math.sqrt(variance) ) )\n const pdf = x.map((x) => ({ x: x, pdf: skellampdf(x, mu1, mu2) }));\n return pdf;\n}","pinned":false,"hidden":false,"mode":"ojs"},{"id":19,"value":"cdf = skellamcdf(quantile, mu1, mu2)","pinned":false,"hidden":false,"mode":"ojs"},{"id":229,"value":"#### Style","pinned":false,"hidden":false,"mode":"md"},{"id":202,"value":"","pinned":false,"hidden":false,"mode":"html"},{"id":58,"value":"viewof md_colour = Inputs.color({label: \"Highlight color cells\", value: \"#F5F5DC\"})","pinned":false,"hidden":false,"mode":"ojs"}],"files":"$W31","stdlib":"1","titling":"auto"}} 2f:"$2e:body"