Search

The propensity methodology is widely used in medical research to compare different treatments in designs with a nonrandomized treatment allocation. The inverse probability weighted (IPW) estimators are a primary tool for estimating the average treatment effect but the large variance of these estimators is often a significant concern for their reliable use in practice. Inspired by Rao-Blackwellization, this paper proposes a method to smooth an IPW estimator by replacing the weights in the original estimator by their mean over a distribution of the potential treatment assignment. In our simulation study, the smoothed IPW estimator achieves a substantial variance reduction over its original version with only a small increased bias, for example two-to-sevenfold variance reduction for the three IPW estimators in Lunceford and Davidian [Statistics in Medicine, 23(19), 2937-2960]. In addition, our proposed smoothing can also be applied to the locally efficient and doubly robust estimator for added protection against model misspecification. An implementation in R is provided.

Electronic health records (EHRs) have become a platform for data-driven granular-level surveillance in recent years. In this paper, we make use of EHRs for early prevention of childhood obesity. The proposed method simultaneously provides smooth disease mapping and outlier information for obesity prevalence, which are useful for raising public awareness and facilitating targeted intervention. More precisely, we consider a penalized multilevel generalized linear model. We decompose regional contribution into smooth and sparse signals, which are automatically identified by a combination of fusion and sparse penalties imposed on the likelihood function. In addition, we weigh the proposed likelihood to account for the missingness and potential non-representativeness arising from the EHR data. We develop a novel alternating minimization algorithm, which is computationally efficient, easy to implement, and guarantees convergence. Simulation studies demonstrate superior performance of the proposed method. Finally, we apply our method to the University of Wisconsin Population Health Information Exchange database. This article is protected by copyright. All rights reserved.

The United Nations' Sustainable Development Goal 3.2 aims to reduce under-5 child mortality to 25 deaths per 1,000 live births by 2030. Child mortality tends to be concentrated in developing regions where much of the information needed to assess achievement of this goal comes from surveys and censuses. In both, women are asked about their birth histories, but with varying degrees of detail. Full birth history (FBH) data contain the reported dates of births and deaths of every surveyed mother's children. In contrast, summary birth history (SBH) data contain only the total number of children born and total number of children who died for each mother. Specialized methods are needed to accommodate this type of data into analyses of child mortality trends. We develop a data augmentation scheme within a Bayesian framework where for SBH data, birth and death dates are introduced as auxiliary variables. Since we specify a full probability model for the data, many of the well-known biases that exist in this data can be accommodated, along with space-time smoothing on the underlying mortality rates. We illustrate our approach in a simulation, showing that uncertainty is reduced when incorporating SBH data over simply analyzing all available FBH data. We also apply our approach to data in the Central region of Malawi. We compare with the well-known Brass method., 50 pages, 18 figures

We consider the proportional hazards model in which the covariates include the discretized categories of a continuous time-dependent exposure variable measured with error. Naively ignoring the measurement error in the analysis may cause biased estimation and erroneous inference. Although various approaches have been proposed to deal with measurement error when the hazard depends linearly on the time-dependent variable, it has not yet been investigated how to correct when the hazard depends on the discretized categories of the time-dependent variable. To fill this gap in the literature, we propose a smoothed corrected score approach based on approximation of the discretized categories after smoothing the indicator function. The consistency and asymptotic normality of the proposed estimator are established. The observation times of the time-dependent variable are allowed to be informative. For comparison, we also extend to this setting two approximate approaches, the regression calibration and the risk-set regression calibration. The methods are assessed by simulation studies and by application to data from an HIV clinical trial.

Drawing inferences for high-dimensional models is challenging as regular asymptotic theories are not applicable. This article proposes a new framework of simultaneous estimation and inferences for high-dimensional linear models. By smoothing over partial regression estimates based on a given variable selection scheme, we reduce the problem to low-dimensional least squares estimations. The procedure, termed as Selection-assisted Partial Regression and Smoothing (SPARES), utilizes data splitting along with variable selection and partial regression. We show that the SPARES estimator is asymptotically unbiased and normal, and derive its variance via a nonparametric delta method. The utility of the procedure is evaluated under various simulation scenarios and via comparisons with the de-biased LASSO estimators, a major competitor. We apply the method to analyze two genomic datasets and obtain biologically meaningful results.

