Closed form for $\sum\limits_{n=1}^{\infty}\frac{n^{4k-1}}{e^{n\pi}-1}-16^k\sum\limits_{n=1}^{\infty}\frac{n^{4k-1}}{e^{4n\pi}-1}$
Suppose we seek to evaluate
$$F(p) = \sum_{n\ge 1} \frac{n^{4p-1}}{e^{\pi n}-1} - 16^p \sum_{n\ge 1} \frac{n^{4p-1}}{e^{4\pi n}-1}.$$
These sums may be evaluated using harmonic summation techniques.
Introduce the sum $$S(x; p) = \sum_{n\ge 1} \frac{n^{4p-1}}{e^{nx}-1}$$ with $p$ a positive integer and $x\ge 0.$
The sum term is harmonic and may be evaluated by inverting its Mellin transform.
Recall the harmonic sum identity $$\mathfrak{M}\left(\sum_{k\ge 1} \lambda_k g(\mu_k x);s\right) = \left(\sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} \right) g^*(s)$$ where $g^*(s)$ is the Mellin transform of $g(x).$
In the present case we have $$\lambda_k = k^{4p-1}, \quad \mu_k = k \quad \text{and} \quad g(x) = \frac{1}{e^x-1}.$$
We need the Mellin transform $g^*(s)$ of $g(x)$ which is $$\int_0^\infty \frac{1}{e^{x}-1} x^{s-1} dx = \int_0^\infty \frac{e^{-x}}{1-e^{-x}} x^{s-1} dx \\ = \int_0^\infty \sum_{q\ge 1} e^{-q x} x^{s-1} dx = \sum_{q\ge 1} \int_0^\infty e^{-q x} x^{s-1} dx \\= \Gamma(s) \sum_{q\ge 1} \frac{1}{q^s} = \Gamma(s) \zeta(s).$$
It follows that the Mellin transform $Q(s)$ of the harmonic sum $S(x,p)$ is given by
$$Q(s) = \Gamma(s) \zeta(s) \zeta(s-(4p-1)) \\ \text{because}\quad \sum_{k\ge 1} \frac{\lambda_k}{\mu_k^s} = \sum_{k\ge 1} k^{4p-1} \frac{1}{k^s} = \zeta(s-(4p-1))$$ for $\Re(s) > 4p.$
The Mellin inversion integral here is $$\frac{1}{2\pi i} \int_{4p+1/2-i\infty}^{4p+1/2+i\infty} Q(s)/x^s ds$$ which we evaluate by shifting it to the left for an expansion about zero.
The two zeta function terms cancel the poles of the gamma function term and we are left with just
$$\begin{align} \mathrm{Res}(Q(s)/x^s; s=4p) & = \Gamma(4p) \zeta(4p) \quad\text{and}\\ \mathrm{Res}(Q(s)/x^s; s=0) & = \zeta(0) \zeta(-(4p-1)). \end{align}$$
Computing these residues we get
$$- (4p-1)! \frac{B_{4p} (2\pi)^{4p}}{2(4p)!} = - \frac{B_{4p} (2\pi)^{4p}}{2\times 4p} \quad\text{and}\quad - \frac{1}{2} \times -\frac{B_{4p}}{4p}.$$
This shows that $$S(x;p) = - \frac{B_{4p} (2\pi)^{4p}}{8p\times x^{4p}} + \frac{B_{4p}}{8p} + \frac{1}{2\pi i} \int_{-1/2-i\infty}^{-1/2+i\infty} Q(s)/x^s ds.$$
To treat the integral recall the duplication formula of the gamma function: $$\Gamma(s) = \frac{1}{\sqrt\pi} 2^{s-1} \Gamma\left(\frac{s}{2}\right) \Gamma\left(\frac{s+1}{2}\right).$$
which yields for $Q(s)$
$$\frac{1}{\sqrt\pi} 2^{s-1} \Gamma\left(\frac{s}{2}\right) \Gamma\left(\frac{s+1}{2}\right) \zeta(s) \zeta(s-(4p-1))$$
Furthermore observe the following variant of the functional equation of the Riemann zeta function: $$\Gamma\left(\frac{s}{2}\right)\zeta(s) = \pi^{s-1/2} \Gamma\left(\frac{1-s}{2}\right) \zeta(1-s)$$
which gives for $Q(s)$ $$\frac{1}{\sqrt\pi} 2^{s-1} \pi^{s-1/2} \Gamma\left(\frac{s+1}{2}\right) \Gamma\left(\frac{1-s}{2}\right) \zeta(1-s)\zeta(s-(4p-1)) \\ = \frac{1}{\sqrt\pi} 2^{s-1} \pi^{s-1/2} \frac{\pi}{\sin(\pi(s+1)/2)} \zeta(1-s)\zeta(s-(4p-1)) \\ = 2^{s-1} \frac{\pi^s}{\sin(\pi(s+1)/2)} \zeta(1-s)\zeta(s-(4p-1)).$$
Now put $s=4p-u$ in the remainder integral to get
$$\frac{1}{x^{4p}} \frac{1}{2\pi i} \int_{4p+1/2-i\infty}^{4p+1/2+i\infty} 2^{4p-1-u} \\ \times \frac{\pi^{4p-u}}{\sin(\pi(4p+1-u)/2)} \zeta(u-(4p-1))\zeta(1-u) x^u du \\ = \frac{2^{4p} \pi^{4p}}{x^{4p}} \frac{1}{2\pi i} \int_{4p+1/2-i\infty}^{4p+1/2+i\infty} 2^{u-1} \\ \times \frac{\pi^{u}}{\sin(\pi(4p+1-u)/2)} \zeta(u-(4p-1))\zeta(1-u) (x/\pi^2/2^2)^u du.$$
Now $$\sin(\pi(4p+1-u)/2) = \sin(\pi(1-u)/2+2\pi p) \\ = \sin(\pi(1-u)/2) = - \sin(\pi(-1-u)/2) = \sin(\pi(u+1)/2).$$
We have shown that $$\bbox[5px,border:2px solid #00A000] {S(x;p) = - \frac{B_{4p} (2\pi)^{4p}}{8p\times x^{4p}} + \frac{B_{4p}}{8p} + \frac{2^{4p} \pi^{4p}}{x^{4p}} S(4\pi^2/x;p)}.$$
In particular we get $$S(\pi;p) = - \frac{B_{4p} 2^{4p}}{8p} + \frac{B_{4p}}{8p} + 2^{4p} S(4\pi;p)$$
and
$$S(4\pi;p) = - \frac{B_{4p} 2^{-4p}}{8p} + \frac{B_{4p}}{8p} + 2^{-4p} S(\pi;p).$$
Therefore
$$S(\pi;p)-2^{4p}S(4\pi;p) \\ = - \frac{B_{4p} (2^{4p}-1)}{8p} + \frac{B_{4p}}{8p} (1-2^{4p}) + (2^{4p} S(4\pi;p)-S(\pi;p)).$$
We finally conclude that
$$F(p) = \frac{B_{4p}}{4p} (1-2^{4p}) - F(p)$$
or
$${\large\color{green}{F(p) = \frac{B_{4p}}{8p} (1-2^{4p})}}.$$
The first few values are
$$\frac{1}{16},{\frac {17}{32}},{\frac {691}{16}}, {\frac {929569}{64}},{\frac {221930581}{16}},\ldots$$
We can use the functional equation to extract additional sum formulae. For example when we choose the pair
$$\sqrt{2}\pi \quad\text{and}\quad 2\sqrt{2}\pi$$
the scalar is $2^{2p}.$ The calculation starts from
$$S(\sqrt{2}\pi;p) = - \frac{B_{4p} 2^{2p}}{8p} + \frac{B_{4p}}{8p} + 2^{2p} S(2\sqrt{2}\pi;p)$$
and
$$S(2\sqrt{2}\pi;p) = - \frac{B_{4p} 2^{-2p}}{8p} + \frac{B_{4p}}{8p} + 2^{-2p} S(\sqrt{2}\pi;p).$$
Therefore
$$S(\sqrt{2}\pi;p)-2^{2p}S(2\sqrt{2}\pi;p) \\ = - \frac{B_{4p} (2^{2p}-1)}{8p} + \frac{B_{4p}}{8p} (1-2^{2p}) + (2^{2p} S(2\sqrt{2}\pi;p)-S(\sqrt{2}\pi;p)).$$
We obtain the formula
$$\sum_{n\ge 1} \frac{n^{4p-1}}{e^{\sqrt{2}\pi n}-1} - 4^p \sum_{n\ge 1} \frac{n^{4p-1}}{e^{2\sqrt{2}\pi n}-1} = \frac{B_{4p}}{8p} (1-2^{2p}).$$
We get for the initial values (factor is $1+2^{2p}$)
$${\frac {1}{80}},\frac{1}{32},{\frac {691}{1040}}, {\frac {3617}{64}},{\frac {5412941}{400}},\ldots$$
Addendum. To illustrate the creation of identities from the functional equation we present a third example, which is
$$\sqrt{3}\pi\quad\text{and}\quad 4\pi/\sqrt{3}.$$
In this example the scalar is $2^{4p} 3^{-2p}.$ The calculation starts from
$$S(\sqrt{3}\pi;p) = - \frac{B_{4p} 2^{4p}}{8p\times 3^{2p}} + \frac{B_{4p}}{8p} + \frac{2^{4p}}{3^{2p}} S(4\pi/\sqrt{3};p)$$
and
$$S(4\pi/\sqrt{3};p) = - \frac{B_{4p} 3^{2p}}{8p\times 2^{4p}} + \frac{B_{4p}}{8p} + \frac{3^{2p}}{2^{4p}} S(\sqrt{3}\pi;p).$$
Therefore
$$S(\sqrt{3}\pi;p)-2^{4p} 3^{-2p} S(4\pi/\sqrt{3};p) \\ = - \frac{B_{4p} (2^{4p} 3^{-2p}-1)}{8p} + \frac{B_{4p}}{8p} (1-2^{4p} 3^{-2p}) + (2^{4p} 3^{-2p} S(\sqrt{3}\pi;p)-S(4\pi/\sqrt{3};p)).$$
We obtain the formula
$$\sum_{n\ge 1} \frac{n^{4p-1}}{e^{n \sqrt{3}\pi}-1} - 2^{4p} 3^{-2p} \sum_{n\ge 1} \frac{n^{4p-1}}{e^{n 4\pi/\sqrt{3}}-1} = \frac{B_{4p}}{8p} (1-2^{4p} 3^{-2p}).$$
The reader is cordially invited to prove this last result by a different method.
Remark. The general pattern for
$$\beta\quad\text{and}\quad 4\pi^2/\beta$$
is
$$\bbox[5px,border:2px solid #00A000] {\sum_{n\ge 1} \frac{n^{4p-1}}{e^{n \beta}-1} - \frac{2^{4p}\pi^{4p}}{\beta^{4p}} \sum_{n\ge 1} \frac{n^{4p-1}}{e^{n 4\pi^2/\beta}-1} = \frac{B_{4p}}{8p} \left(1-\frac{2^{4p}\pi^{4p}}{\beta^{4p}}\right).}$$
Just for fun here is a way of obtaining a closed form expression for the sum (which turns out to be the correct one) by manipulating formulas with a complete disregard for divergent sums (the more divergent the better).
We start with the generating function for the Bernoulli numbers
$$\frac{x^{4k-1}}{e^{x}-1} = \sum_{j\geq 0} \frac{B_j}{j!}x^{j+4k-2}$$
Taking $x=n$ and $x=4n$ respectively we obtain
$$\frac{n^{4k-1}}{e^{n\pi}-1}-16^{k}\frac{n^{4k-1}}{e^{4n\pi}-1} = \frac{1}{\pi}\sum_{j\geq 0} \frac{B_j}{j!}\left[1-4^{j+2k-1}\right]\pi^j \cdot \frac{1}{n^{2-j-4k}}$$
Now by summing this equation over $n$, using the definition of $\zeta$-function, we get $$\frac{1}{\pi}\sum_{j\geq 0} \frac{B_j}{j!}\left[1-4^{j+2k-1}\right]\pi^j\zeta(2-j-4k)$$
Note that the $\zeta$-argument is negative so the sum is as divergent as it's possible to get (each term is a divergent sum), but by interpretting it using the analytical continuation of $\zeta(s)$ we can try to extract some meaning from it. The analytical continuation satisfy $\zeta(2-j-4k) = 0$ when $j$ is even and since $B_{j} = 0$ when $j$ is odd (apart from the term $j=1$ having $B_1 = - \frac{1}{2}$) the formula above reduces to
$$-\frac{1}{2}\left[1-4^{2k}\right]\zeta(1-4k)$$
Finally, by using the $\zeta$ functional equation and the relation for $\zeta(2n)$ when $n$ is an integer this simplifies to
$$\frac{B_{4k}}{8k}\left(1-4^{2k}\right)$$
Remarkably this calculation gives us the correct result (as it often does).