Package 'clinDR'

Description

The functions fitEmax and fitEmaxB fit an Emax model to binary or continuous data using maximum likelihood or Bayesian estimation. They have several generic supporting functions. Functions to produce plots associated with dose response analyses are (plotD, plotB, plot.fitEmax,plot.fitEmaxB). The functions emaxsim and emaxsimB perform simulations of 4- and 3-parameter Emax ML or Bayesian estimation. The ML estimates are replaced with alternative model fits when the primary estimation fails. Several supporting functions are supplied to analyze the output of emaxsim and emaxsimB, including analyses for specific simulated data sets. All of the data sets from dose response meta analyses are included in metaData.

Details

The function compileStanModels must be executed once after the package is installed to create compiled STAN Emax models before the Bayes functions in the package can be executed. This requires 3-10 minutes to complete on most machines. The compiled code is 32-bit or 64-bit specific, and both must be created if both versions of R are used.

The Bayesian computations use the R package rstan. It can be installed from CRAN. Windows users should check the instructions for rstan at the https://mc-stan.org and https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started. Note that Rtools must be installed, which is a simple, but often overlooked step. Instructions for its installation are given in the second URL.

Author(s)

Neal Thomas [aut, cre], Jing Wu[aut]

See Also

DoseFinding

Description

Extract a simulated data set from the output of emaxsim. Data are re-created using the stored random number seed.

Usage

## S3 method for class 'emaxsim'
x[i, ...]

Arguments

Details

Re-creates the ith simulated data set for subsequent analyses. Also returns all analyses done for the ith data set in emaxsim

Value

A list is returned with class(emaxsimobj) containing:

Note

Extraction from a simulation object requires re-creation of the simulated data set. If the extracted object is to be used more than once, it is more efficient to save the extracted object than reuse [].

Author(s)

Neal Thomas

See Also

emaxsim, print.emaxsimobj, plot.emaxsimobj, update.emaxsimobj

Examples

## Not run: 
## code change random number seed

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen.parm<-FixedMean(n,doselev,meanlev,sdy)  

D1 <- emaxsim(nsim,gen.parm,modType=3)
e49<-D1[49]                  #### extract 49th simulation


## End(Not run)

Description

Extract a simulated data set from the output of emaxsimB. Data are re-created using the stored random number seed.

Usage

## S3 method for class 'emaxsimB'
x[i, ...]

Arguments

Details

Re-creates the ith simulated data set for subsequent analyses. Also returns all analyses done for the ith data set in emaxsimB

Value

A list is returned with class(emaxsimBobj) containing:

Note

Extraction from a simulation object requires re-creation of the simulated data set. If the extracted object is to be used more than once, it is more efficient to save the extracted object than reuse [].

Author(s)

Neal Thomas

See Also

emaxsimB, print.emaxsimBobj, plot.emaxsimBobj

Examples

## Not run: 

save.seed<-.Random.seed
set.seed(12357)

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

prior<-emaxPrior.control(epmu=0,epsca=30,difTargetmu=0,
		difTargetsca=30,dTarget=100,p50=50,sigmalow=0.1,
		sigmaup=30,parmDF=5)
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,propInit=0.15,adapt_delta = 0.95)

D1 <- emaxsimB(nsim,gen, prior, modType=3,mcmc=mcmc,check=FALSE)

out<-D1[2]


.Random.seed<-save.seed

## End(Not run)

Description

Bayes posterior predictive test for an Emax (monotone) model fit comparing the best response from lower doses to the response from the highest dose. checkMonoEmax is deprecated. See bpchkMonoEmax.

Usage

bpchkMonoEmax(x, trend='positive', protSel=1)

Arguments

Details

The Bayesian predictive p-value is the posterior probability that a dose group sample mean in a new study with the same sample sizes would yield a higher (or lower for negative trend) difference for one of the lower doses versus the highest dose than was actually obtained from the real sample. There must be at least two non-placebo dose groups (NA returned otherwise). Placebo response is excluded from the comparisons.

The function generates random numbers, so the random number generator/seed must be set before the function is called for exact reproducibility.

When fitEmaxB is applied to first-stage fitted model output with a non-diagonal variance-covariance matrix, the predictive draws are selected from a multivariate model with means computed from the MCMC-generated parameters and input asymptotic variance-covariance matrix vcest. If the fitted model was applied to binary data, the GOF statistic is computed based on the logit rather than observed dose group sample proportion scale. This differs from the setting with patient-level data input to fitEmaxB.

Value

Returns a scalar Bayesian predictive p-value.

Author(s)

Neal Thomas

References

Thomas, N., Sweeney, K., and Somayaji, V. (2014). Meta-analysis of clinical dose response in a large drug development portfolio, Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317. <doi:10.1080/19466315.2014.924876>

Thomas, N., and Roy, D. (2016). Analysis of clinical dose-response in small-molecule drug development: 2009-2014. Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317 <doi:10.1080/19466315.2016.1256229>

Wu, J., Banerjee, A., Jin, B. Menon, M. S., Martin, S. and Heatherington, A. (2017). Clinical dose response for a broad set of biological products: A model-based meta-analysis. Statistical Methods in Medical Research. <doi:10.1177/0962280216684528>

See Also

plot.plotB, plotD, plot.fitEmax

Examples

## Not run: 

data("metaData")
exdat<-metaData[metaData$taid==6 & metaData$poptype==1,]

prior<-emaxPrior.control(epmu=0,epsca=10,difTargetmu=0,difTargetsca=10,dTarget=80.0,
        p50=3.75,sigmalow=0.01,sigmaup=20)
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
parms<-coef(fitout)[,1:4]  #use first intercept

checkMonoEmax(fitout, trend='negative')
      

## End(Not run)

Description

Bayes posterior predictive test for an Emax (monotone) model fit comparing the best response from lower doses to the response from the highest dose. checkMonoEmax is deprecated. See bpchkMonoEmax.

Usage

checkMonoEmax(y, 
      dose, 
      parm, 
      sigma2,
      nvec=rep(1,length(dose)),
      xbase=NULL,
      modelFun=emaxfun,
      trend='positive',
      binary= FALSE,logit=binary)

Arguments

Details

A sample of parameters from the joint posterior distribution must be supplied (typically produced by an MCMC program). The Bayesian predictive p-value is the posterior probability that a dose group sample mean in a new study with the same sample sizes would yield a higher (or lower for negative trend) difference for one of the lower doses versus the highest dose than was actually obtained from the real sample. There must be at least two non-placebo dose groups (NA returned otherwise). Placebo response is excluded from the comparisons.

