Suppose that we have the following independent observations and we know that they come from the same probability density function

k<-c(32,38,42,33,27) #our observations

library('ggplot2')
dat<-data.frame(k,y=0) #we plotted our observations in the x-axis
p<-ggplot(data=dat,aes(x=k,y=y))+geom_point(size=4)
p

We don’t know the exact probability density function, but we know that the shape is the shape of the poisson distribution $f(x|\theta)=f(k|\lambda)=\frac{\lambda^k}{k!}e^{-\lambda}$ Let’s plot two arbitrary probability density functions fixing the parameter of the parametric family

kseq<-seq(20,45,1)
f1<-dpois(kseq,33) #lambda=33
f2<-dpois(kseq,27) #lambda=27
curve1<-data.frame(kseq,f1)
curve2<-data.frame(kseq,f2)

p<-p+geom_segment(data=curve1,aes(x=kseq,xend=kseq,y=0,yend=f1))+
geom_segment(data=curve2,aes(x=kseq+.2,xend=kseq+.2,y=0,yend=f2)) #we introduce a small jitter in the distribution to avoid overplotting
p

It is clear that the probability density function $$f$$ for $$\lambda=33$$ describes the observations better than the other one. We are interested in the $$f$$ that maximizes $$L$$.

The joint density function is $f(k|\lambda)=f(k_1|\lambda)f(k_2|\lambda)...f(k_5|\lambda)=$ $=\frac{\lambda^{k_1}}{k1!}e^{-\lambda}\frac{\lambda^{k_2}}{k2!}e^{-\lambda}...\frac{\lambda^{k_5}}{k5!}e^{-\lambda}=$

that when considered as a function of the parameter is $L(\lambda|k)=L(\lambda|32,38,42,33,27)=f(k|\lambda)=f(32,38,42,33,27|\lambda)=$ $=f(32|\lambda)f(38|\lambda)...f(27|\lambda)=$ $=\frac{\lambda^{32}}{32!}e^{-\lambda}\frac{\lambda^{38}}{38!}e^{-\lambda}...\frac{\lambda^{27}}{27!}e^{-\lambda}$

and $log(L)=log(f(k|\lambda))=log(f(k_1|\lambda))+log(f(k_2|\lambda))+...+log(f(k_5|\lambda))=$ $=log(f(32|\lambda))+log(f(38|\lambda))+...+log(f(27|\lambda))$

logL<-function(p) sum(dpois(k,p,log=T)) # p=lambda

We can calculate the $$log(L)$$ for the two previous examples to verify that $$log(L)$$ is larger for $$\lambda=33$$

logL(33)
## [1] -15.51074
logL(27)
## [1] -20.0261
logL(33)>logL(27)
## [1] TRUE

Let’s plot the $$log(L)$$ for some values of $$\lambda$$

library('plyr')
lambdaSeq<-seq(20,45,.1)
logLike<-laply(lambdaSeq,function(z) logL(z))
dLogL<-data.frame(lambdaSeq,logLike)
ggplot(data=dLogL,aes(x=lambdaSeq,y=logLike))+geom_line()

So it seems that values of $$\lambda$$ around 34 maximizes log(L). Let’s maximize it properly using optimize.

estPar<-optimize(logL,c(10,50),maximum=T) #with optimize we can specify that we search for a max
MLEparameter<-estPar\$maximum
MLEparameter
## [1] 34.4

So, we found that from the parametric family, the probability density function that better characterizes the observations according to MLE is the one described by the parameter $$\lambda=$$ 34.3999994. Let’s plot the distribution in red in the previous graph

fMLE<-dpois(kseq,MLEparameter)
curveMLE<-data.frame(kseq,fMLE)

p<-p+geom_segment(data=curveMLE,aes(x=kseq-.2,xend=kseq-.2,y=0,yend=fMLE),color='red') #we introduce a small jitter in the distribution to avoid overplotting
p

It could be demonstrated that, actually when the parametric family is the poisson distribution, then the MLE of $$\lambda$$ is the mean of the observations

mean(k)
## [1] 34.4