Let
$$
\phi_{r,x}(y)=-\frac{x^2}{2y^{2r}}-\frac{y^2}{2}\tag{1}
$$
Taking the first and second derivatives of $\phi_{r,x}(y)$ yields
$$
\phi_{r,x}^\prime(y)=r\frac{x^2}{y^{2r+1}}-y\tag{2}
$$
and
$$
\phi_{r,x}^{\prime\prime}(y)=-(2r+1)r\frac{x^2}{y^{2r+2}}-1\tag{3}
$$
Using $(2)$, $\phi_{r,x}(y)$ reaches a maximum at $y_0=(rx^2)^{\frac{1}{2r+2}}$. At that point,
$$
\phi_{r,x}(y_0)=-\frac{r+1}{2r}(rx^2)^{\frac{1}{r+1}}\tag{4}
$$
Furthermore, $(3)$ gives that
$$
\frac12\phi_{r,x}^{\prime\prime}(y_0)=-(r+1)\tag{5}
$$
Standard stationary phase methods yield
$$
\begin{align}
&\int_0^\infty\exp\left(-\frac{x^2}{2y^{2r}}-\frac{y^2}{2}\right)\frac{\mathrm{d}y}{y^s}\\
&\sim\exp\left(\phi_{r,x}(y_0)\right)\int_0^\infty\exp\left(-(r+1)(y-y_0)^2\right)\frac{\mathrm{d}y}{y^s}\\
&\sim y_0^{-s}\exp\left(\phi_{r,x}(y_0)\right)\int_{-\infty}^\infty\exp\left(-(r+1)y^2\right)\mathrm{d}y\\
&=(rx^2)^{\frac{-s}{2r+2}}\exp\left(-\frac{r+1}{2r}(rx^2)^{\frac{1}{r+1}}\right)\sqrt{\frac{\pi}{r+1}}\tag{6}
\end{align}
$$
Where $f(x)\sim g(x)$ means that $\lim\limits_{x\to\infty}f(x)/g(x)=1$.
Estimate $(5)$ says that the kind of estimate sought above can be achieved only when $s\le r+1$.
Taking the derivative of $(2)$ yields
$$
\phi_{r,x}^{\prime\prime\prime}(y)=(2r+1)(2r+2)r\frac{x^2}{y^{2r+4}}\tag{7}
$$
which says that the second derivative of the exponent increases monotonically, whereas the second derivative of the quadratic approximation is constant. Since $\phi_{r,x}$ and its first and second derivatives match the quadratic approximation at $y_0$, we get that for $y\ge y_0$,
$$
\phi_{r,x}(y)\ge\phi_{r,x}(y_0)-(r+1)(y-y_0)^2\tag{8}
$$
Furthermore, since $(1+t)^{-s}\ge1-st$ for $t\ge0$, we get
$$
\begin{align}
&\int_0^\infty\exp\left(-\frac{x^2}{2y^{2r}}-\frac{y^2}{2}\right)\frac{\mathrm{d}y}{y^s}\\
&\ge\int_{y_0}^\infty\exp\left(\phi_{r,x}(y_0)-(r+1)(y-y_0)^2\right)y_0^{-s}\left(1-s\frac{y-y_0}{y_0}\right)\mathrm{d}y\\
&=y_0^{-s}\exp(\phi_{r,x}(y_0))\int_0^\infty\exp\left(-(r+1)t^2\right)\left(1-\frac{st}{y_0}\right)\mathrm{d}t\\
&=y_0^{-s}\exp(\phi_{r,x}(y_0))\left(\frac12\sqrt{\frac{\pi}{r+1}}-\frac{s}{2y_0\sqrt{r+1}}\right)\\
&=(rx^2)^{\frac{-s}{2r+2}}\exp\left(-\frac{r+1}{2r}(rx^2)^{\frac{1}{r+1}}\right)\left(\frac12\sqrt{\frac{\pi}{r+1}}-\frac{s}{2y_0\sqrt{r+1}}\right)\tag{9}
\end{align}
$$
For $x\ge x_0$, we get $(9)$ with $y_0=(rx_0^2)^\frac{1}{2r+2}$. This is the bound required as long as $s\le r+1$.
For example, if we set $\displaystyle x_0=\max\left(\frac{s^{r+1}}{\sqrt{r}},1\right)$, for $x\ge x_0$, we get
$$
\begin{align}
&\int_0^\infty\exp\left(-\frac{x^2}{2y^{2r}}-\frac{y^2}{2}\right)\frac{\mathrm{d}y}{y^s}\\
&\ge\left(\frac{\sqrt{\pi}-1}{2\sqrt{r+1}}r^{\frac{-s}{2r+2}}\right)x^{\frac{-s}{r+1}}\exp\left(-\frac{r+1}{2r}(rx^2)^{\frac{1}{r+1}}\right)\\
&\ge\left(\frac{\sqrt{\pi}-1}{2\sqrt{r+1}}r^{\frac{-s}{2r+2}}\right)\frac1x\exp\left(-\frac{r+1}{2r}(rx^2)^{\frac{1}{r+1}}\right)\tag{10}
\end{align}
$$
as long as $s\le r+1$.