Aquí hay otra solución: Vamos a $(I_n)$ por
$$ I_n = \int_{0}^{\infty} \frac{\arctan x}{x(1+x^2)^n} \, dx = \int_{0}^{\frac{\pi}{2}} \frac{\theta}{\sin\theta} \cos^{2n-1}\theta \, d\theta. $$
A continuación, por un simple cálculo,
$$ I_n - I_{n+1} = \int_{0}^{\frac{\pi}{2}} \theta \sin\theta \cos^{2n-1}\theta \, d\theta = \frac{1}{2n} \int_{0}^{\frac{\pi}{2}} \cos^{2n}\theta \, d\theta. $$
Desde $I_n \to 0$$n \to \infty$, nos encontramos con que
$$ I_n = \sum_{k=n}^{\infty} \frac{1}{2k} \int_{0}^{\frac{\pi}{2}} \cos^{2k}\theta \, d\theta. $$
La división de la sumatoria como $\sum_{k=n}^{\infty} = \sum_{k=1}^{\infty} - \sum_{k=1}^{n-1}$, nos encontramos con que
$$ I_n = \frac{\pi}{2}\left( \log 2 - \sum_{k=1}^{n-1} \frac{1}{2k} \frac{(2k-1)!!}{(2k)!!} \right), $$
donde $n!!$ denota el doble factorial.
Edición 1. En general, hemos
$$ \int_{0}^{\infty} \frac{\arctan^s x}{x(1+x^2)^{n+1}} \, dx = \int_{0}^{\frac{\pi}{2}} \theta^s \cot \theta \, d\theta - \sum_{k=1}^{n} \int_{0}^{\frac{\pi}{2}} \theta^s \sin\theta \cos^{2k-1}\theta \, d\theta. \tag{1} $$
Actualmente no tengo idea de cómo obtener una fórmula simple para la siguiente integral
$$ \int_{0}^{\frac{\pi}{2}} \theta^s \sin\theta \cos^{2k-1}\theta \, d\theta, \tag{2} $$
incluso cuando $s = 2$. Por otro lado, para cualquier $s > 0$ $N \geq \lfloor s/2 \rfloor$ hemos
\begin{align*}
\int_{0}^{\frac{\pi}{2}} \theta^s \cot \theta \, d\theta
&= 2^{-s}\cos\left(\frac{\pi s}{2}\right)\Gamma(1+s)\zeta(1+s) \\
&\quad + \left(\frac{\pi}{2}\right)^s \sum_{k=0}^{N} (-1)^k \pi^{-2k} \frac{\Gamma(2k-s)}{\Gamma(-s)} \eta(2k+1) \\
&\quad + \frac{(-1)^{N+1}}{2^s \Gamma(-s)} \int_{0}^{\infty} \frac{t^{2N+1-s}}{1+t^2} \left( \sum_{n=1}^{\infty} \frac{(-1)^{n-1}}{n^{1+s}} e^{-\pi n t} \right) \, dt,
\end{align*}
donde $\eta(s)$ denota la Dirichlet eta función. (Mi solución es un poco complicado, por lo que voy a publicar más tarde si me parece útil para nuestro problema.) En particular, cuando se $s$ es un entero positivo, entonces la integral de la parte se desvanece y la fórmula se vuelve mucho más sencillo. Así, la fórmula (*) da una forma cerrada mientras podemos calcular la integral (2).
Ejemplo 1. Por ejemplo, cuando se $s = 2$ entonces podemos usar $N = 1$ y, a continuación,
\begin{align*}
\int_{0}^{\frac{\pi}{2}} \theta^2 \cot \theta \, d\theta
&= -\frac{1}{2}\zeta(3) + \frac{\pi^2}{4} \log 2 - \frac{1}{2}\eta(3) \\
&= \frac{\pi^2}{4}\log 2 - \frac{7}{8}\zeta(3).
\end{align*}
Ya podemos calcular la integral (2) por $s = 2$$k = 1, \cdots, 4$, podemos obtener fácilmente OP última identidad.
Aquí está otro ejemplo:
Ejemplo 2. Utilizando la fórmula con $s = 6$, se puede comprobar que
\begin{align*}
\int_{0}^{\infty} \frac{\arctan^6 x}{x(1+x^2)^3} \, dx
&= \frac{\pi^6}{64} \log 2 -\frac{45 \pi^4}{128} \zeta(3) + \frac{675 \pi^2}{128} \zeta(5) -\frac{5715}{256} \zeta(7) \\
&\quad - \frac{11 \pi^6}{2048} + \frac{705 \pi^4}{4096} - \frac{8595 \pi^2}{4096} + \frac{135}{16} \\
&\approx 0.0349464822054751922142122595622\cdots.
\end{align*}