The function generates random numbers, so the random number generator/seed must be set before the function is called for exact reproducibility.

Value

Returns a scalar Bayesian predictive p-value.

Author(s)

Neal Thomas

See Also

plot.plotB, plotD, plot.fitEmax

Examples

## Not run: 

data("metaData")
exdat<-metaData[metaData$taid==6 & metaData$poptype==1,]

prior<-emaxPrior.control(epmu=0,epsca=10,difTargetmu=0,difTargetsca=10,dTarget=80.0,
        p50=3.75,sigmalow=0.01,sigmaup=20)
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
parms<-coef(fitout)[,1:4]  #use first intercept

checkMonoEmax(y=exdat$rslt, dose=exdat$dose, parm=parms, sigma2=(sigma(fitout))^2,
      nvec=exdat$sampsize, trend='negative')
      

## End(Not run)

Description

Extract Emax model parameter estimates. MLE for fitEmax. Matrix of MCMC generated parameters for fitEmaxB.

Usage

## S3 method for class 'fitEmax'
coef(object, ...)
	## S3 method for class 'fitEmaxB'
coef(object, local=FALSE, ...)
	## S3 method for class 'emaxsim'
coef(object, ...)
	## S3 method for class 'emaxsimB'
coef(object, local=FALSE, ...)

Arguments

Value

Vector of MLE estimates of model parameter from fitEmax. Matrix of MCMC generated parameters for fitEmaxB. Matrix with posterior median parameter estimates for each emaxsimB simulation: (led50,lambda,emax,e0) or (led50,emax,e0). For emaxsim, a list is returned with the model type fit for each simulation, and a matrix with the corresponding model coefficients. The order of the parameters is given in the emaxsim documentation.

Author(s)

Neal Thomas

See Also

sigma, fitEmax, fitEmaxB, emaxsim, emaxsimB

Examples

doselev<-c(0,5,25,50,100,350)
n<-c(78,81,81,81,77,80)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-8.0
pop<-c(log(ed50),emax,e0)    
dose<-rep(doselev,n)
meanlev<-emaxfun(dose,pop)  

y<-rnorm(sum(n),meanlev,sdy)

testout<-fitEmax(y,dose,modType=4)
coef(testout)

Description

Compile rstan code for Emax models used by fitEmaxB and emaxsimB. This function must be executed once after the clinDR package is installed.

Usage

compileStanModels()

Details

The compiled models are stored in the models sub-directory of the installed clinDR package. The user must have write-access to the package directory. The package can be installed in a user-specified directory if the user does not have write privileges for the default package directory. Execution requires several minutes. The compiled models are 32- or 64- bit specific. Both sets must be compiled if the compiled R type is changed (they are stored in sub-directories comp32 or comp64). It is recommended to execute the function again if the package rstan is updated.

Package rstan must be functional for CompileStanModels to be successful. See https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started. Note especially the instructions for installing Rtools, which is required for execution on a Windows machine.

Value

'basemodel.rds' and 'mrmodel.rds' should be created in the package directory in the sub-directory 'models'.

Author(s)

Neal Thomas

Description

Density plot for distributions conditional on a variable. A grid of values are specified for the conditioning variable, which is plotted on the horizontal axis. The conditioning variable is typically dose or time

Usage

DRDensityPlot(x,qL,qH,qlevL=c(0.025,0.05,0.10,0.25),
		xlim,ylim,xlab='x',ylab='y')

Arguments

Details

The function takes as input percentiles defining confidence intervals or Bayesian probability intervals at different levels (e.g. 5,95, 25,75) for distributions conditional on a variable that is typically dose or time. Regions defined by different confidence/probability levels are represented by different levels of shading. The input parameter, qlevL, is used only to define the input in the matrices qL and qH. The qlevL is not used for any numerical calculations, which must be done before executing the function.

Value

Plotted output only.

Author(s)

Neal Thomas

See Also

plotBdensity

Examples

## Not run: 
data('metaData')
exdat<-metaData[metaData$taid==32,]

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmax(exdat$rslt,exdat$dose,modType=3,count=exdat$sampsize,
								msSat=msSat)

dgrid<-seq(0,100,length=100)
seout95<-predict(fitout,dgrid,clev=0.95)
seout90<-predict(fitout,dgrid,clev=0.9)
seout80<-predict(fitout,dgrid,clev=0.8)
seout50<-predict(fitout,dgrid,clev=0.5)

qlev<-c(0.025,0.05,0.10,0.25)

qL<-cbind(seout95$ubdif,seout90$ubdif,seout80$ubdif,seout50$ubdif)
qH<-cbind(seout95$lbdif,seout90$lbdif,seout80$lbdif,seout50$lbdif)

DRDensityPlot(dgrid,qL,qH,qlevL=qlev,xlab='Dose',ylab='Diff with PBO')


## End(Not run)

Description

ML estimation for 4- and 3-parmeter Emax model. If the 4-parameter model is requested, it is estimated and the 3-parameter model is fit only if the 4-parameter estimation fails. If 3-parameter estimation fails, the linear, log-linear, or exponential model producing the smallest residual SS is substituted. For binary data, the model is fit on the logit scale and then back-transformed.

Usage

emaxalt(y, dose, modType=3,binary=FALSE, 
	iparm=NA,ed50cutoff=2.5*max(doselev), 
	ed50lowcutoff=doselev[2]/1000,switchMod= TRUE,
	truncLambda=6)

Arguments

Details

The partial linear method is used in nls. If it fails, gauss-newton is attempted. If both methods fail, the next simpler model is attempted. For the 4-parameter model, the next step is the 3-parameter model. For the 3-parameter model, a linear, log-linear log(dose+1.0), and exp(dose/max(dose)) are fit using lm, and the 2-parm fit with the smallest residual SS is selected.

Value

A list assigned class "emaxalt" with the following elements:

Author(s)

Neal Thomas

See Also

emaxsim, nls

Examples

save.seed<-.Random.seed
set.seed(12357)

doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
dose<-rep(doselev,n)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanresp<-emaxfun(dose,pop)  
y<-rnorm(sum(n),meanresp,sdy)

simout<-emaxalt(y,dose)

simout2<-emaxalt(y,dose,modType=4)

.Random.seed<-save.seed

Description

Evaluate Emax models for a vector of dose levels for multiple sets of parameters.

