Diffusion magnetic resonance imaging is an imaging technology designed to probe anatomical architectures of biological samples in an in vivo and noninvasive manner through measuring water diffusion. The contribution of this paper is threefold. First, it proposes a new method to identify and estimate multiple diffusion directions within a voxel through a new and identifiable parametrization of the widely used multi-tensor model. Unlike many existing methods, this method focuses on the estimation of diffusion directions rather than the diffusion tensors. Second, this paper proposes a novel direction smoothing method which greatly improves direction estimation in regions with crossing fibers. This smoothing method is shown to have excellent theoretical and empirical properties. Last, this paper develops a fiber tracking algorithm that can handle multiple directions within a voxel. The overall methodology is illustrated with simulated data and a data set collected for the study of Alzheimer’s disease by the Alzheimer’s Disease Neuroimaging Initiative (ADNI).

We consider policy evaluations when the Stable Unit Treatment Value Assumption (SUTVA) is violated due to the presence of interference among units. We propose to explicitly model interference as a function of units’ characteristics. Our approach is applied to the evaluation of a policy implemented in Tuscany (a region in Italy) on small handicraft firms. Results show that the benefits from the policy are reduced when treated firms are subject to high levels of interference. Moreover, the average causal effect is slightly underestimated when interference is ignored. We stress the importance of considering possible interference among units when evaluating and planning policy interventions.

RNA-Seq has emerged as the method of choice for profiling the transcriptomes of organisms. In particular, it aims to quantify the expression levels of transcripts using short nucleotide sequences or short reads generated from RNA-Seq experiments. In real experiments, the label of the transcript, from which each short read is generated, is missing, and short reads are mapped to the genome rather than the transcriptome. Therefore, the quantification of transcript expression levels is an indirect statistical inference problem.

In this article, we propose to use individual exonic base pairs as observation units and, further, to model nonzero as well as zero counts at all base pairs at both the transcript and gene levels. At the transcript level, two-component Poisson mixture distributions are postulated, which gives rise to the Convolution of Poisson mixture (CPM) distribution model at the gene level. The maximum likelihood estimation method equipped with the EM algorithm is used to estimate model parameters and quantify transcript expression levels. We refer to the proposed method as CPM-Seq. Both simulation studies and real data demonstrate the effectiveness of CPM-Seq, showing that CPM-Seq produces more accurate and consistent quantification results than Cufflinks.

Motivated by the important problem of detecting association between genetic markers and binary traits in genome-wide association studies, we present a novel Bayesian model that establishes a hierarchy between markers and genes by defining weights according to gene lengths and distances from genes to markers. The proposed hierarchical model uses these weights to define unique prior probabilities of association for markers based on their proximities to genes that are believed to be relevant to the trait of interest. We use an expectation-maximization algorithm in a filtering step to first reduce the dimensionality of the data and then sample from the posterior distribution of the model parameters to estimate posterior probabilities of association for the markers. We offer practical and meaningful guidelines for the selection of the model tuning parameters and propose a pipeline that exploits a singular value decomposition on the raw data to make our model run efficiently on large data sets. We demonstrate the performance of the model in simulation studies and conclude by discussing the results of a case study using a real-world data set provided by the Wellcome Trust Case Control Consortium.

Early childhood education research often compares a group of children who receive the intervention of interest to a group of children who receive care in a range of different care settings. In this paper, we estimate differential impacts of an early childhood intervention by alternative care type, using data from the Head Start Impact Study, a large-scale randomized evaluation. To do so, we utilize a Bayesian principal stratification framework to estimate separate impacts for two types of Compliers: those children who would otherwise be in other center-based care when assigned to control and those who would otherwise be in home-based care. We find strong, positive short-term effects of Head Start on receptive vocabulary for those Compliers who would otherwise be in home-based care. By contrast, we find no meaningful impact of Head Start on vocabulary for those Compliers who would otherwise be in other center-based care. Our findings suggest that alternative care type is a potentially important source of variation in early childhood education interventions.

