How to numerically evaluate the CDF of this random variable? I have a discrete random variable $X=0,1,2,\ldots$ with the following probability mass function:
$$P(X=x)=\sum_{t=0}^x\binom{n}{t}p^t(1-p)^{n-t}\binom{m-t}{x-t}q^{x-t}(1-q)^{m-x}\tag{1}$$
where $0<p,q<1$, $n\gg0$ and $m=an$ with $a>1$.  The values I am using are $p=2\times10^{-3}$, $q=10^{-4}$, $n=10^5$ and $a=30$ (so $m=3\times10^6$).  I am interested in making $p$ even smaller (as low as $2\times 10^6$ and maybe lower), and increasing $n$ to as much as $10^6$ (with corresponding increase in $m$).  
Essentially, $X$ is a "binomial mixture" of binomial random variables.
I would like to calculate its cumulative distribution $F_X(s)=\sum_{x=0}^{s}P(X=x)$ for relatively small values of $s$ (where most of this distribution should be concentrated).  I would happy with tight upper and lower bounds or approximations of the CDF of $X$.
Right now, I am having a hard time even evaluating (1) at the values above.  Does anyone have ideas?  What am I doing wrong?  What else can I try?
WHAT I'VE DONE


*

*Tried plugging the values above into Wolfram Mathematica, and just
evaluating the sum as is.  Result: overflow (due to binomials).

*Tried sum $\sum_{x=0}^s\hat{p}(x)$ over the following Gaussian
approximation of the p.m.f.:
$$\hat{p}(x)=\int_{0}^{x}\frac{\exp\left[-\frac{(t-np)^2}{2np(1-p)}\right]}{\sqrt{2\pi
    np(1-p)}}\frac{\exp\left[-\frac{(x-t-(m-t)q)^2}{2(m-t)q(1-q)}\right]}{\sqrt{2\pi
    (m-t)q(1-q)}}dt$$
The approximation is due to the local limit theorem for binomial
random variables, however, I am not sure if it applies here.  I
got $\approx2.85\times 10^{-94}$ evaluating the sum for $s=40$ and
the parameters above in Wolfram Mathematica, which isn't right.
 A: 
Let $r=p+q-pq$, then $X=Y+Z$ where $Y$ is binomial $(n,r)$, $Z$ is binomial $(m-n,q)$, and $Y$ and $Z$ are independent. 

Thus, the gaussian approximations are $$Y=nr+\sqrt{nr(1-r)}\,U_n,
$$
and
$$Z=(m-n)q+\sqrt{(m-n)q(1-q)}\,V_{m-n},$$ 
where $U_n\to U$ and $V_{m-n}\to V$ in distribution, for some independent standard normal random variables $U$ and $V$.
Thus, $$X=np(1-q)+mq+\varrho W_{m,n},$$ where $W_{m,n}\to W$ in distribution, for some standard normal random variable $W$, with
$$
\varrho^2=nr(1-r)+(m-n)q(1-q).
$$
Edit: To prove the first assertion, one can identify the generating function of some random variable $X$ with the distribution in the question. This is
$$
E(s^X)=\sum_xP(X=x)s^x=\sum_t\binom{n}{t}(sp)^t(1-p)^{n-t}A_x,
$$
where
$$
A_x=\sum_x\binom{m-t}{x-t}(sq)^{x-t}(1-q)^{m-x}=(sq+1-q)^{m-t},
$$
hence
$$
E(s^X)=(sq+1-q)^{m-n}\sum_t\binom{n}{t}(sp)^t(1-p)^{n-t}(sq+1-q)^{n-t},
$$
and the sum on the RHS is
$$
(sp+(1-p)(sq+1-q))^n=(1-r+rs)^n,
$$
hence $E(s^X)=E(s^Y)E(s^Z)$, as desired.