Prediction modeling for clinical decision making is of great importance and needed to be updated frequently with the changes of patient population and clinical practice. Existing methods are either done in an ad hoc fashion, such as model recalibration or focus on studying the relationship between predictors and outcome and less so for the purpose of prediction. In this article, we propose a dynamic logistic state space model to continuously update the parameters whenever new information becomes available. The proposed model allows for both time-varying and time-invariant coefficients. The varying coefficients are modeled using smoothing splines to account for their smooth trends over time. The smoothing parameters are objectively chosen by maximum likelihood. The model is updated using batch data accumulated at prespecified time intervals, which allows for better approximation of the underlying binomial density function. In the simulation, we show that the new model has significantly higher prediction accuracy compared to existing methods. We apply the method to predict 1 year survival after lung transplantation using the United Network for Organ Sharing data.

Author(s): Bose, Maitreyee; Hodges, James S; Banerjee, Sudipto | Abstract: Gaussian processes (GPs) are widely used as distributions of random effects in linear mixed models, which are fitusing the restricted likelihood or the closely related Bayesian analysis. This article addresses two problems. First, we proposetools for understanding how data determine estimates in these models, using a spectral basis approximation to the GP underwhich the restricted likelihood is formally identical to the likelihood for a gamma-errors GLM with identity link. Second,to examine the data’s support for a covariate and to understand how adding that covariate moves variation in the outcomey out of the GP and error parts of the fit, we apply a linear-model diagnostic, the added variable plot (AVP), both to theoriginal observations and to projections of the data onto the spectral basis functions. The spectral- and observation-domainAVPs estimate the same coefficient for a covariate but emphasize low- and high-frequency data features respectively and thushighlight the covariate’s effect on the GP and error parts of the fit, respectively. The spectral approximation applies to dataobserved on a regular grid; for data observed at irregular locations, we propose smoothing the data to a grid before applyingour methods. The methods are illustrated using the forest-biomass data of Finley et al. (2008).

Multivariate current-status data are frequently encountered in biomedical and public health studies. Semiparametric regression models have been extensively studied for univariate current-status data, but most existing estimation procedures are computationally intensive, involving either penalization or smoothing techniques. It becomes more challenging for the analysis of multivariate current-status data. In this article, we study the maximum likelihood estimations for univariate and bivariate current-status data under the semiparametric probit regression models. We present a simple computational procedure combining the expectation-maximization algorithm with the pool-adjacent-violators algorithm for solving the monotone constraint on the baseline function. Asymptotic properties of the maximum likelihood estimators are investigated, including the calculation of the explicit information bound for univariate current-status data, as well as the asymptotic consistency and convergence rate for bivariate current-status data. Extensive simulation studies showed that the proposed computational procedures performed well under small or moderate sample sizes. We demonstrate the estimation procedure with two real data examples in the areas of diabetic and HIV research.

We consider the optimization of smoothing parameters and variance components in models with a regular log likelihood subject to quadratic penalization of the model coefficients, via a generalization of the method of Fellner (1986) and Schall (1991). In particular: (i) we generalize the original method to the case of penalties that are linear in several smoothing parameters, thereby covering the important cases of tensor product and adaptive smoothers; (ii) we show why the method's steps increase the restricted marginal likelihood of the model, that it tends to converge faster than the EM algorithm, or obvious accelerations of this, and investigate its relation to Newton optimization; (iii) we generalize the method to any Fisher regular likelihood. The method represents a considerable simplification over existing methods of estimating smoothing parameters in the context of regular likelihoods, without sacrificing generality: for example, it is only necessary to compute with the same first and second derivatives of the log-likelihood required for coefficient estimation, and not with the third or fourth order derivatives required by alternative approaches. Examples are provided which would have been impossible or impractical with pre-existing Fellner-Schall methods, along with an example of a Tweedie location, scale and shape model which would be a challenge for alternative methods, and a sparse additive modeling example where the method facilitates computational efficiency gains of several orders of magnitude. This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