Particulate matter (PM) is a class of malicious environmental pollutants known to be detrimental to human health. Regulatory efforts aimed at curbing PM levels in different countries often require high resolution space–time maps that can identify red-flag regions exceeding statutory concentration limits. Continuous spatio-temporal Gaussian Process (GP) models can deliver maps depicting predicted PM levels and quantify predictive uncertainty. However, GP-based approaches are usually thwarted by computational challenges posed by large datasets. We construct a novel class of scalable Dynamic Nearest Neighbor Gaussian Process (DNNGP) models that can provide a sparse approximation to any spatio-temporal GP (e.g., with nonseparable covariance structures). The DNNGP we develop here can be used as a sparsity-inducing prior for spatio-temporal random effects in any Bayesian hierarchical model to deliver full posterior inference. Storage and memory requirements for a DNNGP model are linear in the size of the dataset, thereby delivering massive scalability without sacrificing inferential richness. Extensive numerical studies reveal that the DNNGP provides substantially superior approximations to the underlying process than low-rank approximations. Finally, we use the DNNGP to analyze a massive air quality dataset to substantially improve predictions of PM levels across Europe in conjunction with the LOTOS-EUROS chemistry transport models (CTMs).

We consider the problem of emulating (approximating) computer models (simulators) that produce massive output. The specific simulator we study is a computer model of volcanic pyroclastic flow, a single run of which produces up to $10^{9}$ outputs over a space–time grid of coordinates. An emulator (essentially a statistical model of the simulator—we use a Gaussian Process) that is computationally suitable for such massive output is developed and studied from practical and theoretical perspectives. On the practical side, the emulator does unexpectedly well in predicting what the simulator would produce, even better than much more flexible and computationally intensive alternatives. This allows the attainment of the scientific goal of this work, accurate assessment of the hazards from pyroclastic flows over wide spatial domains. Theoretical results are also developed that provide insight into the unexpected success of the massive emulator. Generalizations of the emulator are introduced that allow for a nugget, which is useful for the application to hazard assessment.

Integrative analysis of multiple experimental datasets measured over a large number of observational units is the focus of large numbers of contemporary genomic and epigenomic studies. The key objectives of such studies include not only inferring a hidden state of activity for each unit over individual experiments, but also detecting highly associated clusters of units based on their inferred states. Although there are a number of methods tailored for specific datasets, there is currently no state-of-the-art modeling framework for this general class of problems. In this paper, we develop the MBASIC (Matrix Based Analysis for State-space Inference and Clustering) framework. MBASIC consists of two parts: state-space mapping and state-space clustering. In state-space mapping, it maps observations onto a finite state-space, representing the activation states of units across conditions. In state-space clustering, MBASIC incorporates a finite mixture model to cluster the units based on their inferred state-space profiles across all conditions. Both the state-space mapping and clustering can be simultaneously estimated through an Expectation-Maximization algorithm. MBASIC flexibly adapts to a large number of parametric distributions for the observed data, as well as the heterogeneity in replicate experiments. It allows for imposing structural assumptions on each cluster, and enables model selection using information criterion. In our data-driven simulation studies, MBASIC showed significant accuracy in recovering both the underlying state-space variables and clustering structures. We applied MBASIC to two genome research problems using large numbers of datasets from the ENCODE project. The first application grouped genes based on transcription factor occupancy profiles of their promoter regions in two different cell types. The second application focused on identifying groups of loci that are similar to a GATA2 binding site that is functional at its endogenous locus by utilizing transcription factor occupancy data and illustrated applicability of MBASIC in a wide variety of problems. In both studies, MBASIC showed higher levels of raw data fidelity than analyzing these data with a two-step approach using ENCODE results on transcription factor occupancy data.

