As an operator on functions, one intuitive way to think about the Laplacian seems to be as an operator that returns the average difference between a function's value at a point and the values of its neighbouring points.
For example, in $\mathbb{R}^n$, this answer explains how the Laplacian $\Delta f(\mathbf{0})$ at $\mathbf{0}$ of a function $f \in C^2(\mathbb{R}^n)$ is limit of the (scaled) average difference $$ \Delta f(\mathbf{0}) = \lim_{r \to 0^+} \dfrac{2nr}{r^2} \dfrac{1}{\omega(S_r)} \int_{S_r} (f(\mathbf{x}) - f(\mathbf{0}) ) d\omega(\mathbf{x}) $$ in a sphere around $\mathbf{0}$.
Suppose now we have a Riemannian manifold $M$, and we construct the De Rham complex $$ \Omega^0(M) \xrightarrow{d^0} \Omega^1(M) \xrightarrow{d^1} \Omega^2(M) \xrightarrow{d^2} \ldots $$ where $d$ is the exterior derivative. The Riemannian structure allows us to define an adjoint to the exterior derivative and the Laplace-Hodge operator $$ \Delta_n = (d^n)^\dagger d^n + d^{n-1} (d^{n-1})^\dagger $$ for any $n \geq 0$.
My Question: For $n > 0$, is there a way to think of $\Delta_n$ as comparing local averages of $n$-forms, similar to the way that $\Delta_0$ does for smooth functions?
In particular, any ideas on whether this works when $n =1$ would be much appreciated!