Simple question regarding the Green's function for the diffusion equation

The differential operator for diffusion in three dimensions is given by $\partial_t - k \nabla^2$ where $k$ is a constant. The Green's function is (according to Wikipedia) $$\theta(t)\left( \frac{1}{4\pi k t} \right)^{3/2} e^{-r^2/4kt}$$ where $\theta$ is the Heaviside step function and $r = |\mathbf{r}|$. Hence, we should have $$(\partial_t - k \nabla^2)\theta(t)\left( \frac{1}{4\pi k t} \right)^{3/2} e^{-r^2/4kt} = \delta(t)\delta^3(\mathbf{r}).$$ I've played around with this for a while, but I simply can't see how evaluating the left-hand side in the above expression gives the right-hand side. Where do the delta functions even come from?

Let's consider the expression we obtain by removing the Heaviside step function from the Green's function. It is a solution to the diffusion equation, viz., $$(\partial_t - k\nabla^{2}) \left(\frac{1}{4\pi k t}\right)^{3/2} e^{-r^2/4kt} = 0$$ Furthermore, one can show that \begin{equation} \lim_{t\rightarrow 0^{+}} \left(\frac{1}{4\pi k t}\right)^{3/2} e^{-r^2/4kt} = \delta^{(3)}(\textbf{r}). \end{equation}
Now we can see that the differential operator $\partial_t - k\nabla^{2}$ almost annihilates the Green's function except when the time derivative acts on $\theta(t)$.
\begin{equation} \begin{split} &(\partial_t - k\nabla^{2})\left[ \theta(t) \left(\frac{1}{4\pi k t}\right)^{3/2} e^{-r^2/4kt} \right]\\ &=\delta(t) \left(\frac{1}{4\pi k t}\right)^{3/2} e^{-r^2/4kt}. \end{split} \end{equation} Then, as $\delta(t)$ demands that $t=0$ in the rest of the expression, what we have obtained is actually just $\delta(t)\delta^{(3)}(\textbf{r})$.