Usage

emaxfun(dose, parm)

Arguments

Details

The Hill parameter is omitted from parm for the hyperbolic model

Value

Returns a matrix of Emax function evaluations. The rows correspond to the parameter replications, and the columns correspond to the dose levels.

Note

The ordering of the parameters was selected to faciliate use of the 'plinear' algorithm in function nls.

Author(s)

Neal Thomas

See Also

dlogis

Examples

doselev<-c(0,5,25,50,100)
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
lambda=2
emax<-solveEmax(diftarget,dtarget,log(ed50),lambda,e0)

parm<-c(log(ed50),lambda,emax,e0)
plot(doselev,emaxfun(doselev,parm))

Description

Set the parameters of the prior distribution for the Emax model implemented in fitEmaxB..

Usage

emaxPrior.control(epmu=NULL,epsca=NULL,
	difTargetmu=NULL,difTargetsca=NULL,
	dTarget=NULL,p50=NULL,
	sigmalow=NULL,sigmaup=NULL,
	effDF=parmDF,parmDF=5,
	loged50mu=0.0,loged50sca=1.73,
	loglammu=0.0,loglamsca=0.425,parmCor=-0.45,
	lowled50=log(0.001),highled50=log(1000),
	lowllam=log(0.3),highllam=log(4.0),
	basemu=NULL,basevar=NULL,binary=FALSE)

Arguments

Details

The prior distribution is based on meta-analyses of dose response described in the references. The E0 and difTarget parameters have independent t-distribution prior distributions. For binary data, these parameters are computed on the logistic scale. The prior means and scales of these parameters must be assigned compound-specific values. The predicted ED50 at the study design stage must must also be specified as 'P50'. For continuous data, the prior distribution for the residual SD is uniform on a user-specifed scale.

The prior distribution of the log(ED50) has a t-distribution centered at log(P50), with scale, degrees of freedom (parmDF), and offset to the P50, defaulting to values given in the references (these can be changed, but they are difficult to interpret outside the context of the meta-analyses). If modType=4, the prior distribution for the Hill parameter is also t-distribution with parmDF degrees of freedom and corParm correlation with the log(ED50).

Value

List of class emaxPrior of prior parameter values for use in fitEmaxB. default is a derived variable set to TRUE when the default values are used for loged50 and loglambda.

Author(s)

Neal Thomas

References

Thomas, N., Sweeney, K., and Somayaji, V. (2014). Meta-analysis of clinical dose response in a large drug development portfolio, Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317. <doi:10.1080/19466315.2014.924876>

Thomas, N., and Roy, D. (2016). Analysis of clinical dose-response in small-molecule drug development: 2009-2014. Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317 <doi:10.1080/19466315.2016.1256229>

Wu, J., Banerjee, A., Jin, B., Menon, S., Martin, S., and Heatherington, A. (2017). Clinical dose-response for a broad set of biological products: A model-based meta-analysis. Vol. 9, 2694-2721. <doi:10.1177/0962280216684528?>

See Also

fitEmaxB

Description

Simulate dose response data and apply 4- or 3- parameter Emax MLE estimation. For binary data, the model is fit on the logit scale and then back-transformed. When MLE estimation fails, models with fewer parameters (including models linear in their parameters) are substituted. Summaries of estimation performance are returned for further analyses. An MCP-MOD test is also performed for each simulated data set.

Usage

emaxsim(
	nsim, 
	genObj, 
	modType=3,
	binary=FALSE,
	seed=12357,
	nproc = parallel::detectCores(),
	negEmax=FALSE,
	ed50contr=NULL,
	lambdacontr=NULL,
	testMods=NULL,
	idmax=length(doselev), 
	iparm=NA,
	ed50cutoff=2.5*max(doselev), 
	ed50lowcutoff=doselev[2]/1000,
	switchMod= TRUE,
	truncLambda=6,
	description="")

Arguments

Details

Continuous data can be simulated from any dose response curve with homogeneous normally distributed residuals. The estimation procedure starts with ML estimation of a 4- or 3- parameter Emax model depending on modType. If modType=3 or 4-parameter estimation fails, a 3 parameter Emax model is fit by maximum likelihood non-linear least squares. If 1) nls fails to converge for a 3 parameter Emax model, 2) the ED50 estimate is <=0, or 3) the ED50 estimate exceeds ed50cutoff, a linear, log-linear (offset of 1.0), or scaled exponental (exp(dose/max(dose))), is fit using simple linear least squares estimation. The model selected has the smallest residual SS.

Binary data are handled similarly using maximum likelihood implemented with the nlm function. The models are fit on the logit scale and then back-transformed for estimation of dose response. Reduced linear models are selected based on the corresponding likelihood deviance.

MCP-MOD tests are created from contrasts based on the Emax function using the DoseFinding package. Different ED50 and lambda (Hill) parameters can be specified to form the contrasts. A contrast matrix output from the DoseFinding package can be specified instead, allowing for other contrast choices.

Value

A list is returned with class(emaxsim) containing:

Author(s)

Neal Thomas

See Also

print.emaxsim, summary.emaxsim, plot.emaxsim, coef.emaxsim, sigma.emaxsim, vcov.emaxsim, predict.emaxsim, emaxfun

Examples

## Not run: 
## emaxsim changes the random number seed
nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

D1 <- emaxsim(nsim,gen,modType=3)
summary(D1,testalph=0.05)

D4 <- emaxsim(nsim,gen,modType=4)
summary(D4,testalph=0.05)

## End(Not run)

Description

Simulate dose response data and apply 4- or 3- parameter sigmoidal or hyperbolic Bayesian estimation. The prior distribution is input by the user with default values for some parameters based on the empirical distribution estimated from dose response meta-analyses. For binary response data, the Emax model is fit on the logit scale, and then back-transformed

Usage

emaxsimB(nsim, genObj, prior, modType = 4, 
	binary = FALSE, seed=12357,
	check = FALSE, nproc=parallel::detectCores(), 
	negEmax = FALSE, ed50contr = NULL, 
	lambdacontr = NULL, testMods = NULL, 
	idmax = length(doselev), 
	mcmc = mcmc.control(), 
	customCode=NULL, customParms=NULL,
	description = "")

Arguments

Details

The Bayesian model fits are implemented in rstan using function fitEmaxB. The function compileStanModels must be executed once to create compiled STAN code before emaxsimB can be used.

Continuous data can be simulated from any dose response curve with homogeneous normally distributed residuals.

