We use a discrete-time proportional hazards model of time to involuntary employment termination. This model enables us to examine both the continuous effect of the age of an employee and whether that effect has varied over time, generalizing earlier work [Kadane and Woodworth J. Bus. Econom. Statist. 22 (2004) 182–193]. We model the log hazard surface (over age and time) as a thin-plate spline, a Bayesian smoothness-prior implementation of penalized likelihood methods of surface-fitting [Wahba (1990) Spline Models for Observational Data. SIAM]. The nonlinear component of the surface has only two parameters, smoothness and anisotropy. The first, a scale parameter, governs the overall smoothness of the surface, and the second, anisotropy, controls the relative smoothness over time and over age. For any fixed value of the anisotropy parameter, the prior is equivalent to a Gaussian process with linear drift over the time–age plane with easily computed eigenvectors and eigenvalues that depend only on the configuration of data in the time–age plane and the anisotropy parameter. This model has application to legal cases in which a company is charged with disproportionately disadvantaging older workers when deciding whom to terminate. We illustrate the application of the modeling approach using data from an actual discrimination case.

Recently-developed genotype imputation methods are a powerful tool for detecting untyped genetic variants that affect disease susceptibility in genetic association studies. However, existing imputation methods require individual-level genotype data, whereas, in practice, it is often the case that only summary data are available. For example, this may occur because, for reasons of privacy or politics, only summary data are made available to the research community at large; or because only summary data are collected, as in DNA pooling experiments. In this article we introduce a new statistical method that can accurately infer the frequencies of untyped genetic variants in these settings, and indeed substantially improve frequency estimates at typed variants in pooling experiments where observations are noisy. Our approach, which predicts each allele frequency using a linear combination of observed frequencies, is statistically straightforward, and related to a long history of the use of linear methods for estimating missing values (e.g., Kriging). The main statistical novelty is our approach to regularizing the covariance matrix estimates, and the resulting linear predictors, which is based on methods from population genetics. We find that, besides being both fast and flexible—allowing new problems to be tackled that cannot be handled by existing imputation approaches purpose-built for the genetic context—these linear methods are also very accurate. Indeed, imputation accuracy using this approach is similar to that obtained by state-of-the-art imputation methods that use individual-level data, but at a fraction of the computational cost.

This paper compares two methods built on a hierarchical Bayesian foundation and designed for modeling hourly ozone concentrations over the eastern United States. One, a dynamic linear state space model (DLM) that has been proposed earlier, lies in a very contemporary setting where two historical paths to temporal process models, the Kalman filter and the dynamic system with random perturbations, converge. The other, which we call the Bayesian spatial predictor (BSP), is a Bayesian alternative to the purely spatial method of kriging. The DLM as a dynamic system model has parameters that are states of the process which generate the ozone and change with time. More specifically, the model includes a time-varying site invariant mean field as well as time-varying coefficients for 24 and 12 hour diurnal cyclic components. The resulting model’s great flexibility comes at the cost of complexity, forcing the use of an MCMC approach and very time-consuming computations. Thus, the size of the DLM’s spatial domain of applicability has to be restricted and the number of monitoring sites that can be treated limited. The paper’s assessment of the DLM reveals other difficulties that point to the need to consider a less flexible competitor, a Bayesian spatial predictor (BSP). The two methods are compared in a variety of ways and overall conclusions given. In particular, the conclusions point to the BSP as the more practical alternative for spatial prediction.

Genomic imprinting has been thought to play an important role in seed development in flowering plants. Seed in a flowering plant normally contains diploid embryo and triploid endosperm. Empirical studies have shown that some economically important endosperm traits are genetically controlled by imprinted genes. However, the exact number and location of the imprinted genes are largely unknown due to the lack of efficient statistical mapping methods. Here we propose a general statistical variance components framework by utilizing the natural information of sex-specific allelic sharing among sibpairs in line crosses, to map imprinted quantitative trait loci (iQTL) underlying endosperm traits. We propose a new variance components partition method considering the unique characteristic of the triploid endosperm genome, and develop a restricted maximum likelihood estimation method in an interval scan for estimating and testing genome-wide iQTL effects. Cytoplasmic maternal effect which is thought to have primary influences on yield and grain quality is also considered when testing for genomic imprinting. Extension to multiple iQTL analysis is proposed. Asymptotic distribution of the likelihood ratio test for testing the variance components under irregular conditions are studied. Both simulation study and real data analysis indicate good performance and powerfulness of the developed approach.

