Function minimizers

Generating some random Poisson data:

myData<-rpois(10,1)
myData

##  [1] 2 1 0 1 0 0 1 1 0 1

mean(myData)

## [1] 0.7

Recapturing the lambda paramter of the Poisson data.

First, writing a negative log likelihood function for the Poisson:

nllPois<-function(parVec,data){
  nllik<- -sum(dpois(x=data,lambda=parVec,log=TRUE))
  # cat("nllik= ",nllik,sep=" ",fill=T);cat(" ",sep=" ",fill=T)
  return(nllik)
}

Then minimizing the negative log likelihood (equivalent to maximizing the log likelihood):

library(Rvmmin)
# ?Rvmmin
# Rvmmin(par, fn, gr, lower, upper, bdmsk, control = list(), ...)
parVec<-0.1
outPoisR<-Rvmmin(par=parVec,fn=nllPois,
               lower=c(0.01),upper=c(Inf),
               data=myData)
outPoisR$value

## [1] 10.18987

outPoisR$par

## [1] 0.7

The outPois value is the negative log likelihood associated with the mle. The outPois par is the mle estimate. The mle lambda parameter is equal to the mean of the data.

Minimizing using another minimizer (optim) just to demonstrate that there are other minimizer functions:

parVec<-0.1
outPois<-optim(par=parVec,fn=nllPois,method="L-BFGS-B",
               lower=c(0.01),upper=c(Inf),
               data=myData)
outPois$value

## [1] 10.18987

outPois$par

## [1] 0.7000005

Calculating the associated AIC with the model fit from Rvmmin minimizer. AIC = 2k-2log(likelihood) where k is number of fitted parameters.

2*1-2*outPoisR$value

## [1] -18.37974

Using AIC to decide whether a Poisson or Normal distribution is a better fit to these data.

The likelihood function for the Normal distribution:

nllNorm<-function(parVec,data){
  myMean<-parVec[1]
  mySd<-parVec[2]
  nllik<- -sum(dnorm(x=data,mean=myMean,sd=mySd,log=TRUE))
  return(nllik)
}

Now estimating the mean and sd paramters:

require(Rvmmin)
parVec<-c(0.1,2.0)
outNormR<-Rvmmin(par=parVec,fn=nllNorm,
                 lower=c(-Inf,0.0),upper=c(Inf,Inf),
                 data=myData)
outNormR$value

## [1] 9.731395

outNormR$par

## [1] 0.7000000 0.6403124

Calculating AIC:

2*2-2*outNormR$value

## [1] -15.46279

The lowest AIC is the favored model. Thus the Poisson is a better model for these data as expected.