$\int_{0}^{\infty} \frac{1}{1 + x^r}\:dx = \frac{1}{r}\Gamma\left( \frac{r - 1}{r}\right)\Gamma\left( \frac{1}{r}\right)$ As part of a recent question I posted, I decided to try and generalise for a power of $2$ to any $r \in \mathbb{R}$. As part of the method I took, I had to solve the following integral:
\begin{equation}
 I = \int_{0}^{\infty} \frac{1}{1 + x^r}\:dx
\end{equation}
I believe what I've done is correct, but I'm concerned that I may missed something (in particular whether it holds for all $r \neq 0$). So, here I have two questions (1) Is what I've done correct? and (2) What other methods can be employed that doesn't rely on complex analysis?
Here is the method I took:
First make the substitution $u = x^{\frac{1}{r}}$ to arrive at
\begin{equation}
 I = \frac{1}{n} \int_{0}^{\infty} \frac{1}{1 + u} \cdot u^{1 -\frac{1}{r}}\:du
\end{equation}
We now substitute $t = \frac{1}{1 + u}$ to arrive at:
\begin{align}
 I &= \frac{1}{r} \int_{1}^{0} t \cdot \left(\frac{1 - t}{t}\right)^{\frac{1}{r} -1}\frac{1}{t^2}\:dt = \frac{1}{r}\int_{0}^{1}t^{-\frac{1}{r}}\left(1 - t\right)^{ \frac{1}{r} - 1}\:dt \\
&= \frac{1}{r}B\left(1 - \frac{1}{n}, 1 + \frac{1}{r} - 1\right) = \frac{1}{r} B\left(\frac{r - 1}{r},\frac{1}{r}\right) \\
&= \frac{1}{r} B\left(\frac{r - 1}{r},\frac{1}{r}\right)
\end{align}
Wheer $B(a,b)$ is the Beta function. Using the relationship between the Beta and Gamma function we arrive at:
\begin{equation}
I = \frac{1}{r} \frac{\Gamma\left( \frac{r - 1}{r}\right)\Gamma\left( \frac{1}{r}\right)}{\Gamma\left(\frac{r - 1}{r} + \frac{1}{r}\right)} = \frac{1}{r}\Gamma\left( \frac{r - 1}{r}\right)\Gamma\left( \frac{1}{r}\right)
\end{equation}
And so, we arrive at:
\begin{equation}
 I = \int_{0}^{\infty} \frac{1}{1 + x^r}\:dx = \frac{1}{r}\Gamma\left( \frac{r - 1}{r}\right)\Gamma\left( \frac{1}{r}\right)
\end{equation}
for $r > 1$
As per KemonoChen's comment and others, we can employ Euler's Reflection Formula to position this result for $\frac{1}{r} \not \in \mathbb{Z}$ Here, as $r \in \mathbb{R}, r > 1 \rightarrow \frac{1}{r} \not \in \mathbb{Z}$ and so our formula holds.
\begin{equation}
 I = \int_{0}^{\infty} \frac{1}{1 + x^r}\:dx = \frac{1}{r}\Gamma\left( \frac{r - 1}{r}\right)\Gamma\left( \frac{1}{r}\right) = \frac{\pi}{r\sin\left(\frac{\pi}{r} \right)}
\end{equation}
Thank you also to Winther, Jjagmath, and MrTaurho's for their comments and corrections/clarifications.

NOT A FULL SOLUTION:
I've been working with special cases of the integral. 
Here we will consider $r = 2m$ where $m \in \mathbb{N}$. In doing so, we observe that the roots of the denominator are $m$ pairs of complex roots $(z, c(z))$ where $c(z)$ is the conjugate of $z$. To verify this:
\begin{align}
 x^{2m} + 1 = 0 \rightarrow x^{2m} = e^{\pi i}
\end{align}
By De Moivre's formula, we observe that:
\begin{align}
 x = \exp\left({\frac{\pi + 2\pi j}{2m} i} \right) \mbox{ for } j = 0\dots 2m - 1
