You can find a counterexample in this link, however I would like to note that for every $0<\beta<\alpha$ we do have convergence. Indeed, define for $x\neq y$ $$T(\epsilon,x,y)=\int_{B_1}\left|\varphi(z)\frac{g(x-\epsilon z)-g(x)-(g(y-\epsilon z)-g(y))}{|x-y|^\beta}\right|.$$
Case 1: $|x-y|\ge \delta_1$
In this case we have that for any $\eta>0$, there is $\epsilon_1>0$ such that if $\epsilon<\epsilon_1$ then
$$T(\epsilon,x,y)\leq 2\frac{\|\varphi\|_\infty}{\delta_1^{\alpha-\beta}}\left|B_1\right|\eta.\tag{1}$$
To prove $(1)$, use the uniform continuity of $g$ and the fact that $\frac{1}{|x-y|}\leq \frac{1}{\delta_1}$.
Case 2: $|x-y|<\delta_1$
In this case we use the $\alpha$-Hölder continuity of $g$ to conclude that for any $\epsilon>0$, $$T(\epsilon,x,y)\leq 2\|\varphi\|_\infty \left|B_1\right|\left|x-y\right|^{\alpha-\beta}\leq 2\|\varphi\|_\infty \left|B_1\right| \delta_1^{\alpha-\beta}.\tag{2}$$
We choose a suitable $\delta_1$ and combine $(1)$ and $(2)$ to conclude that $$g\star\varphi_\epsilon\to g\ \mbox{in}\ C^{\beta}(\overline{B_1}).$$
I know a different approach, and I think it is more "standard." Basically it is sufficient to prove that if $\Omega \subset \mathbb{R}^n$ is limited open, then $(-\Delta)^{-1}$ is a compact and injective self-adjoint operator on $L^2(\Omega)$ and on $H^1_{0}(\Omega)$ and then after applies Hilbert-Schmidt spectral theorem. Recalling a little proof of the Hilbert-Schmidt spectral theorem, I think that what you ask lies precisely in the proof of this theorem. More precisely there is the following theorem
"If $H$ is a real (it is true also in the complex case) Hilbert spaces and $K:H \longrightarrow H$ is a compact adjoint operator, then exists an orthonormal basis of eigenvectors $\lbrace u_n \rbrace$ of $K$ with eigenvalues $\lbrace \lambda_n \rbrace$ and it has the representation
$\displaystyle Kx = \sum_{n \geq 1} \lambda_n (x, u_n)_H u_n$ $(x \in H)$"
Now, the facts you say apply generally to elliptic operators in divergence form, ie type
$\displaystyle Lu:= - \sum_{i,j=1}^n (a_{ij} u_{x_i})_{x_j} + \sum_{i=1}^n (b_i u)_{x_i} +c u$
where $a_{ij}, b_i, c \in L^\infty(\Omega)$, and it assumes that $L$ is uniformly elliptic. Basically the case of the Laplace operator is a particular case of $L$. After we introduce the weak solutions, assume $b_i=c=0$, there is the following theorem
"If $a_{ij}=a_{ji} \in L^\infty(\Omega)$, and considering $L^{-1} : L^2(\Omega) \longrightarrow H^1_{0}(\Omega) \subset L^2(\Omega)$, then $L^{-1}$ is a compact and injective self-adjoint operator. In addition there is an orthonormal basis $\lbrace \phi_k : k \in \mathbb{N} \rbrace$ of $L^2(\Omega)$ of eigenfunctions associated to the eigenvalues of $L^{-1}$ and
$\displaystyle L^{-1}f =\sum_{k=1}^\infty \lambda_k (f, \phi_k)_{L^2} \phi_k$
where $Lu=f$ wih $f \in L^2(\Omega)$. In particular $\lim_{k \rightarrow \infty} \lambda_k =0$".
This whole theory is very well explained in the book "Lecture Notes on Functional Analysis: With Applications to Linear Partial Differential Equations" by A. Bressan.
Best Answer
Let $\psi$ be a standard non-negative positive mollifier, $C^\infty$, supported on the unit ball in $\mathbb{R}^n$, with $\int_{\mathbb{R}^n} \psi = 1$. Define $\psi_k(x) = k^n \psi(kx)$ as usual. Extend $f$ trivially outside the unit ball and define $\tilde f_k(x) = f(kx/(k-1))$ for all $x$, for $K > 1$. This function is now supported on the ball with radius $1 - k^{-1}$ and is non-negative and in $H^1_0$. Show that $\tilde f_k \to f$ in $H^1_0$. (It is enough to show weak convergence in $L^2$ and convergence of the norm, since we are working in a Hilbert space.)
Then define
$$ \varphi_k(x) = \int_{\mathbb{R}^n} \psi_k(x-y)\tilde f_k(y) dy $$ as usual. Then $\varphi_k$ is $C^\infty$, supported in the unit ball, and non-negative. Apply a triangle inequality argument to show that $\varphi_k \to f$ in $H^1_0$.