In this paper we propose network methodology to infer prognostic cancer biomarkers based on the epigenetic pattern DNA methylation. Epigenetic processes such as DNA methylation reflect environmental risk factors, and are increasingly recognised for their fundamental role in diseases such as cancer. DNA methylation is a gene-regulatory pattern, and hence provides a means by which to assess genomic regulatory interactions. Network models are a natural way to represent and analyse groups of such interactions. The utility of network models also increases as the quantity of data and number of variables increase, making them increasingly relevant to large-scale genomic studies. We propose methodology to infer prognostic genomic networks from a DNA methylation-based measure of genomic interaction and association. We then show how to identify prognostic biomarkers from such networks, which we term “network community oncomarkers”. We illustrate the power of our proposed methodology in the context of a large publicly available breast cancer dataset.

We propose a unified framework for conducting inference on complex aggregated data in high-dimensional settings. We assume the data are a collection of multiple non-Gaussian realizations with underlying undirected graphical structures. Using the concept of median graphs in summarizing the commonality across these graphical structures, we provide a novel semiparametric approach to modeling such complex aggregated data, along with robust estimation of the median graph, which is assumed to be sparse. We prove the estimator is consistent in graph recovery and give an upper bound on the rate of convergence. We further provide thorough numerical analysis on both synthetic and real datasets to illustrate the empirical usefulness of the proposed models and methods.

Maternal smoking is well known to adversely affect birth outcomes, and there is considerable spatial variation in the rates of maternal smoking in the city of Glasgow, Scotland. This spatial variation is a partial driver of health inequalities between rich and poor communities, and it is of interest to determine the extent to which these inequalities have changed over time. Therefore in this paper we develop a Bayesian hierarchical model for estimating the spatio-temporal pattern in smoking incidence across Glasgow between 2000 and 2013, which can identify the changing geographical extent of clusters of areas exhibiting elevated maternal smoking incidences that partially drive health inequalities. Additionally, we provide freely available software via the R package CARBayesST to allow others to implement the model we have developed. The study period includes the introduction of a ban on smoking in public places in 2006, and the results show an average decline of around 11% in maternal smoking rates over the study period.

In an observational study of the effects caused by treatments, a sensitivity analysis asks about the magnitude of bias from unmeasured covariates that would need to be present to alter the conclusions of a naive analysis that presumes adjustments for measured covariates remove all biases. When there are two or more outcomes in an observational study, these outcomes may be unequally sensitive to unmeasured biases, and the least sensitive finding may concern a combination of several outcomes. A method of sensitivity analysis is proposed using Scheffé projections that permits the investigator to consider all linear contrasts in two or more scored outcomes while controlling the family-wise error rate. In sufficiently large samples, the method will exhibit insensitivity to bias that is greater than or equal to methods, such as the Bonferroni–Holm procedure, that focus on individual outcomes; that is, Scheffé projections have larger design sensitivities. More precisely, if the least sensitive linear combination is a single one of the several outcomes, then the design sensitivity using Scheffé projections equals that using a Bonferroni correction, but if the least sensitive combination is a nontrivial combination of two or more outcomes, then Scheffé projections have larger design sensitivities. This asymptotic property is examined in terms of finite sample power of sensitivity analyses using simulation. The method is applied to a replication with recent data of a well-known study of the effects of smoking on periodontal disease. In the example, the comparison that is least sensitive to bias from unmeasured covariates combines results for lower and upper teeth, but emphasizes lower teeth. This pattern would be difficult to anticipate prior to examining the data, but Scheffé’s method permits use of this unanticipated pattern without fear of capitalizing on chance.