Summary This article investigates a generalized semiparametric varying-coefficient model for longitudinal data that can flexibly model three types of covariate effects: time-constant effects, time-varying effects, and covariate-varying effects. Different link functions can be selected to provide a rich family of models for longitudinal data. The model assumes that the time-varying effects are unspecified functions of time and the covariate-varying effects are parametric functions of an exposure variable specified up to a finite number of unknown parameters. The estimation procedure is developed using local linear smoothing and profile weighted least squares estimation techniques. Hypothesis testing procedures are developed to test the parametric functions of the covariate-varying effects. The asymptotic distributions of the proposed estimators are established. A working formula for bandwidth selection is discussed and examined through simulations. Our simulation study shows that the proposed methods have satisfactory finite sample performance. The proposed methods are applied to the ACTG 244 clinical trial of HIV infected patients being treated with Zidovudine to examine the effects of antiretroviral treatment switching before and after HIV develops the T215Y/F drug resistance mutation. Our analysis shows benefits of treatment switching to the combination therapies as compared to continuing with ZDV monotherapy before and after developing the 215-mutation.

Summary Many diseases arise due to exposure to one of multiple possible pathogens. We consider the situation in which disease counts are available over time from a study region, along with a measure of clinical disease severity, for example, mild or severe. In addition, we suppose a subset of the cases are lab tested in order to determine the pathogen responsible for disease. In such a context, we focus interest on modeling the probabilities of disease incidence given pathogen type. The time course of these probabilities is of great interest as is the association with time-varying covariates such as meteorological variables. In this set up, a natural Bayesian approach would be based on imputation of the unsampled pathogen information using Markov Chain Monte Carlo but this is computationally challenging. We describe a practical approach to inference that is easy to implement. We use an empirical Bayes procedure in a first step to estimate summary statistics. We then treat these summary statistics as the observed data and develop a Bayesian generalized additive model. We analyze data on hand, foot, and mouth disease (HFMD) in China in which there are two pathogens of primary interest, enterovirus 71 (EV71) and Coxackie A16 (CA16). We find that both EV71 and CA16 are associated with temperature, relative humidity, and wind speed, with reasonably similar functional forms for both pathogens. The important issue of confounding by time is modeled using a penalized B-spline model with a random effects representation. The level of smoothing is addressed by a careful choice of the prior on the tuning variance.

Array-based CGH experiments are designed to detect genomic aberrations or regions of DNA copy-number variation that are associated with an outcome, typically a state of disease. Most of the existing statistical methods target on detecting DNA copy number variations in a single sample or array. We focus on the detection of group effect variation, through simultaneous study of multiple samples from multiple groups. Rather than using direct segmentation or smoothing techniques, as commonly seen in existing detection methods, we develop a sequential model selection procedure that is guided by a modified Bayesian information criterion. This approach improves detection accuracy by accumulatively utilizing information across contiguous clones, and has computational advantage over the existing popular detection methods. Our empirical investigation suggests that the performance of the proposed method is superior to that of the existing detection methods, in particular, in detecting small segments or separating neighboring segments with differential degrees of copy-number variation.

Motivated by ultrahigh-dimensional biomarkers screening studies, we propose a model-free screening approach tailored to censored lifetime outcomes. Our proposal is built upon the introduction of a new measure, survival impact index (SII). By its design, SII sensibly captures the overall influence of a covariate on the outcome distribution, and can be estimated with familiar nonparametric procedures that do not require smoothing and are readily adaptable to handle lifetime outcomes under various censoring and truncation mechanisms. We provide large sample distributional results that facilitate the inference on SII in classical multivariate settings. More importantly, we investigate SII as an effective screener for ultrahigh-dimensional data, not relying on rigid regression model assumptions for real applications. We establish the sure screening property of the proposed SII-based screener. Extensive numerical studies are carried out to assess the performance of our method compared with other existing screening methods. A lung cancer microarray data is analyzed to demonstrate the practical utility of our proposals.