Binary data are handled similarly. The models are fit on the logit scale and then back-transformed for estimation of dose response. Reduced linear models are selected based on the corresponding likelihood deviance.

MCP-MOD tests are created from contrasts based on the Emax function using the DoseFinding package. Different ED50 and lambda (Hill) parameters can be specified to form the contrasts. A contrast matrix output from the DoseFinding package can be specified instead, allowing for other contrast choices.

Customized code:

For binary data, the inputs to the function customCode for each simulated data set will be (parms,pVal,dose,y), where parms is the matrix of parameters generated from the posterior distribution with columns in the order given in function emaxfun, pVal is the MCP-MOD p-value, dose and y are the patient-level simulated data. For continuous data, the inputs are (parms,residSD,pVal,dose,y), where residSD are the variance parameters generated from their posterior distribution. The customParms supply other user-inputs such as a target efficacy level. When it is not null, the customCode inputs must be (parms,pVal,dose,y,customParms) or (parms,residSD,pVal,dose,y,customParms).

Value

A list is returned with class(emaxsim) containing:

Note

The default modType was changed from 3 to 4 for clinDR version >2.0

Author(s)

Neal Thomas

References

Thomas, N., Sweeney, K., and Somayaji, V. (2014). Meta-analysis of clinical dose response in a large drug development portfolio, Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317. <doi:10.1080/19466315.2014.924876>

Thomas, N., and Roy, D. (2016). Analysis of clinical dose-response in small-molecule drug development: 2009-2014. Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317 <doi:10.1080/19466315.2016.1256229>

See Also

print.emaxsimB, summary.emaxsimB, plot.emaxsimB, coef.emaxsimB, sigma.emaxsimB, emaxfun

Examples

## Not run: 

### emaxsimB changes the random number seed

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

prior<-emaxPrior.control(epmu=0,epsca=30,difTargetmu=0,
		difTargetsca=30,dTarget=100,p50=50,sigmalow=0.1,
		sigmaup=30,parmDF=5)
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,
		propInit=0.15,adapt_delta = 0.95)

### custom code to compute the distribution of the dose yielding
### a target diff with pbo
customCode<-function(parms,residSD,pVal,dose,y,customParms){
	target<-customParms
	ed50<-exp(parms[,1])
	emax<-parms[,2]
	td<-ifelse(emax-target>0,ed50*(target/(emax-target)),Inf)
	tdest<-median(td)
	lb<-quantile(td,0.1)
	ub<-quantile(td,0.9)
	return(c(td=tdest,lb=lb,ub=ub))
}

D1 <- emaxsimB(nsim,gen, prior, modType=4,seed=12357,mcmc=mcmc,check=FALSE,
				customCode=customCode,customParms=1.0)
D1

## End(Not run)

Description

Solve the Emax function for dose or Emax to yield a specified response.

Usage

solveEmax(target,dose,led50,lambda,e0,pboadj=TRUE)
solveDose(target,led50,lambda,emax,e0,pboadj=TRUE)

Arguments

Author(s)

Neal Thomas

See Also

fitEmax, fitEmaxB, emaxsim, emaxsimB

Examples

e0<-10
dose<-1
led50<-log(0.5)
lambda<-2
target<- -1.5
emax<-solveEmax(target,dose,led50,lambda,e0)
emax

dose1<-solveDose(target,led50,lambda,emax,e0)
dose1

emaxfun(dose=dose1,parm=c(led50,lambda,emax,e0)) - e0

Description

Calls Newton-Raphson optimizers, nls and nlm, for a hyperbolic or sigmoidal Emax model. Different intercepts for multiple protocol-data are supported. For binary data, the Emax model is on the logit scale.

Usage

fitEmax(y,dose,iparm,xparm,modType=4,
	prot=rep(1,length(y)),count=rep(1,length(y)),xbase=NULL,
	binary=FALSE,diagnostics=TRUE,msSat=NULL,
	pboAdj=rep(FALSE,max(prot)),optObj=TRUE)

Arguments

Details

Fits the 3- or 4- Emax model using nls. A newton-raphson algorithm is tried first followed by a partial linear optimatization if needed. Binary data are fit using nlm.

Value

A list assigned class "fitEmax" with:

Author(s)

Neal Thomas

See Also

nls, nlm, nllogis, predict.fitEmax, plot.fitEmax, coef.fitEmax

Examples

## the example changes the random number seed

doselev<-c(0,5,25,50,100,350)
n<-c(78,81,81,81,77,80)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)


sdy<-8.0
pop<-c(log(ed50),emax,e0)    
dose<-rep(doselev,n)
meanlev<-emaxfun(dose,pop)  

y<-rnorm(sum(n),meanlev,sdy)

testout<-fitEmax(y,dose,modType=4)

Description

Uses Rpackage rstan to fit a Bayesian hyperbolic or sigmoidal Emax model. Different intercepts for multiple protocol-data are supported. For binary data, the Emax model is on the logit scale.

Usage

fitEmaxB(y, dose, prior, modType = 4, prot = rep(1, length(y)), 
	count = rep(1, length(y)), xbase=NULL,
	binary = FALSE, msSat = NULL,vcest=NULL,
	pboAdj = FALSE, mcmc = mcmc.control(), estan = NULL, 
	diagnostics = TRUE, nproc = getOption("mc.cores", 1L))

Arguments

Details

The function compileStanModels must be executed once to create compiled STAN code before fitEmaxB can be used.

MCMC fit of a Bayesian hyperbolic or sigmoidal Emax model. The prior distributions available are based on the publication Thomas, Sweeney, and Somayaji (2014), Thomas and Roy (2016), and Wu, et al (2017).

The posterior distributions are complex because the distributions of the Emax and ED50 parameters change substantially as a function of the lambda, often creating 'funnel' type conditions. Small numbers of divergences are common with the 4-parameter model and do not appear easily avoided. Extensive simulation using evaluations with emaxsimB support the utility of the resulting approximate posterior distributions. The number of divergences can be viewed using diagnostics=TRUE. The usual convergence diagnostics should always be checked.

Value

A list assigned class "fitEmaxB" with:

Note

The default modType was changed from 3 to 4 for clinDR version >2.0

Author(s)

Neal Thomas

References

Thomas, N., Sweeney, K., and Somayaji, V. (2014). Meta-analysis of clinical dose response in a large drug development portfolio, Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317. <doi:10.1080/19466315.2014.924876>

