{'fraction_non_alphanumeric': 0.0333442838370565, 'fraction_numerical': 0.01540735873850197, 'mean_word_length': 4.6143489487274065, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Fast Accelerated Failure Time Modeling for Case-Cohort Data

Semiparametric accelerated failure time (AFT) models directly relate the predicted failure times to covariates and are a useful alternative to models that work on the hazard function or the survival function. For case-cohort data, much less development has been done with AFT models. In addition to the missing covariates outside of the sub-cohort in controls, challenges from AFT model inferences with full cohort are retained. The regression parameter estimator is hard to compute because the most widely used rank-based estimating equations are not smooth. Further, its variance depends on the unspecified error distribution, and most methods rely on computationally intensive bootstrap to estimate it. We propose fast rank-based inference procedures for AFT models, applying recent methodological advances to the context of case-cohort data. Parameters are estimated with an induced smoothing approach that smooths the estimating functions and facilitates the numerical solution. Variance estimators are obtained through efficient resampling methods for nonsmooth estimating functions that avoids full blown bootstrap. Simulation studies suggest that the recommended procedure provides fast and valid inferences among several competing procedures. Application to a tumor study demonstrates the utility of the proposed method in routine data analysis.


Introduction
A case-cohort design (Prentice, 1986) is an effective and economical design which reduces the effort and cost of a full-scale cohort study. Such design originated to allow efficient analysis of studies where it is too expensive and time consuming to collect and analyze data on all subjects. Cases and controls refer to subjects who have and have not, respectively, developed the disease of interest by the end of the study period. A case-cohort design is typically composed of two steps. First, a subset called sub-cohort is randomly selected from the whole cohort regardless of their disease status. Second, the remaining cases in the cohort are added to the sub-cohort. Cases and controls refer to subjects who have and have not, respectively, developed the disease of interest by the end of the study period. Measurement on the main risk factors are taken only on subjects in the sub-cohort and the remaining cases outside of the sub-cohort. This leads to substantial reduction in the effort and cost of conducting large scale cohort studies, especially when the disease of interest is rare or the main risk factors are expensive to measure.
A semiparametric accelerated failure time (AFT) model is a log-linear model for the failure times with unspecified error distribution. It directly relates the failure time to covariates such that the effect of a covariate is to multiply the predicted failure time by a constant. For failure time data from case-cohort studies, most statistical methods have focused on semiparametric models that work on either the hazard function (Barlow, 1994;Kang and Cai, 2009;Kulich and Lin, 2000;Lin and Ying, 1993;Prentice, 1986;Self and Prentice, 1988;Sun et al., 2004;Therneau and Li, 1999), or the survival function (Chen, 2001a,b;Kong et al., 2004;Lu and Tsiatis, 2006). Parametric AFT models were considered by Kalbfleisch and Lawless (1988). Inferences about semiparametric AFT models for case-cohort data are much less developed, with only a few recent works (Kong and Cai, 2009;Nan et al., 2006;Yu, 2011;Yu et al., 2007).
Inferences for semiparametric AFT models have been difficult for not only case-cohort data but also for complete data. The most important estimator is the rank-based estimator motivated from inverting the weighted log-rank test (Prentice, 1978), with asymptotic properties rigorously studied (Tsiatis, 1990;Ying, 1993). Nevertheless, the estimator has not been as widely used as it should be due to lack of efficient and reliable computing algorithm to obtain both parameter estimates and their standard errors.
The parameter estimates are hard to compute because the most widely used rank-based estimating equations are not smooth. Recent works shed light on bringing AFT models into routine data analysis practice, including case-cohort studies. Jin et al. (2003) exploited that the rank-based estimating equation with Gehan's weight is the gradient of an objective function and obtained estimates by solving it with linear programming. This approach was adapted to case-cohort data by Kong and Cai (2009). Nevertheless, the optimization with linear programming is still computationally very demanding, especially for larger sample sizes. A more computing efficient approach for rank-based inference is the induced smoothing procedure of Brown and Wang (2007). This approach is an application of the general induced smoothing method of Brown and Wang (2005), where the discontinuous estimating equations are replaced with a smoothed version, whose solutions are asymptotically equivalent to those of the former. The smoothed estimating equations are differentiable, thus facilitates rapid numerical solution.
Direct estimation of the variance is difficult because it involves nonparametric estimation of the unspecified error distribution. Most existing methods rely on bootstrap which is very computing intensive. Jin et al. (2003) estimated the variance through a multiplier resampling method, which requires a large bootstrapping sample in order to obtain a reliable variance estimate. For case-cohort data, Kong and Cai (2009) adopted a specially designed bootstrap procedure (Wacholder et al., 1989). The demanding computing task in linear programming is amplified because it requires solving estimating equations for each bootstrap sample. Huang (2002) proposed an easy-to-compute variance estimator based on the asymptotic linearity property of the estimating equations. A decomposition matrix of the variance matrix is estimated by solving estimating equations, but the number of the estimating equations to solve is much smaller; it is just the dimension of the parameters. For general nonsmooth estimating functions, Zeng and Lin (2008) proposed a resampling strategy that does not require solving estimating equations or minimizing objective functions. Instead, it only involves evaluations of estimating functions and simple linear regression in estimating the slope matrix. The resulting variance estimators are computationally more efficient and stable than those from existing resampling methods.
In this article, we propose a fast rank-based inference procedure for semiparametric AFT models in the context of case-cohort studies. The parameters are estimated with an induced smoothing approach. Variance estimators are obtained through an efficient resampling methods for nonsmooth estimating functions that avoids full blown bootstrap. Of course, the methods also apply to full cohort data.
The rest of this article is organized as follows. Point estimation procedures based on smoothed estimating equations for case-cohort data when the sub-cohort is a simple random sample from the full cohort are proposed in Section 2. Four variance estimation procedures, one based on full multiplier bootstrap and three based on possibly multiplier bootstrapaided sandwich variance estimator, are proposed in Section 3. A large scale simulation study is reported in Section 4, comparing the performances of the variance estimator and their timings. The methods are applied to a tumor study with both case-cohort data and full cohort data in Section 5. A discussion concludes in Section 6.

Point Estimation
Let {T i , C i , X i }, i = 1, . . . , n, be n independent copies of {T, C, X}, where T i and C i are logtransformed failure time and log-transformed censoring time, X i is a p × 1 covariate vector, and given X, C and T are assumed to be independent. A semiparametric AFT model has the form where β is an unknown p × 1 vector of regression parameters, ǫ i 's are independent and identically distributed random variables with an unspecified distribution. It is also assumed that ǫ i 's are independent of X i . In a full cohort study, due to censoring, the observed data are (Y i , ∆ i , X i ), i = 1, . . . , n, where Y i = min(T i , C i ), ∆ i = I[T i < C i ], and I[·] is the indicator function. A rank based estimating equation with Gehan's weight is where e i (β) = Y i − X ⊤ i β. The root of (1) is consistent to the true parameter β 0 , and is asymptotically normal (Tsiatis, 1990). Despite these nice properties, even for the most promising method to date that solves it via linear programming (Jin et al., 2003), the computing burden increases drastically when bootstrapping is used to estimate the variance of the estimator.
For a case-cohort study, the covariate vector X i 's are not completely available for each individual. Measurement of some covariates is taken only on the subjects in the sub-cohort and cases outside the sub-cohort, and, thus, estimating function (1) cannot be evaluated. Using the observed data naively in (1) would lead to misleading results because the casecohort sample is biased -it includes all cases but only a fraction of controls. It is possible, however, to adjust the biases by incorporating a weight that depends on the selection scheme of case-cohort samples. Suppose we select a sub-cohort of sizeñ by simple random sampling without replacement from the whole cohort. Let ξ i be the sub-cohort indicator; ξ i = 1 if the ith observation is in the sub-cohort and ξ i = 0 otherwise. Let p = lim n→∞ p n , where p n =ñ/n is the sub-cohort inclusion probability. Under these assumptions, the desired case-cohort weight is The solution to (2),β n , remains to be consistent and asymptotically normal (Kong and Cai, 2009). For full cohort data, a computationally more efficient approach for rank-based inference with Gehan's weight is the induced smoothing procedure of Brown and Wang (2007). Such smoothing method leads to continuously differentiable estimating equations that can be solved with standard numerical methods. Let Z be a p-dimensional standard normal random vector. The estimating function U n (β) in (1) is replaced with E[U n (β + n −1/2 Z)], where the expectation is taken with respect to Z. This lead tõ where r 2 ij = n −1 (X i − X j ) ⊤ (X i − X j ) and Φ(·) denotes the standard normal cumulative distribution function. The solution to (3) is consistent to β 0 and has the same asymptotic distribution as the solution to (1) (Johnson and Strawderman, 2009).
For case-cohort data, we propose a smoothed version of (2) by adapting the idea of Brown and Wang (2007). Specifically, we replace U c n (β) with E[U c n (β + n −1/2 Z)] to obtain the induced smooth version of (2), The solutionβ n to (4) is a consistent estimator to β 0 and is asymptotically normal. Furthermore, the asymptotic distribution ofβ n is also the same as that ofβ n . These arguments can be justified similarly as those in Johnson and Strawderman (2009).

Variance Estimation
The asymptotic variance ofβ n is even harder to estimate for case-cohort data than for full cohort data because of the extra complexity caused by the data structure. The terms in the summation inŨ n (β) are not independent since the sub-cohort is drawn from the full cohort without replacement. We propose four variance estimators; one is fully resampling based while the other three use resampling to a component of the sandwich variance estimator.

Multiplier Bootstrap
The multiplier bootstrap estimator of Jin et al. (2003) is adapted to case-cohort data by inserting proper case-cohort weights, h i 's, in the multiplier bootstrap estimating equations. Let η i , i = 1, . . . , n, be independent and identically distributed positive random variables For a realization of (η 1 , . . . , η n ), the solution to (5) provides one draw ofβ n from its asymptotic distribution. By repeating this process a large number B times, the variance matrix of β n can be estimated directly by the sampling variance matrix of the bootstrap sample ofβ n . Since the asymptotic variance ofβ n is the same as that ofβ n , the covariance matrix ofβ n can also be estimated by (2) through multiplier bootstrap. This is, however, not recommended because it would need to solve a large number B nonsmooth estimating equations. As will be seen in our simulation study, even with the computationally more efficient smoothing estimating equations, the multiplier bootstrap approach can still be very time consuming, especially for larger sample sizes or more covariates.

Sandwich Estimator
To improve the computational efficiency, we consider alternative variance estimation procedures based on the sandwich form that avoid solving estimating equations repetitively. The asymptotic variances ofβ n andβ n are the same, both having a sandwich form. Under some regularity conditions (Zeng and Lin, 2008), uniformly in a neighborhood of β 0 , equation (2) can be expressed as where S i (β 0 ) is a zero-mean random vector, and A is asymptotic slope matrix of n −1/2Ũ c n (β 0 ). The analytical details of S i (β 0 ) for case-cohort data is presented in the Appendix. The asymptotic variance matrix of . Estimation of Σ involves estimating V and A by estimator V n and A n , respectively. The variance estimator then has the sandwich formΣ n = A −1 n V n (A −1 n ) ⊤ .

Estimation of V
Matrix V can be estimated either through a closed-form estimator or through bootstrapping the estimating equations. For case-cohort data, due to the correlated feature of ξ i 's in h i 's, V is different from its full cohort counterpart. There are two sources of variations contributing to V : variation due to the sampling of a full cohort (V 1 ) and variation due to the sampling of a sub-cohort within the full cohort (V 2 ). In particular, we have where V 2 vanishes if full cohort data are available.

Closed-form With explicit expressions for
andŜ i (β n ) is obtained by replacing unknown quantities in S i (β) with their sample counterparts.
Multiplier Bootstrap WhenŜ i (β n ) have complicated expressions, it is more convenient and perhaps more accurate to estimate V via bootstrap (Zeng and Lin, 2008). Because U c n andŨ c n have the same asymptotic distribution, we apply the multiplier bootstrap approach toŨ c n . Evaluation of (5) atβ n with each realization of (η 1 , . . . , η n ) provides one bootstrap replicate ofŨ c * n (β n ). With B replicates, we estimate V by the sample variance of the bootstrap sample ofŨ c * n (β n ). The bootstrap here is much less demanding than the full multiplier bootstrap above, because it only involves evaluations of estimating equations instead of solving them to obtain each bootstrap replicate.

Estimation of A
With V estimated by V n , we next propose three approaches to estimate the slope matrix A. Depending whether V n is based on closed-form or multiplier bootstrap, we will have two versions of estimator of Σ for each approach of slope matrix estimation.
Induced Smoothing WithŨ c n , the smoothed version of U c n , the slope matrix A can be estimated directly by The close-form expression of A n can be evaluated easily. The variance estimator then has the sandwich formΣ Smoothed Huang's (2002) Approach Huang (2002) avoided the difficulty in estimating the slope matrix of nonsmooth estimating equations by exploiting the asymptotic linearity of the estimating equations. Nevertheless, this approach still requires solving p nonsmooth estimating equations, whose convergence may be a problem. We adapt Huang's approach by replacing the p nonsmooth estimating equations with their smoothed versions. Let V n = L ⊤ n L n be the Cholesky decomposition of V n . Let q nj be the solution to the following estimating equations for γ, j = 1, . . . , p, n −1Ũ c n (γ) = n −1/2 l j , where l j is the jth column of L n . The solutions can be obtained with from general purpose nonlinear equation solvers; in our implementation we used R packages nleqslv (Hasselman, 2012) and BB (Varadhan and Gilbert, 2009). Let Q n be the matrix whose jth column is q nj −β n . Then Q ⊤ n Q n is an estimate of Σ. With the adaptation to smooth estimating equations, this approach has an advantage compared to the induced smoothing approach in that the closed-form derivative matrix is not required, and, hence, can be applied to more general nonsmooth estimating equations. Zeng and Lin (2008) proposed to estimate the slope matrix by regressing the perturbed estimating functions on the perturbations. Let Z b , b = 1, . . . , B, be B realizations of a p-dimensional standard normal random vector. For casecohort data, let U c nj be the jth component of U c n . We estimate the jth row of A, j = 1, . . . , p, by A nj , the least squares estimate of the regression coefficients when regressing n −1/2 U nj (β n + n −1/2 Z b ) on Z b , i = 1, . . . , n. The variance estimator also has the sandwich formΣ n = A −1 n V n (A −1 n ) ⊤ . This approach differs from the induced smoothing approach in that the slope matrix A is estimated via a resampling procedure that involves p least squares regressions, instead of taking the derivatives of a smooth function. It can be viewed as an empirical version of the induced smoothing approach.

Simulation
We conducted an extensive simulation study to assess the performance of the our point and variance estimators. Failure time T was generated from AFT model log(T ) = 2 + X 1 + X 2 + X 3 + ǫ, where X 1 was Bernoulli with rate 0.5, X 2 and X 3 were uncorrelated standard normal variables. Censoring time C was generated from unif(0, τ ) where τ was tuned to achieve desired censoring rate C p . The distribution of ǫ had three types: standard normal, standard logistic, or standard Gumbel, abbreviated by N, L, and G, respectively. The censoring rate C p had two levels, 90% and 97%, representing a mildly rare disease and a very rare disease, respectively. For the mildly rare disease, the full cohort size was set to be 1500 and the case-cohort size was set tom = 300 on average. For the very rare disease, the full cohort sizes were set to be 1500 and 3000, each with case-cohort sizes averaged atm ∈ {150, 300}. The sub-cohort sampling proportion p n was set to yield the desired average case-cohort size given censoring rate and full cohort size. For each viable combination, we generated 1000 datasets.
Given a dataset, point estimates of regression coefficients were obtained from both nonsmooth and smoothed estimating equations. The estimator from the nonsmooth version was obtained using linear programming (Jin et al., 2003), denoted by LP. The estimator from the induced smoothing approach with estimating equations (4) was obtained using R package nleqslv (Hasselman, 2012), denoted by IS. The two estimators are expected to be asymptotically the same, but with the IS estimator obtained much faster. Eight variance estimates were computed for the point estimate. The first two were full multiplier bootstrap estimates, denoted by MB, one based on the LP approach and the other based on the IS approach. The rest six were sandwich estimates constructed by combinations of three approaches to estimate A and two approaches to estimate V . We use abbreviations IS, SH, and ZL to denote the induced smoothing, smoothed Huang's, and Zeng and Lin's approach for A, respectively. We use abbreviations CF and MB to denote the closed-form estimate approach and the multiplier bootstrap approach for V , respectively.
Results for the mildly rare disease case with censoring percentage C p = 90%, full cohort size 1500, and average case-cohort sizem = 300 are summarized in Table 1. Both the LP and the IS estimators appear to be virtually unbiased. In fact, they agreed with each other closely on a 45 degree line (not shown). Consequently, their empirical standard errors agreed with each other, and their bootstrap based standard errors agreed with each other. The bootstrap standard errors and the empirical standard errors match closely, suggesting that the bootstrap variance estimators provide good estimation of the empirical variantion. The other six standard errors based on sandwich variance estimators agreed quite well with the empirical standard errors too. The associated 95% confidence intervals based on all eight standard errors had empirical coverage percentages reasonably close to the nominal level. These observations were invariant to the error distributions. Table 2 summarizes the results for the very rare disease case with censoring rate 97% and full cohort size 3000. The results for full cohort size 1500 were similar and not reported. The two point estimates, their empirical standard errors, and their average bootstrap standard errors still agree with each other. The bootstrap standard errors for case-cohort size 150, however, are underestimating the true variation, and as a result, the 95% confidence intervals had coverage percentage smaller than the nominal level. Not surprisingly, the six sandwich variance estimators performed no better than the two multiplier bootstrap variance estimators. When the case-cohort size was increased to 300, all variance estimators performed reasonably well in estimating the true variation and the coverage percentage was reasonably close to the nominal level. Among all the sandwich variance estimators, the IS-MB and ZL-MB approaches seem to provide confidence intervals with the best coverage percentage. The SH-MB approach is slight inferior, which might be explained by the fact that this approach has two layers of approximation -one from asymptotic linear approximation and the other     from induced smoothing. Of more interest is Table 3, which summarizes the timing results in seconds averaged from 1000 replicates for both point estimation and variance estimation on a 2GHz linux machine. For point estimation with full cohort size 1500 and censoring percentage C p = 0.90%, the IS approach was up to 4.5 times as fast as the LP approach (with normal error distribution). The multipler bootstrap variance estimation with the IS approach was up to 3.6 times as fast as the LP approach (again with normal error). Nevertheless, the multiplier bootstrap IS approach still needed about 9 minutes on average to obtain a variance estimator. All sandwich variance estimators are strikingly much faster, especially with the closed-form approach: ZL-CF approach took about 2 seconds on average; the IS-CF and SH-CF approaches took about 3 seconds on average. For each sandwich variance estimator, the version with CF estimation of V is over 5 times faster than the version with MB estimation of V . Using the LP approach as benchmark, the IS-CF and SH-CF estimators is 695 times faster and the ZL-CF estimators 1003 times faster. Since the performance of all variance estimators are similar for this setting, the IS-CF, SH-CF and ZL-CF approaches are obviously preferred for this setup with a mildly rare event.
The timing results for full cohort size 3000 and censoring percentage C p = 0.97% follow a similar pattern. Compare to case C p = 0.90%, time for point estimation is shorter because the number of cases decreases in the case C p = 0.97% even when the average case-cohort size were both at 300. Sandwich variance estimators with CF estimation of V is up to 8 times faster than with those with MB estimation of V , at the expense of slightly worse performance in coverage percentage. The IS-MB and ZL-MB approaches yield the most reliable variance estimates but IS-MB is slightly faster than faster than ZL-MB. As the average case-cohort size doubles, the computing time of the sandwich variance estimates with MB estimation for V appear to double accordingly, in contrast to those with CF estimation for V , which do not necessarily double linearly. In summary, based on the performance and speed, our recommended variance estimator is the IS-MB estimator.

National Wilm's Tumor Study
We demonstrate the performance of our proposed methods with an application to the cohort study conducted by the National Wilm's Tumor Study Group (NWTSG) (D'Angio et al., 1989;Green et al., 1998). Wilm's tumor is a rare kidney cancer in young children. The interest of the study was to assess the relationship between the tumor histology and the outcome, time to tumor relapse. Tumor histology can be classified into two categories, favorable or unfavorable, depending on the cell type. The central histological diagnosis was made by an individual pathologist at the central pathology center, which was believed to be more accurate than a local diagnosis yet more expensive to measure and required more efforts to obtain. Although in the full version of the data, the central histology measurement was available for all the cohort members, it was only available for a case-cohort sample in the case-cohort version. We take advantage of the full version in this example. Other covariates that were available for all cohort members were patient age, disease stage and study group. According to the staging system employed by NWTSG, four stages (I -IV) of Wilms' tumors, with Stage IV as the latest stage, indicated the spread of the tumor. Each subject came from one of the two study groups, NWTSG-3 and NWTSG-4. The case-cohort version of the data was analyzed with Cox models (Breslow et al., 2009;Kulich and Lin, 2004) and additive hazards models (Kulich and Lin, 2000), respectively.
There were a total of 4028 subjects in the full cohort. Among them, 571 were cases who experienced the relapse of tumor -a censoring rate of about 86%. We considered an AFT model for the time to relapse with the following covariates: central histology measurement (1 = favorable, 0 = unfavorable), age (measure in year) at diagnosis, three tumor stages indicators (Stage I as reference) and a study group indicator (NWTSG-3 as reference). The case-cohort version of the data had 668 patients selected as sub-cohort sample and the total case-cohort sample size was 1154. To take advantage of availability of full cohort data, we drew 1000 new sub-cohort samples with size 668 and formed a case-cohort by including the remaining cases for each replicate. We then averaged these estimates and estimated standard errors from the 1000 replicates of the case-cohort analysis.
The results of the average from 1000 replicates of case-cohort analyses are summarized in Table 4. Due to its poor timing performance, the LP approach was not considered. Since the MB standard error is considered to reflect the true variation quite well from the simulation study, we are interested in how close the various sandwich standard errors to the MB standard error. For all three sandwich estimators, the CF versions systematically underestimate noticeably, although the underestimation is less severe in the IS and ZL approach than in the SH approach. The MB versions of the sandwich estimates appear to agree with the MB standard error closely, and again the agreement appears to be better for the IS and ZL approach than for the SH approach. In particular, the SH standard error for the age effect is Based on the IS-MB standard errors, the coefficients of central histological diagnosis, age, and all three stage indicators were found to be significantly different from zero with p-values 0.000, 0.009, 0.019, 0.011 and 0.002, respectively. No significant difference was found between the two study groups. In terms of timing, the MB-IS standard error took over a hour whereas the MB-based sandwich estimates only took 14-15 seconds on average. For comparison purpose, we also analyzed the full cohort data with the same approaches and reported the results in Table 4. Point estimates are close to these in case-cohort analysis, with their standard errors taken into consideration. All the standard errors decrease compared to the case-cohort analyses, which is expected as full information became available for all covariates. The best sandwich variance estimators are still IS-MB and ZL-MB, both closely approximates the full blown MB standard error. With full cohort size 4028 and censoring rate 86%, the IS point estimates took 4 and a half minutes, the MB variance estimation took 35.36 hours, while IS-MB only took 7 and a half minutes.