In many metropolitan areas efforts are made to count the homeless to ensure proper provision of social services. Some areas are very large, which makes spatial sampling a viable alternative to an enumeration of the entire terrain. Counts are observed in sampled regions but must be imputed in unvisited areas. Along with the imputation process, the costs of underestimating and overestimating may be different. For example, if precise estimation in areas with large homeless c ounts is critical, then underestimation should be penalized more than overestimation in the loss function. We analyze data from the 2004–2005 Los Angeles County homeless study using an augmentation of L₁ stochastic gradient boosting that can weight overestimates and underestimates asymmetrically. We discuss our choice to utilize stochastic gradient boosting over other function estimation procedures. In-sample fitted and out-of-sample imputed values, as well as relationships between the response and predictors, are analyzed for various cost functions. Practical usage and policy implications of these results are discussed briefly.

Given the potential for illicit nuclear material being used for terrorism, most ports now inspect a large number of goods entering national borders for radioactive cargo. The U.S. Department of Homeland Security is moving toward one hundred percent inspection of all containers entering the U.S. at various ports of entry for nuclear material. We propose a Bayesian classification approach for the real-time data collected by the inline Polyvinyl Toluene radiation portal monitors. We study the computational and asymptotic properties of the proposed method and demonstrate its efficacy in simulations. Given data available to the authorities, it should be feasible to implement this approach in practice.

Neuroimaging produces data that are continuous in one or more dimensions. This calls for an inference framework that can handle data that approximate functions of space, for example, anatomical images, time–frequency maps and distributed source reconstructions of electromagnetic recordings over time. Statistical parametric mapping (SPM) is the standard framework for whole-brain inference in neuroimaging: SPM uses random field theory to furnish p-values that are adjusted to control family-wise error or false discovery rates, when making topological inferences over large volumes of space. Random field theory regards data as realizations of a continuous process in one or more dimensions. This contrasts with classical approaches like the Bonferroni correction, which consider images as collections of discrete samples with no continuity properties (i.e., the probabilistic behavior at one point in the image does not depend on other points). Here, we illustrate how random field theory can be applied to data that vary as a function of time, space or frequency. We emphasize how topological inference of this sort is invariant to the geometry of the manifolds on which data are sampled. This is particularly useful in electromagnetic studies that often deal with very smooth data on scalp or cortical meshes. This application illustrates the versatility and simplicity of random field theory and the seminal contributions of Keith Worsley (1951–2009), a key architect of topological inference.

Government statistical agencies often apply statistical disclosure limitation techniques to survey microdata to protect the confidentiality of respondents. There is a need for valid and practical ways to assess the protection provided. This paper develops some simple methods for disclosure limitation techniques which perturb the values of categorical identifying variables. The methods are applied in numerical experiments based upon census data from the United Kingdom which are subject to two perturbation techniques: data swapping (random and targeted) and the post randomization method. Some simplifying approximations to the measure of risk are found to work well in capturing the impacts of these techniques. These approximations provide simple extensions of existing risk assessment methods based upon Poisson log-linear models. A numerical experiment is also undertaken to assess the impact of multivariate misclassification with an increasing number of identifying variables. It is found that the misclassification dominates the usual monotone increasing relationship between this number and risk so that the risk eventually declines, implying less sensitivity of risk to choice of identifying variables. The methods developed in this paper may also be used to obtain more realistic assessments of risk which take account of the kinds of measurement and other nonsampling errors commonly arising in surveys.

The generation of multi-step density forecasts for non-Gaussian data mostly relies on Monte Carlo simulations which are computationally intensive. Using aggregated wind power in Ireland, we study two approaches of multi-step density forecasts which can be obtained from simple iterations so that intensive computations are avoided. In the first approach, we apply a logistic transformation to normalize the data approximately and describe the transformed data using ARIMA–GARCH models so that multi-step forecasts can be iterated easily. In the second approach, we describe the forecast densities by truncated normal distributions which are governed by two parameters, namely, the conditional mean and conditional variance. We apply exponential smoothing methods to forecast the two parameters simultaneously. Since the underlying model of exponential smoothing is Gaussian, we are able to obtain multi-step forecasts of the parameters by simple iterations and thus generate forecast densities as truncated normal distributions. We generate forecasts for wind power from 15 minutes to 24 hours ahead. Results show that the first approach generates superior forecasts and slightly outperforms the second approach under various proper scores. Nevertheless, the second approach is computationally more efficient and gives more robust results under different lengths of training data. It also provides an attractive alternative approach since one is allowed to choose a particular parametric density for the forecasts, and is valuable when there are no obvious transformations to normalize the data.

