A relationship between matrices, bernoulli polynomials, and binomial coefficients

We define the following polynomials, for $n≥0$: $$P_n(x)=(x+1)^{n+1}-x^{n+1}=\sum_{k=0}^{n}{\binom{n+1}{k}x^k}$$ For $n=0,1,2,3$ this gives us, $$P_0(x)=1\enspace P_1(x)=2x+1\enspace P_2(x)=3x^2+3x+1\enspace P_3(x)=4x^3+6x^2+4x+1$$

We then define the set $P_{(3)}=\{P_0,P_1,P_2,P_3\}$. It can be easily shown that this set is a basis over the vector space of polynomials of degree $3$ and lower. We take $3$ for the sake of brevity.

Taking the coefficients of these polynomials and turning them into column vectors, we can construct the matrix (coefficients from the lowest term to the highest term) $$\large{M_{P_{(3)}}}=\begin{pmatrix} 1 & 1 & 1 & 1 \\ 0 & 2 & 3 & 4 \\ 0 & 0 & 3 & 6 \\ 0 & 0 & 0 & 4 \end{pmatrix}$$ We'll call this matrix the pascal in the context of this post, and the above polynomials as pascal polynomials. The inverse of this matrix is the matrix, $$M_{P_{(3)}}^{-1}=\begin{pmatrix} 1 & -\frac{1}{2} & \frac{1}{6} & 0 \\ 0 & \frac{1}{2} & -\frac{1}{2} & \frac{1}{4} \\ 0 & 0 & \frac{1}{3} & -\frac{1}{2} \\ 0 & 0 & 0 & \frac{1}{4} \end{pmatrix}$$

We'll factor this matrix into two matrices as follows:

$$M_{P_{(3)}}^{-1}=\begin{pmatrix} 1 & 0 & 0 & 0 \\ 0 & \frac{1}{2} & 0 & 0 \\ 0 & 0 & \frac{1}{3} & 0 \\ 0 & 0 & 0 & \frac{1}{4} \end{pmatrix}×\begin{pmatrix} 1 & -\frac{1}{2} & \frac{1}{6} & 0 \\ 0 & 1 & -1 & \frac{1}{2} \\ 0 & 0 & 1 & -\frac{3}{2} \\ 0 & 0 & 0 & 1 \end{pmatrix}$$ We can see the Bernoulli numbers in the first row of the matrix. Every column is a coefficient vector of a Bernoulli polynomial

The following are extended versions of these matrices:

Are there accepted names for these matrices and polynomials? What is the meaning of these relationships?

In particular, is there some treatment of using these matrices as change of basis transformations between representations of polynomials? E.g. from a linear combination of pascal polynomials to a linear combination of monomial terms.

-
+1 for a great question! –  Alex Nelson Dec 25 '12 at 2:49
Many thanks to whoever fixed the formatting. I'm not very good with latex :( And thanks for the compliment :P –  GreĝRos Dec 25 '12 at 2:57
Oh, I hate to burst the magic: your matrix factorization is incorrect. If you carry out the matrix multiplication, you don't recover the correct matrix :( –  Alex Nelson Dec 25 '12 at 4:09
This looks like it is related. –  Mike Spivey Dec 25 '12 at 5:28
I've done a long deal with that matrix-representation myself; the discussions in go.helms-net.de/math/pascal/bernoulli_en.pdf might be of interest. One step further one can use that matrices to define "zeta"-polynomials for the summing of like powers: go.helms-net.de/math/binomial_new/04_3_SummingOfLikePowers.pdf . And even the Euler/McLaurin-formula can be found in that matrices... go.helms-net.de/math/divers/EulerMacLaurin.pdf –  Gottfried Helms Dec 25 '12 at 10:44

Here's another way to look at this. The Bernoulli polynomials can be defined by the property

$$\int_x^{x+1} B_n(u) \, du = x^n.$$