Due to recent advances in technology, medical diagnosis data are becoming increasingly complex and, nowadays, applications where measurements are curves or images are ubiquitous. Motivated by the need of modeling a functional covariate on a metabolic syndrome case study, we develop a nonparametric functional regression model for the area under the specificity receiver operating characteristic curve. This partial area is a meaningful summary measure of diagnostic accuracy for cases in which misdiagnosis of diseased subjects may lead to serious clinical consequences, and hence it is critical to maintain a high sensitivity. Its normalized value can be interpreted as the average specificity over the interval of sensitivities considered, thus summarizing the trade-off between sensitivity and specificity. Our methods are motivated by, and applied to, a metabolic syndrome study that investigates how restricting the sensitivity of the gamma-glutamyl-transferase, a metabolic syndrome marker, to certain clinical meaningful values, affects its corresponding specificity and how it might change for different curves of arterial oxygen saturation. Application of our methods suggests that oxygen saturation is key to gamma-glutamyl transferase’s performance and that some of the different intervals of sensitivities considered offer a good trade-off between sensitivity and specificity. The simulation study shows that the estimator associated with our model is able to recover successfully the true overall shape of the functional covariate-adjusted partial area under the curve in different complex scenarios.

We introduce a dependent Bayesian nonparametric model for the probabilistic modeling of membership of subgroups in a community based on partially replicated data. The focus here is on species-by-site data, that is, community data where observations at different sites are classified in distinct species. Our aim is to study the impact of additional covariates, for instance, environmental variables, on the data structure, and in particular on the community diversity. To this end, we introduce dependence a priori across the covariates and show that it improves posterior inference. We use a dependent version of the Griffiths–Engen–McCloskey distribution defined via the stick-breaking construction. This distribution is obtained by transforming a Gaussian process whose covariance function controls the desired dependence. The resulting posterior distribution is sampled by Markov chain Monte Carlo. We illustrate the application of our model to a soil microbial data set acquired across a hydrocarbon contamination gradient at the site of a fuel spill in Antarctica. This method allows for inference on a number of quantities of interest in ecotoxicology, such as diversity or effective concentrations, and is broadly applicable to the general problem of community response to environmental variables.

Bayesian Melding (BM) and downscaling are two Bayesian approaches commonly used to combine data from different sources for statistical inference. We extend these two approaches to combine accurate but sparse direct observations with another set of high-resolution but biased calculated observations. We use our methods to estimate the path of a moving or evolving object and apply them in a case study of tracking northern fur seals. To make the BM approach computationally feasible for high-dimensional (big) data, we exploit the properties of the processes along with approximations to the likelihood to break the high-dimensional problem into a series of lower dimensional problems. To implement the alternative, downscaling approach, we use R-INLA to connect the two sources of observations via a linear mixed effect model. We compare the predictions of the two approaches by cross-validation as well as simulations. Our results show that both approaches yield similar results—both provide accurate, high resolution estimates of the at-sea locations of the northern fur seals, as well as Bayesian credible intervals to characterize the uncertainty about the estimated movement paths.

Imaging genetics has rapidly emerged as a promising approach for investigating the genetic determinants of brain mechanisms that underlie an individual’s behavior or psychiatric condition. In particular, for early detection and targeted treatment of schizophrenia, it is of high clinical relevance to identify genetic variants and imaging-based biomarkers that can be used as diagnostic markers, in addition to commonly used symptom-based assessments. By combining single-nucleotide polymorphism (SNP) arrays and functional magnetic resonance imaging (fMRI), we propose an integrative Bayesian risk prediction model that allows us to discriminate between individuals with schizophrenia and healthy controls, based on a sparse set of discriminatory regions of interest (ROIs) and SNPs. Inference on a regulatory network between SNPs and ROI intensities (ROI–SNP network) is used in a single modeling framework to inform the selection of the discriminatory ROIs and SNPs. We use simulation studies to assess the performance of our method and apply it to data collected from individuals with schizophrenia and healthy controls. We found our approach to outperform competing methods that do not link the ROI–SNP network to the selection of discriminatory markers.

Individuals of protected species, such as amphibians and reptiles, often need to be removed from sites before development commences. Usually, the population is considered to be closed. All individuals are assumed to (i) be present and available for detection at the start of the study period and (ii) remain at the site until the end of the study, unless they are detected. However, the assumption of population closure is not always valid. We present new removal models which allow for population renewal through birth and/or immigration, and population depletion through sampling as well as through death/emigration. When appropriate, productivity may be estimated and a Bayesian approach allows the estimation of the probability of total population depletion. We demonstrate the performance of the models using data on common lizards, Zootoca vivipara, and great crested newts, Triturus cristatus.

