0
$\begingroup$

The classic example of the sum of product binomial coefficients is Vandermonde's Identity:

$$ \sum_{a_{1,1}+a_{1,2}+a_{1,3}=x}\binom{n}{a_{1,1}}\binom{n}{a_{1,2}}\binom{n}{a_{1,3}}=\binom{3n}{x} $$

However, I'm interested in the following twist. To motivate the question, consider the simple case in which I'd like to compute the following sum. Assume that $0\leq c_1\leq n$, $0\leq c_2\leq 2c_1$, $n$ and $x\leq 3c_1$ are given.

$$ \sum_{\bigcap\left(\substack{a_{1,1}\leq c_1,a_{1,2}\leq c_1,a_{1,3}\leq c_1\\ a_{1,1}+a_{1,2}=a_{2,1}\leq c_2,a_{1,2}+a_{1,3}=a_{2,2}\leq c_2\\a_{1,1}+a_{1,2}+a_{1,3}=x}\right)}\binom{n}{a_{1,1}}\binom{n}{a_{1,2}}\binom{n}{a_{1,3}} $$ After some algebra, this can be simplified to:

$$ \sum_{\bigcap\left(\substack{x-c_2\leq a_{1,1}\leq c_1\\x-c_1\leq a_{2,1}\leq c_2\\a_{2,1}-a_{1,1}\leq c_1}\right)}\binom{n}{a_{1,1}}\binom{n}{a_{2,1}-a_{1,1}}\binom{n}{x-a_{2,1}} $$

which only involves two terms $a_{1,1}$ and $a_{2,1}$ that we can sum over.

This is fairly simple to code and corresponds to the case $L=3$. However, I will eventually run into the complicated case where I have the following constrained sum for general $L$:

$$ \sum_{\bigcap\left(\substack{\bigcup_{k=1}^{L-1}\bigcup_{r=1}^{L-k+1}\{a_{k,r}\leq c_k \}\\ a_{L,1}=x}\right)}\binom{n}{a_{1,1}}\binom{n}{a_{1,2}}\binom{n}{a_{1,3}}\cdots\binom{n}{a_{1,L}} $$

or alternatively, like above,

$$ \sum_{\bigcap\left(\substack{ \bigcup_{k=1}^{L-1}\{x-c_{L-k}\leq a_{k,1}\leq c_k \} \\ \bigcup_{k=1}^{L-2}\bigcup_{r=1}^{L-k-1}\{a_{k+r,1}-a_{r,1}\leq c_k \}}\right)}\binom{n}{a_{1,1}}\binom{n}{a_{2,1}-a_{1,1}}\binom{n}{a_{3,1}-a_{2,1}}\cdots\binom{n}{a_{L-1,1}-a_{L-2,1}}\binom{n}{x-a_{L-1,1}} $$

where $a_{k,r}=\sum_{s=r}^{r+k-1} a_{1,s}$

I'm interested in exact solutions to the problem. Is there an efficient way to calculate this sum for general $L$ using R (or similar languages as MATLAB and Python)?

$\endgroup$
0
$\begingroup$

I believe I have a working solution that works for general $L>1$:

binom.coeff <- function(x,l,size,cvec){
 sum(apply(valid.counts(l=l,cvec=cvec,x=x),1,function(j) prod(choose(size,j))))
}
valid.counts <- function(l,cvec,x){
  vec = seq(from=0,to=cvec[1])
  lst = lapply(numeric(l), function(i) vec)
  mat = as.matrix(expand.grid(lst))

  if(l==2){
    return(mat[which(rowSums(mat)==x),])
  }
  else{
    #k=1 is redundant since all must be less than or equal to cvec[1]
    for (k in seq(from=2, to=l-1)){
      partial = NULL
      for (r in seq(l-k+1)){
        partial = which(rowSums(mat[,r:(r+k-1),drop=FALSE])<=cvec[k])
        mat = mat[unique(partial),]
      }
    }
    return(mat[which(rowSums(mat)==x),])  
  }
}
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.