# Particle approximation to probability density functions: Dirac delta function representation

**Scientific Memo**, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

In the previous post, I have briefly shown the idea of using dirac delta function for discrete data representation. In the second example there, a histogram locations for a given set of points are presented as spike trains, where as heights are somehow given in a second sum. This is hard to follow and visualise, of course if you are not that good in reading formulation with different indexes. Pedagocial reasons, an easier representation of arbitrary probability density function (PDF), $p(x)$, one would simply need to couple each discrete points with a corresponding weight.

Hence, a set $\{x_{i}, \omega^{i}\}_{i=1}^{N}$ would be an estimation of PDF, $\hat{p}(x)$ . At this point we can invoke dirac delta function,

$ \hat{p}(x) = \sum_{i=1}^{N} \omega^{i} \delta(x-x_{i})$

Let’s revisit the R code given there, this time let’s draw uniform numbers between $[-2, 2]$ to get 100 x_{i} values. Simply these numbers will indicate the locations on the x-axis, a spike train. For simplicity, let’s use Gausian distribution for target PDF, $ \mathcal{N}(0, 1)$. Than, for weights we need to draw numbers using the spike locations. This approach is easier to understand compare to my previous double index notation.

**R Example code**

Above explained procedure is trivial to implement in R.

1 2 3 4 5 6 7 8 9 10 |
# Generate locations 100 x locations # out out 1000 points in [-2.0, 2.0] set.seed(42) # Domain where Dirac comb operates Xj = seq(-2,2,0.002) Xi = sample(Xj, 100) # Now generate weights from N(0,1) at those given locations Wi = dnorm(Xi) # Now visualise plot(Xi, Wi, type="h",xlim=c(-2.0,2.0),ylim=c(0,0.6),lwd=2,col="blue",ylab="p") |

**Conclusion**

Above notation introduces second abuse of notation while actually there must be a secondary regular grid that pics $x_{i}$ values using dirac delta in practice. Because, the argument of $\hat{p}(x)$, is in the discrete domain. So a little better notation that reflects the above code would be

$ \hat{p}(x_j) = \sum_{i=1}^{N} \omega^{i} \delta(x_{j}-x_{i})$

The set $x_j$ is simply defined in a certain domain, for example regularly. Hence I only recommend not to introduce dirac delta for explaining a particle approximation to PDFs for novice students in the class. It will only confuse them even more.

Figure: Spike trains with weights $\hat{p}(x) = \sum_{i=1}^{N} \omega^{i} \delta(x-x_{i})$ |

**leave a comment**for the author, please follow the link and comment on their blog:

**Scientific Memo**.

R-bloggers.com offers

**daily e-mail updates**about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.

Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.