Gene regulatory networks represent the regulatory relationships between genes and their products and are important for exploring and defining the underlying biological processes of cellular systems. We develop a novel framework to recover the structure of nonlinear gene regulatory networks using semiparametric spline-based directed acyclic graphical models. Our use of splines allows the model to have both flexibility in capturing nonlinear dependencies as well as control of overfitting via shrinkage, using mixed model representations of penalized splines. We propose a novel discrete mixture prior on the smoothing parameter of the splines that allows for simultaneous selection of both linear and nonlinear functional relationships as well as inducing sparsity in the edge selection. Using simulation studies, we demonstrate the superior performance of our methods in comparison with several existing approaches in terms of network reconstruction and functional selection. We apply our methods to a gene expression dataset in glioblastoma multiforme, which reveals several interesting and biologically relevant nonlinear relationships.

Summary Linear regressions are commonly used to calibrate the signal measurements in proteomic analysis by mass spectrometry. However, with or without a monotone (e.g., log) transformation, data from such functional proteomic experiments are not necessarily linear or even monotone functions of protein (or peptide) concentration except over a very restricted range. A computationally efficient spline procedure improves upon linear regression. However, mass spectrometry data are not necessarily homoscedastic; more often the variation of measured concentrations increases disproportionately near the boundaries of the instruments measurement capability (dynamic range), that is, the upper and lower limits of quantitation. These calibration difficulties exist with other applications of mass spectrometry as well as with other broad-scale calibrations. Therefore the method proposed here uses a functional data approach to define the calibration curve and also the limits of quantitation under the two assumptions: (i) that the variance is a bounded, convex function of concentration; and (ii) that the calibration curve itself is monotone at least between the limits of quantitation, but not necessarily outside these limits. Within this paradigm, the limit of detection, where the signal is definitely present but not measurable with any accuracy, is also defined. An iterative approach draws on existing smoothing methods to account simultaneously for both restrictions and is shown to achieve the global optimal convergence rate under weak conditions. This approach can also be implemented when convexity is replaced by other (bounded) restrictions. Examples from Addona et al. (2009, Nature Biotechnology 27, 663–641) both motivate and illustrate the effectiveness of this functional data methodology when compared with the simpler linear regressions and spline techniques.

Summary. We propose a method to test the correlation of two random fields when they are both spatially autocorrelated. In this scenario, the assumption of independence for the pair of observations in the standard test does not hold, and as a result we reject in many cases where there is no effect (the precision of the null distribution is overestimated). Our method recovers the null distribution taking into account the autocorrelation. It uses Monte-Carlo methods, and focuses on permuting, and then smoothing and scaling one of the variables to destroy the correlation with the other, while maintaining at the same time the initial autocorrelation. With this simulation model, any test based on the independence of two (or more) random fields can be constructed. This research was motivated by a project in biodiversity and conservation in the Biology Department at Stanford University.

Functional data are defined as realizations of random functions (mostly smooth functions) varying over a continuum, which are usually collected with measurement errors on discretized grids. In order to accurately smooth noisy functional observations and deal with the issue of high-dimensional observation grids, we propose a novel Bayesian method based on the Bayesian hierarchical model with a Gaussian-Wishart process prior and basis function representations. We first derive an induced model for the basis-function coefficients of the functional data, and then use this model to conduct posterior inference through Markov chain Monte Carlo. Compared to the standard Bayesian inference that suffers serious computational burden and unstableness for analyzing high-dimensional functional data, our method greatly improves the computational scalability and stability, while inheriting the advantage of simultaneously smoothing raw observations and estimating the mean-covariance functions in a nonparametric way. In addition, our method can naturally handle functional data observed on random or uncommon grids. Simulation and real studies demonstrate that our method produces similar results as the standard Bayesian inference with low-dimensional common grids, while efficiently smoothing and estimating functional data with random and high-dimensional observation grids where the standard Bayesian inference fails. In conclusion, our method can efficiently smooth and estimate high-dimensional functional data, providing one way to resolve the curse of dimensionality for Bayesian functional data analysis with Gaussian-Wishart processes., Comment: Under review