Soil erosion and sediment transport into waterways and the ocean can adversely affect water clarity, leading to the deterioration of marine ecosystems such as the iconic Great Barrier Reef (GBR) in Australia. Quantifying a sediment load and its associated uncertainty is an important task in delineating how changes in management practices can contribute to improvements in water quality, and therefore continued sustainability of the GBR. However, monitoring data are spatially (and often temporally) sparse, making load estimation complicated, particularly when there are lengthy periods between sampling or during peak flow periods of major events when samples cannot be safely taken.

We develop a spatio-temporal statistical model that is mechanistically motivated by a process-based deterministic model called Dynamic SedNet. The model is developed within a Bayesian hierarchical modelling framework that uses dimension reduction to accommodate seasonal and spatial patterns to assimilate monitored sediment concentration and flow data with output from Dynamic SedNet. The approach is applied in the Upper Burdekin catchment in Queensland, Australia, where we obtain daily estimates of sediment concentrations, stream discharge volumes and sediment loads at 411 spatial locations across 20 years. Our approach provides a method for assimilating both monitoring data and modelled output, providing a statistically rigorous method for quantifying uncertainty through space and time that was previously unavailable through process-based models.

Mapping molecular QTLs has emerged as an important tool for understanding the genetic basis of cell functions. With the increasing availability of functional genomic data, it is natural to incorporate genomic annotations into QTL discovery. Discovering molecular QTLs is typically framed as a multiple hypothesis testing problem and solved using false discovery rate (FDR) control procedures. Currently, most existing statistical approaches rely on obtaining $p$-values for each candidate locus through permutation-based schemes, which are not only inconvenient for incorporating highly informative genomic annotations but also computationally inefficient. In this paper, we discuss a novel statistical approach for integrative QTL discovery based on the theoretical framework of Bayesian FDR control. We use a Bayesian hierarchical model to naturally integrate genomic annotations into molecular QTL mapping and propose an empirical Bayes-based computational procedure to approximate the necessary posterior probabilities to achieve high computational efficiency. Through theoretical arguments and simulation studies, we demonstrate that the proposed approach rigorously controls the desired type I error rate and greatly improves the power of QTL discovery when incorporating informative annotations. Finally, we demonstrate our approach by analyzing the expression-genotype data from 44 human tissues generated by the GTEx project. By integrating the simple annotation of SNP distance to transcription start sites, we discover more genes that harbor expression-associated SNPs in all 44 tissues, with an average increase of 1485 genes per tissue.

A multilevel functional data method is adapted for forecasting age-specific mortality for two or more populations in developed countries with high-quality vital registration systems. It uses multilevel functional principal component analysis of aggregate and population-specific data to extract the common trend and population-specific residual trend among populations. If the forecasts of population-specific residual trends do not show a long-term trend, then convergence in forecasts may be achieved. This method is first applied to age- and sex-specific data for the United Kingdom, and its forecast accuracy is then further compared with several existing methods, including independent functional data and product-ratio methods, through a multi-country comparison. The proposed method is also demonstrated by age-, sex- and state-specific data in Australia, where the convergence in forecasts can possibly be achieved by sex and state. For forecasting age-specific mortality, the multilevel functional data method is more accurate than the other coherent methods considered. For forecasting female life expectancy at birth, the multilevel functional data method is outperformed by the Bayesian method of Raftery, Lalic and Gerland [Demogr. Res. 30 (2014) 795–822]. For forecasting male life expectancy at birth, the multilevel functional data method performs better than the Bayesian methods in terms of point forecasts, but less well in terms of interval forecasts. Supplementary materials for this article are available online.

