Given the density function of the Gaussian distribution $\mathcal{N}(0, \sigma)$:
$$ \frac{1}{\sigma\sqrt{2\pi}} \exp\left(-0.5\,\frac{x^2}{\sigma^2}\right)$$
I have noted that emprically, summing discretized values, with $x\in\mathbb{Z}$, leads to 1 irrespectively of the value of $\sigma$ (as Michael has pointed out in the comments, for certain values of $\sigma<1$ this does not hold), i.e.:
$$
\begin{align}
1 &= \frac{1}{\sigma\sqrt{2\pi}} \sum_{x=-\infty}^\infty \exp\left(-0.5\,\frac{x^2} {\sigma^2}\right), \qquad x\in\mathbb{Z} \\[10pt]
%
&= \frac{1}{\sigma\sqrt{2\pi}} \sum \left( \ldots,\,\exp\left(-0.5\,\frac{(-1)^2} {\sigma^2}\right),\, 1,\, \exp\left(-0.5\,\frac{1^2} {\sigma^2}\right),\,\ldots \right)
\end{align} \tag{1}\label{1}
$$
Is this a known fact? My concern is that I am not sure if this this is true since I was not able to prove it / find a proof.
Motivation: In image processing, when applying a Gaussian filter to an image, the filter/kernel used must properly average the neighbouring pixel values (e.g. this great post explains its importance). Thereby the sum of the weights (filter/kernel values) must sum to $1$.
As a result, by sampling enough values ($\geq ceil(3\sigma)+1$, given the increasing decay of the function), and if we followed eq. $\eqref{1}$, we would end-up always with a proper filter without the need of any normalization.
This is the result of the empirical test that I've mentioned above:
where the number of samples (x-axis) to extract from the normal distribution is determined by $n = 2\text{ceil}(\text{cutoff}\cdot\sigma) + 1$.
As can be seen, by growing the number of samples (increasing the cutoff), we end up with a sum that tends to $1$.
The python code for reproducing it is the following:
from math import sqrt, ceil
import numpy as np
import matplotlib.pyplot as plt
def sample_normal_mean0(sigma, cutoff=3.0):
"""
Create an array, centered at 0, of length determined by the cutoff freqency
"""
# length of the array (number of samples)
n = 2 * ceil(cutoff * sigma) + 1
# discrete values
x = np.linspace(-(n//2), n//2, n)
# density function at these values
f = np.exp(-0.5 * (x**2) / (sigma**2)) / (sigma * sqrt(2. * np.pi))
return f
def emprical_test(cutoffs, sigma, ax):
""" Plot the sum of the discretized values as a function of \sigma """
# XY coordinates:
x = 2 * np.ceil(cutoffs * sigma) + 1
y = [sample_normal_mean0(sigma, cutoff_i).sum() for cutoff_i in cutoffs]
# plot
ax.plot(x, y, '-o', lw=2, label=f'$\sigma$={sigma:.1f}')
if __name__ == "__main__":
# standard deviations
sigmas_to_test = np.linspace(start=1, stop=20, num=20)
# cutoff frequencies
cutoffs_to_test = np.arange(start=1, stop=10, step=1.0)
# create plot
fig, ax = plt.subplots()
for sigma_i in sigmas_to_test:
emprical_test(cutoffs_to_test, sigma_i, ax)
ax.set_xlabel('#samples = 2 * ceil(cutoff * sigma) + 1', fontsize=12)
ax.set_ylabel('Discretized sum', fontsize=12)
ax.legend(loc='center right', ncol=2)
ax.grid()
plt.show()
Best Answer
You are looking at a consequence of the remarkable Poisson Summation Formula. In your context, with $f_\sigma(x)$ denoting the above Gaussian distribution density, its Fourier transform (using the normalization as in the link) is given by $$ \hat{f}_\sigma(w) = \exp\left(-\frac12 (2\pi \sigma)^2 w^2\right) ,$$ and the P-S formula simply states that $$ \sum_{k\in {\Bbb Z}} f_\sigma(k) = \sum_{w\in {\Bbb Z}} \hat{f}_\sigma(w).$$
As an example, for $\sigma=1$ you get $$ \sum_{k\in {\Bbb Z}} f_1(k) = \sum_{w\in {\Bbb Z}} \exp(-2\pi^2 w^2) = 1 + 2 \exp(-2\pi^2) + 2 \exp(-8\pi^2) + \cdots.$$ The second term on the RHS is of order $5.35 \times 10^{-9}$ and the third of order $10^{-34}$. Approximating by $1$ is good for $\sigma$ large, but bad when $\sigma$ is small.