The claimed expression for the determinant indeed holds. In fact, one may explicitly determine the eigenvalue decomposition of the infinite matrix $M_\infty$ viewed as a self-adjoint operator on a Hilbert space. Its eigenvalues are
$$ \left(\frac{2}{q_x^{p+\frac12}+q_x^{-p-\frac12}}\right)_{p=0}^\infty,\qquad q_k = e^{-\pi\frac{K\left(\sqrt{1-x^2}\right)}{K(x)}}, \tag{1}
$$ where $q_k$ is the elliptic nome of modulus $x$ and $K(x)$ the complete elliptic integral of the first kind.
Once this is known one may compute
$$
\begin{aligned}
\det(I - M_\infty^2) &= \prod_{p=0}^\infty\left(1-\left(\frac{2}{q_x^{p+\frac12}+q_x^{-p-\frac12}}\right)^2\right) \\
&=\prod_{p=0}^\infty\left(\frac{1-q_x^{2p+1}}{1+q_x^{2p+1}}\right)^2 = \frac{\theta_4(0,q_x)}{\theta_3(0,q_x)}=\sqrt[4]{1-x^2},
\end{aligned}
$$
where the last two equalities follow from standard properties of the Jacobi theta functions $\theta_i$.
To see how the eigenvalue decomposition can be obtained, it is convenient to look at the Dirichlet space $\mathcal{D}$ of complex analytic functions $f$ on the open unit disk that vanish at $0$ and that have finite norm with respect to the Dirichlet inner product
$$\langle f,g\rangle_{\mathcal{D}} = \frac{1}{\pi} \int_{|z|<1} \overline{f'(z)}g'(z) \mathrm{d}^2z = \sum_{n=1}^\infty n\,\overline{[z^n]f(z)}\,[z^n]g(z).$$
A basis is given by the monomials $(e_n:=z^n)_{n\geq 1}$, which satisfy $\langle e_n,e_m\rangle_{\mathcal{D}}=n\, \delta_{n,m}$. Let us also introduce the bounded linear operator on $\mathcal{D}$ given by
$$\mathbf{R}_x f = f \circ \psi_x, \qquad \psi_x(z) = \frac{1-\sqrt{1-x z^2}}{\sqrt{x}\,z}.$$
By power series expansion one may check that
$$ \mathbf{R}_x e_p = \sum_{\ell=p}^\infty \left(\frac{x}{4}\right)^{\ell/2} \frac{p}{\ell}\binom{\ell}{(\ell+p)/2}\, 1_{\{\ell+p\text{ even}\}} e_\ell$$
and the adjoint with respect to $\langle \cdot,\cdot\rangle_{\mathcal{D}}$ is determined by
$$ \mathbf{R}_x^\dagger e_p = \sum_{\ell=1}^p \left(\frac{x}{4}\right)^{p/2} \binom{p}{(\ell+p)/2}\, 1_{\{\ell+p\text{ even}\}} e_\ell.$$
In particular, $\mathbf{R}_x\mathbf{R}_x^\dagger$ is a self-adjoint operator that preserves the even and odd functions in $\mathcal{D}$.
Up to a factor of $2$ its matrix elements on the odd monomials are precisely $[M_\infty]_{i,j}$:
$$
\begin{aligned}
2\mathbf{R}_x\mathbf{R}_x^\dagger e_{2i+1} &= \sum_{j=0}^\infty \sum_{\ell=0}^{\min(i,j)} 2\left(\frac{x}{4}\right)^{i+j+1} \binom{2i+1}{i+\ell+1}\frac{2\ell+1}{2j+1} \binom{2j+1}{j+\ell+1} e_{2j+1}\\
&= \sum_{j=0}^\infty\frac{2i+1}{2(i+j+1)}\binom{i-1/2}i\binom{j-1/2}jx^{i+j+1}e_{2j+1}\\
&= \sum_{j=0}^\infty [M_\infty]_{i,j} e_{2j+1}.
\end{aligned}
$$
In T. Budd, Winding of simple walks on the square lattice, arXiv:1709.04042, Section 2.1, I examined a closely related operator $\mathbf{J}_x = \mathbf{R}_x^\dagger\mathbf{R}_x$, whose matrix elements count diagonal walks on $\mathbb{Z}^2$ starting on the positive $x$-axis and hitting the $y$-axis at prescribed height. By a convenient elliptic parametrization of the disc one can determine all eigenvectors of $\mathbf{J}_x$. Since $\mathbf{R}_x$ is injective, $\mathbf{R}_x\mathbf{R}_x^\dagger$ has the same eigenvalue decomposition (after applying $\mathbf{R}_x$ to the eigenvectors of $\mathbf{J}_x$) given by Proposition 9, of which the eigenvectors with odd label $m=1,3,5,\ldots$
span the odd functions in $\mathcal{D}$. The corresponding eigenvalues are precisely the ones given above in (1).
Added (22 Jan): One may also check the formula for $\det(I\pm M_\infty)$ proposed in the post by Hucht. Using the theta function product formulas we find
$$
\det(I-M_\infty) = \prod_{p=0}^\infty \frac{(1-q_x^{p+\frac12})^2}{1+q^{2p+1}} = \left(\frac{\theta_2(0,\sqrt{q_x})}{\theta_1'(0,\sqrt{q_x})}\right)^{\frac14} \frac{\theta_4(0,\sqrt{q_x})}{\theta_3(0,q_x)^{\frac12}}.
$$
With the help of $\theta_1'(0,\sqrt{q_x}) = \theta_2(0,\sqrt{q_x})\theta_3(0,\sqrt{q_x})\theta_4(0,\sqrt{q_x})$ and
$$
\theta_3(0,\sqrt{q_x})^2 = (1+x)\theta_3(0,q_x)^2,\quad \theta_4(0,\sqrt{q_x})^2 = (1-x)\theta_3(0,q_x)^2
$$
this yields
$$
\det(I- M_\infty) = \left(\frac{(1-x)^3}{1+x}\right)^{1/8}.
$$
Similarly
$$
\det(I-M_\infty) = \prod_{p=0}^\infty \frac{(1+q_x^{p+\frac12})^2}{1+q^{2p+1}} = \left(\frac{\theta_2(0,\sqrt{q_x})}{\theta_1'(0,\sqrt{q_x})}\right)^{\frac14} \frac{\theta_3(0,\sqrt{q_x})}{\theta_3(0,q_x)^{\frac12}}= \left(\frac{(1+x)^3}{1-x}\right)^{1/8}.
$$