This project aims to explore which combinations of meteorological conditions are associated with extreme ground level ozone conditions. Our approach focuses only on the tail by optimizing the tail dependence between the ozone response and functions of meteorological covariates. Since there is a long list of possible meteorological covariates, the space of possible models cannot be explored completely. Consequently, we perform data mining within the model selection context, employing an automated model search procedure. Our study is unique among extremes applications, as optimizing tail dependence has not previously been attempted, and it presents new challenges, such as requiring a smooth threshold. We present a simulation study which shows that the method can detect complicated conditions leading to extreme responses and resists overfitting. We apply the method to ozone data for Atlanta and Charlotte and find similar meteorological drivers for these two Southeastern US cities. We identify several covariates which help to differentiate the meteorological conditions which lead to extreme ozone levels from those which lead to merely high levels.

Large-scale multiple testing tasks often exhibit dependence. Leveraging the dependence between individual tests is still one challenging and important problem in statistics. With recent advances in graphical models, it is feasible to use them to capture the dependence among multiple hypotheses. We propose a multiple testing procedure which is based on a Markov-random-field-coupled mixture model. The underlying true states of hypotheses are represented by a latent binary Markov random field, and the observed test statistics appear as the coupled mixture variables. The model can be learned by a novel EM algorithm. The next step is to infer the posterior probability that each hypothesis is null (termed local index of significance), and the false discovery rate can be controlled accordingly. We also provide a semiparametric variation of the graphical model which is useful in the situation where $f_{1}$ (the density function of the test statistic under the alternative hypothesis) is heterogeneous among multiple hypotheses. This semiparametric approach exactly generalizes the local FDR procedure [J. Amer. Statist. Assoc. 96 (2001) 1151–1160] and connects with the BH procedure [J. Roy. Statist. Soc. Ser. B 57 (1995) 289–300]. Simulations show that the numerical performance of multiple testing can be improved substantially by using our procedure. We apply the procedure to a real-world genome-wide association study on breast cancer, and we identify several SNPs with strong association evidence.

Space–time Hawkes point process models for the conditional rate of earthquake occurrences traditionally make many parametric assumptions about the form of the triggering function for the rate of aftershocks following an earthquake. As an alternative, Marsan and Lengliné [Science 319 (2008) 1076–1079] developed a completely nonparametric method that provides an estimate of a homogeneous background rate for mainshocks, and a histogram estimate of the triggering function. At each step of the procedure the model estimates rely on computing the probability each earthquake is a mainshock or aftershock of a previous event. The focus of this paper is the improvement and assessment of Marsan and Lengliné’s method in the following ways: (a) the proposal of novel ways to incorporate a spatially inhomogeneous background rate; (b) adding error bars to the histogram estimates which quantify the sampling variability in the estimation of the underlying seismic process. A simulation study is designed to evaluate and validate the ability of our methods to recover the triggering function and spatially varying background rate. An application to earthquake data from the Tohoku District in Japan is discussed at the end, and the results are compared to a well-established parametric model of seismicity for this region.

In this paper, we consider an estimation problem arising in the measurement of bifurcated flow of the Teesta, a trans-boundary river flowing through India and Bangladesh. The location of measurement is an Indian Barrage, where a part of the flow is diverted from the main stream to a canal. The flows through the two channels are regulated by different control structures and are measured indirectly from the height of the water level and the dimensions of the control structures. The computational formula for the measurement involves a hydrological constant used as a multiplier. Empirical findings indicate that incorrect multipliers are currently used in the computational formula for the two channels. For implementing any water sharing treaty between the two countries, the measurements need to be brought to a common scale. For this purpose, we present a model with carefully considered assumptions to estimate the correction factor. The model permits diagnostic tests for validation of the assumptions. We provide a nonparametric and consistent estimator of the desired factor.

Analysis of historical flow data shows that a main stream flow measured as 100 cumec would be measured as 76 cumec if it is diverted through the canal. Adjustment of emerging measurements through this finding would help the governments of India and Bangladesh to effectively implement and monitor any water sharing agreement.