Learn R Programming

surveillance (version 1.5-4)

twinstim_profile: Profile Likelihood Computation and Confidence Intervals for twinstim objects

Description

Function to compute estimated and profile likelihood based confidence intervals for twinstim objects. Computations might be cumbersome!

Usage

## S3 method for class 'twinstim':
profile(fitted, profile, alpha = 0.05,
        control = list(fnscale = -1, factr = 10, maxit = 100),
        do.ltildeprofile=FALSE, ...)

Arguments

fitted
an object of class "twinstim".
profile
a list with elements being numeric vectors of length 4. These vectors must have the form c(index, lower, upper, gridsize). [object Object],[object Object],[object Object]
alpha
$(1-\alpha)%$ profile likelihood based confidence intervals are computed. If alpha
control
control object to use in optim for the profile log-likelihood computations. It might be necessary to control maxit or reltol in order to obtain results in finite time.
do.ltildeprofile
If TRUE calculate profile likelihood as well. This might take a while, since an optimisation for all other parameters has to be performed. Useful for likelihood based confidence intervals. Default: FALSE.
...
unused (argument of the generic).

Value

  • list with profile log-likelihood evaluations on the grid and highest likelihood and wald confidence intervals. The argument profile is also returned.

encoding

latin1

Examples

Run this code
# the following call takes a while
#Load the twinstim model fitted to the IMD data
data("imdepifit")

#Generate profiling object for a list of parameters for the new model
names <- c("h.(Intercept)","e.typeC")
coefList <- lapply(names, function(name) {
  c(pmatch(name,names(coef(imdepi.fit))),NA,NA,11)
})

#Profile object (necessary to specify a more loose convergence
#criterion). Speed things up by using do.ltildeprofile=FALSE (the default)
prof <- profile(imdepi.fit, coefList,control=list(fnscale=-1,maxit=50,
   reltol=0.1,REPORT=1,trace=5),do.ltildeprofile=TRUE)

#Plot result for one variable
par(mfrow=c(1,2))
for (name in names) {
  with(as.data.frame(prof$lp[[name]]),matplot(grid,cbind(profile,estimated,wald),
    type="l",xlab=name,ylab="loglik"))
  legend(x="bottomleft",c("profile","estimated","wald"),lty=1:3,col=1:3)
}

Run the code above in your browser using DataLab