Description Usage Arguments Details Value Author(s) References See Also Examples
View source: R/jointSurroCopPenal.R
Joint FrailtyCopula model for Surrogacy definition
Fit the onestep Joint frailtycopula surrogate model for the evaluation of a canditate surrogate endpoint,
with different integration methods on the random effects, using a semiparametric penalized
likelihood estimation. This approach extends that of Burzykowski et al.
(2001) by
including in the bivariate copula model the random effects treatmentbytrial interaction.
Assume S_{ij} and T_{ij} the failure times associated respectively
with the surrogate and the true endpoints, for subject j(j = 1,..., n
_{i})
belonging to
the trial i (i = 1,..., G)
.
Let v_{i} = (u_{i}, v_{Si}, v_{Ti}) be the vector of trial level random effects; Z_{S,ij} = (Z_{Sij1}, ..., Z_{Sijp})^{'} and Z_{T,ij} = (Z_{Tij1}, ..., Z_{Tijp})^{'} be covariates associated with S_{ij} and T_{ij}. The joint frailtycopula model is defined as follows:
where
and the conditional survival functions are given by
in which
u_{i} ~ N(0,γ), u_{i} ⊥ v_{Si}, u_{i} ⊥ v_{Ti}; (v_{Si}, v_{Ti})^{T} ~ N(0,Σ_{v})
with
In this model, λ_{0s}(x) is the baseline hazard function associated with the surrogate endpoint and β_{S} the fixed effects (or loghazard ratio) corresponding to the covariates Z_{S,ij}; λ_{0T}(x) is the baseline hazard function associated with the true endpoint and β_{T} the fixed treatment effects corresponding to the covariates Z_{T,ij}. The copula model serves to consider dependence between the surrogate and true endpoints at the individual level. In the copula model, θ is the copula parameter used to quantify the strength of association. u_{i} is a shared frailty effect associated with the baseline hazard function that serve to take into account the heterogeneity between trials of the baseline hazard function, associated with the fact that we have several trials in this metaanalytical design. The power parameter α distinguishes triallevel heterogeneity between the surrogate and the true endpoint. v_{Si} and v_{Ti} are two correlated random effects treatmentbytrial interactions. Z_{Sij1} or Z_{Tij1} represents the treatment arm to which the patient has been randomized.
For simplicity, we focus on the Clayton and GumbelHougaard copula functions. In Clayton's model, the copula function has the form
and in Gumbel's model, the copula function has the form
Surrogacy evaluation
We propose to base validation of a candidate surrogate endpoint on Kendall's τ at the individual level and
coefficient of determination at the trial level, as in the classical approach (Burzykowski et al.
, 2001).
The formulations are given below.
Individuallevel surrogacy
From the proposed model, according to the copula function, it can be shown that Kendall's τ is defined as :
where θ is the copula parameter. Kendall's τ is the difference between the probability of concordance and the probability of discordance of two realizations of S_{ij} and T_{ij}. It belongs to the interval [1,1] and assumes a zero value when S_{ij} and T_{ij} are independent.
Triallevel surrogacy
The key motivation for validating a surrogate endpoint is to be able to predict the effect of treatment on the true endpoint, based on the observed effect of treatment on the surrogate endpoint. As shown by Buyse et al. (2000), the coefficenient of determination obtains from the covariance matrix Σ_{v} of the random effects treatmentbytrial interaction can be used to evaluate underlined prediction, and therefore as surrogacy evaluation measurement at triallevel. It is defined by:
The SEs of R_{trial}^{2} and τ are calculated using the Deltamethod. We also propose R_{trial}^{2} and 95% CI computed using the parametric bootstrap. The use of deltamethod can lead to confidence limits violating the [0,1], as noted by (Burzykowski et al., 2001). However, using other methods would not significantly alter the findings of the surrogacy assessment
1 2 3 4 5 6 7 8 9 10  jointSurroCopPenal(data, maxit = 40, indicator.alpha = 1,
frail.base = 1, n.knots = 6, LIMparam = 0.001, LIMlogl = 0.001,
LIMderiv = 0.001, nb.mc = 1000, nb.gh = 20, nb.gh2 = 32,
adaptatif = 0, int.method = 0, nb.iterPGH = 5, true.init.val = 0,
thetacopula.init = 1, sigma.ss.init = 0.5, sigma.tt.init = 0.5,
sigma.st.init = 0.48, gamma.init = 0.5, alpha.init = 1,
betas.init = 0.5, betat.init = 0.5, scale = 1,
random.generator = 1, kappa.use = 4, random = 0,
random.nb.sim = 0, seed = 0, init.kappa = NULL, ckappa = c(0,0),
typecopula = 1, nb.decimal = 4, print.times = TRUE, print.iter = FALSE)

data 
A

