I am going to add the restrictions $b>a>0$ and $x$ instead of $x^4$
$$\int^b_a \frac{x e^{x}}{e^{x}-1}\,dx = \int^b_a x \,dx+\int^b_a \frac{x}{e^{x}-1}\,dx = \frac{b^2-a^{2}}{2}+ I(a,b)$$
Now for that integral we can rewrite it as
$$I(a,b)=\int^b_0 \frac{x}{e^x-1}\,dx-\int^a_0\frac{x}{e^x-1}\,dx$$
Now let us look at the following
$$I(a) = \int^a_0 \frac{x}{e^x-1}\,dx= \int^a_0xe^{-x} \sum_{n\geq 0}e^{-nx} = \sum_{n\geq 0}\int^a_0 xe^{-x(n+1)}\,dx$$
$$=\sum_{n\geq 0}\frac{1}{(n+1)^2}-\frac{e^{-(n+1) a}}{(n+1)^2} -\frac{ae^{-(n+1) a}}{(n+1)} = \zeta(2)+a\log(1-e^{-a})-\mathrm{Li}_2(e^{-a})$$
Let us take $a \to \infty$
$$\int^\infty_0 \frac{x}{e^{x}-1}= \zeta(2)+\lim_{a \to \infty}a\log(1-e^{-a})-\mathrm{Li}_2(0)=\frac{\pi^2}{6}$$
Hence
$$I(a,b)= b\log(1-e^{-b})-a\log(1-e^{-a})+\mathrm{Li}_2(e^{-a}) -\mathrm{Li}_2(e^{-b})$$
The method could be generalized for
$$\int^b_a \frac{x^{n}e^{x}}{e^{x}-1}\,dx$$