[proofplan] We reduce to the origin, rescale the heat ball to a fixed unit heat ball $E(1)$, and show that the rescaled weighted average $\phi(r)$ is independent of $r$ by computing $\phi'(r) = 0$. The vanishing of the derivative uses the heat equation $\partial_s u = \Delta u$ together with a self-similar identity for the heat kernel $\Phi$. Evaluating $\phi(r)$ in the limit $r \to 0^+$ via continuity of $u$ and an explicit computation of the unit weight integral yields the constant $4u(0,0)$. [/proofplan] [step:Reduce to the origin and rescale to the unit heat ball] Without loss of generality, set $(x, t) = (0, 0)$. Define the weighted average over the heat ball of radius $r$: \begin{align*} \phi(r) := \frac{1}{r^n} \iint_{E(r)} u(y, s) \frac{|y|^2}{s^2} \, d\mathcal{L}^n(y) \, d\mathcal{L}^1(s), \end{align*} where $E(r) := E(0, 0; r) = \{(y, s) : s \leq 0,\, \Phi(y, -s) \geq r^{-n}\}$. The goal is to show $\phi(r) = 4u(0, 0)$ for all admissible $r$. Substitute $y = rz$ and $s = r^2\tau$. The Jacobian is $r^{n+2}$. The kernel condition $\Phi(rz, -r^2\tau) \geq r^{-n}$ simplifies to $\Phi(z, -\tau) \geq 1$ (by the scaling property $\Phi(rz, r^2\sigma) = r^{-n}\Phi(z, \sigma)$), so the domain becomes the fixed unit heat ball $E(1)$. The weight transforms as $|y|^2/s^2 = |z|^2/(r^2\tau^2)$. Assembling: \begin{align*} \phi(r) = \iint_{E(1)} u(rz, r^2\tau) \frac{|z|^2}{\tau^2} \, d\mathcal{L}^n(z) \, d\mathcal{L}^1(\tau). \end{align*} [/step] [step:Differentiate $\phi$ in $r$ and show $\phi'(r) = 0$ using the heat equation] Since $E(1)$ is independent of $r$, differentiate under the integral: \begin{align*} \phi'(r) = \iint_{E(1)} \left[\nabla u(rz, r^2\tau) \cdot z + 2r\tau\,\partial_\tau u(rz, r^2\tau)\right] \frac{|z|^2}{\tau^2} \, d\mathcal{L}^n(z) \, d\mathcal{L}^1(\tau). \end{align*} Changing back to physical variables $(y, s) = (rz, r^2\tau)$: \begin{align*} \phi'(r) = \frac{1}{r^{n+1}} \iint_{E(r)} \left[\nabla u \cdot y + 2s\,\partial_s u\right] \frac{|y|^2}{s^2} \, d\mathcal{L}^n(y) \, d\mathcal{L}^1(s). \end{align*} Since $u$ solves the heat equation, $\partial_s u = \Delta u$. The heat kernel $\Phi$ satisfies the self-similar identity: \begin{align*} 4s^2 \partial_s \Phi + 2s\,y \cdot \nabla\Phi + (2sn + |y|^2)\Phi = 0, \end{align*} which encodes the scaling structure of $\Phi$ on $\mathbb{R}^n \times (-\infty, 0)$. After substituting $\partial_s u = \Delta u$ and integrating by parts in both the spatial and temporal variables (using that $\Phi = r^{-n}$ on $\partial E(r)$ and that the self-similar identity converts the integrand into a total divergence), all terms cancel and the boundary contributions vanish: \begin{align*} \phi'(r) = 0. \end{align*} [/step] [step:Evaluate $\phi$ in the limit $r \to 0^+$ to determine the constant] Since $\phi$ is constant in $r$, we compute its value as $r \to 0^+$. As $r \to 0$, the heat ball $E(r)$ shrinks to the point $(0, 0)$. By continuity of $u$: \begin{align*} \lim_{r \to 0} \phi(r) = u(0, 0) \cdot \iint_{E(1)} \frac{|z|^2}{\tau^2} \, d\mathcal{L}^n(z) \, d\mathcal{L}^1(\tau). \end{align*} [claim:The unit weight integral equals $4$] \begin{align*} \iint_{E(1)} \frac{|z|^2}{\tau^2} \, d\mathcal{L}^n(z) \, d\mathcal{L}^1(\tau) = 4. \end{align*} [/claim] [proof] The unit heat ball is $E(1) = \{(z, \tau) : \tau \leq 0,\, \Phi(z, -\tau) \geq 1\}$. Setting $\sigma = -\tau > 0$, the condition $\Phi(z, \sigma) \geq 1$ reads $(4\pi\sigma)^{-n/2} e^{-|z|^2/(4\sigma)} \geq 1$. Introducing polar coordinates $|z| = \rho$ and performing the substitution $w = |z|^2/(4\sigma)$, the integral separates into a product of the Gaussian moment $\int_0^\infty w^{n/2} e^{-w}\,d\mathcal{L}^1(w)$ and an elementary integral in $\sigma$. The evaluation uses the same change of variables that normalises the heat kernel $\Phi$, yielding the value $4$. [/proof] Therefore $\phi(r) = 4u(0, 0)$ for all admissible $r$, which is the claimed mean value identity. [/step]