Pinheiro, J., Bornkamp, B., Glimm, E., and Bretz, F. (2014). Model-based dose finding under model uncertainty using general parametric models, Vol. 33, No. 10, 1646-1661 <doi:/10.1002/sim.6052>

Thomas, N., and Roy, D. (2016). Analysis of clinical dose-response in small-molecule drug development: 2009-2014. Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317 <doi:10.1080/19466315.2016.1256229>

Wu, J., Banerjee, A., Jin, B. Menon, M. S., Martin, S. and Heatherington, A. (2017). Clinical dose response for a broad set of biological products: A model-based meta-analysis. Statistical Methods in Medical Research. <doi:10.1177/0962280216684528>

See Also

fitEmax, predict.fitEmaxB, plot.fitEmaxB, coef.fitEmaxB

Examples

## Not run: 

data("metaData")
exdat<-metaData[metaData$taid==1,]

prior<-emaxPrior.control(epmu=0,epsca=4,difTargetmu=0,difTargetsca=4,dTarget=20,
				p50=(2+5)/2,
				sigmalow=0.01,sigmaup=3)
										
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,prot=exdat$protid,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
plot(fitout)

## End(Not run)

Description

Creates a list object that contains inputs and a function to create simulated data sets with a common mean (proportion) for use in emaxsim with normal or continuous data

Usage

FixedMean(n, doselev, meanlev, resSD, parm = NULL, binary=FALSE)

Arguments

Value

A list of length 2. The first element is itself a list named genP that contains named elments n, resSD, doselev, dose, parm, binary, and the element meanlev, which is specific to FixedMean. The second element is a function named genFun that takes genP as input and returns a list with named elements meanlev, parm, resSD, y.

Author(s)

Neal Thomas

See Also

emaxsim, RandEmax

Examples

## Not run: 
##  example changes the random number seed

doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)  
  
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
genp<-FixedMean(n,doselev,meanlev,sdy,pop)  

### binary example
n<-rep(500,5)
doselev<-c(0,5,25,50,1000)
dose<-rep(doselev,n)

e0<- qlogis(0.2)
ed50<-20
diftarget<-qlogis(0.6)-qlogis(0.2)
lambda<-2
dtarget<-100
emax<-solveEmax(diftarget,dtarget,log(ed50),lambda,e0)

pop<-c(log(ed50),lambda,emax,e0)
meanlev<-plogis(emaxfun(doselev,pop))

genp<-FixedMean(n,doselev,meanlev,sdy,pop,binary=TRUE)  

tapply(genp$genFun(genp$genP)$y,dose,mean)
meanlev


## End(Not run)

Description

Set MCMC controls. Also control spread of initial parameter values.

Usage

mcmc.control(chains = 1, thin = 1, 
	warmup = 1000, iter = 3333* thin+warmup, 
	propInit = 0.50, seed = 12357, adapt_delta = 0.9)

Arguments

Note

Some defaults were changed with version>=2.0. For earlier versions, warmup = 500, iter = 5000* thin, and adapt_delta=0.8

Description

Dose response data from over 200 compounds included in published meta-analyses. The data are aggregated in a single data frame in a common format.

Usage

data('metaData')

Format

The data frame has one row for each compound, protocol within compound, and dose group within protocol. Compound and protocol level descriptors are repeated on each row of the data frame.

drugid

A numerical ID identifying each drug

taid

A drug can be studied in more than one therapeutic area. The taid ID identifies each TA/drug combination.

protid

Numerical (1,2,3,...) ID for protocols specific to each TAID.

gname

Generic drug name

bname

Branded(USA) drug name

drugtype

Drug classified as SMALL MOLECULE, BIOLOGIC, OTHER

route

Route of administration, e.g., oral, subcutaneous,...

routeShort

Abbreviated format for route

oralForm

Formulation (e.g., TABLET, POWDER,...) for drugs with oral administration.

fdaapproved

NA if status was not yet determined

metasource

Meta-analysis contributing compounds. BIO14: biological compounds through 2014; FDA914: FDA approved small molecules and 'other' 2009-2014; FDA1417: FDA approved compounds 2014-2017; Pfizer P2 compounds 1998-2009; PFIZERUPDATE18: Pfizer compounds 2009-2018

protno

Sponsor assigned protocol name/number

nctno

Clintrial.gov protocol ID

protyear

When available, year of first patient/first visit. In some cases, date of journal publication

design

PARELLEL, CROSSOVER,...

actcomp

Indicator if an active comparator was included in the protocol

etype

etype=1 for the designated primary endpoint. For completeness, where there was ambiguity in the selection of the endpoint, additional endpoint data was included on separate rows and indicated by etype=2,3,... Most analyses subset on etype=1

poptype

For a compound and TA, there can be distinctly different populations with anticipated response differences, e.g., treatment-naive and pre-treated patients. The population with the most studied doses has poptype=1. For completeness, additional populations are included and identified by poptype=2,3,.... Most analyses subset on poptype=1

primsource

IRO/PRO investigator/patient reported outcome; L lab, V vitals

primtype

Primary endpoint is BINARY, CONTINOUS, TIMETOEVEN

primtime

time units to primary endpoint from randomization

timeunit

DAY, HR, MIN, MONTH, WK for primary endpoint

indication

Disease description

broadta

Broad TA classification of the indication

endpointLong,endpointShort

Endpoint name and an abbreviated form using for example, cfb and pcfb for change from baseline and percent change from baseline

dose

Total daily dose for small molecules, total weekly dose for biologics in mg or mg/kg for weight-based dosing.

tload

Amount of any loading dose

nload

Number of visits with a loading dose

regimen

Dosing frequency

primregimen

primregimen=1 for most doses/regimens, but primregimen=2 for a few regimens that clearly differed from the most common regimen for the same total dose. Most analyses subset on primregimen=1

rslt

The sample dose group mean (continuous) or proportion (binary) of the primary endpoint. Analyses of the time-to-event endpoints was compound specific (either a mean or a proportion was estimated).

se

Standard error of rslt

sd

Dose group sample standard deviation for continuous data

lcl, ucl, alpha

alpha-level interval (lcl,ucl) when confidence intervals were extracted from the original data source because se were not reported

sampsize

Sample size reported for rslt. The handling of missing data by the protocol sponsors varied, but 'completers' was most common.

ittsize

The number randomized. The counts are usually available, except for internal data before 2009, where it was not collected.

pmiss

Percent of missing data.

Details

Compound sampling plans and other details are given in the publications:

