optim Function

optim is the first function in concern. Documentation of the optim function is useful and it indicates that this function is only used on one input situations only. This means our EstMLEBetaBin function has to be modified. Reason for this is only the parameters that should be estimated need to be input values but our EstMLEBetaBin function has four parameters which are a,b,x(Binomial Random Variable) and freq(corresponding frequency values).

While using optim function first index refers to shape parameter a and second index refers to shape parameter b. Further, we have to input the observations or in our case the Binomial random variable values and their respective frequencies. I think it is inconvenient to modify the EstMLEBetaBin function, because if we want to estimate parameters for different data-sets it would become tedious. After modification we have a new function foroptim which can be used for demonstration and comparison.

Below is the code to estimation and going through the outputs. It should be noted that we have to provide initial parameter values as an input to the optim function, and it is best to provide values in the domain of shape parameter values which we want to estimate.

Here the shape parameters a and b are in the region of greater than zero but less than positive infinity (\(+\infty >a,b>0\)). So for the initial parameters of a=0.1 and b=0.2 we are finding parameters which would minimize the Negative Log Likelihood function of Beta-Binomial distribution with the Alcohol Consumption data.

# new function to facilitate optim criteria
# only one input but has two elements
foroptim<-function(a)
  {
  EstMLEBetaBin(x=Alcohol_data$Days, freq=Alcohol_data$week1,a=a[1],b=a[2])
  }

# optimizing values for a,b using default mathematical method
optim_answer<-optim(par=c(0.1,0.2),fn=foroptim)

# obtaining class of output
class(optim_answer)

#length of output
length(optim_answer)

# the outputs
optim_answer$par # estimated values for a, b
optim_answer$value # minimized function value 
optim_answer$counts  # see the documentation to understand
optim_answer$convergence # indicates successful completion
optim_answer$message # additional information

# fitting the Beta-Binomial distribution with estimated shape parameter  values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           optim_answer$par[1],optim_answer$par[2])

So the foroptim function can be used as above and parameters are estimated for \(\alpha\) and \(\beta\) (or a, b) for the Alcohol Consumption data week 1. Further we have scrutinized the function as below.

• package : stats
• No of Inputs: 7
• Minimum required Inputs : 2
• Class of output : list
• No of outputs: 5
• No of Analytical Methods : 6
• Default Method : Nelder-Mead

nlm Function

nlm function is also similar to optim but only one analytical method will be used, which is a Newton-type Algorithm. Here also there needs to be changes made to our EstMLEBetaBin function as previously. After making those changes we have called the new Negative Log-likelihood function of Beta-Binomial distribution as fornlm. Then we can use the nlm function and estimate a and b for the initial shape parameter values of 0.1 and 0.2 respectively. Documentation of nlm function is very useful so that we can understand how it works.

Below is the code for using nlm function appropriately and fiddling with the results.

#new function to facilitate nlm criteria
fornlm<-function(a)
  {
  EstMLEBetaBin(x=Alcohol_data$Days, freq=Alcohol_data$week1,a=a[1],a[2])
  }

#optimizing values for a,b using the only analytical method
nlm_answer<-nlm(f=fornlm,p=c(0.1,0.2))

#obtaining class of output
class(nlm_answer)

#length of output
length(nlm_answer)

# the outputs
nlm_answer$estimate # estimated values for a, b
nlm_answer$minimum # minimized function value 
nlm_answer$gradient  # gradient at the estimated minimum of given funciton
nlm_answer$code # indicates successful completion
nlm_answer$iterations # number of iterations performed

# fitting the Beta-Binomial distribution with estimated shape parameter  values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           nlm_answer$estimate[1],nlm_answer$estimate[2])

Similarly for the fornlm function we have estimated values for a and b which would fit the Alcohol consumption data of week 1. Below is a point form summary of nlm function.

• package : stats
• No of Inputs: 12
• Minimum required Inputs : 2
• Class of output : list
• No of outputs: 5
• No of Analytical Methods : 1
• Default Method : Newton-type Algorithm

nlminb Function

nlminb is also similar and requires the EstMLEBetaBin function to be restructured as similar to previous situations. After this task we now have a function called fornlminb. nlminb function is based on analytical method of unconstrained and box-constrained optimization using PORT routines.

After choosing initial parameter values for a and b, which are respectively 0.1 and 0.2 the estimation was done following the process below

# new function to facilitate nlminb criteria
fornlminb<-function(a)
  {
  EstMLEBetaBin(x=Alcohol_data$Days, freq=Alcohol_data$week1,a=a[1],a[2])
}

# optimizing values for a,b using default analytical method
nlminb_answer<-nlminb(start=c(0.1,0.2), objective=fornlminb)

# obtaining class of output
class(nlminb_answer)

# length of output
length(nlminb_answer)

# the outputs
nlminb_answer$par # estimated values for a, b
nlminb_answer$objective # minimized function value 
nlminb_answer$evaluations  # see the documentation to understand
nlminb_answer$convergence # indicates successful completion
nlminb_answer$message # additional information
nlminb_answer$iterations # number of iterations performed

# fitting the Beta-Binomial distribution with estimated shape parameter values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           nlminb_answer$par[1],nlminb_answer$par[2])

Documentation includes the information related to the function therefore referring it will be useful. After estimating values for a and b using nlminb function these were noticed regarding the function in concern

• package : stats
• No of Inputs: 8
• Minimum required Inputs : 2
• Class of output : list
• No of outputs: 6
• No of Analytical Methods : 1
• Default Method : Unconstrained and box-constrained optimization using PORT routines

ucminf Function

Package ucminf produces the ucminf function, as previously mentioned functions here also we have to change the EstMLEBetaBin function. After making the changes we will be using the forucminf function to the estimation process of shape parameters a and b.

When initial parameter values are set to a=0.1 and b=0.2 we will obtain results from ucminf function, which will minimize the Negative Log-likelihood value of Beta-Binomial distribution. Below is the code for estimation and using the results to understand the function ucminf.

library(ucminf)

# new function to facilitate ucminf criteria
forucminf<-function(a)
  {
  EstMLEBetaBin(x=Alcohol_data$Days,freq = Alcohol_data$week1,a=a[1],a[2])
}

# optimizing values for a,b using default analytical method
ucminf_answer<-ucminf(par=c(0.1,0.2),fn=forucminf)

#obtaining class
class(ucminf_answer)

# length of output
length(ucminf_answer)

# the outputs
ucminf_answer$par # estimated values for a, b
ucminf_answer$value # minimized function value 
ucminf_answer$invhessian.lt  # see the documentation understand
ucminf_answer$convergence # indicates successful completion
ucminf_answer$message # additional information
ucminf_answer$info

# fitting the Beta-Binomial distribution with estimated shape parameter values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           ucminf_answer$par[1],ucminf_answer$par[2])

With the help of R Documentation the estimation was done for shape parameter values a and b using the ucminf function. Below is the initial understanding of the ucminf function

• package : ucminf
• No of Inputs: 5
• Minimum required Inputs : 2
• Class of output : list
• No of outputs: 6
• No of Analytical Methods : 1
• Default Method : Quasi-Newton Algorithm type with BFGS updating

maxLik Function

maxLik function is from the maxLik package, which only maximizes the Log Likelihood function. Therefore we have to restructure EstMLEBetaBin as previously mentioned, but as an addition a negative sign is added for the output. This new function will be called as formaxLik.

For the initial parameter values where a=0.1 and b=0.2 the maxLik function will be used and results will be evaluated as below.

library(maxLik)

# new function to facilitate maxLik criteria
formaxLik<-function(a)
  {
  -EstMLEBetaBin(x=Alcohol_data$Days,freq = Alcohol_data$week1,a=a[1],a[2])
  }

# optimizing values for a,b using default analytical method
maxLik_answer<-maxLik(logLik =formaxLik, start = c(0.1,0.2))

# obtaining class of output
class(maxLik_answer)

# length of output
length(maxLik_answer)

# the outputs
maxLik_answer$estimate # estimated values for a, b
maxLik_answer$maximum # minimized function value 
maxLik_answer$iterations  # no of iterations to succeed
maxLik_answer$gradient # last gradient value which was calculated
maxLik_answer$message # additional information
maxLik_answer$hessian # hessian matrix
maxLik_answer$code # indicates successful completion
maxLik_answer$fixed # logical vector indicating which parameters are constants
maxLik_answer$type # type of maximization
maxLik_answer$last.step # list describing the last unsuccessful step
maxLik_answer$control # see the documentation understand

# fitting the Beta-Binomial distribution with estimated shape parameter values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           maxLik_answer$estimate[1],maxLik_answer$estimate[2])

Using the Documentation of maxLik and the Documentation of maxNR function the above analysis was done for the maxLik function. According to the above code, below are the findings from the maxLik function in point form

• package : maxLik
• No of Inputs: 6
• Minimum required Inputs : 2
• Class of output : list or class of maxim or class of maxLik
• No of outputs: 11
• No of Analytical Methods : 7
• Default Method : Automatically chosen

mle Function