So if we let $T$ be the operator from the set of polynomials to itself given by $(Tf)(x) = \int_x^{x+1} f(u) \, du$, then we have $(TB_n)(x) = x^n$. The operator $T$ sends $x^n$ to $$\int_x^{x+1} u^n \, du = \frac{1}{n+1}\left((x+1)^{n+1} - x^{n+1}\right) = \sum_{k=0}^{n} \frac{1}{n+1} {n+1 \choose k}x^k.$$

Writing $T$ as an infinite matrix with respect to the basis $1,x,x^2, \ldots$, gives

$$T = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{rrrr} 1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{4} & \ldots\\ 0 & 1 & 1 & 1 & \ldots\\ 0 & 0 & 1 & \frac{3}{2} & \ldots\\ 0 & 0 & 0 & 1& \ldots\\ \vdots & \vdots & \vdots & \vdots & \ddots \end{array}\right).$$

We can find approximate the inverse by taking the inverse of this truncation, giving the matrix form of $T^{-1}$:

$$T^{-1} = \newcommand{\Bold}[1]{\mathbf{#1}}\left(\begin{array}{rrrr} 1 & -\frac{1}{2} & \frac{1}{6} & 0 & \ldots\\ 0 & 1 & -1 & \frac{1}{2} & \ldots\\ 0 & 0 & 1 & -\frac{3}{2} & \ldots\\ 0 & 0 & 0 & 1& \ldots\\ \vdots & \vdots & \vdots & \vdots & \ddots \end{array}\right)$$

This operator sends $x^n$ to $B_n(x)$, so the columns are the coefficients of Bernoulli polynomials.


The matrix on the right is your $M_P$ (infinitely extended). Calling the diagonal matrix on the left $D$, we have $T = DM_p$, or $DT^{-1} = M_P^{-1}$ which is your last equation.

To my mind, this is the most natural way to compute Bernoulli polynomials. A few years ago I was playing around with this idea before I knew what a Bernoulli polynomial was. I was slightly disappointed, although not too surprised, to hear that someone else had discovered this first. I was beaten by some three hundred years, no less. It ended up as part of my undergraduate thesis. :)

-
+1. What was the main topic of your undergraduate thesis, if you don't mind me asking? –  Mike Spivey Dec 25 '12 at 18:18
Well, that was the main topic - Bernoulli numbers. I did a little bit with related topics/applications like Euler-Maclaurin summation and the umbral calculus. –  Jair Taylor Dec 25 '12 at 18:45

This is a very fascinating problem, but in the $P_{(4)}$ case, we have $$\tag{1} M_{P_{(4)}}^{-1}=\begin{pmatrix} 1 & -1/2 & 1/6 & 0 & -1/30\\ & 1/2 & -1/2 & 1/4 & 0\\ & & 1/3 & -1/2 & 1/3\\ & & & 1/4 & -1/2\\ & & & & 1/5 \end{pmatrix}$$ This factors out as $$M_{P_{(4)}}^{-1}=\begin{pmatrix} 1 & -1 & 1 & 0 & -1/6\\ & 1 & -3/2 & 1 & 0\\ & & 1 & -2 & 5/3\\ & & & 1 & -5/2\\ & & & & 1 \end{pmatrix} \begin{pmatrix} 1 & & & & \\ & 1/2 & & & \\ & & 1/3 & & \\ & & & 1/4 & \\ & & & & 1/5 \end{pmatrix}$$ Unfortunately, as we can read off, the columns in the matrix on the left do not give us Bernoulli polynomials :(

Although, please note the first row of the matrix in (1) gives us Bernoulli numbers. Perhaps the first row of the inverse matrix does give you the Bernoulli numbers, but you do not obtain the Bernoulli polynomials.

Again, read this with caution and suspicion: the good reader should compute this and verify it for him or herself!

Write out

$$M_{P_{(n)}} = DU$$

where $D$ is a diagonal matrix, and $U=I+X$ is a unit upper triangular matrix (here $X$ is the strictly upper triangular part of $U$). The claim is that

$$(DU)^{-1}=(I+X)^{-1}D^{-1}=(I-X+X^{2}-X^{3}+\dots+(-1)^{n}X^{n})D^{-1}$$

produces the coefficients of Bernoulli polynomials in the columns, and the first row is the Bernoulli numbers. I claim the first (Bernoulli polynomial data) implies the latter (the first row consists of Bernoulli numbers), and claim this is obvious.

Conjecture

It seems that the OP found the matrix version for the Bernoulli polynomials described on Wikipedia using forward differences.

I would have urged the OP to change notation, and use lower triangular matrices. Why? Because then you could write $$M_{P_{(n)}}\begin{bmatrix}0\\ 0\\\vdots\\0\\ x^{n}\end{bmatrix} \cong (1+x)^{n+1}-x^{n+1}$$ in vector form. Current notation demands we use row vectors for polynomials.

I'm about to go to bed, so I do not believe I have time to prove my conjecture. If no one has proven it by tomorrow, I'll try writing up a proof.

-
You're absolutely right. I'm very sorry, I made a silly error. The factorization needs to be reversed. –  GreĝRos Dec 25 '12 at 4:21
I just ran everything through Mathematica. Alex, your calcs are correct, and Greg's factorization does work if the order is reversed. The presence of Bernoulli polynomials also seems to hold for every n I've checked (up to 100). –  AndrewG Dec 25 '12 at 4:34
@AndrewGibson: Many thanks for double checking my work! I really appreciate it :) –  Alex Nelson Dec 25 '12 at 4:42
@AlexNelson no problem! :) I have no idea how to go about investigating this rigorously, so at least plugging in some code gives me something to do. I'm really fascinated by this though and would love to understand what's going on. –  AndrewG Dec 25 '12 at 5:21