Given a large number of t-statistics, we consider the problem of approximating the distribution of noncentrality parameters (NCPs) by a continuous density. This problem is closely related to the control of false discovery rates (FDR) in massive hypothesis testing applications, e.g., microarray gene expression analysis. Our methodology is similar to, but improves upon, the existing approach by Ruppert, Nettleton, and Hwang (2007, Biometrics, 63, 483-495). We provide parametric, nonparametric, and semiparametric estimators for the distribution of NCPs, as well as estimates of the FDR and local FDR. In the parametric situation, we assume that the NCPs follow a distribution that leads to an analytically available marginal distribution for the test statistics. In the nonparametric situation, we use convex combinations of basis density functions to estimate the density of the NCPs. A sequential quadratic programming procedure is developed to maximize the penalized likelihood. The smoothing parameter is selected with the approximate network information criterion. A semiparametric estimator is also developed to combine both parametric and nonparametric fits. Simulations show that, under a variety of situations, our density estimates are closer to the underlying truth and our FDR estimates are improved compared with alternative methods. Data-based simulations and the analyses of two microarray datasets are used to evaluate the performance in realistic situations.

Summary We consider frailty models with additive semiparametric covariate effects for clustered failure time data. We propose a doubly penalized partial likelihood (DPPL) procedure to estimate the nonparametric functions using smoothing splines. We show that the DPPL estimators could be obtained from fitting an augmented working frailty model with parametric covariate effects, whereas the nonparametric functions being estimated as linear combinations of fixed and random effects, and the smoothing parameters being estimated as extra variance components. This approach allows us to conveniently estimate all model components within a unified frailty model framework. We evaluate the finite sample performance of the proposed method via a simulation study, and apply the method to analyze data from a study of sexually transmitted infections (STI).

Summary Applied scientists often like to use ordinary differential equations (ODEs) to model complex dynamic processes that arise in biology, engineering, medicine, and many other areas. It is interesting but challenging to estimate ODE parameters from noisy data, especially when the data have some outliers. We propose a robust method to address this problem. The dynamic process is represented with a nonparametric function, which is a linear combination of basis functions. The nonparametric function is estimated by a robust penalized smoothing method. The penalty term is defined with the parametric ODE model, which controls the roughness of the nonparametric function and maintains the fidelity of the nonparametric function to the ODE model. The basis coefficients and ODE parameters are estimated in two nested levels of optimization. The coefficient estimates are treated as an implicit function of ODE parameters, which enables one to derive the analytic gradients for optimization using the implicit function theorem. Simulation studies show that the robust method gives satisfactory estimates for the ODE parameters from noisy data with outliers. The robust method is demonstrated by estimating a predator–prey ODE model from real ecological data.

In this work we propose a fully Bayesian semiparametric method to estimate the intensity of an inhomogeneous spatial point process. The basic idea is to first convert intensity estimation into a Poisson regression setting via binning data points on a regular grid, and then model the log intensity semiparametrically using an adaptive version of Gaussian Markov random fields to smooth the corresponding counts. The inference is carried by an efficient Markov chain Monte Carlo simulation algorithm. Compared to existing methods for intensity estimation, for example, parametric modeling and kernel smoothing, the proposed estimator not only provides inference regarding the dependence of the intensity function on possible covariates, but also uses information from the data to adaptively determine the amount of smoothing at the local level. The effectiveness of using our method is demonstrated through simulation studies and an application to a rainforest dataset.

Summary Oysters from the Pacific Northwest coast of British Columbia, Canada, contain high levels of cadmium, in some cases exceeding some international food safety guidelines. A primary goal of this article is the investigation of the spatial and temporal variation in cadmium concentrations for oysters sampled from coastal British Columbia. Such information is important so that recommendations can be made as to where and when oysters can be cultured such that accumulation of cadmium within these oysters is minimized. Some modern statistical methods are applied to achieve this goal, including monotone spline smoothing, functional principal component analysis, and semi-parametric additive modeling. Oyster growth rates are estimated as the first derivatives of the monotone smoothing growth curves. Some important patterns in cadmium accumulation by oysters are observed. For example, most inland regions tend to have a higher level of cadmium concentration than most coastal regions, so more caution needs to be taken for shellfish aquaculture practices occurring in the inland regions. The semi-parametric additive modeling shows that oyster cadmium concentration decreases with oyster length, and oysters sampled at 7 m have higher average cadmium concentration than those sampled at 1 m.