maxit 
maximum number of iterations for the Marquardt algorithm.
The default being 
indicator.alpha 
A binary, indicating whether the power's parameter α should
be estimated (1) or not (0). If 
frail.base 
A binary, indicating whether the heterogeneity between trial on the baseline
risk is considered ( 
n.knots 
integer giving the number of knots to use. Value required in
the penalized likelihood estimation. It corresponds to the (n.knots+2)
splines functions for the approximation of the hazard or the survival
functions. We estimate I or Msplines of order 4. When the user set a
number of knots equals to k (n.knots=k) then the number of interior knots
is (k2) and the number of splines is (k2)+order. Number of knots must be
between 4 and 20. (See 
LIMparam 
Convergence threshold of the Marquardt algorithm for the
parameters, 10^{3} by default (See 
LIMlogl 
Convergence threshold of the Marquardt algorithm for the
loglikelihood, 10^{3} by default (See 
LIMderiv 
Convergence threshold of the Marquardt algorithm for the gradient,
10^{3} by default
(See 
nb.mc 
Number of samples considered in the MonteCarlo integration. Required in the event

nb.gh 
Number of nodes for the GaussianHermite quadrature. It can
be chosen among 5, 7, 9, 12, 15, 20 and 32. The default is 
nb.gh2 
Number of nodes for the GaussHermite quadrature used to reestimate the model,
in the event of nonconvergence, defined as previously. The default is 
adaptatif 
A binary, indicates whether the pseudo adaptive GaussianHermite quadrature 
int.method 
A numeric, indicates the integration method: 
nb.iterPGH 
Number of iterations before the reestimation of the posterior random effects,
in the event of the twosteps pseudoadaptive Gaussianhermite quadrature. If set to 
true.init.val 
Numerical value. Indicates if the given initial values to parameters 
thetacopula.init 
Initial values for the copula parameter (θ), required if 
sigma.ss.init 
Initial values for
σ^{2}_{vS}, required if 
sigma.tt.init 
Initial values for
σ^{2}_{vT}, required if 
sigma.st.init 
Initial values for
σ_{vST}, required if 
gamma.init 
Initial values for γ, required if 
alpha.init 
Initial values for α, required if 
betas.init 
Initial values for β_{S}, required if 
betat.init 
Initial values for β_{T}, required if 
scale 
A numeric that allows to rescale (by multiplication) the survival times, to avoid numerical
problems in the event of some convergence issues. If no change is needed the argument is set to 1, the default value.
eg: 
random.generator 
Random number generator used by the Fortran compiler,

kappa.use 
A numeric, that indicates how to manage the smoothing parameters
k_{1}
and k_{2} in the event of convergence issues. If it is set to 
random 
A binary that says if we reset the random number generation with a different environment
at each call 
random.nb.sim 
If 
seed 
The seed to use for data (or samples) generation. required if 
init.kappa 
smoothing parameter used to penalized the loglikelihood. By default (init.kappa = NULL) the values used are obtain by crossvalidation. 
ckappa 
Vector of two fixed values to add to the smoothing parameters. By default it is set to (0,0). this argument allows to well manage the smoothing parameters in the event of convergence issues. 
typecopula 
The copula function used, can be 1 for clayton or 2 for GumbelHougaard. The default is 
nb.decimal 
Number of decimal required for results presentation. 
print.times 
a logical parameter to print estimation time. Default is TRUE. 
print.iter 
a logical parameter to print iteration process. Default is FALSE. 
The estimated parameter are obtained using the robust Marquardt algorithm (Marquardt, 1963) which is a combination between a NewtonRaphson algorithm and a steepest descent algorithm. The iterations are stopped when the difference between two consecutive loglikelihoods was small (< 10^{3} ), the estimated coefficients were stable (consecutive values (< 10^{3} )), and the gradient small enough (< 10^{3} ), by default. Cubic Msplines of order 4 are used for the hazard function, and Isplines (integrated Msplines) are used for the cumulative hazard function.
The inverse of the Hessian matrix is the variance estimator and to deal with the positivity constraint of the variance component and the spline coefficients, a squared transformation is used and the standard errors are computed by the Δmethod (Knight & Xekalaki, 2000). The smooth parameter can be chosen by maximizing a likelihood cross validation criterion (Joly and other, 1998).
We proposed based on the joint surrogate model a new definition of the Kendall's τ. Moreover, distinct numerical integration methods are available to approximate the integrals in the marginal loglikelihood.
Nonconvergence case management procedure
Special attention must be given to initializing model parameters, the choice of the number of
spline knots, the smoothing parameters and the number of quadrature points to solve convergence
issues. We first initialized parameters using the user's desired strategy, as specified
by the option true.init.val
. When numerical or convergence problems are encountered,
with kappa.use
set to 4
, the model is fitted again using a combination of the following strategies:
vary the number of quadrature point (nb.gh
to nb.gh2
or nb.gh2
to nb.gh
)
in the event of the use of the Gaussian Hermite quadrature integration (see int.method
);
divided or multiplied the smoothing parameters ( k_{1},
k_{2}) by 10 or 100 according to
their preceding values, or used parameter vectors obtained during the last iteration (with a
modification of the number of quadrature points and smoothing parameters). Using this strategy,
we usually obtained during simulation the rejection rate less than 3%. A sensitivity analysis
was conducted without this strategy, and similar results were obtained on the converged samples,
with about a 23% rejection rate.
This function return an object of class jointSurroPenal with elements :
EPS 
A vector containing the obtained convergence thresholds with the Marquardt algorithm, for the parameters, the loglikelihood and for the gradient; 
b 
A vector containing estimates for the splines parameter's; elements of the
lower triangular matrix (L) from the Cholesky decomposition such that Σ = LL^{T}, with Σ
the covariance of the random effects (v_{Si},v_{Ti});
the coefficient α (if 
varH 
The variance matrix of all parameters in 
varHIH 
The robust estimation of the variance matrix of all parameters in 
loglikPenal 
The complete marginal penalized loglikelihood; 
LCV 
the approximated likelihood crossvalidation criterion in the semiparametric case (with 
xS 
vector of times for surrogate endpoint where both survival and hazard function are estimated. By default seq(0,max(time),length=99), where time is the vector of survival times; 
lamS 
array (dim = 3) of hazard estimates and confidence bands, for surrogate endpoint; 
survS 
array (dim = 3) of baseline survival estimates and confidence bands, for surrogate endpoint; 
xT 
vector of times for true endpoint where both survival and hazard function are estimated. By default seq(0, max(time), length = 99), where time is the vector of survival times; 
lamT 
array (dim = 3) of hazard estimates and confidence bands, for true endpoint; 
survT 
array (dim = 3) of baseline survival estimates and confidence bands, for true endpoint; 
n.iter 
number of iterations needed to converge; 
theta 
Estimate for θ; 
gamma 
Estimate for γ; 
alpha 
Estimate for α; 
zeta 
A value equals to 
sigma.s 
Estimate for σ^{2}_{vS}; 
sigma.t 
Estimate for σ^{2}_{vT}; 
sigma.st 
Estimate for σ_{vST}; 
beta.s 
Estimate for β_{S}; 
beta.t 
Estimate for β_{T}; 
ui 
A binary, that indicates if the heterogeneity between trial on the baseline risk
has been Considered ( 
ktau 
The Kendall's τ with the correspondant 95 \% CI obtained from the deltamethod; 
R2.boot 
The

Coefficients 
The estimates with the corresponding standard errors and the 95 \% CI 
kappa 
Positive smoothing parameters used for convergence. These values could be different to initial
values if 
scale 
The value used to rescale the survival times 
data 
The dataset used in the model 
varcov.Sigma 
Covariance matrix of the estimates of the estimates of (σ^{2}_{vS},σ^{2}_{vT}, σ_{vST}) obtained from the deltamethod 
parameter 
List of all arguments used in the model 
type.joint 
A code 
Casimir Ledoux Sofeu casimir.sofeu@ubordeaux.fr, scl.ledoux@gmail.com and Virginie Rondeau virginie.rondeau@inserm.fr
Burzykowski, T., Molenberghs, G., Buyse, M., Geys, H., and Renard, D. (2001). Validation of surrogate end points in multiple randomized clinical trials with failure time end points. Journal of the Royal Statistical Society: Series C (Applied Statistics) 50, 405422.
Buyse, M., Molenberghs, G., Burzykowski, T., Renard, D., and Geys, H. (2000). The validation of surrogate endpoints in metaanalyses of randomized experiments. Biostatistics 1, 4967
Sofeu, C. L., Emura, T., and Rondeau, V. (2019). Onestep validation method for surrogate endpoints using data from multiple randomized cancer clinical trials with failuretime endpoints. Statistics in Medicine 38, 29282942.
R. B. Nelsen. An introduction to Copulas. Springer, 2006
Prenen, L., Braekers, R., and Duchateau, L. (2017). Extending the archimedean copula methodology to model multivariate survival data grouped in clusters of variable size. Journal of the Royal Statistical Society: Series B (Statistical Methodology) 79, 483505.
Sofeu, C. L., Emura, T., and Rondeau, V. (2020). A joint frailtycopula model for metaanalytic
validation of failure time surrogate endpoints in clinical trials. Under review
jointSurrCopSimul
, summary.jointSurroPenal
, jointSurroPenal
, jointSurroPenalSimul
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16  ## Not run:
# Data from the advanced ovarian cancer randomized clinical trials.
data(dataOvarian)
joint.surro.Gumbel < jointSurroCopPenal(data = dataOvarian, int.method = 0,
n.knots = 8, maxit = 50, kappa.use = 4, nb.mc = 1000, typecopula = 2,
print.iter = F, scale = 1/365)
print(joint.surro.Gumbel)
joint.surro.Clayton < jointSurroCopPenal(data = dataOvarian, int.method = 0,
n.knots = 8, maxit = 50, kappa.use = 4, nb.mc = 1000, typecopula = 1,
print.iter = F, scale = 1/365)
print(joint.surro.Clayton)
## End(Not run)

Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.