It is not hard to show, by induction on $m\in\mathbb N$, that there exist a sequence $(B_n)_{n\geq0}$ of rational numbers such that

$$\sum_{k=1}^nk^m=\frac1{m+1}\sum_{k=0}^m\binom{m+1}kB_k\,n^{m+1-k},\ \style{font-family:inherit;}{\text{for all}}\ n\geq1\,.$$

These are the Bernoulli numbers of the second kind, that is, $B_1=+1/2$. A corollary of the proof (by induction) of the fact above is a recurrence formula for such numbers $B_n$, which are known as Bernoulli numbers:

$$\sum_{k=0}^{n-1}\binom nkB_k=0\,.\tag{$\ast$}$$

On the other hand, there are a number of explicit formulae for $B_n$, which are obtained using the equality $\frac x{e^x-1}=\sum_{n=0}^\infty B_n\,\frac{x^n}{n!}$; see here and here. For example:

$$B_n=\sum_{k=0}^n\frac1{k+1}\sum_{r=0}^k(-1)^r\binom krr^n\,.\tag{$\ast\ast$}$$

Is there some clever way to manipulate the recurrence $(\ast)$ in order to obtain $(\ast\ast)$?

  • $\begingroup$ ??? How is (*) true for $n=1$ ? We have then $ \binom10 \cdot B_0 = 0 $ but it must be $ \binom10 \cdot B_0 = 1 $ $\endgroup$ Commented Jul 25, 2013 at 21:23
  • 1
    $\begingroup$ @GottfriedHelms No, $B_0=1$. $\endgroup$
    – Pedro
    Commented Jul 25, 2013 at 22:02
  • $\begingroup$ It seems, I'll never understand, why in many openly available places this formula is not stated as any other mathematical formula with an explicite statement for the range of its validity. $\endgroup$ Commented Jul 25, 2013 at 22:07
  • 2
    $\begingroup$ A simple correction of (*) to make it a fully valid formula is $ \displaystyle \\ \sum_{k=0}^{n-1} \binom nk B_k = \delta_{n,1} \\ $ . I'll try to generate an answer based on this... $\endgroup$ Commented Jul 25, 2013 at 22:24
  • $\begingroup$ @GottfriedHelms Ah, the formula is true for $n>1$; yes. Do you think you can work out something out of what I did? $\endgroup$
    – Pedro
    Commented Jul 25, 2013 at 22:55

2 Answers 2


[20-08] I will add something that sheds light on why the latter works out nicely. In what follows we provide a proof of Faulhaber's Formula using the operation of binomial convolution on sequences of real numbers (in fact, we needn't worry about the whole real numbers, but merely the rational numbers). This operation arises naturally when multiplying (formal) exponential generating series of the form $$\sum_{n\geqslant 0}\gamma_n \frac{z^n}{n!}$$ and collecting powers of $z$. The problem of finding $1^n+2^n+\cdots+m^n$ motivates us to look into each of the sequences $1,j,j^2,j^3,\ldots$ that have the EGF $e^{jz}$. Thus the sum has EGF $$\sum_{j=1}^m e^{jz}=e^z\frac{e^{mz}-1}{e^{z}-1}=\frac{e^{mz}-1}z \frac{(-z)}{e^{-z}-1}$$

It is easy to obtain explicit generating sequences for $\dfrac{e^{mz}-1}z$; indeed $$\frac{e^{mz}-1}z=\sum_{n\geqslant 0}\frac{m^{n+1}}{n+1}\frac{z^n}{n!}$$

In particular when $m=1$ we get the Harmonic numbers generate $$H(z)=\frac{e^{z}-1}z$$

The Bernoulli numbers arise naturally as the inverses of the Harmonic numbers, in the sense that $$B(z)=\frac{z}{e^z-1}=\sum_{n\geqslant 0}B_n\frac{z^n}{n!}$$

The first equation then gives Faulhaber's Formula immediately: we have $$\sum_{j=1}^m e^{jz}=e^z\frac{e^{mz}-1}{e^{z}-1}=B(-z)\frac{e^{mz}-1}z $$

Equating coefficients and using the convolution gives $$\sum_{j=1}^m j^n =\sum_{k=0}^n\binom nk (-1)^k B_k\frac{m^{n-k+1}}{n-k+1}$$

This is the shortest proof I know of. We work things formally, without worrying about convergence. We use the fact that a sequence with non-zero first term $a_0$ admits a unique inverse (as is the case of the Harmonic numbers). This is proven by observing the convolution gives a recurrence for the terms of this inverse.

The above is very similar to what will follow. The first lemma gives $$H(z)e^{jz}=\frac{e^{(j+1)z}-1}z-\frac{e^{jz}-1}z$$ and the second gives $$e^z\frac{e^{mz}-1}z=\frac{e^{(m+1)z}-1}z-\frac{e^{z}-1}z$$ Telescoping gives $$H(z)\sum\limits_{j = 1}^m {{e^{jz}}} = \frac{{{e^{(m + 1)z}} - {e^z}}}{z} = {e^z}\frac{{{e^{mz}} - 1}}{z}$$