Technological advances in genotyping have given rise to hypothesis-based association studies of increasing scope. As a result, the scientific hypotheses addressed by these studies have become more complex and more difficult to address using existing analytic methodologies. Obstacles to analysis include inference in the face of multiple comparisons, complications arising from correlations among the SNPs (single nucleotide polymorphisms), choice of their genetic parametrization and missing data. In this paper we present an efficient Bayesian model search strategy that searches over the space of genetic markers and their genetic parametrization. The resulting method for Multilevel Inference of SNP Associations, MISA, allows computation of multilevel posterior probabilities and Bayes factors at the global, gene and SNP level, with the prior distribution on SNP inclusion in the model providing an intrinsic multiplicity correction. We use simulated data sets to characterize MISA’s statistical power, and show that MISA has higher power to detect association than standard procedures. Using data from the North Carolina Ovarian Cancer Study (NCOCS), MISA identifies variants that were not identified by standard methods and have been externally “validated” in independent studies. We examine sensitivity of the NCOCS results to prior choice and method for imputing missing data. MISA is available in an R package on CRAN.

Boundary distance (BD) plotting is a technique for making orientation invariant comparisons of the spatial distribution of biochemical markers within and across cells/nuclei. Marker expression is aggregated over points with the same distance from the boundary. We present a suite of tools for improved data analysis and statistical inference using BD plotting. BD is computed using the Euclidean distance transform after presmoothing and oversampling of nuclear boundaries. Marker distribution profiles are averaged using smoothing with linearly decreasing bandwidth. Average expression curves are scaled and registered by x-axis dilation to compensate for uneven lighting and errors in nuclear boundary marking. Penalized discriminant analysis is used to characterize the quality of separation between average marker distributions. An adaptive piecewise linear model is used to compare expression gradients in intra, peri and extra nuclear zones. The techniques are illustrated by the following: (a) a two sample problem involving a pair of voltage gated calcium channels (Cav1.2 and AB70) marked in different cells; (b) a paired sample problem of calcium channels (Y1F4 and RyR1) marked in the same cell.

Presence-only data, point locations where a species has been recorded as being present, are often used in modeling the distribution of a species as a function of a set of explanatory variables—whether to map species occurrence, to understand its association with the environment, or to predict its response to environmental change. Currently, ecologists most commonly analyze presence-only data by adding randomly chosen “pseudo-absences” to the data such that it can be analyzed using logistic regression, an approach which has weaknesses in model specification, in interpretation, and in implementation. To address these issues, we propose Poisson point process modeling of the intensity of presences. We also derive a link between the proposed approach and logistic regression—specifically, we show that as the number of pseudo-absences increases (in a regular or uniform random arrangement), logistic regression slope parameters and their standard errors converge to those of the corresponding Poisson point process model. We discuss the practical implications of these results. In particular, point process modeling offers a framework for choice of the number and location of pseudo-absences, both of which are currently chosen by ad hoc and sometimes ineffective methods in ecology, a point which we illustrate by example.

Modeling species abundance patterns using local environmental features is an important, current problem in ecology. The Cape Floristic Region (CFR) in South Africa is a global hot spot of diversity and endemism, and provides a rich class of species abundance data for such modeling. Here, we propose a multi-stage Bayesian hierarchical model for explaining species abundance over this region. Our model is specified at areal level, where the CFR is divided into roughly 37, 000 one minute grid cells; species abundance is observed at some locations within some cells. The abundance values are ordinally categorized. Environmental and soil-type factors, likely to influence the abundance pattern, are included in the model. We formulate the empirical abundance pattern as a degraded version of the potential pattern, with the degradation effect accomplished in two stages. First, we adjust for land use transformation and then we adjust for measurement error, hence misclassification error, to yield the observed abundance classifications. An important point in this analysis is that only 28% of the grid cells have been sampled and that, for sampled grid cells, the number of sampled locations ranges from one to more than one hundred. Still, we are able to develop potential and transformed abundance surfaces over the entire region.

In the hierarchical framework, categorical abundance classifications are induced by continuous latent surfaces. The degradation model above is built on the latent scale. On this scale, an areal level spatial regression model was used for modeling the dependence of species abundance on the environmental factors. To capture anticipated similarity in abundance pattern among neighboring regions, spatial random effects with a conditionally autoregressive prior (CAR) were specified. Model fitting is through familiar Markov chain Monte Carlo methods. While models with CAR priors are usually efficiently fitted, even with large data sets, with our modeling and the large number of cells, run times became very long. So a novel parallelized computing strategy was developed to expedite fitting. The model was run for six different species. With categorical data, display of the resultant abundance patterns is a challenge and we offer several different views. The patterns are of importance on their own, comparatively across the region and across species, with implications for species competition and, more generally, for planning and conservation.

We consider the problem of constructing optimal designs for population pharmacokinetics which use random effect models. It is common practice in the design of experiments in such studies to assume uncorrelated errors for each subject. In the present paper a new approach is introduced to determine efficient designs for nonlinear least squares estimation which addresses the problem of correlation between observations corresponding to the same subject. We use asymptotic arguments to derive optimal design densities, and the designs for finite sample sizes are constructed from the quantiles of the corresponding optimal distribution function. It is demonstrated that compared to the optimal exact designs, whose determination is a hard numerical problem, these designs are very efficient. Alternatively, the designs derived from asymptotic theory could be used as starting designs for the numerical computation of exact optimal designs. Several examples of linear and nonlinear models are presented in order to illustrate the methodology. In particular, it is demonstrated that naively chosen equally spaced designs may lead to less accurate estimation.

Individual-level health data are often not publicly available due to confidentiality; masked data are released instead. Therefore, it is important to evaluate the utility of using the masked data in statistical analyses such as regression. In this paper we propose a data masking method which is based on spatial smoothing techniques. The proposed method allows for selecting both the form and the degree of masking, thus resulting in a large degree of flexibility. We investigate the utility of the masked data sets in terms of the mean square error (MSE) of regression parameter estimates when fitting a Generalized Linear Model (GLM) to the masked data. We also show that incorporating prior knowledge on the spatial pattern of the exposure into the data masking may reduce the bias and MSE of the parameter estimates. By evaluating both utility and disclosure risk as functions of the form and the degree of masking, our method produces a risk-utility profile which can facilitate the selection of masking parameters. We apply the method to a study of racial disparities in mortality rates using data on more than 4 million Medicare enrollees residing in 2095 zip codes in the Northeast region of the United States.

Assessment of circulating CD4 count change over time in HIV-infected subjects on antiretroviral therapy (ART) is a central component of disease monitoring. The increasing number of HIV-infected subjects starting therapy and the limited capacity to support CD4 count testing within resource-limited settings have fueled interest in identifying correlates of CD4 count change such as total lymphocyte count, among others. The application of modeling techniques will be essential to this endeavor due to the typically nonlinear CD4 trajectory over time and the multiple input variables necessary for capturing CD4 variability. We propose a prediction-based classification approach that involves first stage modeling and subsequent classification based on clinically meaningful thresholds. This approach draws on existing analytical methods described in the receiver operating characteristic curve literature while presenting an extension for handling a continuous outcome. Application of this method to an independent test sample results in greater than 98% positive predictive value for CD4 count change. The prediction algorithm is derived based on a cohort of n = 270 HIV-1 infected individuals from the Royal Free Hospital, London who were followed for up to three years from initiation of ART. A test sample comprised of n = 72 individuals from Philadelphia and followed for a similar length of time is used for validation. Results suggest that this approach may be a useful tool for prioritizing limited laboratory resources for CD4 testing after subjects start antiretroviral therapy.

Graphs and networks are common ways of depicting biological information. In biology, many different biological processes are represented by graphs, such as regulatory networks, metabolic pathways and protein–protein interaction networks. This kind of a priori use of graphs is a useful supplement to the standard numerical data such as microarray gene expression data. In this paper we consider the problem of regression analysis and variable selection when the covariates are linked on a graph. We study a graph-constrained regularization procedure and its theoretical properties for regression analysis to take into account the neighborhood information of the variables measured on a graph. This procedure involves a smoothness penalty on the coefficients that is defined as a quadratic form of the Laplacian matrix associated with the graph. We establish estimation and model selection consistency results and provide estimation bounds for both fixed and diverging numbers of parameters in regression models. We demonstrate by simulations and a real data set that the proposed procedure can lead to better variable selection and prediction than existing methods that ignore the graph information associated with the covariates.

