Computing the sums of powers of the integers
useful in statistics and data reduction
Note that for
and
Link to Tabulation
The brute force method for computing the polynomial coefficients
1a  definition of  
1b  declaration of b_{j}.There is no b_{0} because S_{r}(0) = 0and r+1 is highest power  
2a  =  substitute n+1 for n in 1a  
2b  =  extract final term of sum  
2c  =  apply 1a in reverse  
2d  =  apply 1b  
2e  =  binomial expansion of 1st term  
=  extract first term of k sum and last of j 

=  independent indexes can be merged into a single summation 

3a  =  substitute n+1 for n in 1b  
3b  =  binomial expansion  
3c  =  b_{j} does not depend on k 
Since eqn 3_{*} equals eqn 2_{*} (1a=1b) we can write r+2 equations for the coefficient
of where k=0..r+1
k  from 2e  from 3c  equating coefficients of  
4a  0  1  =  the 0th term of the inner loop of 3c is included for each j 

4b  =  for all x.  
5  r+1  =  *\binom{r+1}{r+1}  in 3c only j=R+1 and k=j generate an R+1 term in the sum 

=  not a useful equation!but we had one more than we needed so we still have enough to solve for all coefficients 

6a  p=1..R  C_{R,p+}b_{p}  =  Sum[j=p..R+1]{b_{j}*C_{j,p}}  the inner sum of 3c has a term with n^{p} when j is at least p, and has only one such term for each j 

6b  C_{R,p+}b_{p}  =  b_{p}*C_{p,p+}Sum[j=p+1..R+1]{b_{j}*C_{j,p}}  extract first term of summation  
C_{R,p}  =  Sum[j=p+1..R+1]{bj*C_{j},_{p}}  Cp,p=1 so b_{p}*C_{p,p}=b_{p} so b_{p}can be subracted from both sides 

6c  1  =  Sum[j=p+1..R+1]{b_{j}*C_{j,p}/C_{R,p}}  divide both sides by C_{R,p}  
6d  1  =  Sum[j=p+1..R+1]{b_{j}*[j!(Rp)!/R!(jp)!]}  resolve C_{x,y}‘s into factorials, cancelling p! term from each’s denominator 

Solving these equations is fairly straightforward in that b_{R+1} appears by itself in an equation (p=R), and each lower order b_{j} appears in an equation with just it and higher order terms.Therefore we can solve the equations from highest order descending with simple algebra at each step, no simultaneous solving is required. 

7  p=R  1  =  Sum[j=R+1..R+1]{b_{j}*[j!(RR)!/R!(jR)!]}  
=  b_{R+1}*[(R+1)!(RR)!/R!(R+1R)!]  only term to sum is j=R+1  
=  b_{R+1}*(R+1)  resolve factorials  
b_{R+1}  =  1/(R+1)  which is what you would expect from integral calculus 

8  p=1..R1  1  =  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+[(R+1)!(Rp)!/R!(R+1p)!]/(R+1)  substitute 7 into 6d  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+(R+1)!(Rp)!/(R+1)!(R+1p)!  merge (R+1)! =(R+1)*R!  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+[(Rp)!/(Rp+1)!]  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+1/(Rp+1) 
Each remaining equation will be of the form 1= bx * (factorials in R
and x) +
expression created by substituting all known bx’s solutions into the
definition eqn of 6d adjusting the summation range.
Expression = Sum[k=x+1..R+1]{b_{k}*[k!(Rp)!/R!(kp)!]} with
each b_{k} symbol replaced with the known expression for it.
Rearranging terms gives: bx= (1 expression)/ (4 factorials in R
and x)
9  b_{p+1}  =  (1expression)  /  (p+1)!(Rp)!/R!(p+1p)!  term from 6d exposed 
=  (1expression)  /  (p+1)p!(Rp)!/R!  (n+1)! = (n+1) n! from definition of ! and p+1p = 1 and 1!=1 

=  (1expression)  /  (p+1)/C_{R,p}  
=  (C_{R,p}/(p+1))_{ * }(1expression) 
Let’s work some examples to check our algebra so far.
The q below is a convenient variable for condensing the algebra.
q  b_{Rq}  (C_{R,p}/(p+1))  *(1  sums of j!(Rp)!/R!(jp)! (see following table) 
)  
0  b_{R}  =  (C_{R,R1}/(R))  *(1(  1/(0+2)  )  ) 
=  (R/R*1!)  *(1(  1/2  )  )  
=  1  *  1/2  
b_{R}  =  1/2  2nd coefficient is always 1/2 regardless of power 

1  b_{R1}  =  (C_{R,R2}/(R1))  *(1(  1/(1+2) + 1/2  )  ) 
=  (R(R1)/(R1)2!)  *(1(  1/3 + 1/2  )  )  
=  R/2  *(1  5/6  )  
=  R/12  3rd coefficient is simple 

