Data Engineering and BI courses are free this week!
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.
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))A list assigned class "fitEmaxB" with:
The rstan object with the model fit.
emaxPrior specified]
Input values.
Outcome for each patient. Missing Y values
are are not permitted. Dose/protocol group means for grouped continuous data. For binary
data, y must be 0/1 and counts must be supplied for each 0/1 value.
Dose for each patient.
Prior specification through an object of type 'emaxPrior' or 'prior'.
See emaxPrior.control and prior.control for details.
The 'emaxPrior' specifies
the magnitude of the potential effect for a specified dose (typically the highest
anticipated dose and/or the dose in a POC study), while the 'prior' specifies the
theoretical maximum effect (the emax parameter). The 'prior' specification is
deprecated and will be removed.
modType=3 (default) for the 3-parameter hyperbolic Emax model. modType=4 for the 4-parameter sigmoidal Emax model.
Protocol (group) membership used to create multiple intercepts. The default is a single protocol. The prior disribution for the placebo response is re-used independently for each intercept.
Counts for the number of patients when the Y are dose
continuous group means or binary 0/1 values.
Default is 1 (ungrouped data).
A matrix of baseline covariates with rows corresponding to
y that enter as linear additive predictors. The baseline covariates must
be centered about
their (protocol-specific) means. xbase does not include an intercept or
protocol indicators. Covariates cannnot be specified with PBO adjusted or
aggregated input.
When TRUE, the y are assumed to be coded 0/1, and
the means reported are proportions. The Emax model is specified on the logit
scale, and proportions are estimated from the model by back-transformation.
If continuous Y are dose/protocol group means rather than
individual measurements, the within group variance, msSat, should be
supplied. This variance is the mean square from the model saturated in
dose and protocol. It is used to
improve the residual variance estimate for the Emax model. If it is not
supplied, the
residual SD (and associated SE) will have low degrees of freedom.
The input, Y, can be estimates of dose group
responses from a first-stage
model. The vcest is the variance-covariance matrix of the
model-based estimates. The most common usage is when a saturated model
is fit using maximum likelihood estimation to longitudinal data
to produce dose group estimates that are valid under the MAR assumption
for missing values. Other applications are possible, see the Pinheiro,
et al reference. The count and msSat are ignored when
vcest is specified. Covariates xbase cannot be specified
with vcest, but covariates can be included in the first stage
modeling.
For published data with only pbo-adjusted dose group means and
SEs, the model is fit without an intercept(s). If initial parameters
are supplied, the intercept (E0) should be assigned 0.
A zero for the placebo mean should not be included in Y.
This option is not available for binary data. Potential correlation between between
placebo-adjusted means is ignored.
Inputs controlling rstan execution. See mcmc.control
for details.
The compiled rstan Emax model is usually loaded
automatically. It can be load to an object using the function
selEstan and passed to fitEmaxB for repeated executions to
improve efficiency and stability.
Printed output from rstan. See Details for more information.
The number of processor requested for STAN MCMC
computations. Defaults to the value set by the
rstan installation. When set explicitly, nproc
is usually 1 or the number of MCMC chains. If greater than the
number of chains, it is set to the number of chains.
Neal Thomas
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.
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>
fitEmax, predict.fitEmaxB, plot.fitEmaxB, coef.fitEmaxB
if (FALSE) {
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)
}
Run the code above in your browser using DataLab