<s>[INST] Combinatorial identity: $\sum_{i,j \ge 0} \binom{i+j}{i}^2 \binom{(a-i)+(b-j)}{a-i}^2=\frac{1}{2} \binom{(2a+1)+(2b+1)}{2a+1}$ In my research, I found this identity and as I experienced, it's surely right. But I can't give a proof for it.
Could someone help me?
This is the identity:
let $a$ and $b$ be two positive integers; then:
$\sum_{i,j \ge 0} \binom{i+j}{i}^2 \binom{(a-i)+(b-j)}{a-i}^2=\frac{1}{2} \binom{(2a+1)+(2b+1)}{2a+1}$.
 [/INST] Denote $h(x,y)=\sum_{i,j\geqslant 0} \binom{i+j}i x^iy^j=\frac1{1-(x+y)}$, $f(x,y)=\sum_{i,j\geqslant 0} \binom{i+j}i^2 x^iy^j$. We want to prove that $2xyf^2(x^2,y^2)$ is an odd (both in $x$ and in $y$) part of the function $h(x,y)$. In other words, we want to prove that $$2xyf^2(x^2,y^2)=\frac14\left(h(x,y)+h(-x,-y)-h(x,-y)-h(-x,y)\right)=\frac{2xy}{1-2(x^2+y^2)+(x^2-y^2)^2}.$$
So, our identity rewrites as $$f(x,y)=(1-2(x+y)+(x-y)^2)^{-1/2}=:f_0(x,y)$$
This is true for $x=0$, both parts become equal to $1/(1-y)$. Next, we find a differential equation in $x$ satisfied by the function $f_0$. It is not a big deal: $$\left(f_0(1-2(x+y)+(x-y)^2)\right)'_x=(x-y-1)f_0.$$
Since the initial value $f_0(0,y)$ and this relation uniquely determine the function $f_0$, it remains to check that this holds for $f(x,y)$, which is a straightforward identity with several binomials. Namely, comparing the coefficients of $x^{i-1}y^j$ we get
$$ i\left(\binom{i+j}j^2-2\binom{i+j-1}j^2-2\binom{i+j-1}i^2+\binom{i+j-2}i^2+\binom{i+j-2}j^2-2\binom{i+j-2}{i-1}^2\right) 
$$ 
for $(f(1-2(x+y)+(x-y)^2))'_x$ and $$\binom{i+j-2}j^2-\binom{i+j-1}j^2-\binom{i+j-2}{j-1}^2$$
for $(x-y-1)f$. Both guys are equal to $$-2\frac{j}{i+j-1}\binom{i+j-1}{j}^2.$$
</s>