Regarding the interpretation of the matrices: Let $$X = \left(\begin{array}{c} 1 \\ x \\ \vdots \\ x^n \end{array}\right) \qquad \textrm{and}\qquad A = \left(\begin{array}{c} a_0 \\ a_1 \\ \vdots \\ a_n \end{array}\right)$$ so $$X^T A = \sum_{k=0}^n a_k x^k.$$ (I've written $X^T A$ instead of the preferred $A^T X$ to agree with your notation and avoid things like $(M^{-1})^T$.) Then $$\begin{eqnarray*} X^T A &=& X^T M M^{-1} A \\ &=& {X'}^T A' \end{eqnarray*}$$ where $X' = M^T X$ is the new basis and $A' = M^{-1} A$ are the coefficients in the new basis. This interpretation holds for any invertible transformation $M$. The formalism can be generalized in a straightforward way to transformations between any two bases of polynomials.

For your transformation we find the coefficients of $\sum_{k=0}^3 a_k x^k$ in the "Pascal basis" are $$M_{P_{(3)}}^{-1}A = \left( \begin{array}{c} a_0 - \frac{1}{2} a_1 + \frac{1}{6} a_2 \\ \frac{1}{2}a_1 - \frac{1}{2} a_2 + \frac{1}{4} a_3 \\ \frac{1}{3} a_2 - \frac{1}{2} a_3 \\ \frac{1}{4} a_3 \end{array} \right).$$ Thus, $$\begin{eqnarray*} \sum_{k=0}^3 a_k x^k &=& \left(a_0 - \frac{1}{2} a_1 + \frac{1}{6} a_2\right) P_0(x) \\ &&+ \left(\frac{1}{3} a_2 - \frac{1}{2} a_3\right) P_1(x) \\ &&+ \left(\frac{1}{3} a_2 - \frac{1}{2} a_3\right) P_2(x) \\ &&+ \frac{1}{4} a_3 P_3(x). \end{eqnarray*}$$

-