since the 2nd feedback term is always 1/2 we can combine that with the (1 so we start now with (1/2 

2  b_{R2}  =  (C_{R,R3}/(R2))  *(1/2(  1/(2+2) +(2+1)/12  )  ) 
=  (C_{R,R3}/(R2))  *(1/2(  1/2  )  )  
=  0  4th coefficient is always 0 ! 

3  b_{R3}  =  (C_{R,R4}/(R3))  *(1/2  1/(3+2) + (3+1)/12  )  ) 
=  (R(R1)(R2)/4!  *(15/30(  6/30 + 10/30  )  )  
=  R(R1)(R2)/6!  6!= 720. 

we can skip use of the Cx,y symbol as we are always expanding it into its factorials,and we will directly resolve the summed terms instead of showing the fine detail of the substitutions 

4  b_{R4}  =  (R!/(R4)!5!)  *(1/2(  1/6+5/121/12  )  ) 
=  (R!/(R4)!5!)  *(6/12  (2+51)/12  )  
=  0  6th coefficient is always 0 ! 

5  b_{R5}  =  (R!/(R5)!6!)  *(1/2(  1/7+ 1/21/6  )  ) 
=  R!/(R5)!7!*6  7!*6 = 30240  
6  b_{R6}  =  (R!/(R6)!7!)  *(1/2(  1/8+7/127/24+1/12  )  ) 
=  0  8th coefficient is always 0! a pattern seems to be emerging 
feedback components
calculation
formula  *  j!  (Rp)!  /  R!  (jp)!  net  q=Rp1  
b_{R+1}  1/(R+1)  *  (R+1)!  (Rp)!  /  R!  (R+1p)!  1/(R+1p)  1/(q+2)  
b_{R}  1/2  *  R!  (Rp)!  /  R!  (Rp)!  1/2  1/2  1/Prod[i=2..2] 
b_{R1}  R/12  *  (R1)!  (Rp)!  /  R!  (R1p)!  (Rp)/12  (q+1)/12  2* Prod[i=q+1..q+1]/Prod[i=2..4] 
b_{R3}  R(R1)(R2)/720  *  (R3)!  (Rp)!  /  R!  (R3p)!  (Rp)(Rp1)(Rp2)(1/720)  (q+1)(q)(q1)/6!  Prod[i=q1..q+1]/Prod[i=2..6] 
b_{R5}  R!/(R5)!7!6  (R5)!  (Rp)!  /  R!  (R5p)!  (Rp)!/(R5p)!(1/6*7!)  (q+1)!/(q4)!6*7!  (4/3)*Prod[i=q3..q+1]/Prod[i=2..8] 
The ‘expression’ of eqn 9 is a function of the step that can be reduced
to a rational number. That means the form of b_{p} is a
combinatorial in R times a rational number.
That ‘expression’ also is zero for even values of our ‘q’ index.
Which means, that_{ }b_{R2k }= 0 for
k>=1.
Direct
proof of 4b: Sum[j=1..R+1]
{b_{j}} =1
A:  S_{R}(1)  =  Sum[j=1..1]{j^{R}}  =  1^{R}  =  1  
B:  S_{R}(1)  =  Sum[j=1..R+1]{b_{j}1^{j}}  =  Sum[j=1..R+1] {b_{j}}  
A=B so  Sum[j=1..R+1] {b_{j}}  =  1 
Checking the formulae
The column headers are derived from the algebra, the inline comments at
the ends of rows are the sum of the coefficients, which by 4c are equal
to 1.
power  1/  1/  R/  R(R1)(R2)/  R(R1)(R2)(R3)(R4)/  R(R1)(R2)(R3)(R4)(R5)(R6)/  R(R1)(R2)(R3)(R4)(R5)(R6)(R7)(R8)/  
(R+1)  2  12 sep 1  5*6*4*3!  6*7*6*5!  5*6*8*7!  10*9!*66/5  B_{C }is “Bernoulli number” 

=  C_{R,0}/2  C_{R,1}/2*2*3  C_{R,3}/3*5*6  C_{R,5}/6*7*6  C_{R,7}/8*(5*6)  C_{R,9}/10*11*6/5  C_{R,C}/(C+1)B_{C}  
0  n  
1  1/2n^{2}  +1/2n  11/2 = 1/2 

2  1/3n^{3}  +1/2n^{2}  + 1/6n  11/31/2 = 1/6 = 2/x=>x=12 

3  1/4n^{4}  +1/2n^{3}  + 1/4n^{2}  
4  1/5n^{5}  +1/2n^{4}  + 1/3n^{3} 
− 1/30n  ^{ }1(1/5+1/2+1/3)=(30(6+15+10))/30 =1/30=(4*3*2)/x =>x = 24*30 = 720 = 6! 

5  1/6n^{6}  +1/2n^{5}  + 5/12n^{4}  − 1/12n^{2}  
6  1/7n^{7}  +1/2n^{6}  + 1/2n^{5}  − 1/6 n^{3}  + 1/42 n  1(1/7+1/2+1/21/6)=(4262424+7)/42=1/42=6*5*4*3*2/x =>x=42*6!=6*7!  
7  1/8n^{8}  +1/2n^{7}  + 7/12n^{6}  − 7/24n^{4}  + 1/12 n^{2}  
8  1/9n^{9}  +1/2n^{8}  + 2/3n^{7}  – 7/15n^{5}  + 2/9 n^{3}  − 1/30 n 
1(1/9+1/2+2/37/15+2/9)=(90104560+4220)/90=3/90=1/30=8!/x => x=30*8!  
9  1/10n^{10}  +1/2n^{9}  + 3/4n^{8}  – 7/10n^{6}  + 1/2n^{4}  – 3/20n^{2}  
10  1/11n^{11}  +1/2n^{10}  + 5/6n^{9} 
– n^{7}  + n^{5}  – 1/2n^{3}  + 5/66n  5/66 = 10!/x =>x=10!*66/5 
11  1/12n^{12}  +1/2n^{11}  +11/12n^{10}  – 11/8n^{8}  + 11/6n^{6}  – 11/8n^{4 }  + 5/12 n^{2} 
Incremental computation of all coefficients up to power of interest:
 coefficient(row, col>=row)=0
 coefficient(row,0) is always 1/row. (could be folded into teh
subsequent rules, but is worth optimizing)  coefficient(row,1) is always 1/2. since we will be computing
1sum(…) where … includes this 1/2  coefficient(row, col>1 &&odd)=0
 for col<row
coefficient(row,col>2)=coefficient(row,col)*row/(rowcol)  if row is even then
coefficient(row,row1)=1/21/(row+1)sum(col=1..row3){coefficient(row,col)}
Once the first nonzero quantity in a column is calculated the
remaining terms in that column can be calculated.
The first nonzero quantity in a column is calculated from the already
calculated members of its row. So another algorithm is:
 coefficient(row, col>=row)=0
 coefficient(row,0) is always 1/row. (could be folded into teh
subsequent rules, but is worth optimizing)  coefficient(row,1) is always 1/2. since we will be computing
1sum(…) where … includes this 1/2  coefficient(row, col>1 &&odd)=0
 for col<row
coefficient(row,col>2)=coefficient(row,col)*row/(rowcol)  if row is even then
coefficient(row,row1)=1/21/(row+1)sum(col=1..row3){coefficient(row,col)}
(=
Bournoulli
number (row) )
If you don’t need all the lower
order Sk’s then use this method, derived by T.J. Pfaff from
combinatorial reasoning.
Pfaff’s method counts combinations via two different approaches and
shows that they are equal.
S_{k}(n)  =  Sum [j=1 .. k]{ Sum[i=0..j] {(1)^{ji} i^{k}Cj,iC_{n+1,j+1}}} 
TJ Pfaff  
if k>0  =  Sum [j=1 .. k]{ C_{n+1,j+1 * }Sum[i=1..j] {(1)^{ji} i^{k}Cj,i}} 
can drop useless 0 term of inner sumand move term that doesn’t depend upon sum outside of the sum 

S_{1}(n)  =  Sum [j=1 .. 1]{ C_{n+1,j+1 * }Sum[i=1..j] {(1)^{ji} i^{1}Cj,i}} 
k>1  
=  C_{n+1,1+1 * }Sum[i=1..1] {(1)^{1i} i^{1}C1,i} 
expand sum over j  
=  C_{n+1,2 * }(1)^{11} 1C_{1,1} 
expand sums over i  
=  C_{n+1,2} = (n+1)(n)/2 
correct!  
S_{2}(n)  =  Sum [j=1 .. 2]{ Sum[i=1..j] {(1)^{ji} i^{2}Cj,iC_{n+1,j+1}}} 
R>2  
=  C_{n+1,1+1 }* Sum[i=1..1] {(1)^{1i} i^{2}C_{1},i}+C_{n+1,2+1 }* Sum[i=1..2] {(1)^{2i} i^{2}C2,i} 
expand sum over j  
=  Cn+1,2* (1)^{11} 1^{2}C_{1},_{1}_{
}+C_{n+1,3} * ((1)^{21} 1^{2}C_{2,1} +(1)^{22} 2^{2}C_{2,2} ) 
expand sums over i  
C_{n+1,2 }+ C_{n+1,3 }*(2+ 4)  simple reductions  
C_{n+1,2 }+2*C_{n+1,2} * (n1)/3 
Cx,r+1 = Cx,r (xr)/(r+1)  
C_{n+1,2 }*(1+(2*(n1)/3)  
C_{n+1,2 }*(3+2n2)/3)  
C_{n+1,2 }*(2n+1)/3  
n(n+1)*(2n+1)/(2*3)  correct! 
Pfaff’s approach readily
generates the expressions factored.
power  integer factorization  convenient relationships  nested monomial form 
1  n(n+1)/2  S1  
2  n(n+1)(2n+1)/6  S1*(2n+1)/3  
3  n^{2}(n+1)^{2}/4  S1^{2}  
4  n(n+1)(2n+1)(3n^{2}+3n1)/(6*5)  S2*(3n^{2}+3n1)/5  S2*((n+1)3n1)/5 
5  n^{2}(n+1)^{2 }(2n^{2}+2n1)/(4*3)  S3*(2n^{2}+2n1)/3  S3*((n+1)2n1)/3 
6  n(n+1)(2n+1)(3n^{4}+6n^{3}3n+1)/(6*7)  S2*(3n^{4}+6n^{3}3n+1)/7  S2*((n+2)n^{2}1)3n+1)/7 
7  n^{2}(n+1)^{2}(3n^{4}+6n^{3}n^{2}4n+2)/(4*6)  S3*(3n^{4}+6n^{3}n^{2}4n+2)/6  S3*((((n+2)3n1)n4)n+2)/6 
8  n(n+1)(2n+1)(5n^{6}+15n^{5}+5n^{4}15n^{3}n^{2}+9n3)/(6*15)  S2*(5n^{6}+15n^{5}+5n^{4}15n^{3}n^{2}+9n3)/15  S2*((((((n+3)n+1)n3)5n1)n+9)n3)/15 
9  n^{2}(n+1)^{2}(n^{2}+n1)(2n4+4n3n23n+3)/(4*5)  S3*(n^{2}+n1)(2n^{4}+4n^{3}n^{2}3n+3)/5  S3*((n+1)n1)*((((n+2)2n1)n3)n+3)/5 
10  n(n+1)(2n+1)(n^{2}+n1)(3n^{6}+9n^{5}+2n^{4}11n^{3}+3n^{2}+10n5)/(6*11)  S2*(n^{2}+n1)(3n^{6}+9n^{5}+2n^{4}11n^{3}+3n^{2}+10n5)/11  S3*((n+1)n1)*((((((n+3)3n+2)n11)n+3)n+10)n5)/11 
Observations
 all sums have a factor of n as the ‘0’ term of the polynomial has
a
coefficient of 0.  all sums for R>0 are also divisible by (n+1).
 all sums for Even R>0 are divisible by S_{2}: S_{2k}(n)=S_{2}(n)R_{k}(n)
where
R_{k}(n) is a polynomial of degree 2*(k1).  all sums for Odd R>3 are divisible by S_{3}: S_{2k+1}(n)=S_{3}(n)*T_{k}(n) where T_{k}(n) is a
polynomial of degree 2*k.
S_{R}(n)=Sum[i=1..R]{
Sum[j=0..i]{1^{j} (ij)^{R}
C_{n+pi+1,ni} C_{R+1,j}}}
sumsOfPowersOfIntegers
Computing the sums of powers of the integers
useful in statistics and data reduction
Note that for
and
first, a tabulation
power  as sum  factored 
0  n  n 
1  1/2n^{2} + 1/2n  n(n+1)/2 
2  1/3n^{3} + 1/2n^{2} + 1/6n 
n(n+1)(2n+1)/6 
3  1/4n^{4} + 1/2n^{3} + 1/4n^{2} 
n^{2}(n+1)^{2}/4 
4  1/5n^{5} + 1/2n^{4
}+ 1/3n^{3} − 1/30n 
n(n+1)(2n+1)(3n^{2}+3n1)/(6*5) 
5  1/6n^{6} + 1/2n^{5} + 5/12n^{4}− 1/12 n^{2} 
n^{2}(n+1)^{2 }(2n^{2}+2n1)/(4*3) 
6  1/7n^{7} + 1/2n^{6} + 1/2n^{5} − 1/6 n^{3} + 1/42 n 
n(n+1)(2n+1)(3n^{4}+6n^{3}3n+1)/(6*7) 
7  1/8n^{8} + 1/2n^{7} + 7/12n^{6}− 7/24 n^{4} + 1/12 n^{2} 
n^{2}(n+1)^{2}(3n^{4}+6n^{3}n^{2}4n+2)/(4*6) 
8  1/9n^{9} + 1/2n^{8} + 2/3n^{7} – 7/15 n^{5} + 2/9 n^{3} − 1/30 n 
n(n+1)(2n+1)(5n^{6}+15n^{5}+5n^{4}15n^{3}n^{2}+9n3)/(6*15) 
9  1/10n^{10}+1/2n^{9} + 3/4n^{8} – 7/10 n^{6} + 1/2n^{4} – 3/20n^{2} 
n^{2}(n+1)^{2}(n^{2}+n1)(2n^{4}+4n^{3}n^{2}3n+3)/(4*5) 
10  1/11n^{11}+1/2n^{10}+ 5/6n^{9} – n^{7} + n^{5} – 1/2n^{3} + 5/66 n 
n(n+1)(2n+1)(n^{2}+n1)(3n^{6}+9n^{5}+2n^{4}11n^{3}+3n^{2}+10n5)/(6*11) 
Some convenient relationships exist:
For odd r>=3: $latex S_r = S_3*R_{r2}$ where is a polynomial of degree k.
For even r>=2 where is a
polynomial of degree k.
The brute force method
1a  definition of  
1b  declaration of b_{j}.There is no b_{0} because S_{r}(0) = 0and r+1 is highest power  
2a  =  substitute n+1 for n in 1a  
2b  =  extract final term of sum  
2c  =  apply 1a in reverse  
2d  =  apply 1b  
2e  =  binomial expansion of 1st term  
=  extract first term of k sum and last of j 

=  independent indexes can be merged into a single summation 

3a  =  substitute n+1 for n in 1b  
3b  =  $latex \sum_{j=1}^{r+1} {b_{j}( \sum_{k=0}^j{\binom{j}{k}n^k})}$ 
binomial expansion  
3c  =  b_{j} does not depend on k 
Since eqn 3_{*} equals eqn 2_{*} (1a=1b) we can write
R+2
equations for the coefficient
of n^{k} where k=0..R+1
k  from 2e  from 3c  equating coefficients of n^{j}  
4a  0  1  =  Sum[j=1..R+1]{ b_{j}*C_{j,0}}  the 0th term of the inner loop of 3c is included for each j 

4b  =  Sum[j=1..R+1] {b_{j}}  C_{x,0}=1 for all x.  
5  R+1  b_{R+1}  =  b_{R+1}*C_{R+1,R+1}  in 3c only j=R+1 and k=j generate an R+1 term in the sum 

b_{R+1}  =  b_{R+1}*1  not a useful equation!but we had one more than we needed so we still have enough to solve for all coefficients 

6a  p=1..R  C_{R,p+}b_{p}  =  Sum[j=p..R+1]{b_{j}*C_{j,p}}  the inner sum of 3c has a term with n^{p} when j is at least p, and has only one such term for each j 

6b  C_{R,p+}b_{p}  =  b_{p}*C_{p,p+}Sum[j=p+1..R+1]{b_{j}*C_{j,p}}  extract first term of summation  
C_{R,p}  =  Sum[j=p+1..R+1]{bj*C_{j},_{p}}  Cp,p=1 so b_{p}*C_{p,p}=b_{p} so b_{p}can be subracted from both sides 

6c  1  =  Sum[j=p+1..R+1]{b_{j}*C_{j,p}/C_{R,p}}  divide both sides by C_{R,p}  
6d  1  =  Sum[j=p+1..R+1]{b_{j}*[j!(Rp)!/R!(jp)!]}  resolve C_{x,y}‘s into factorials, cancelling p! term from each’s denominator 

Solving these equations is fairly straightforward in that b_{R+1} appears by itself in an equation (p=R), and each lower order b_{j} appears in an equation with just it and higher order terms.Therefore we can solve the equations from highest order descending with simple algebra at each step, no simultaneous solving is required. 

7  p=R  1  =  Sum[j=R+1..R+1]{b_{j}*[j!(RR)!/R!(jR)!]}  
=  b_{R+1}*[(R+1)!(RR)!/R!(R+1R)!]  only term to sum is j=R+1  
=  b_{R+1}*(R+1)  resolve factorials  
b_{R+1}  =  1/(R+1)  which is what you would expect from integral calculus 

8  p=1..R1  1  =  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+[(R+1)!(Rp)!/R!(R+1p)!]/(R+1)  substitute 7 into 6d  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+(R+1)!(Rp)!/(R+1)!(R+1p)!  merge (R+1)! =(R+1)*R!  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+[(Rp)!/(Rp+1)!]  
=  Sum[j=p+1..R]{b_{j}*[j!(Rp)!/R!(jp)!]}+1/(Rp+1) 
Each remaining equation will be of the form 1= bx * (factorials in R
and x) +
expression created by substituting all known bx’s solutions into the
definition eqn of 6d adjusting the summation range.
Expression = Sum[k=x+1..R+1]{b_{k}*[k!(Rp)!/R!(kp)!]} with
each b_{k} symbol replaced with the known expression for it.
Rearranging terms gives: bx= (1 expression)/ (4 factorials in R
and x)
9  b_{p+1}  =  (1expression)  /  (p+1)!(Rp)!/R!(p+1p)!  term from 6d exposed 
=  (1expression)  /  (p+1)p!(Rp)!/R!  (n+1)! = (n+1) n! from definition of ! and p+1p = 1 and 1!=1 

=  (1expression)  /  (p+1)/C_{R,p}  
=  (C_{R,p}/(p+1))_{ * }(1expression) 
Let’s work some examples to check our algebra so far.
The q below is a convenient variable for condensing the algebra.
q  b_{Rq}  (C_{R,p}/(p+1))  *(1  sums of j!(Rp)!/R!(jp)! (see following table) 
)  
0  b_{R}  =  (C_{R,R1}/(R))  *(1(  1/(0+2)  )  ) 
=  (R/R*1!)  *(1(  1/2  )  )  
=  1  *  1/2  
b_{R}  =  1/2  2nd coefficient is always 1/2 regardless of power 

1  b_{R1}  =  (C_{R,R2}/(R1))  *(1(  1/(1+2) + 1/2  )  ) 
=  (R(R1)/(R1)2!)  *(1(  1/3 + 1/2  )  )  
=  R/2  *(1  5/6  )  
=  R/12  3rd coefficient is simple 

since the 2nd feedback term is always 1/2 we can combine that with the (1 so we start now with (1/2 

2  b_{R2}  =  (C_{R,R3}/(R2))  *(1/2(  1/(2+2) +(2+1)/12  )  ) 
=  (C_{R,R3}/(R2))  *(1/2(  1/2  )  )  
=  0  4th coefficient is always 0 ! 

3  b_{R3}  =  (C_{R,R4}/(R3))  *(1/2  1/(3+2) + (3+1)/12  )  ) 
=  (R(R1)(R2)/4!  *(15/30(  6/30 + 10/30  )  )  
=  R(R1)(R2)/6!  6!= 720. 

we can skip use of the Cx,y symbol as we are always expanding it into its factorials,and we will directly resolve the summed terms instead of showing the fine detail of the substitutions 

4  b_{R4}  =  (R!/(R4)!5!)  *(1/2(  1/6+5/121/12  )  ) 
=  (R!/(R4)!5!)  *(6/12  (2+51)/12  )  
=  0  6th coefficient is always 0 ! 

5  b_{R5}  =  (R!/(R5)!6!)  *(1/2(  1/7+ 1/21/6  )  ) 
=  R!/(R5)!7!*6  7!*6 = 30240  
6  b_{R6}  =  (R!/(R6)!7!)  *(1/2(  1/8+7/127/24+1/12  )  ) 
=  0  8th coefficient is always 0! a pattern seems to be emerging 
feedback components
calculation
formula  *  j!  (Rp)!  /  R!  (jp)!  net  q=Rp1  
b_{R+1}  1/(R+1)  *  (R+1)!  (Rp)!  /  R!  (R+1p)!  1/(R+1p)  1/(q+2)  
b_{R}  1/2  *  R!  (Rp)!  /  R!  (Rp)!  1/2  1/2  1/Prod[i=2..2] 
b_{R1}  R/12  *  (R1)!  (Rp)!  /  R!  (R1p)!  (Rp)/12  (q+1)/12  2* Prod[i=q+1..q+1]/Prod[i=2..4] 
b_{R3}  R(R1)(R2)/720  *  (R3)!  (Rp)!  /  R!  (R3p)!  (Rp)(Rp1)(Rp2)(1/720)  (q+1)(q)(q1)/6!  Prod[i=q1..q+1]/Prod[i=2..6] 
b_{R5}  R!/(R5)!7!6  (R5)!  (Rp)!  /  R!  (R5p)!  (Rp)!/(R5p)!(1/6*7!)  (q+1)!/(q4)!6*7!  (4/3)*Prod[i=q3..q+1]/Prod[i=2..8] 
The ‘expression’ of eqn 9 is a function of the step that can be reduced
to a rational number. That means the form of b_{p} is a
combinatorial in R times a rational number.
That ‘expression’ also is zero for even values of our ‘q’ index.
Which means, that_{ }b_{R2k }= 0 for
k>=1.
Direct
proof of 4b: Sum[j=1..R+1]
{b_{j}} =1
A:  S_{R}(1)  =  Sum[j=1..1]{j^{R}}  =  1^{R}  =  1  
B:  S_{R}(1)  =  Sum[j=1..R+1]{b_{j}1^{j}}  =  Sum[j=1..R+1] {b_{j}}  
A=B so  Sum[j=1..R+1] {b_{j}}  =  1 
Checking the formulae
The column headers are derived from the algebra, the inline comments at
the ends of rows are the sum of the coefficients, which by 4c are equal
to 1.
power  1/  1/  R/  R(R1)(R2)/  R(R1)(R2)(R3)(R4)/  R(R1)(R2)(R3)(R4)(R5)(R6)/  R(R1)(R2)(R3)(R4)(R5)(R6)(R7)(R8)/  
(R+1)  2  12 sep 1  5*6*4*3!  6*7*6*5!  5*6*8*7!  10*9!*66/5  B_{C }is “Bernoulli number” 

=  C_{R,0}/2  C_{R,1}/2*2*3  C_{R,3}/3*5*6  C_{R,5}/6*7*6  C_{R,7}/8*(5*6)  C_{R,9}/10*11*6/5  C_{R,C}/(C+1)B_{C}  
0  n  
1  1/2n^{2}  +1/2n  11/2 = 1/2 

2  1/3n^{3}  +1/2n^{2}  + 1/6n  11/31/2 = 1/6 = 2/x=>x=12 

3  1/4n^{4}  +1/2n^{3}  + 1/4n^{2}  
4  1/5n^{5}  +1/2n^{4}  + 1/3n^{3} 
− 1/30n  ^{ }1(1/5+1/2+1/3)=(30(6+15+10))/30 =1/30=(4*3*2)/x =>x = 24*30 = 720 = 6! 

5  1/6n^{6}  +1/2n^{5}  + 5/12n^{4}  − 1/12n^{2}  
6  1/7n^{7}  +1/2n^{6}  + 1/2n^{5}  − 1/6 n^{3}  + 1/42 n  1(1/7+1/2+1/21/6)=(4262424+7)/42=1/42=6*5*4*3*2/x =>x=42*6!=6*7!  
7  1/8n^{8}  +1/2n^{7}  + 7/12n^{6}  − 7/24n^{4}  + 1/12 n^{2}  
8  1/9n^{9}  +1/2n^{8}  + 2/3n^{7}  – 7/15n^{5}  + 2/9 n^{3}  − 1/30 n 
1(1/9+1/2+2/37/15+2/9)=(90104560+4220)/90=3/90=1/30=8!/x => x=30*8!  
9  1/10n^{10}  +1/2n^{9}  + 3/4n^{8}  – 7/10n^{6}  + 1/2n^{4}  – 3/20n^{2}  
10  1/11n^{11}  +1/2n^{10}  + 5/6n^{9} 
– n^{7}  + n^{5}  – 1/2n^{3}  + 5/66n  5/66 = 10!/x =>x=10!*66/5 
11  1/12n^{12}  +1/2n^{11}  +11/12n^{10}  – 11/8n^{8}  + 11/6n^{6}  – 11/8n^{4 }  + 5/12 n^{2} 
Incremental computation of all coefficients up to power of interest:
 coefficient(row, col>=row)=0
 coefficient(row,0) is always 1/row. (could be folded into teh
subsequent rules, but is worth optimizing)  coefficient(row,1) is always 1/2. since we will be computing
1sum(…) where … includes this 1/2  coefficient(row, col>1 &&odd)=0
 for col<row
coefficient(row,col>2)=coefficient(row,col)*row/(rowcol)  if row is even then
coefficient(row,row1)=1/21/(row+1)sum(col=1..row3){coefficient(row,col)}
Once the first nonzero quantity in a column is calculated the
remaining terms in that column can be calculated.
The first nonzero quantity in a column is calculated from the already
calculated members of its row. So another algorithm is:
 coefficient(row, col>=row)=0
 coefficient(row,0) is always 1/row. (could be folded into teh
subsequent rules, but is worth optimizing)  coefficient(row,1) is always 1/2. since we will be computing
1sum(…) where … includes this 1/2  coefficient(row, col>1 &&odd)=0
 for col<row
coefficient(row,col>2)=coefficient(row,col)*row/(rowcol)  if row is even then
coefficient(row,row1)=1/21/(row+1)sum(col=1..row3){coefficient(row,col)}
(=
Bournoulli
number (row) )
If you don’t need all the lower
order Sk’s then use this method, derived by T.J. Pfaff from
combinatorial reasoning.
Pfaff’s method counts combinations via two different approaches and
shows that they are equal.
S_{k}(n)  =  Sum [j=1 .. k]{ Sum[i=0..j] {(1)^{ji} i^{k}Cj,iC_{n+1,j+1}}} 
TJ Pfaff  
if k>0  =  Sum [j=1 .. k]{ C_{n+1,j+1 * }Sum[i=1..j] {(1)^{ji} i^{k}Cj,i}} 
can drop useless 0 term of inner sumand move term that doesn’t depend upon sum outside of the sum 

S_{1}(n)  =  Sum [j=1 .. 1]{ C_{n+1,j+1 * }Sum[i=1..j] {(1)^{ji} i^{1}Cj,i}} 
k>1  
=  C_{n+1,1+1 * }Sum[i=1..1] {(1)^{1i} i^{1}C1,i} 
expand sum over j  
=  C_{n+1,2 * }(1)^{11} 1C_{1,1} 
expand sums over i  
=  C_{n+1,2} = (n+1)(n)/2 
correct!  
S_{2}(n)  =  Sum [j=1 .. 2]{ Sum[i=1..j] {(1)^{ji} i^{2}Cj,iC_{n+1,j+1}}} 
R>2  
=  C_{n+1,1+1 }* Sum[i=1..1] {(1)^{1i} i^{2}C_{1},i}+C_{n+1,2+1 }* Sum[i=1..2] {(1)^{2i} i^{2}C2,i} 
expand sum over j  
=  Cn+1,2* (1)^{11} 1^{2}C_{1},_{1}_{
}+C_{n+1,3} * ((1)^{21} 1^{2}C_{2,1} +(1)^{22} 2^{2}C_{2,2} ) 
expand sums over i  
C_{n+1,2 }+ C_{n+1,3 }*(2+ 4)  simple reductions  
C_{n+1,2 }+2*C_{n+1,2} * (n1)/3 
Cx,r+1 = Cx,r (xr)/(r+1)  
C_{n+1,2 }*(1+(2*(n1)/3)  
C_{n+1,2 }*(3+2n2)/3)  
C_{n+1,2 }*(2n+1)/3  
n(n+1)*(2n+1)/(2*3)  correct! 
Pfaff’s approach readily
generates the expressions factored.
power  integer factorization  convenient relationships  nested monomial form 
1  n(n+1)/2  S1  
2  n(n+1)(2n+1)/6  S1*(2n+1)/3  
3  n^{2}(n+1)^{2}/4  S1^{2}  
4  n(n+1)(2n+1)(3n^{2}+3n1)/(6*5)  S2*(3n^{2}+3n1)/5  S2*((n+1)3n1)/5 
5  n^{2}(n+1)^{2 }(2n^{2}+2n1)/(4*3)  S3*(2n^{2}+2n1)/3  S3*((n+1)2n1)/3 
6  n(n+1)(2n+1)(3n^{4}+6n^{3}3n+1)/(6*7)  S2*(3n^{4}+6n^{3}3n+1)/7  S2*((n+2)n^{2}1)3n+1)/7 
7  n^{2}(n+1)^{2}(3n^{4}+6n^{3}n^{2}4n+2)/(4*6)  S3*(3n^{4}+6n^{3}n^{2}4n+2)/6  S3*((((n+2)3n1)n4)n+2)/6 
8  n(n+1)(2n+1)(5n^{6}+15n^{5}+5n^{4}15n^{3}n^{2}+9n3)/(6*15)  S2*(5n^{6}+15n^{5}+5n^{4}15n^{3}n^{2}+9n3)/15  S2*((((((n+3)n+1)n3)5n1)n+9)n3)/15 
9  n^{2}(n+1)^{2}(n^{2}+n1)(2n4+4n3n23n+3)/(4*5)  S3*(n^{2}+n1)(2n^{4}+4n^{3}n^{2}3n+3)/5  S3*((n+1)n1)*((((n+2)2n1)n3)n+3)/5 
10  n(n+1)(2n+1)(n^{2}+n1)(3n^{6}+9n^{5}+2n^{4}11n^{3}+3n^{2}+10n5)/(6*11)  S2*(n^{2}+n1)(3n^{6}+9n^{5}+2n^{4}11n^{3}+3n^{2}+10n5)/11  S3*((n+1)n1)*((((((n+3)3n+2)n11)n+3)n+10)n5)/11 
Observations
 all sums have a factor of n as the ‘0’ term of the polynomial has
a
coefficient of 0.  all sums for R>0 are also divisible by (n+1).
 all sums for Even R>0 are divisible by S_{2}: S_{2k}(n)=S_{2}(n)R_{k}(n)
where
R_{k}(n) is a polynomial of degree 2*(k1).  all sums for Odd R>3 are divisible by S_{3}: S_{2k+1}(n)=S_{3}(n)*T_{k}(n) where T_{k}(n) is a
polynomial of degree 2*k.
S_{R}(n)=Sum[i=1..R]{
Sum[j=0..i]{1^{j} (ij)^{R}
C_{n+pi+1,ni} C_{R+1,j}}}