Thus $$\sum\limits_{j = 1}^m {{e^{jz}}} = \frac{{{e^{(m + 1)z}} - {e^z}}}{z} = {e^z}B\left( z \right)\frac{{{e^{mz}} - 1}}{z}$$ and using $B(-z)=e^zB(z)$ $-$ this $1\star B=\mu B$ which was quite ugly to prove without EGFs $-$ finishes things to give Faulhaber's Formula $$\sum\limits_{j = 1}^m {{e^{jz}}} = B\left( { - z} \right)\frac{{{e^{mz}} - 1}}{z}$$

DEF Let ${\bf x} =(x_0,x_1,x_2,\ldots,x_n,\ldots)$, ${\bf y}=(y_0,y_1,y_2,\ldots,y_n,\ldots)$ be sequences. We define a new sequence, their binomial convolution, by the formula $$\left({\bf x}\star {\bf y}\right)_n:=\sum_{k=0}^n \binom{n}{k}x_ky_{n-k}$$

OBS The binomial convolution is associative and commutative.

DEF We define the special sequences $\mu=(1,-1,1,-1,1,-1,1,\ldots)$, ${\bf 1 }=(1,1,1,1,1,\ldots)$, ${\rm id}=(1,0,0,0,0,0,\ldots)$. Note then that $$\mu\star {\bf 1} =\rm{id}$$ $${\rm id}\star {\bf x}={\bf x}$$

We're ready to prove the

THM (Inversion formula) Let ${\bf x} =(x_0,x_1,x_2,\ldots,x_n,\ldots)$, ${\bf y}=(y_0,y_1,y_2,\ldots,y_n,\ldots)$ be sequences. Then ${\bf x}\star {\bf 1}={\bf y}$ if, and only if ${\bf x}={\bf y}\star \mu$.

P By the above $$\begin{align}{\bf x}\star {\bf 1}&={\bf y}\\({\bf x}\star {\bf 1})\star \mu&={\bf y}\star \mu\\{\bf x}\star ({\bf 1}\star\mu)&={\bf y}\star \mu\\{\bf x}\star {\rm id }&={\bf y}\star \mu\\{\bf x}&={\bf y}\star \mu\end{align}$$

and the other direction is analogous.

[7-28] I have decided to define the Bernoulli numbers as the inverse of the Harmonic numbers, since this puts in evidence their importance when proving Faulhaber's formula. I will keep looking for a proof of the closed form formula.

DEF Denote by $\bf H$ the Harmonic numbers $\left(1,\frac 1 2,\frac 13,\frac 14 \ldots\right)$. We define the Bernoulli numbers by $${\bf H}\star {\bf B}=\rm id$$

This is well-defined since inverses are unique. A few values are $$\left(1,-\dfrac 1 2,\dfrac 1 6,0,-\dfrac 1 {30},0,\dfrac1{42},0,-\dfrac{1}{30},\ldots\right)$$

OBS The definition says that $\sum_{k=0}^n\binom{n+1}kB_k=[n=0]$ and ${\bf B}\star {\bf 1}={\bf B}+[n=1]$

PROP For $n\geqslant 1$ we have $B_{2n+1}=0$.

P (Credits to Rob) From ${\bf B}\star {\bf 1}={\bf B}+[n=1]$ we obtain by inversion that ${\bf B} ={\bf B}\star \mu+[n=1]\star \mu$. But $[n=1]\star \mu=-(-1)^nn$; so ${\bf B}\star \mu=B+\mu{\bf N}$. Since ${\bf x}\star \mu {\bf y}=\mu(\mu{\bf x}\star {\bf y})$ we obtain $\mu{\bf B}\star {\bf 1}=\mu{\bf B}+{\bf N}$. Now we write things explicitly $$\sum_{k=0}^n\binom{n+1}{k}B_k=[n=0]\;\;,\;\;\sum_{k=0}^{m-1}\binom mk(-1)^kB_k=m$$

Thus by $m\mapsto n+1$ $$\sum_{k=0}^n\binom{n+1}{k}B_k=[n=0]\;\;,\;\;\sum_{k=0}^{n}\binom {n+1}k(-1)^kB_k=n+1$$

Thus $$n+1-[n=0]=\sum_{k=0}^n\left((-1)^k-1\right)\binom{n+1}kB_k$$ and $n\mapsto 2n+1$ gives $$2n+2=\sum_{k=0}^{2n+1}\left((-1)^k-1\right)\binom{2n+2}kB_k$$ but the $k=1$ term is $-2\binom{2n+2}1\left(-\frac1 2\right)=2n+2$ so, since even terms vanish we get $$0=-2\sum_{k=1}^{n}\binom{2n+2}{2k+1}B_{2k+1}$$ and induction does the rest. $\blacktriangle$

COR The formula ${\bf B}\star{\bf 1}=\mu{\bf B}$ holds.

First, we have:

LEMMA1 Fix $j\geqslant 0$. Then, $$\left({\bf{H}} \star {j^k}\right)_n = \frac{{{{\left( {j + 1} \right)}^{n + 1}} - {j^{n + 1}}}}{{n + 1}}$$

P This follows from direct calculation and use of the binomial theorem.

LEMMA2 Fix $m$. Then $$\left( {1 \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right)_n = \frac{{{{\left( {m + 1} \right)}^{n + 1}} - 1}}{{n + 1}}$$

P This should also be fairly obvious by integrating the binomial expansion and finding the appropriate constant of integration.

COR By the two previous lemmas: $$\left( {1 \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right) = \sum\limits_{j = 1}^m {{{\left( {{\bf{H}} \star {j^k}} \right)}_n}} $$

Now we can prove

THM (Faulhaber's Formula) $$\left( {\mu {\bf{B}} \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right) = \sum\limits_{j = 1}^m {{j^n}} $$

P We have that $$\sum\limits_{j = 1}^m {{{\left( {{\bf{H}} \star {j^k}} \right)}_n}} = \left( {1 \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right)$$

By the distributive property $$\left( {1 \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right) = {\bf{H}} \star \sum\limits_{j = 1}^m {{j^n}} $$ and after multiplication by $\bf B$, this is $$\eqalign{ & {\bf{B}} \star \left( {1 \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right) = \sum\limits_{j = 1}^m {{j^n}} \cr & \left( {\mu {\bf{B}} \star \frac{{{m^{k + 1}}}}{{k + 1}}} \right) = \sum\limits_{j = 1}^m {{j^n}} \cr} $$

as desired. $\blacktriangle$

  • $\begingroup$ After learning about dirichlet convolutions, I thought it might be fun to see if it could be useful to try and do this with power series, but never got around to it lol $\endgroup$ Commented Jul 24, 2013 at 6:02
  • $\begingroup$ The relation between $ B_k$ and $ B^*_k$ could be expressed as $ B_k = (-1)^k B^*_k $ but I do not know how this would be expressed in the convolution-notation. But perhaps this gives another impulse ... $\endgroup$ Commented Jul 26, 2013 at 5:27
  • $\begingroup$ @GottfriedHelms The correction has help me get back on track. In fact, I think I am close to a proof of Faulhaber's formula using only what I posted. $\endgroup$
    – Pedro
    Commented Jul 26, 2013 at 17:16
  • $\begingroup$ That would be really good! Just go on! $\endgroup$ Commented Jul 26, 2013 at 19:35
  • $\begingroup$ @GottfriedHelms I succeeded in proving Faulhaber's formula! =) $\endgroup$
    – Pedro
    Commented Jul 27, 2013 at 0:30

To avoid further editing to the other answer, I give you some details on the closed form formula of the Bernoulli numbers. Note that the Stirling numbers of the second kind which ounts the number of nonempty partitions in sets of size $k$ of a size $n$ set can be given explicitly by an inclusion exclusion arguement $-$ by counting of surjective functions (this is the proof I know, and it is not much complicated) $-$ as $$\left\{\begin{matrix}n\\k\end{matrix}\right\}=\frac 1 {k!}\sum_{j=0}^k\binom kj(-1)^{j-k}j^n$$

Thus the formula is $$B_n=\sum_{k=0}^n \frac{(-1)^k}{k+1}k!\left\{\begin{matrix}n\\k\end{matrix}\right\}$$

Moreover, if we define the Bernoulli polynomials as $$B_n(x)=\sum_{k=0}^n\binom nk B_k x^{n-k}$$ so that $B_n(0)=B_n$ we have to see if it is possible to prove things like $$\Delta {B_n}\left( x \right) = n{x^{n - 1}}$$ $${B_n}\left( x \right) = \sum\limits_{k = 0}^n {\frac{{{{\left( { - 1} \right)}^k}}}{{k + 1}}{\Delta ^n}{x^k}} $$ where $${\Delta ^n}{x^k} = \sum\limits_{j = 0}^n\binom nj {{{\left( { - 1} \right)}^{n - j}}{{\left( {x + j} \right)}^k}} $$

Note then that the special case $x=0$ gives the Stirling Number of the Second Kind, so $$B_n=B_n(0)=\sum_{k=0}^n \frac{(-1)^k}{k+1}k!\left\{\begin{matrix}n\\k\end{matrix}\right\}$$

is what you're after. Since EGFs have turned up, recall that the Bernoulli polynomials $B_n(t)$ have the exponential generating function $$\frac{ze^{tz}}{e^z-1}=\sum_{n\geqslant 0}B_n(t)\frac{z^n}{n!}$$

that is, using convolution, $$B_n(t)=\sum_{k=0}^n\binom nk B_k t^{n-k}$$

On the other hand, the Stirling numbers arise as $$\frac{(e^z-1)^k}{k!}=\sum_{n\geqslant 0}\left\{\begin{matrix}n\\k\end{matrix}\right\} \frac{z^n}{n!}$$

Surely this can clean up a road for a proof.


You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .