saemix (version 0.96.1)

conddist.saemix: Estimate conditional mean and variance of individual parameters using the MCMC algorithm

Description

When the parameters of the model have been estimated, we can estimate the individual parameters (psi_i).

Let hat{theta} be the estimated value of theta computed with the SAEM algorithm and let p(phi_i |y_i; hat{theta}) be the conditional distribution of phi_i for 1

Usage

conddist.saemix(saemixObject,nsamp=1,max.iter=NULL,...)

Arguments

saemixObject

an object returned by the saemix function

nsamp

Number of samples to be drawn in the conditional distribution for each subject. Defaults to 1

max.iter

Maximum number of iterations for the computation of the conditional estimates. Defaults to twice the total number of iterations (sum(saemixObject["options"]$nbiter.saemix)*2)

...

optional arguments passed to the plots. Plots will appear if the option displayProgress in the saemixObject object is TRUE

Details

See PDF documentation for details of the computation. Briefly, the MCMC algorithm is used to obtain samples from the individual conditional distributions of the parameters. The algorithm is initialised for each subject to the conditional estimate of the individual parameters obtained at the end of the SAEMIX fit. A convergence criterion is used to ensure convergence of the mean and variance of the conditional distributions. When nsamp>1, several chains of the MCMC algorithm are run in parallel to obtain samples from the conditional distributions, and the convergence criterion must be achieved for all chains. When nsamp>1, the estimate of the conditional mean is obtained by averaging over the different samples.

The shrinkage for any given parameter for the conditional estimate is obtained as Sh=1-var(eta_i)/omega(eta)

where var(eta_i) is the empirical variance of the estimates of the individual random effects, and omega(eta) is the estimated variance.

The function adds or modifies the following elements in the results: [object Object],[object Object],[object Object],[object Object],[object Object],[object Object] A warning is output if the maximum number of iterations is reached without convergence (the maximum number of iterations is saemix.options$nbiter.saemix[2]).

References

Kuhn, E., and Lavielle, M. Maximum likelihood estimation in nonlinear mixed effects models. Computational Statistics and Data Analysis 49, 4 (2005), 1020-1038.

Monolix32_UsersGuide.pdf (http://software.monolix.org/sdoms/software/)

See Also

SaemixData,SaemixModel, SaemixObject,saemixControl,saemix

Examples

data(theo.saemix)

saemix.data<-saemixData(name.data=theo.saemix,header=TRUE,sep="",na=NA,
  name.group=c("Id"),name.predictors=c("Dose","Time"),
  name.response=c("Concentration"),name.covariates=c("Weight","Sex"),
  units=list(x="hr",y="mg/L",covariates=c("kg","-")), name.X="Time")

model1cpt<-function(psi,id,xidep) { 
	  dose<-xidep[,1]
	  tim<-xidep[,2]  
	  ka<-psi[id,1]
	  V<-psi[id,2]
	  CL<-psi[id,3]
	  k<-CL/V
	  ypred<-dose*ka/(V*(ka-k))*(exp(-k*tim)-exp(-ka*tim))
	  return(ypred)
}

saemix.model<-saemixModel(model=model1cpt,
  description="One-compartment model with first-order absorption", 
  psi0=matrix(c(1.,20,0.5,0.1,0,-0.01),ncol=3, byrow=TRUE,dimnames=list(NULL, 
  c("ka","V","CL"))),transform.par=c(1,1,1), 
  covariate.model=matrix(c(0,1,0,0,0,0),ncol=3,byrow=TRUE),fixed.estim=c(1,1,1),
  covariance.model=matrix(c(1,0,0,0,1,0,0,0,1),ncol=3,byrow=TRUE), 
  omega.init=matrix(c(1,0,0,0,1,0,0,0,1),ncol=3,byrow=TRUE), error.model="constant")

saemix.options<-list(seed=632545,save=FALSE,save.graphs=FALSE)

saemix.fit<-saemix(saemix.model,saemix.data,saemix.options)

saemix.fit<-conddist.saemix(saemix.fit,nsamp=3)
# First sample from the conditional distribution 
# (a N (nb of subject) by nb.etas (nb of parameters) matrix)
saemix.fit["results"]["phi.samp"][,,1]

# Second sample
saemix.fit["results"]["phi.samp"][,,2]