x <- c(2,3,5)
p <- c(.3,.2,.5)

multinomial <- function(x,p) {
	n <- sum(x)
	if (length(p) != length(x))print("x and p are different lengths")
	if (sum(p) != 1) print("p does not sum to 1")
	multi.fn <- prod(p^x)*factorial(n)/prod(factorial(x))
}
p
prob <- multinomial(x,p)
prob

x <- c(2,2,2,2,2)
p <- rep(1/5,5)
prob <- multinomial(x,p)
prob

x <- c(1,3,2,2,2)
prob <- multinomial(x,p)
prob

x <- c(1,3,1,3,2)
prob <- multinomial(x,p)
prob

f <- c(44,43,9,4)
p <- f/sum(f)
x <- c(4,4,1,1)
prob <- multinomial(x,p)
prob

x <- 0:15
pois.cdf <- ppois(lambda=5,x)
plot(pois.cdf,type="b",pch=16,xlab="x",ylab="P(X less than or equal to x)")

n <- 500
p <- .01

x <- 0:20
bin.probs <- dbinom(x,n,p)
pois.probs <- dpois(x,n*p,log = F)
cbind(bin.probs,pois.probs)
plot(bin.probs,pois.probs,pch=16)
abline(a=0,b=1)


x <- 5:20
bin.probs <- dbinom(x,n,p)
pois.probs <- dpois(x,n*p,log = F)
cbind(bin.probs,pois.probs)
plot(bin.probs,pois.probs,pch=16)
abline(a=0,b=1)