Thomas, N., Sweeney, K., and Somayaji, V. (2014). Meta-analysis of clinical dose response in a large drug development portfolio, Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317. <doi:10.1080/19466315.2014.924876>

Thomas, N., and Roy, D. (2016). Analysis of clinical dose-response in small-molecule drug development: 2009-2014. Statistics in Biopharmaceutical Research, Vol. 6, No.4, 302-317 <doi:10.1080/19466315.2016.1256229>

Wu, J., Banerjee, A., Jin, B. Menon, M. S., Martin, S. and Heatherington, A. (2017). Clinical dose response for a broad set of biological products: A model-based meta-analysis. Statistical Methods in Medical Research. <doi:10.1177/0962280216684528>

Examples

data('metaData')
names(metaData)

Description

The negative log likelihood function evaluated with a single input set of parameters for the binary Emax model on the logistic scale. For use with function fitEmax

Usage

nllogis(parms,y,dose,
          prot=rep(1,length(y)),
          count=rep(1,length(y)),
          xbase=NULL)

Arguments

Details

The negative log likelihood for the 3- or 4- Emax model on the logit scale for binary data. Note the ordering of the parameters and their transformations. A 3 vs 4 parameter model is deterimined by the length of parms.

Value

Negative log likelihood value is returned.

Author(s)

Neal Thomas

See Also

nlm, fitEmax

Examples

data('metaData')
exdat<-metaData[metaData$taid==8,]

cy<-round(exdat$sampsize*exdat$rslt)
y<-c(rep(1,length(cy)),rep(0,length(cy)))
cy<-c(cy,exdat$sampsize-cy)
drep<-c(exdat$dose,exdat$dose)
plotD(exdat$rslt,exdat$dose,se=FALSE)
nllogis(parms=c(log(2.5),-3.26,-0.15), y, drep,count=cy)

Description

A Q-Q plot of the dose response estimate of the mean at a specified dose minus the population value divided by the standard error of the estimator (computed using the delta method). Estimates based on alternative models when the Emax estimation fails are highlighted in red.

Usage

## S3 method for class 'emaxsim'
plot(x, id = x$idmax, plotDif= TRUE, ...)

Arguments

Value

No output is returned.

Author(s)

Neal Thomas

See Also

emaxsim, print.emaxsim, summary.emaxsim

Examples

## Not run: 
nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop.parm<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop.parm)  

###FixedMean is specialized constructor function for emaxsim
gen.parm<-FixedMean(n,doselev,meanlev,sdy)  
D1 <- emaxsim(nsim,gen.parm)

plot(D1,id=3)

## End(Not run)

Description

A Q-Q plot of the posterior mean of the mean dose response at a specified dose minus the population value divided by the posterior SD of the mean difference.

Usage

## S3 method for class 'emaxsimB'
plot(x, id = x$idmax, plotDif= TRUE, ...)

Arguments

Value

ggplot object is returned

Author(s)

Neal Thomas

See Also

emaxsimB, print.emaxsimB, summary.emaxsimB

Examples

## Not run: 
## emaxsimB changes the random number seeds
nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

prior<-emaxPrior.control(epmu=0,epsca=30,difTargetmu=0,
		difTargetsca=30,dTarget=100,p50=50,sigmalow=0.1,
		sigmaup=30,parmDF=5)
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,propInit=0.15,adapt_delta = 0.95)

D1 <- emaxsimB(nsim,gen, prior, modType=3,mcmc=mcmc,check=FALSE)

plot(D1,id=3)

## End(Not run)

Description

Plot of population dose response curve, sample dose group means, posterior and posterior predictive intervals, and the model-based estimated (posterior means) dose response curve.

Usage

## S3 method for class 'emaxsimBobj'
plot(
	x, clev=0.9, plotDif=FALSE, 
	plotPop=c('m','3','4'), 
	logScale=FALSE, plotResid=FALSE, 
	plot=TRUE, ... )

Arguments

Note

The estimated curve is the posterior mean evaluated along a grid of dose values.

Examples

## Not run: 

## emaxsimB changes the random number seed

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 
dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

prior<-emaxPrior.control(epmu=0,epsca=30,difTargetmu=0,
		difTargetsca=30,dTarget=100,p50=50,sigmalow=0.1,
		sigmaup=30,parmDF=5)
		
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,propInit=0.15,adapt_delta = 0.95)

D1 <- emaxsimB(nsim,gen, prior, modType=3,mcmc=mcmc,check=FALSE)

plot(D1,id=3)
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,propInit=0.15,adapt_delta = 0.95)

D1 <- emaxsimB(nsim,gen, prior, modType=3,mcmc=mcmc,check=FALSE)

plot(D1[2])

## End(Not run)

Description

Plot of population dose response curve, dose group means with CIs, predictive intervals, and the model-based estimated dose response curve.

Usage

## S3 method for class 'emaxsimobj'
plot(
    x, xlim, xat=NULL, ylim, xlab, ylab,
    plotDif=FALSE,
    plotResid=FALSE,
    clev = 0.9, 
    plotPop=c('m','3','4'),     
    negC = FALSE, 
    logScale=FALSE,
    predict=TRUE, 
    plot=TRUE, ...)

Arguments

Value

ggplot object is returned

Author(s)

Neal Thomas

See Also

emaxsim, print.emaxsimobj, summary.emaxsimobj, update.emaxsimobj

Examples

## Not run: 
## emaxsim changes the random number seed

nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen.parm<-FixedMean(n,doselev,meanlev,sdy)  
D1 <- emaxsim(nsim,gen.parm)
e49<-D1[49]

plot(e49,clev=0.8)  

## End(Not run)

Description

Plot an Emax model stored in an object created by function fitEmax.

Usage

## S3 method for class 'fitEmax'
plot(
          x,int=0,plotResid=FALSE,clev=0.9,
          predict=TRUE,plotci=TRUE,plotDif=FALSE,
          xlab='Dose',
          ylab=ifelse(plotResid,'Residuals',ifelse(plotDif,
                'Difference With Placebo','Response')),
          ncol=NULL,
          symbol=NULL,symbolLabel='Group',symbolShape=8,
          symbolColor='red',symbolSize=4,
          bwidth=NULL,
          xlim=NULL,
          xat=NULL,
          ylim=NULL,
          logScale=FALSE,       
          ngrid=200,
          plot=TRUE, ...)

Arguments

Details

Model estimates, standard errors, and confidence bounds are computed using function SeEmax.