\end{align}
Which we can express as the set
\begin{align}
 S &= \Bigg\{  \exp\left({\frac{\pi + 2\pi \cdot 0}{2m} i} \right) , \:\exp\left({\frac{\pi + 2\pi \cdot 1}{2m} i} \right),\dots,\:\exp\left({\frac{\pi + 2\pi \cdot (2m - 2)}{2m} i} \right)\\
&\qquad\:\exp\left({\frac{\pi + 2\pi \cdot (2m - 1)}{2m} i} \right)\Bigg\} \\
\end{align}
Which can be expressed as the set of 2-tuples
\begin{align}
 S &= \left\{ \left( \exp\left({\frac{\pi + 2\pi j}{2m} i} \right) , \:\exp\left({\frac{\pi + 2\pi(2m - 1 - j )}{2m} i} \right)\right)\: \bigg|\: j = 0 \dots m - 1\right\}\\
& = \left\{ (z_j, c\left(z_j\right)\:|\: j = 0 \dots m - 1 \right\}
\end{align}
From here, we can factor $x^{2m} + 1$ into the form 
\begin{align}
x^{2m} + 1 &= \prod_{r \in S} \left(x + r_j\right)\left(x + c(r_j)\right) \\
 &= \prod_{r \in S} \left(x^2 + \left(r_j + c(r_j)\right)x + r_j c(r_j)\right) \\
 &= \prod_{r \in S} \left(x^2 + 2\Re\left(r_j\right)x + \left|r_j \right|^2\right)
\end{align}
For our case here $\left|r_j \right|^2 = 1$ and $\Re\left(r_j\right) = \cos\left({\frac{\pi + 2\pi j}{2m} } \right)$. Hence, 
\begin{align}
\frac{1}{x^{2m} + 1} =  \prod_{j = 0}^{m - 1}\frac{1}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}
\end{align}
From here, to evaluate the integral we must employ Partial Fraction Decomposition:
\begin{align}
\frac{1}{x^{2m} + 1} =  \prod_{j = 0}^{m - 1}\frac{1}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1} = \sum_{j = 0}^{m - 1}\frac{\alpha_j + \beta_jx}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}
\end{align}
And solve for $\alpha_j$ and $\beta_j$. Putting the coefficents to the side we can find general expressions for the integral:
\begin{align}
\frac{1}{x^{2m} + 1} =  \prod_{j = 0}^{m - 1}\frac{1}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}
\end{align}
From here, to evaluate the integral we must employ Partial Fraction Decomposition:
\begin{align}
\int_{0}^{\infty}\frac{1}{x^{2m} + 1}\:dx &= \int_{0}^{\infty}\sum_{j = 0}^{m - 1}\frac{\alpha_j + \beta_jx}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}\:dx \\
&= \sum_{j = 0}^{m - 1}\left[ \int_{0}^{\infty}\frac{\alpha_j}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}\:dx + \int_{0}^{\infty}\frac{\beta_jx}{ x^2 + 2\cos\left({\frac{\pi + 2\pi j}{2m} } \right)x + 1}\:dx\right] \\
&=\sum_{j = 0}^{m - 1}\left[ I_1 + I_2\right]
\end{align}
Evaluating each individually:
\begin{align}
\int_{0}^{\infty} \frac{\alpha_j}{ x^2 + 2\cos\left(\frac{\pi + 2\pi j}{2m}  \right)x + 1}\:dx &= \left[ \csc\left(\frac{\pi + 2\pi j}{2m} \right)\arctan\left(\frac{(x - 1)\tan\left(\frac{\pi + 2\pi j}{4m}  \right)}{x + 1} \right)\right]_{0}^{\infty} \\
&= \csc\left(\frac{\pi + 2\pi j}{2m} \right)\left( \frac{\pi + 2\pi j}{2m}\right)
\end{align}
Now if we consider the second integral on it's own we find that as a general expression the integral doesn't converge. This doesn't take away from this method, but to proceed we must solve for the unknown coefficients to proceed.