stats4 package is installed so that mle function can be operated for the purpose of estimating a and b parameters. Here also as previously we need to make some changes as below and create a new Negative Log Likelihood function called formle.

For the initial parameter values of a=0.1 and b=0.2 the Negative Log Likelihood value of Beta-Binomial distribution has been minimized using mle function. Below is the code for estimation and investigation from the outputs of mle after estimation.

library(stats4)

# new function to facilitate mle criteria
formle<-function(a,b)
  {
  EstMLEBetaBin(x=Alcohol_data$Days,freq = Alcohol_data$week1,a,b)
  }

# optimizing values for a,b using default analytical method
mle_answer<-mle(minuslogl=formle,start = list(a=0.1,b=0.2))

# obtaining class
class(mle_answer)

# length of output
length(mle_answer)

# the outputs
mle_answer@call # inputs i have used 
mle_answer@coef # estimated values for a,b
mle_answer@fullcoef # all values, even the fixed values we did not want to estimate
mle_answer@vcov # variance covariance matrix for a,b
mle_answer@min # minimized function value
mle_answer@details # details after estimation process
mle_answer@nobs # number of observations to be used for computing only if given 
mle_answer@method # optimization methods used

# Methods used
confint(mle_answer) # confidence intervals for estimated values
logLik(mle_answer) # Negative loglikelihood value for estimated values 
profile(mle_answer) # Likelihood profile generation.
nobs(mle_answer) # number of observations to be used for computing only if given
show(mle_answer) # display object briefly
summary(mle_answer) # generate a summary
#update()   # updating if we have new data and need to estimate new values
vcov(mle_answer) # variance covariance matrix

# fitting the Beta-Binomial distribution with estimated shape parameter values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           mle_answer@coef[1],mle_answer@coef[2])

It should be noted that in all 5 previous functions we had only outputs in the form of lists, but through mle we are seeing a class of mle output.The Documention explains the inputs and outputs of mle. Documentation for mle-class explains further about the methods that can be used. Below is a list of findings based on the outputs of the mle function.

• package : stats4
• No of Inputs: 5
• Minimum required Inputs : 2
• Class of output : class of mle
• No of outputs: 9
• Methods for output: 8
• No of Analytical Methods :6
• Default Method : Nelder and Mead

mle2 Function

mle2 function is advanced than mle function and it is from the package bbmle. Here, there is no need to modify the EstMLEBetabin function from the fitODBOD package to estimation as all previous situations .

For the initial parameter values of a=0.1 and b=0.2 Negative Log Likelihood function of Beta-Binomial distribution will be minimized where outputs will be investigated and methods related to output of class mle2 will be used.

Below is the code for using mle2 function and scrutinizing the output and methods related to it.

library(bbmle)

# optimizing values for a,b using default analytical method
mle2_answer<-mle2(minuslogl= EstMLEBetaBin,start = list(a=0.1,b=0.2),
                  data = list(x=Alcohol_data$Days,freq=Alcohol_data$week1))

# obtaining class
class(mle2_answer)

# length of output
length(mle2_answer)

# the outputs
mle2_answer@call # inputs generally considered 
mle2_answer@call.orig # inputs i have given
mle2_answer@coef # estimated values for a,b
mle2_answer@fullcoef # all values, even the fixed values we did not want to estimate
mle2_answer@vcov # variance covariance matrix for a,b
mle2_answer@min # minimized function value
mle2_answer@details # details after estimation process
mle2_answer@method # optimization methods used
mle2_answer@data # data used for estimation 
mle2_answer@formula # if a formula was specified in the input 
mle2_answer@optimizer # function used for optimizing

# Methods used
coef(mle2_answer) # extrat the estimated values
confint(mle2_answer) # confidence intervals for estimated values
show(mle2_answer) # display object briefly
summary(mle2_answer) # generate a summary
#update()   #updating if we have new data and need to estimate new values
vcov(mle2_answer) # variance covariance matrix
#formula(mle2_answer) # if a formula was specified in the input 
#plot(mle2_answer) # plot the profile
logLik(mle2_answer) # Negative loglikelihood value for estimated values 
profile(mle2_answer) # profile of estimated values

# fitting the Beta-Binomial distribution with estimated shape parameter values
fitBetaBin(Alcohol_data$Days,Alcohol_data$week1,
           mle2_answer@coef[1],mle2_answer@coef[2])

Documention of mle2 function and Documentation of mle2-class
provide a few findings as mentioned below in point form.

• package : bbmle
• No of Inputs: 22
• Minimum required Inputs : 3
• Class of output : class of mle2
• No of outputs: 12
• Methods for output: 9
• No of Analytical Methods : 6
• Default Method : Nelder and Mead