We in fact have a finite bound for all $\lambda>0$ and even just assuming $\phi\in L^{1}$. We have the L2 of the convolution
$$\int_{0}^{T}\int |P_{u}\ast\phi_{\lambda}(y)|^{2}dy du.$$
We will use The Dirichlet heat semigroup, $L^1_\delta$, and the dimension shift phenomenon we have nice Lp estimates
For $\varphi \in L^p(U)$ and $1\le p \le q < \infty$, we have $$\| e^{-t A} \varphi\|_{q} \lesssim \|\varphi\|_p t^{-\frac{n}{2}(\frac{1}{p}-\frac{1}{q})}, \quad t >0.$$
and also found here
Let $N_t:\mathbb{R}^N\to \mathbb{R}$, $t>0$, be the function defined by
$$N_t(x)=(4\pi t)^{-N/2}e^{-|x|^2/4t}.$$
Since
$$\int_{\mathbb{R}^N} e^{-a|x|^2}dx=\left(\frac{\pi}{a}\right)^{N/2},\tag{1}\label{1}$$
we can see that $N_t\in L^1(\mathbb{R}^N)$ and $\|N_t\|_{ L^1(\mathbb{R}^N)}=1$.
We know that $S(t)v=N_t\ast v$. From Young's Inequality, we have
$$\|S(t)v\|_{ L^p(\Omega)}\leq \|N_t\ast v\|_{ L^p(\Omega)}\leq \|N_t\|_{ L^m(\Omega)}\|v\|_{ L^q(\Omega)},$$
where $1+\frac{1}{p}=\frac{1}{m}+\frac{1}{q}.$
Now, we just have to estimate $\|N_t\|_{ L^m(\Omega)}$. From \eqref{1}, we can see that
$$\|N_t\|_{ L^m(\Omega)}=(4\pi t)^{-N/2}\left(\int_{\mathbb{R}^N} e^{-\frac{m}{4t}|x|^2}dx\right)^{1/m}=(4\pi t)^{-N/2}\left(\frac{\pi}{\frac{m}{4t}}\right)^{N/2m}=C_{m,N}t^{-\frac{N}{2}\left(1-\frac{1}{m}\right)}=C_{m,N}t^{-\frac{N}{2}\left(\frac{1}{q}-\frac{1}{p}\right)}.$$
Hence, we have the result.
and "Parabolic L p − L q estimates" by Dietmar A. Salamon. We set $q=2$ and so we get
$$\int^{T}\int |P_{u}\ast\phi_{\lambda}(y)|^{2}dydu\leq \lambda^{-\frac{p-1}{p}2}\int^{T} u^{-(\frac{2-p}{2p})}du=\lambda^{-\frac{p-1}{p}2}T^{1-(\frac{2-p}{2p})}.$$
To bring in $\beta$, we set $2\beta=1+\frac{p-1}{p}2$ and $2\epsilon:=(3p-2)/2p$. So to get positive $\epsilon>0$ we need $2\geq p>2/3$ (which is true since $p\geq 1$). This restricts $\beta\in [\frac{1}{2},1]$.
And if we take $p=1$, the $\lambda$ disappears and so we have a finite bound and so it works for all $\beta$.
Another way to see it is that for larger $\beta>1$, since it shows up as negative and $\lambda\in (0,1)$, we can just bound
$$\lambda^{1/2-\tilde{\beta}}<\lambda^{1/2-\beta}$$
where $\tilde{\beta}\in [1/2,1]$