I have a data set $X$ which can be broken up into, say, $n$ bins and an expected value for each bin, $a_i$. I am modeling the probability within each bin with a Poisson distribution. I would like to be able to calculate the pdf for a given $Y\subset X$ with $y_i$ events in each bin. The approach I have been using is to just write, $$P(Y)=\prod_{i=1}^nf(a_i; y_i)$$ where $f$ is the Poisson pdf, $$f(\lambda; k)=\frac{\lambda^k\exp(-\lambda)}{k!}$$ My concern is that I have tried this numerically, and have calculated maximum of $P(Y)$ for a fixed size of $Y$ over different configurations with the same size. Then, I compared this for different sizes of $Y$ to maximize the probability. I expected to get a maximum when the size of $Y$ was equal to $\sum a_i$, but that doesn't seem to be the case.

What would be the correct way to combine Poisson pdfs in this fashion?