\documentclass[notitlepage]{article}
\usepackage{../Math323}
\usepackage{listings}
\usepackage{numprint}
\usepackage{enumerate}
\usepackage{bbm}

\usepackage{amsfonts,amsmath}

\setlength{\parindent}{0pt}

\def\E{\Expect}


\lstloadlanguages{R}

\definecolor{keywordcolor}{rgb}{0,0.6,0.6}
\definecolor{delimcolor}{rgb}{0.461,0.039,0.102}
\definecolor{Rcommentcolor}{rgb}{0.101,0.043,0.432}

\lstdefinestyle{Rsettings}{
  basicstyle=\ttfamily,
  breaklines=true,
  showstringspaces=false,
  keywords={if, else, function, theFunction, tmp}, % Write as many keywords
  otherkeywords={},
  commentstyle=\itshape\color{Rcommentcolor},
  keywordstyle=\color{keywordcolor},
  moredelim=[s][\color{delimcolor}]{"}{"},
}

\lstset{basicstyle=\ttfamily, numbers=none, literate={~} {$\sim$}{2}}

\begin{document}


<<setup, cache=FALSE, echo=FALSE>>=
library(knitr)
# global chunk options
opts_chunk$set(cache=TRUE, autodep=TRUE)
options(scipen=999)
options(repos=c(CRAN="https://cloud.r-project.org/"))
@
\begin{center}
{\textsclarge{The Poisson distribution}}
\end{center}

<<pois00,comment='+', fig.width=10, fig.height=6.666, fig.align='center',echo=FALSE>>=
lam<-2
@
The Poisson distribution arises as the limiting case of the Binomial distribution whose parameter $p$ varies with $n$.  Specifically, we consider a fixed constant $\lambda$, and consider the $Binomial(n,p)$ distribution, and observe what happens when $n$ increases with $p=\lambda/n$ decreasing so that $n p = \lambda$ stays fixed.    We have that as $n \longrightarrow \infty$
\[
\dbinom{n}{y} \left( \frac{\lambda}{n} \right)^y \left(1- \frac{\lambda}{n} \right)^{n-y} \longrightarrow e^{-\lambda} \frac{\lambda^y}{y!}
\]
In the following plots, we generate binary sequences repeated a number of times (repeats are in the rows of the plots), with 1s or 'successes' marked as the red squares, and 0s or `failures' marked as white squares.  In the simulations, $\lambda=\Sexpr{lam}$ is fixed, and $n=30,60,\ldots,300$.

<<pois01,comment='+', fig.width=10, fig.height=6.666, fig.align='center'>>=
lam<-2
 par(mar=c(5.1, 4.1, 4.1, 8.1), xpd=TRUE)
for(n in 1:10){
	xn<-30*n
	yn<-20*n
	p<-lam/xn
	xmat<-matrix(rbinom(xn*yn,1,p),nrow=xn,ncol=yn)
	image((1:xn/xn)-0.5/xn,(1:yn)/yn-0.5/yn,xmat,xlab=substitute(paste(nval,' repeats'),
            list(nval = yn)),ylab='',axes=FALSE,col=rev(rgb(c(1,1),c(0,1),c(0,1))))
	for(i in 0:yn){lines(c(0,1),c(i,i)/yn)}
	for(i in 0:xn){lines(c(i,i)/xn,c(0,1))}
	title(substitute(paste(n==nval,' ',p==pval), list(nval = xn, pval=p)))
	legend(1.1,1.0,c('1','0'),fill=rgb(c(1,1),c(0,1),c(0,1)),cex=1.2,bty='n')
}
@
We remove the gridlines to make the plot easier to see.

<<pois02,comment='+', fig.width=10, fig.height=6.666, fig.align='center',echo=FALSE>>=
par(mar=c(5.1, 4.1, 4.1, 8.1), xpd=TRUE)
image((1:xn/xn)-0.5/xn,(1:yn)/yn-0.5/yn,xmat,xlab=substitute(paste(nval,' repeats'), list(nval = yn)),
        ylab='',axes=FALSE,col=rev(rgb(c(1,1),c(0,1),c(0,1))))
for(i in c(0,yn)){lines(c(0,1),c(i,i)/yn)}
for(i in c(0,xn)){lines(c(i,i)/xn,c(0,1))}

title(substitute(paste(n==nval,' ',p==pval), list(nval = xn, pval=p)))
legend(1.1,1.0,c('1','0'),fill=rgb(c(1,1),c(0,1),c(0,1)),cex=1.2,bty='n')
@

We can investigate how well the Poisson distribution approximates the binomial distribution: in each case, $p=\lambda/n$, with $\lambda = \Sexpr{lam}$. Entries in the table are
\[
p(y) = P(Y=y) \qquad y=0,1,2,\ldots,9
\]
for $n$ up to 600.

<<pois03,comment='+', fig.width=10, fig.height=6.666, fig.align='center',echo=FALSE>>=
x<-c(0:9)
bvec<-matrix(0,nrow=length(x),ncol=20)
for(n in 1:20){
    p<-lam/(30*n)
    bvec[,n]<-dbinom(x,30*n,p)
}
bvec<-data.frame(cbind(bvec,dpois(x,lam)))
rownames(bvec)<-paste('y=',as.character(x),sep='')
colnames(bvec)<-c(paste('n=',seq(1,20)*30,sep=''),'Poisson')
round(t(bvec),4)
@
The probabilities match to the third decimal place.

\end{document}