Discussion
In AFT modeling of case-cohort data, both point estimation and variance estimation are challenging with the nonsmooth estimating equations. Resampling methods are commonly used to estimate the variance, which are time consuming even with a computationally efficient point estimator such as our induced smoothing approach with rank-based estimating equations. We have proposed six sandwich variance estimators and compared their performances with the bootstrap variance estimator in numerical studies. The IS-MB and ZL-MB approaches were found to provide good approximation to the true variation and are computationally very efficient. All the methods are implemented in an R package aftgee (Chiou et al., 2012). The package had the potential to bring AFT modeling of case-cohort data into routine analysis.
The IS approach was built on Gehan's weight for rank-based estimating equations, in which case closed-form expectations of the perturbed estimating equations are available. Alternative weights such as the logrank weight are possible, though the computation is less straightforward than that for Gehan's weight. Incorporating a general, possibly optimized weight in the IS approach merits further investigation for both full cohort and case-cohort data. The estimates from Gehan's weight always serve as a good initial value in numerical equations solving.
Some extensions of the proposed methods are worth considering. Auxiliary covariates are often available for the entire cohort, which can be used to construct strata for subcohort members selection. The resulting estimators from the stratified case-cohort design has been shown to be more efficient than their traditional case-cohort counterpart for the Cox model (Kulich and Lin, 2004) and the additive hazards model (Kulich and Lin, 2000). An extension of the proposed methods with the AFT model to a stratified case-cohort design may lead to efficiency improvement too. When more than one diseases are considered in a casecohort design, a multivariate extension will be needed, and a possible dependence among the multivariate failure times needs to be taken into account. For the Cox model, Kang and Cai