Here is a simple proof of the mean-value theorem for harmonic functions of two variables. We start with Green's Third Identity
$$u(\vec x_0)=\oint_C \left(u(\vec x')\frac{\partial G(\vec x_0,\vec x')}{\partial n'}-G(\vec x_0,\vec x')\frac{\partial u(\vec x')}{\partial n'}\right)\,d \ell' \tag 1$$
where $u(\vec x)$ is harmonic in a region $S$ bounded by the smooth contour $C$ and $G(\vec x_0,\vec x')-\frac{1}{2\pi}\log|\vec x_0-\vec x'|$ is the Green's Function for the Laplacian.
If we choose $S$ to be the region $|\vec x_0-\vec x'|\le r$, and thus $C$ to be the circle of radius $r$ and center $\vec x_0$, then on $C$ we have
$$\begin{align}
G(\vec x_0,\vec x')&=\frac{1}{2\pi}\log(r) \tag 2\\\\
\frac{\partial G(\vec x_0,\vec x')}{\partial n'}&=\frac{1}{2\pi r} \tag 3
\end{align}$$
Using $(2)$ and $(3)$ in $(1)$ yields
$$\begin{align}
u(\vec x_0)&=\frac{1}{2\pi r}\int_0^{2\pi} u(\vec x')\,d\ell'-\frac{1}{2\pi}\log(r)\oint_C \frac{\partial u(\vec x')}{\partial n'}d \ell' \tag 4\\\\
&=\frac{1}{2\pi r}\int_0^{2\pi} u(\vec x')\,d\ell'-\frac{1}{2\pi}\log(r)\int_{S} \nabla '\cdot \nabla 'u(\vec x')\,dS' \tag 5\\\\
&=\frac{1}{2\pi r}\int_0^{2\pi} u(\vec x')\,d\ell'-\frac{1}{2\pi}\log(r)\int_{S} \nabla '^2 u(\vec x')\,dS' \\\\
&=\frac{1}{2\pi r}\int_0^{2\pi} u(\vec x')\,d\ell' \tag 6
\end{align}$$
In going from $(4)$ to $(5)$, we invoked the Divergence Theorem, while in arriving at $(6)$ we used the fact that $u$ is harmonic.
Finally, if we multiply both sides of $(6)$ by $2\pi r$ and integrate with respect to $r$ we find
$$\pi r^2\,u(\vec x_0)=\int_S u(\vec x')\,dS'$$
whereupon dividing by $\pi r^2$ yields the coveted mean-value theorem
$$u(\vec x_0)=\frac{1}{\pi r^2}\int_{|\vec x_0-\vec x'|\le r} u(\vec x')\,dS'$$