The function generates random numbers when predict=TRUE, so the random number generator/seed must be set before the function is called for exact reproducibility.

Value

A list with ggplot object, and a matrix with the confidence and prediction interval limits.

Author(s)

Neal Thomas

See Also

nls

Examples

### example changes the random number seed

doselev<-c(0,5,25,50,100,350)
n<-c(78,81,81,81,77,80)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-8.0
pop.parm<-c(log(ed50),emax,e0)    
dose<-rep(doselev,n)
meanlev<-emaxfun(dose,pop.parm)  

y<-rnorm(sum(n),meanlev,sdy)

testout<-fitEmax(y,dose,modType=4)

plot(testout)

Description

Plot an Emax model stored in an object created by function fitEmaxB.

Usage

## S3 method for class 'fitEmaxB'
plot(
      x,int=0,plotResid=FALSE,clev=0.9,
      predict=TRUE,plotci=TRUE,plotDif=FALSE,
      xlab='Dose',
      ylab=ifelse(plotResid,'Residuals',ifelse(plotDif,
            'Difference With Placebo','Response')),
      ncol=NULL,
      symbol=NULL,symbolLabel='Group',symbolShape=8,
      symbolColor='red',symbolSize=4,
      bwidth=NULL,
      xlim=NULL,
      xat=NULL,
      ylim=NULL,
      logScale=FALSE,       
      ngrid=200,
      plot=TRUE, ...)

Arguments

Details

Model-based medians, standard deviations, and interval bounds for the dose groups means/proportions based on the MCMC parameters evaluated in the Emax function.

The function generates random numbers when predict=TRUE, so the random number generator/seed must be set before the function is called for exact reproducibility.

If baseline covaraties were included in the fit, then the mean of the predictions for the protocol given by int is plotted. This can be computationally intensive when the dosing grid is dense, the MCMC sample size is large, and the input sample size is large. Consider reducing ngrid in this situation. Note that the protocol must be specified, or the prediction defaults to patients from the first protocol.

Value

A list with ggplot object, and posterior and prediction interval limits.

Author(s)

Neal Thomas

See Also

fitEmaxB

Examples

## Not run: 

data("metaData")
exdat<-metaData[metaData$taid==1,]

prior<-emaxPrior.control(epmu=0,epsca=4,difTargetmu=0,difTargetsca=4,dTarget=20,
				p50=(2+5)/2,
				sigmalow=0.01,sigmaup=3)
										
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,prot=exdat$protid,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
plot(fitout)

## End(Not run)

Description

Plot a dose response curve fit by Bayes MCMC methods (with optional posterior interval bars). Also plot dose group means (with optional CI bars)

Usage

## S3 method for class 'plotB'
plot( x,
	plotDif= FALSE, plotMed= FALSE, 
	plotResid=FALSE, predict= TRUE,
	logScale=FALSE,
	xlim,
	xat=NULL,
	ylim,  
	xlab, 
	ylab,labac='Act Comp',shapeac=8,colac='red',  
	symbolLabel='Group',symbolShape=8,
	symbolColor='red',symbolSize=4, ...)

Arguments

Details

Produce additional plots from output of plotB without any re-computing. A plot is produced by default on return from the function. When active control is specified, the plot is 'printed' within the function. If there is a symbol group variable, it must be specified when plotB is executed. The symbol label, shape, color, and size must be re-specified in subsequent plot requests.

Value

ggplot object of the dose response curve, which will be plotted by default unless the output of the plot is assigned. When an active control group is present, the value returned is an invisible list with the ggplot for the dosing data, and a second ggplot for the ac data.

Note

PlotB can also be used with draws from a prior distribution to evaluate the prior dose response curve.

Author(s)

Neal Thomas

See Also

plotB, plotD, plot.fitEmax

Examples

## Not run: 
data("metaData")
exdat<-metaData[metaData$taid==6 & metaData$poptype==1,]

prior<-emaxPrior.control(epmu=0,epsca=100,difTargetmu=0,difTargetsca=100,dTarget=80.0,
        p50=3.75,sigmalow=0.01,sigmaup=20)
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
parms<-coef(fitout)[,1:4]  #use first intercept

outB<-plotB(exdat$rslt,exdat$dose,parms, sigma2=(sigma(fitout))^2,
			ylab="Change in EDD")
			
plot(outB,plotDif=TRUE)

## End(Not run)

Description

Plot a dose response curve fit by Bayes MCMC methods (with optional posterior interval bars). Also plot dose group means (with optional CI bars)

Usage

plotB(y, 
	dose, 
	parm, 
	sigma2,
	count=rep(1,length(y)),
	dgrid=sort(unique(c(seq(0,max(dose),length=50), dose))),
	predict= TRUE,plotDif=FALSE,plotMed=FALSE,
	plotResid=FALSE,clev=0.9,
	binary=c('no','logit','probit','BinRes'),BinResLev,
	BinResDir=c('>','<'),
	activeControl=FALSE,ac,yac,
	countac=rep(1,length(yac)),
	labac='Act Comp',shapeac=8,colac='red',
	symbol,symbolLabel='Group',symbolShape=8,
	symbolColor='red',symbolSize=4,
	xlim,ylim,xat=NULL,xlab="Dose",
	ylab=ifelse(plotDif,"Diff with Comparator","Mean"),
	modelFun=emaxfun,makePlot=TRUE,
	...)

Arguments

Details

A sample of parameters from the joint posterior distribution must be supplied (typically produced by BUGS). The Bayesian dose response curve is the Bayes posterior mean (or median) at each value on dgrid. The bar (interval) is the (clev/2,1-clev/2) Bayes posterior interval (which can differ from the Bayes HPD interval). The intervals are plotted only at the dose levels included in the study. Predictive intervals are formed by adding independent random draws from the sampling distributions of the dose group sample means to the population means.

The function generates random numbers when predict=TRUE, so the random number generator/seed must be set before the function is called for exact reproducibility.

Value

Returns an object of class plotB. Three inputs are saved for later plotting: doses in the original design, dgrid, and clev. The following matrices are saved:

Note

PlotB can also be used with draws from a prior distribution to evaluate the prior dose response curve.

Author(s)

Neal Thomas

References

Spiegelhalter, D., Thomas, A., Best, N., and Lunn, D. (2003), WinBUGS User Manual Version 1.4, Electronic version www.mrc-bsu.cam.ac.uk/bugs

See Also

plot.plotB, plotD, plot.fitEmax