In this article we develop a nonparametric estimation procedure for the varying coefficient model when the within-subject covariance is unknown. Extending the idea of iterative reweighted least squares to the functional setting, we iterate between estimating the coefficients conditional on the covariance and estimating the functional covariance conditional on the coefficients. Smoothing splines for correlated errors are used to estimate the functional coefficients with smoothing parameters selected via the generalized maximum likelihood. The covariance is nonparametrically estimated using a penalized estimator with smoothing parameters chosen via a Kullback-Leibler criterion. Empirical properties of the proposed method are demonstrated in simulations and the method is applied to the data collected from an ovarian tumor study in mice to analyze the effects of different chemotherapy treatments on the volumes of two classes of tumors.

Attachment loss (AL), the distance down a tooth's root that is no longer attached to surrounding bone by periodontal ligament, is a common measure of periodontal disease. In this article, we develop a spatiotemporal model to monitor the progression of AL. Our model is an extension of the conditionally autoregressive (CAR) prior, which spatially smooths estimates toward their neighbors. However, because AL often exhibits a burst of large values in space and time, we develop a nonstationary spatiotemporal CAR model that allows the degree of spatial and temporal smoothing to vary in different regions of the mouth. To do this, we assign each AL measurement site its own set of variance parameters and spatially smooth the variances with spatial priors. We propose a heuristic to measure the complexity of the site-specific variances, and use it to select priors that ensure parameters in the model are well identified. In data from a clinical trial, this model improves the fit compared to the usual dynamic CAR model for 90 of 99 patients' AL measurements.

Smoothing over a domain with irregular boundaries or interior gaps and holes is addressed. Consider the problem of estimating mercury in sediment concentrations in the estuarine waters in New Hampshire. A modified version of low-rank thin plate splines (LTPS) is introduced where the geodesic distance is applied to evaluate dissimilarity of any two data observations: loosely speaking, distances between locations are not measured as the crow flies, but as the fish swims. The method is compared with competing smoothing techniques, LTPS, and finite element L-splines.

AFLP is a DNA fingerprinting technique frequently used in plant and animal sciences. A drawback of the technique is the occurrence of multiple DNA fragments of the same length in a single AFLP lane, which we name a collision. In this article we quantify the problem. The well-known birthday problem plays a role. Calculation of collision probabilities requires a fragment length distribution (fld). We discuss three ways to estimate the fld: based on theoretical considerations, on in-silico determination using DNA sequence data from Arabidopsis thaliana, or on direct estimation from AFLP data. In the latter case we use a generalized linear model with monotone smoothing of the fragment length probabilities. Collision probabilities are calculated from two perspectives, assuming known fragment counts and assuming known band counts. We compare results for a number of fld's, ranging from uniform to highly skewed. The conclusion is that collisions occur often, with higher probabilities for higher numbers of bands, for more skewed distributions, and, to a lesser extent, for smaller scoring ranges. For a typical plant genome an AFLP with 19 bands is likely to contain the first collision. Practical implications of collisions are discussed. AFLP examples from lettuce and chicory are used for illustration.

A linear mixed model with a smooth random effects density is proposed. A similar approach to P-spline smoothing of Eilers and Marx (1996, Statistical Science 11, 89-121) is applied to yield a more flexible estimate of the random effects density. Our approach differs from theirs in that the B-spline basis functions are replaced by approximating Gaussian densities. Fitting the model involves maximizing a penalized marginal likelihood. The best penalty parameters minimize Akaike's Information Criterion employing Gray's (1992, Journal of the American Statistical Association 87, 942-951) results. Although our method is applicable to any dimensions of the random effects structure, in this article the two-dimensional case is explored. Our methodology is conceptually simple, and it is relatively easy to fit in practice and is applied to the cholesterol data first analyzed by Zhang and Davidian (2001, Biometrics 57, 795-802). A simulation study shows that our approach yields almost unbiased estimates of the regression and the smoothing parameters in small sample settings. Consistency of the estimates is shown in a particular case.

