Closed form of $\int_0^\infty \left(\frac{\arctan x}{x}\right)^ndx$
Used formulas and definitions:
$\displaystyle\Re(i^n)=\cos\left(\frac{\pi n}{2}\right)\enspace$ for integer $~n\enspace$ ; $\enspace\displaystyle\tan x = -i\frac{1-e^{-i2x}}{1+e^{-i2x}}$
$\displaystyle \sum\limits_{j=0}^{-1} a_j := \sum\limits_{j=0}^0 a_j - a_0 = 0\enspace$ ; $\enspace\displaystyle \left(\sum\limits_{k=0}^{\infty} a_k\right) \left(\sum\limits_{k=0}^{\infty} b_k\right) = \left(\sum\limits_{k=0}^{\infty} \sum\limits_{v=0}^k a_v b_{k-v}\right) $
$\displaystyle \int\limits_0^a x^m e^{zx} dx = \frac{m!}{(-z)^{m+1}} - e^{az}\sum\limits_{v=0}^m\frac{m!a^{m-v}}{(m-v)!(-z)^{v+1}} \enspace$ for $\enspace m\in\mathbb{N}_0 , \, a\in\mathbb{R} , \, z\in\mathbb{C}\setminus\{0\} $
Stirling numbers of the first kind $\displaystyle \begin{bmatrix}n\\k\end{bmatrix}$ defined by
$\hspace{4cm}\displaystyle\sum\limits_{j=0}^n\begin{bmatrix}n\\j\end{bmatrix}x^j:=\prod\limits_{k=0}^{n-1}(x+k) \enspace$ for $\enspace n\in\mathbb{N}_0 , \, x\in\mathbb{C} $
$\displaystyle \sum\limits_{j=0}^n f(j) \sum\limits_{l=0}^j g(j,l) k^l = \sum\limits_{j=0}^n k^j \sum\limits_{l=j}^n f(l)g(l,j) \enspace\enspace$ formal summation permutation
$\displaystyle \sum\limits_{v=-1}^k {\binom {n+1} {k-v}}{\binom {n+v} v} = \displaystyle \sum\limits_{v=0}^{n+1} {\binom {n+1} v}{\binom {n+k-v} n}\enspace$ for $\enspace k,n\geq 0$
$\displaystyle c_{n,j} := \frac{1}{n!} \sum\limits_{v=0}^{n+1} {\binom {n+1} {v}} \sum\limits_{l=j}^m \begin{bmatrix}n+1\\l+1\end{bmatrix} {\binom l j} (-v)^{l-j} \enspace$ for $\enspace 0\leq j\leq n\enspace$ with $\enspace 0^0:=1$
$\displaystyle b_{n,k} := \sum\limits_{v=0}^k {\binom {n+1} {k-v}} {\binom {n+v} v} = -0^{n+k} +\sum\limits_{j=0}^n k^j c_{n,j} \enspace$ for $\enspace k,n \geq 0 \,~~ ; \enspace b_{n,0 }=1$
It follows:
$\displaystyle\sum\limits_{k=1}^\infty \frac{b_{n,k}~x^k}{k^s} = \sum\limits_{j=0}^n c_{n,j}~\text{Li}_{s-j}(x) \enspace$ for $\enspace s\in\mathbb{C}\enspace$ where $\enspace \text{Li}_s(x)\enspace$ is the Polylogarithm
with the special cases $\enspace \text{Li}_s(1)\equiv \zeta(s)\,$ and
$\hspace{4.1cm}\text{Li}_s(-1)\equiv -\eta(s) = \left(2^{1-s}-1\right)\zeta(s) \enspace$ with $\enspace\eta(1)=\ln 2$
$\underline{\text{Solution:}}$
Let $\enspace\displaystyle |a| \leq \frac{\pi}{2} \, , \enspace 0\leq n\leq m \,$ .
$\displaystyle \int\limits_0^{\tan a} \frac{\arctan^m x}{x^n} dx = \int\limits_0^a \frac{x^m}{\tan^n x} dx + \int\limits_0^a \frac{x^m}{\tan^{n-2} x} dx \enspace$ for $\enspace n\geq 2$
$$\int\limits_0^a \frac{x^m}{\tan^n x}dx = i^n\frac{a^{m+1}}{m+1} + \frac{i^{n-m-1}m!}{2^{m+1}}\sum\limits_{j=0}^{n-1}~c_{n-1,j}~\text{Li}_{m+1-j}~(1) $$
$$\hspace{2.5cm} -\sum\limits_{v=0}^m \frac{i^{n-v-1}m!a^{m-v}}{(m-v)!2^{v+1}} \sum\limits_{j=0}^{n-1}~c_{n-1,j}~\text{Li}_{v+1-j}~(e^{-i2a})$$
$\displaystyle a:=\frac{\pi}{2}\,$ :
$$\int\limits_0^{\pi/2} \frac{x^m}{\tan^n x}dx = $$
$$ \cos\left(\frac{\pi n}{2}\right)\frac{(\pi/2)^{m+1}}{m+1} +\cos\left(\frac{\pi (n-m-1)}{2}\right) \frac{m!}{2^{m+1}}\sum\limits_{j=0}^{n-1}~c_{n-1,j}~\zeta(m+1-j) $$
$$\hspace{2.5cm} +\sum\limits_{v=0}^m \cos\left(\frac{\pi (n-v-1)}{2}\right)\frac{m!(\pi/2)^{m-v}}{(m-v)!2^{v+1}} \sum\limits_{j=0}^{n-1}~c_{n-1,j}~\eta(v+1-j)\hspace{1cm}$$
Finally we get:
$$\int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^n dx = $$
$$\frac{1}{2^{n+1}} \sum\limits_{v=0}^n \left(\cos\frac{\pi(n-v-1)}{2}\right) \frac{n!\pi^{n-v}}{(n-v)!}~\cdot $$
$$\cdot~\left( \sum\limits_{j=0}^{n-1} c_{n-1,j}~\eta(v-j+1) - \sum\limits_{j=0}^{n-3} c_{n-3,j}~\eta(v-j+1) \right)$$
Analytical continuations $(s\in\mathbb{C})$ :
$$\begin{align} \zeta(1-s)&=\dfrac{2}{(2\pi)^s}\cos\bigg(\dfrac{\pi s}{2}\bigg)~\Gamma(s)~\zeta(s) \\\\ \eta(1-s)&=\dfrac{2^s-1}{1-2^{s-1}}~\pi^{-s}\cos\bigg(\dfrac{\pi s}{2}\bigg)~\Gamma(s)~\eta(s) \end{align} $$
Simplifications $(k\in\mathbb{N}_0)$ :
$$\begin{align} \eta(1)=\ln 2 \enspace ; \enspace \eta(1-k)~&=~\frac{2^k-1}{k}~B_k \end{align}$$
$$\begin{align} \eta(2k)~&=~(-1)^{k-1}~\frac{2^{2k-1}-1}{(2k)!}~B_{2k}~\pi^{2k} \end{align}$$
$$\begin{align} \eta(2k+1)~&=~\bigg(1-\frac{1}{2^{2k}}\bigg)\zeta(2k+1) \end{align}$$
Examples:
$\displaystyle c_{0,0} = \frac{2}{0!}(1)$
$\displaystyle (c_{1,0}~;~ c_{1,1}) = \frac{2}{1!} \left(0~;~2\right)$
$\displaystyle (c_{2,0}~;~c_{2,1}~;~c_{2,2}) = \frac{2}{2!} \left(2~;~0~;~4\right)$
$\displaystyle (c_{3,0}~;~c_{3,1}~;~c_{3,2}~;~c_{3,3}) = \frac{2}{3!} \left(0~;~16~;~0~;~8\right)$
$\displaystyle (c_{4,0}~;~c_{4,1}~;~c_{4,2}~;~c_{4,3}~;~c_{4,4}) = \frac{2}{4!} \left(24~;~0~;~80~;~0~;~16\right)$
$\displaystyle (c_{5,0}~;~c_{5,1}~;~c_{5,2}~;~c_{5,3}~;~c_{5,4}~;~c_{5,5}) = \frac{2}{5!} \left(0~;~368~;~0~;~320~;~0~;~32\right)$
$\displaystyle (c_{6,0}~;~c_{6,1}~;~c_{6,2}~;~c_{6,3}~;~c_{6,4}~;~c_{6,5}~;~c_{6,6}) = \frac{2}{6!} \left(720~;~0~;~3136~;~0~;~1120~;~0~;~64\right)$
$\displaystyle (c_{7,0}~;~c_{7,1}~;~c_{7,2}~;~c_{7,3}~;~c_{7,4}~;~c_{7,5}~;~ c_{7,6}~;~ c_{7,7}) $
$\hspace{7cm}\displaystyle =\frac{2}{7!} \left(0~;~16896~;~0~;~19712~;~0~;~3584~;~0~;~128\right)$
...
$\displaystyle \int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^2 dx = \pi\ln 2$
$\displaystyle \int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^3 dx = -\frac{\pi^3}{16} + \frac{3\pi}{2}\ln 2$
$\displaystyle \int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^4 dx = -\frac{\pi^3}{12}(2\ln 2 + 1) + \frac{\pi}{4}(3\zeta(3) + 8\ln 2)$
$\displaystyle \int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^5 dx = \frac{\pi^5}{128} - \frac{5\pi^3}{48}(8\ln 2 + 1) + \frac{5\pi}{4}(3\zeta(3) + 2\ln 2)$
$\displaystyle \int\limits_0^\infty \left(\frac{\arctan x}{x}\right)^6 dx = $
$\displaystyle\hspace{8mm} =\frac{3\pi^5}{320}(4\ln 2 + 3) - \frac{\pi^3}{16}(9\zeta(3) + 40\ln 2 + 2) + \frac{3\pi}{32}(45\zeta(5) + 120\zeta(3) + 32\ln 2)$
...
$$\int_{0}^{+\infty}\left(\frac{\arctan x}{x}\right)^3\,dx \stackrel{x\mapsto\tan\theta}{=}\int_{0}^{\pi/2}\frac{\theta^3\,d\theta}{\tan^3\theta\cos^2\theta}\stackrel{\text{IBP}}{=}-\frac{\pi^3}{16}+\frac{3}{2}\int_{0}^{\pi/2}\frac{\theta^2\,d\theta}{\sin^2\theta}$$ and in general the computation of the wanted integrals boils down to the computation of $\int_{0}^{\pi/2}\left(\frac{\theta}{\sin\theta}\right)^m\,d\theta$ or the computation of $\oint_\gamma \frac{\log^m z}{z\left(z-\frac{1}{z}\right)^m}\,dz$ where $\gamma$ is the quarter circle joining $1$ and $i$. The integration of $\frac{\log^m z}{z\left(z-\frac{1}{z}\right)^m}$ along the segments joining $0$ and $1$ or $0$ and $i$ can be simply performed through Maclaurin series; in particular the wanted integrals can be always expressed in terms of standard or alternating Euler sums.
It wasn't requested, but instead of exact representations the OP might want an asymptotic expression for $n \to \infty.$ This one works well with the technique of Depoissonization. Make an exponential power series and analyze it asymptotically: $$ \sum_{n=0}^\infty \frac{y^n}{n!} C_n = \int_0^\infty \exp{\big(\frac{y}{x} \text{arctan}(x) \big)} dx , \quad C_n=\int_0^\infty \Big(\frac{\text{arctan}(x)}{x}\Big)^n dx$$ Now $\text{arctan}(x)/x = 1-x^2/3+x^4/5+...$ and with $y$ large and keeping on the first term in the asymptotic expansion $$ e^{-y} \sum_{n=0}^\infty \frac{y^n}{n!} C_n \sim \int_0^\infty \exp{\big(-y\,\frac{x^2}{3}\big)} dx = \frac{1}{2} \sqrt{\frac{3\pi}{y}}.$$ By Depoissonization we can conclude that $$ C_n \sim \frac{1}{2} \sqrt{\frac{3\pi}{n}} .$$ For $n=50$ the asymptotic expression is within 2% of the value from a numerical integration.