Fitting data to a Poisson distribution, what are the errors? I have large data set of $\approx 10^6$ points where each point contains the information of a
(year, count)

of a particular event. There are many counts for a each year. Looking at the histograms of counts for a given year, I've conjectured that they follow a Poisson distribution. By generating the Poisson distribution with the maximum likelihood parameter $P(\lambda_{MLE})$, $\lambda_{MLE}=\frac{1}{n}\sum{k_i}$. I think what I need to do next is compare the values of this expected distribution to my normalized histogram and compute a $\chi^2$ statistic. Wiki however warns that:

In cases where the expected value, E, is found to be small (indicating a small underlying population probability, and/or a small number of observations), the normal approximation of the multinomial distribution can fail, and in such cases it is found to be more appropriate to use the G-test, a likelihood ratio-based test statistic.

This will certainly be the case for larger values of $\lambda$ at the more extreme counts. How do I calculate how "good" my data will fit to the Poisson distribution? Ideally what I would like is to show a plot of ($\lambda$ vs year) with an error bar related to this "goodness".
 A: The direct answer to the question - how to deal with small expected counts - is that one might either 
(a) combine ranges of cells at the end (a very common approach), 
(b) use a different (and perhaps better) goodness of fit test, or 
(c) consider dropping the chi-square approximation, and see if one can deal with the discrete distribution of the test statistic more directly, perhaps by simulation.
Approach (a) can be found in many texts. There are many ways you can go about combining cells, but many people simply work from one end or the other, combining cells either to the left, or to the right, until the expected counts are sufficiently high for their purpose.
(However, the most commonly cited rule of thumb for the expected number - that it should be at least 5 for the chi-square approximation to hold - is unnecessarily strict for the sort of approximation most people would require. Many subsequent papers have suggested less stringent rules.)
The other answer by user36381 suggests that with such a large sample size, goodness of fit tests are almost certain to reject; this is true. However, I'm not so sure comparing to other reference distributions will help, since they, too, would almost certainly be rejected by a decent goodness of fit test.
(Why are you testing whether it's Poisson? If you have around a million data points, the sample itself contains a lot of information about distributional shape - do you actually need a name for the distribution?)
A: This doesn't entirely answer your question regarding how to implement the goodness of fit test. But with a sample size that large, your data will almost certainly fail any test comparing it to a Poisson distribution even if your distribution closely resembles a Poisson distribution (unless, of course, you simulated Poisson data).  The power of your test to detect minute deviations from the Poisson distribution would be very high.
It may be better to compare it to other reference distributions to see which distribution it fits best.
A: One way to go is by using the "deviance residual" from glms.  This is a likelihood ratio for each count.  You can easily get these from the glm () function in R.  The glm you would fit is one with no intercept, and dummy variables for each year.  However with a million observstions this may run into memory problems.
You can calculate the deviance residuals for the poisson model as
$$ d_{ti}^2=2 (\hat{\lambda}_t-n_{ti})+2n_{ti}\log\left (\frac {n_{ti}}{\hat {\lambda}_t}\right) $$
I have use "t" to denote the year ("time"), "i" to denote the observations, $ n_{ti} $ is the "count" variable, and $\hat {\lambda}_t $ is the MLE for year "t".  This is approximately equal to the usual "pearson chi-square" residual when lambda is large, but has more "robust" behaviour for small lambda due to the logarithm.  The standard pearson residual is $r_{ti}=\frac {(\hat{\lambda}_t-n_{ti})^2}{\hat{\lambda}_t} $ which diverges like $\frac {1}{\hat{\lambda}_t}$ - much faster than logarithm.
Now if your Poisson model is adequate you should roughly have $\sum_{ti} d_{ti}^2\approx df $ where "df" is the number of observed counts minus the number of lambdas you fitted.
Another way you can test the poisson model is to fit a negative binomial model, and see how close the "concentration" parameter is to infinity.
A: You could try fitting an exponential distribution to the inter-arrival times (if you have the times of arrivals in your dataset). In this case, you can do a simple KS Test for the exponential fit. If it fits, then the Poisson fits your count data. See the following reference:
ON THE KOLMOGOROV-SMIRNOV TEST FOR THE EXPONENTIAL DISTRIBUTION WITH MEAN UNKNOWN -HUBERT W. LILLIEFORS, AMERICAN STATISTICAL ASSOCIATION JOURNAL, MARCH 1969