Summary. The fish abundance index over an ocean region is defined here to be the integral of expected catch per unit effort (CPUE), approximated by the sum of expected CPUE over grid squares. When trawl surveys are done within grid squares selected according to a probability sampling design, several other sources of variation such as the fish population dynamics and the catching process are also involved. In such situations model-assisted methods for estimating abundance, assessed under both design and model perspectives, have some advantages over purely design-based methods such as the Horvitz–Thompson (HT) estimator or purely model-based prediction approaches. This article develops model-assisted empirical likelihood (EL) methods via loglinear regression and nonparametric smoothing. The methods are applied to grid surveys of the Grand Bank region carried out annually by Fishery Products International from 1996 through 2002. The HT and EL methods produce similar point estimates of abundance indices. Simulation results, however, indicate that the EL estimator under local linear smoothing is associated with smaller standard errors.

Cox's regression model is the standard regression tool for survival analysis in most applications. Often, however, the model only provides a rough summary of the effect of some covariates. Therefore, if the aim is to give a detailed description of covariate effects and to consequently calculate predicted probabilities, more flexible models are needed. In another article, Scheike and Zhang (2002, Scandinavian Journal of Statistics 29, 75-88), we suggested a flexible extension of Cox's regression model, which aimed at extending the Cox model only for those covariates where additional flexibility are needed. One important advantage of the suggested approach is that even though covariates are allowed a nonparametric effect, the hassle and difficulty of finding smoothing parameters are not needed. We show how the extended model also leads to simple formulae for predicted probabilities and their standard errors, for example, in the competing risk framework.

This article introduces a new approach for estimating the hazard function for possibly interval- and right-censored survival data. We weakly parameterize the log-hazard function with a piecewise-linear spline and provide a smoothed estimate of the hazard function by maximizing the penalized likelihood through a mixed model-based approach. We also provide a method to estimate the amount of smoothing from the data. We illustrate our approach with two well-known interval-censored data sets. Extensive numerical studies are conducted to evaluate the efficacy of the new procedure.

The model that specifies that cancer incidence, I, is the convolution of the preclinical incidence, g, and the density of time in the preclinical phase, f, has frequently been utilized to model data from cancer screening trials and to estimate such quantities as sojourn time, lead time, and sensitivity. When this model is fit to the above data, the parameters of f as well as the parameter(s) governing screening sensitivity must be estimated. Previously, g was either assumed to be equal to clinical incidence or assumed to be a constant or exponential function that also had to be estimated. Here we assume that the underlying incidence, I, in the study population (in the absence of screening) is known. With I known, g then becomes a function of f, which can be solved for using (numerical) deconvolution, thus eliminating the need to estimate g or make assumptions about it. Since numerical deconvolution procedures may be highly unstable, however, we incorporate a smoothing procedure that produces a realistic g function while still closely reproducing the original incidence function I upon convolution with f. We have also added the concept of competing mortality to the convolution model. This, along with the realistic preclinical incidence function described above, results in more accurate estimates of sojourn time and lead time and allows for estimation of quantities related to overdiagnosis, which we define here.

Kaplan-Meier curves provide an effective means of presenting the distributional pattern in a sample of survival data. However, in order to assess the effect of a covariate, a standard scatterplot is often difficult to interpret because of the presence of censored observations. Several authors have proposed a running median as an effective way of indicating the effect of a covariate. This article proposes a form of kernel estimation, employing double smoothing, that can be applied in a simple and efficient manner to construct an estimator of a percentile of the survival distribution as a function of one or two covariates. Permutations and bootstrap samples can be used to construct reference bands that help identify whether particular features of the estimates indicate real features of the underlying curve or whether this may be due simply to random variation. The techniques are illustrated on data from a study of kidney transplant patients.