In many medical studies, patients are followed longitudinally and interest is on assessing the relationship between longitudinal measurements and time to an event. Recently, various authors have proposed joint modeling approaches for longitudinal and time-to-event data for a single longitudinal variable. These joint modeling approaches become intractable with even a few longitudinal variables. In this paper we propose a regression calibration approach for jointly modeling multiple longitudinal measurements and discrete time-to-event data. Ideally, a two-stage modeling approach could be applied in which the multiple longitudinal measurements are modeled in the first stage and the longitudinal model is related to the time-to-event data in the second stage. Biased parameter estimation due to informative dropout makes this direct two-stage modeling approach problematic. We propose a regression calibration approach which appropriately accounts for informative dropout. We approximate the conditional distribution of the multiple longitudinal measurements given the event time by modeling all pairwise combinations of the longitudinal measurements using a bivariate linear mixed model which conditions on the event time. Complete data are then simulated based on estimates from these pairwise conditional models, and regression calibration is used to estimate the relationship between longitudinal data and time-to-event data using the complete data. We show that this approach performs well in estimating the relationship between multivariate longitudinal measurements and the time-to-event data and in estimating the parameters of the multiple longitudinal process subject to informative dropout. We illustrate this methodology with simulations and with an analysis of primary biliary cirrhosis (PBC) data.

In this article we describe a method for carrying out Bayesian estimation for the double Pareto lognormal (dPlN) distribution which has been proposed as a model for heavy-tailed phenomena. We apply our approach to estimate the dPlN / M / 1 and M / dPlN / 1 queueing systems. These systems cannot be analyzed using standard techniques due to the fact that the dPlN distribution does not possess a Laplace transform in closed form. This difficulty is overcome using some recent approximations for the Laplace transform of the interarrival distribution for the Pareto / M / 1 system. Our procedure is illustrated with applications in internet traffic analysis and risk theory.

We investigate the effect that the choice of measurement scale has upon inference and extrapolation in extreme value analysis. Separate analyses of variables from a single process on scales which are linked by a nonlinear transformation may lead to discrepant conclusions concerning the tail behavior of the process. We propose the use of a Box–Cox power transformation incorporated as part of the inference procedure to account parametrically for the uncertainty surrounding the scale of extrapolation. This has the additional feature of increasing the rate of convergence of the distribution tails to an extreme value form in certain cases and thus reducing bias in the model estimation. Inference without reparameterization is practicably infeasible, so we explore a reparameterization which exploits the asymptotic theory of normalizing constants required for nondegenerate limit distributions. Inference is carried out in a Bayesian setting, an advantage of this being the availability of posterior predictive return levels. The methodology is illustrated on both simulated data and significant wave height data from the North Sea.

We develop a new principal components analysis (PCA) type dimension reduction method for binary data. Different from the standard PCA which is defined on the observed data, the proposed PCA is defined on the logit transform of the success probabilities of the binary observations. Sparsity is introduced to the principal component (PC) loading vectors for enhanced interpretability and more stable extraction of the principal components. Our sparse PCA is formulated as solving an optimization problem with a criterion function motivated from a penalized Bernoulli likelihood. A Majorization–Minimization algorithm is developed to efficiently solve the optimization problem. The effectiveness of the proposed sparse logistic PCA method is illustrated by application to a single nucleotide polymorphism data set and a simulation study.

Stochastic processes often exhibit sudden systematic changes in pattern a short time before certain failure events. Examples include increase in medical costs before death and decrease in CD4 counts before AIDS diagnosis. To study such terminal behavior of stochastic processes, a natural and direct way is to align the processes using failure events as time origins. This paper studies backward stochastic processes counting time backward from failure events, and proposes one-sample nonparametric estimation of the mean of backward processes when follow-up is subject to left truncation and right censoring. We will discuss benefits of including prevalent cohort data to enlarge the identifiable region and large sample properties of the proposed estimator with related extensions. A SEER–Medicare linked data set is used to illustrate the proposed methodologies.