The Jacobi theta function $\theta_4$ is defined by $$\displaystyle \theta_4(q)=\sum_{n \in \mathbb{Z}} (-1)^n q^{n^2} \tag{1}$$ For this question, set $q=\large e^{-\pi x}$ and $\theta_4 \equiv \theta_4(q)$. Define $\theta_3(q)=\theta_4(-q)$.
Using Lambert-Series representation for powers of $\theta_4$ (which I will describe in a moment) and integrating term by term, I have obtained a family of neat identites: $$ \int_0^{\infty} \frac{\theta_4^2}{1+x^2} dx=1 \tag{2}$$ $$ \int_0^{\infty} \frac{\theta_4^4}{1+x^2} dx=\frac{4 \ln2}{\pi} \tag{3}$$ $$ \int_0^{\infty} \frac{\theta_4^6}{1+x^2} dx=\frac{16 G}{\pi^2}-\frac23 \tag{4}$$ $$ \int_0^{\infty} \frac{\theta_4^8}{1+x^2} dx=\frac{20 \zeta(3)}{\pi^3} \tag{5}$$ Here, $G$ is Catalan's constant. These identities check out to a lot of digits, according to Mathematica.
As I mentioned, in order to drive these identities I used Lambert series representations for powers of $\theta_4$ which I have found online. For instance:
$$\theta_4^2 = 1+4\sum_{n=1}^{\infty} \frac{(-1)^n q^n}{1+q^{2n}} \tag{6}$$ (See e.g. [1], or use $r_2(n)=4\sum_{d \mid n} \sin(\frac{\pi}{2} d)$ ([2]), and switch order of summation)
$$\theta_4^6=1+16\sum_{n=0}^{\infty}\frac{(-1)^n n^2 q^n}{1+q^{2n}}+4\sum_{n=0}^{\infty} \frac{(-1)^n (2n+1)^2 q^{2n+1}}{1+q^{2n+1}} \tag{7}$$ (A proof is given in [3], together with proofs for similar formulas for the $4$th power of $\theta_4$, and the $8$th power.)
I am interested in a closed form for the integral $$I(n)= \int_0^{\infty} \frac{\theta_4(e^{-\pi x})^n}{1+x^2} dx.\tag{8}$$
1st Question
I have searched the web for a similar Lambert series combination for $\theta_4^{10}$, but all my efforts remain fruitless. Also, I wasn't able to derive one using the approach taken in [3]. Can we find a formula, possibly in the form of $(6)$ or $(7)$, for $\theta_4^{10}(q)$? Such a formula may be used to find the closed form of $I(10)$.
2nd Question
Can we find closed forms for $I(n)$ for other values of $n$?
In general, can we find a Lambert series representation of $\theta_4^n$ for each even $n$? If not, for which $n$ we can find one and for which we cannot?
References
[1]: Proving $\left(\sum_{n=-\infty}^\infty q^{n^2} \right)^2 = \sum_{n=-\infty}^\infty \frac{1}{\cos(n \pi \tau)}$ ,(Note that $\theta_4(q)=\theta_3(-q)$)
[2]: Eric W. Weisstein, Sum of Squares Function,Mathworld.(25),(The line below references proofs). Link
[3]: George E. Andrews, Richard Lewis and Zhi-Guo Liu, An identity relating a theta function to a sum of Lambert series, (7)-(9).
This question has been edited. For more information about how I obtained these results, check the original version of this post.
$\displaystyle \large \mathbf{Progress}\,\mathbf{Report}$
Thanks to Paramanand Singh (see his answer below) I've been able to make some progress. Again, we set $q=e^{-\pi x}$ and $$\theta_2 \equiv \theta_2(q)= \sum_{n \in \mathbb{Z}} q^{(n+\frac12)^2}\\\theta_3\equiv \theta_3(q)=\sum_{n \in \mathbb{Z}} q^{n^2}\\\theta_4\equiv \theta_4(q)=\sum_{n \in \mathbb{Z}}(-1)^n q^{n^2}$$ We also have the relations $$ \theta_2^4+\theta_4^4=\theta_3^4 \\ \theta_2(e^{-\pi/x})=\frac1{\sqrt{x}}\theta_4(e^{-\pi x}) \\ \theta_3(e^{-\pi/x})=\frac1{\sqrt{x}}\theta_3(e^{-\pi x})\tag{9} $$ I will omit all the intermediate steps and just state what I've found so far, and also some other related identities which I've not proven, but match numerically. So, using Singh's formulas (changing $q$ to $-q$ and rewriting in terms of hyperbolic functions), I get that $$\int_0^{\infty} \frac{\theta_4^{10}}{1+x^2}dx=\frac{768 \beta(4)}{5 \pi^4}-\frac{32}{75}-2\int_0^{\infty} \frac{\theta_2^4\theta_4^6}{1+x^2}dx\tag{10}$$ $$\int_0^{\infty} \frac{\theta_4^{12}}{1+x^2}dx=\frac{450 \zeta(5)}{\pi^5}-\int_0^{\infty} \frac{\theta_2^4\theta_3^4\theta_4^4}{1+x^2}dx\tag{11}$$ Notice that $(10)$ and $(11)$ are a bit different from the expressions I wrote in the comments on Singh's answer. Numerical evidence suggests that: $$\int_0^{\infty} \frac{\theta_2^4\theta_3^4\theta_4^4}{1+x^2}dx=\frac23\int_0^{\infty} \theta_2^4\theta_3^4\theta_4^4 dx\tag{12}$$ $$\int_0^{\infty} \frac{\theta_2^4\theta_3^4\theta_4^2}{1+x^2}dx=\frac83\int_0^{\infty} \theta_2^2\theta_3^4\theta_4^4 dx\tag{13}$$ $$\int_0^{\infty} \theta_2^4\theta_4^2dx=\int_0^{\infty} \theta_2^2\theta_4^2dx=1\tag{14}$$ $$\int_0^{\infty} \frac{\theta_2^4\theta_4^4}{1+x^2}dx=\frac{8\zeta(3)}{\pi^3}\tag{15}$$ $$\int_0^{\infty} \frac{\theta_2^4\theta_4^2}{1+x^2}dx=\frac23\tag{16}$$ All of these just make me more confident that an ultimate colsed form for both $(10)$ ,$(11)$ and higher power combinations of theta functions exist. As a side note, these integrals identities can be translated into identites about lattice sums, by expanding the theta functions into their series definition.