Examples

## Not run: 
data("metaData")
exdat<-metaData[metaData$taid==6 & metaData$poptype==1,]

prior<-emaxPrior.control(epmu=0,epsca=100,difTargetmu=0,difTargetsca=100,dTarget=80.0,
        p50=3.75,sigmalow=0.01,sigmaup=20)
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
parms<-coef(fitout)[,1:4]  #use first intercept

outB<-plotB(exdat$rslt,exdat$dose,parms, sigma2=(sigma(fitout))^2,
			ylab="Change in EDD")
			
plot(outB,plotDif=TRUE)

## End(Not run)

Description

Density plot over a grid of doses displaying the prior or posterior distribution for the mean dose response computed from simulated input model parameters.

Usage

plotBdensity(dgrid, 
      parm, 
      modelFun=emaxfun,
      qlevL=c(0.025,0.05,0.10,0.25),
      plotDif= FALSE,
      logit= FALSE, ...)

Arguments

Details

A sample of parameters from the joint prior or posterior distribution must be supplied (typically produced by BUGS). A density plot with contours corresponding to the perentiles in qlevL created by function DRDensityPlot .

Value

A list containing two matrices with the number of rows equal to the number dose grid points, and columns corresponding to percentiles in qlevL:

Author(s)

Neal Thomas

References

Spiegelhalter, D., Thomas, A., Best, N., and Lunn, D. (2003), WinBUGS User Manual Version 1.4, Electronic version www.mrc-bsu.cam.ac.uk/bugs

See Also

plot.plotB, plotD, plot.fitEmax, DRDensityPlot

Examples

## Not run: 
data("metaData")
exdat<-metaData[metaData$taid==6 & metaData$poptype==1,]

prior<-emaxPrior.control(epmu=0,epsca=10,difTargetmu=0,difTargetsca=10,dTarget=80.0,
        p50=3.75,sigmalow=0.01,sigmaup=20)
mcmc<-mcmc.control(chains=3)

msSat<-sum((exdat$sampsize-1)*(exdat$sd)^2)/(sum(exdat$sampsize)-length(exdat$sampsize))
fitout<-fitEmaxB(exdat$rslt,exdat$dose,prior,modType=4,
				count=exdat$sampsize,msSat=msSat,mcmc=mcmc)
parms<-coef(fitout)[,1:4]  #use first intercept

dgrid<-seq(0,1,length=100)

pout<-plotBdensity(dgrid,parm=parms)

pout2<-plotBdensity(dgrid,parm=parms,plotDif=TRUE,
       xlab='Dose',ylab='Dif with PBO')

## End(Not run)

Description

Plot dose group means vs dose with options to connect points by lines, and include CI about each dose group mean based on within-group SDs

Usage

plotD(y, dose, baseline, se = TRUE, line = TRUE, 
    meansOnly=FALSE,sem=NULL,clev = 0.9, 
    xlab='Dose',ylab='Response', logScale=FALSE)

Arguments

Value

Returns a list with the ggplot object and two vectors with the dose group means and their standard errors.

Author(s)

Neal Thomas

See Also

plot.fitEmax, plotB

Examples

data("metaData")
exdat<-metaData[metaData$taid==2 & metaData$etype==1,]
with(exdat,plotD(rslt,dose,meansOnly=TRUE,se=TRUE,sem=se,ylab=
"Y",xlab="Dose(mg)"))

Description

Estimated mean and standard error for specified doses computed from the output of a model fit by function emaxalt. Also returns mean difference with placebo and their standard errors.

Usage

## S3 method for class 'emaxalt'
predict(object,dose, dref=0, ...)

Arguments

Value

A list containing:

Author(s)

Neal Thomas

See Also

emaxalt, predict.emaxsimobj, predict.emaxsim

Examples

## Not run: 
## random number seed changed by this example

doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
dose<-rep(doselev,n)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop.parm<-c(log(ed50),e0,emax)    
meanresp<-emaxfun(dose,pop.parm)  
y<-rnorm(sum(n),meanresp,sdy)

simout<-emaxalt(y,dose)
predict(simout,c(75,150))

simout2<-emaxalt(y,dose,modType=4)
predict(simout2,c(75,150))

## End(Not run)

Description

Estimated mean/proportion and standard error for each simulated data set in an emaxsim object. Also returns mean difference with placebo and their standard errors.

Usage

## S3 method for class 'emaxsim'
predict(object, 
         dose, dref=0, ...)

Arguments

Value

A list containing:

Author(s)

Neal Thomas

See Also

emaxsim, summary.emaxsim, plot.emaxsim

Examples

## Not run: 
## random number seed changed by this example
nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-9.032497
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop.parm<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop.parm)  

###FixedMean is specialized constructor function for emaxsim
gen.parm<-FixedMean(n,doselev,meanlev,sdy)  
D1 <- emaxsim(nsim,gen.parm)

predout<-predict(D1,c(75,150))

## End(Not run)

Description

Return warning and explanation that only predicted values at doses included in the study are available. The code needed to obtain predicted values at other doses is indicated.

Usage

## S3 method for class 'emaxsimB'
predict(object, 
          dose, dref=0, ...)

Arguments

Value

No output.

Author(s)

Neal Thomas

See Also

emaxsimB, summary.emaxsimB, plot.emaxsimB

Examples

## Not run: 
nsim<-50
idmax<-5
doselev<-c(0,5,25,50,100)
n<-c(78,81,81,81,77)
Ndose<-length(doselev)

### population parameters for simulation
e0<-2.465375 
ed50<-67.481113 

dtarget<-100
diftarget<-2.464592
emax<-solveEmax(diftarget,dtarget,log(ed50),1,e0)

sdy<-7.967897
pop<-c(log(ed50),emax,e0)    
meanlev<-emaxfun(doselev,pop)  

###FixedMean is specialized constructor function for emaxsim
gen<-FixedMean(n,doselev,meanlev,sdy)  

prior<-emaxPrior.control(epmu=0,epsca=30,difTargetmu=0,
		difTargetsca=30,dTarget=100,p50=50,sigmalow=0.1,
		sigmaup=30,parmDF=5)
mcmc<-mcmc.control(chains=1,warmup=500,iter=5000,seed=53453,
		propInit=0.15,adapt_delta = 0.95)


D1 <- emaxsimB(nsim,gen, prior, modType=3,seed=12357,mcmc=mcmc,check=FALSE)

predict(D1,dose=20)

## End(Not run)