Showing that $2 \Gamma(a) \zeta(a) \left(1-\frac{1}{2^{a}} \right) = \int_{0}^{\infty}\left( \frac{x^{a-1}}{\sinh x} - x^{a-2} \right) \mathrm dx$

First step, write the integral as $$ I = \int_{0}^{\infty}\left(-{\frac {{x}^{a-2}{{\rm e}^{2\,x}}}{{{\rm e}^{2\,x}}-1}}+2\,{\frac {{x }^{a-1}{{\rm e}^{x}}}{{{\rm e}^{2\,x}}-1}}+{\frac {{x}^{a-2}}{{{\rm e} ^{2\,x}}-1}}\right) dx. $$

Second step, use the change of variables $ 2x=u $. Third step, use the identity (the Hurwitz zeta function)

$$ \zeta(s,q)=\frac{1}{\Gamma(s)} \int_0^\infty \frac{t^{s-1}e^{-qt}}{1-e^{-t}}dt.$$

Here are related problems: I, II. By the way, I tried Maple and WolframAlpha, but no answer.

$$ \color{blue}{\Gamma(s-N)\zeta(s-N)=\int_{0}^{\infty}x^{s-N-2}\left[\frac{x}{e^x-1}-\left(\sum_{n=0}^{N}B_{n}\frac{x^n}{n!}\right)\right]\,dx} $$ $$ {\small \,0\lt\,Re\{s\}\,\lt1 ,\quad N\in\{\,0,\,1,\,2,\,\cdots\,\} ,\quad B_{n}\,\,{Bernoulli\,Number} ,\quad B_{1}=-1/2} $$

$$ \begin{align} I\,& =\,\int_{0}^{\infty}\left(\frac{x^{s-1}}{\sinh{x}}-x^{s-2}\right)\,dx =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{e^x}{e^{2x}-1}-\frac{1}{2x}\right)\,dx \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{e^x\color{red}{+1-1}}{e^{2x}-1}-\frac{1}{2x}\color{red}{-\frac{1}{2x}+\frac{1}{2x}}\right)\,dx \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}-\frac{1}{e^{2x}-1}+\frac{1}{2x}\right)\,dx \end{align} $$

For ${\small\color{red}{\,0\le\,Re\{s\}\,\le+1\,}}\colon\quad\Gamma(s)\zeta(s)=\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^x-1}-\frac{1}{x}\right)\,dx\,$; Thus: $$ \begin{align} I\,& =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}-\frac{1}{e^{2x}-1}+\frac{1}{2x}\right)\,dx \qquad\color{blue}{\small\{\text{both converge}\}} \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}\right)\,dx\,-\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{2x}-1}-\frac{1}{2x}\right)\,dx \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}\right)\,dx\,-\,\frac{2}{2^{s}}\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}\right)\,dx \\[2mm] & =\,\color{red}{\left(2-2^{1-s}\right)\,\Gamma(s)\zeta(s)} \end{align} $$

For ${\small\color{red}{\,-1\le\,Re\{s\}\,\le0\,}}\colon\quad\Gamma(s)\zeta(s)=\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^x-1}-\frac{1}{x}+\frac{1}{2}\right)\,dx\,$; Thus: $$ \begin{align} I\,& =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}\color{red}{+\frac{1}{2}}-\frac{1}{e^{2x}-1}+\frac{1}{2x}\color{red}{-\frac{1}{2}}\right)\,dx \qquad\color{blue}{\small\{\text{both converge}\}} \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}+\frac{1}{2}\right)\,dx\,-\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{2x}-1}-\frac{1}{2x}+\frac{1}{2}\right)\,dx \\[2mm] & =\,2\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}+\frac{1}{2}\right)\,dx\,-\,\frac{2}{2^{s}}\int_{0}^{\infty}x^{s-1}\left(\frac{1}{e^{x}-1}-\frac{1}{x}+\frac{1}{2}\right)\,dx \\[2mm] & =\,\color{red}{\left(2-2^{1-s}\right)\,\Gamma(s)\zeta(s)} \end{align} $$