A Bayesian multivariate mixture model for high throughput spatial transcriptomics.

High throughput spatial transcriptomics (HST) is a rapidly emerging class of experimental technologies that allow for profiling gene expression in tissue samples at or near single-cell resolution while retaining the spatial location of each sequencing unit within the tissue sample. Through analyzing HST data, we seek to identify sub-populations of cells within a tissue sample that may inform biological phenomena. Existing computational methods either ignore the spatial heterogeneity in gene expression profiles, fail to account for important statistical features such as skewness, or are heuristic-based network clustering methods that lack the inferential benefits of statistical modeling. To address this gap, we develop SPRUCE: a Bayesian spatial multivariate finite mixture model based on multivariate skew-normal distributions, which is capable of identifying distinct cellular sub-populations in HST data. We further implement a novel combination of Pólya-Gamma data augmentation and spatial random effects to infer spatially correlated mixture component membership probabilities without relying on approximate inference techniques. Via a simulation study, we demonstrate the detrimental inferential effects of ignoring skewness or spatial correlation in HST data. Using publicly available human brain HST data, SPRUCE outperforms existing methods in recovering expertly annotated brain layers. Finally, our application of SPRUCE to human breast cancer HST data indicates that SPRUCE can distinguish distinct cell populations within the tumor microenvironment. An R package spruce for fitting the proposed models is available through The Comprehensive R Archive Network.

Contrasting principal stratum and hypothetical strategy estimands in multi-period crossover trials with incomplete data.

Complete case analyses of complete crossover designs provide an opportunity to make comparisons based on patients who can tolerate all treatments. It is argued that this provides a means of estimating a principal stratum strategy estimand, something which is difficult to do in parallel group trials. While some trial users will consider this a relevant aim, others may be interested in hypothetical strategy estimands, that is, the effect that would be found if all patients completed the trial. Whether these estimands differ importantly is a question of interest to the different users of the trial results. This paper derives the difference between principal stratum strategy and hypothetical strategy estimands, where the former is estimated by a complete-case analysis of the crossover design, and a model for the dropout process is assumed. Complete crossover designs, that is, those where all treatments appear in all sequences, and which compare t treatments over p periods with respect to a continuous outcome are considered. Numerical results are presented for Williams designs with four and six periods. Results from a trial of obstructive sleep apnoea-hypopnoea (TOMADO) are also used for illustration. The results demonstrate that the percentage difference between the estimands is modest, exceeding 5% only when the trial has been severely affected by dropouts or if the within-subject correlation is low.

Bayesian additive regression trees for multivariate skewed responses.

This paper introduces a nonparametric regression approach for univariate and multivariate skewed responses using Bayesian additive regression trees (BART). Existing BART methods use ensembles of decision trees to model a mean function, and have become popular recently due to their high prediction accuracy and ease of use. The usual assumption of a univariate Gaussian error distribution, however, is restrictive in many biomedical applications. Motivated by an oral health study, we provide a useful extension of BART, the skewBART model, to address this problem. We then extend skewBART to allow for multivariate responses, with information shared across the decision trees associated with different responses within the same subject. The methodology accommodates within-subject association, and allows varying skewness parameters for the varying multivariate responses. We illustrate the benefits of our multivariate skewBART proposal over existing alternatives via simulation studies and application to the oral health dataset with bivariate highly skewed responses. Our methodology is implementable via the R package skewBART, available on GitHub.

Bayesian Noise Modelling for State Estimation of the Spread of COVID-19 in Saudi Arabia with Extended Kalman Filters.

The epistemic uncertainty in coronavirus disease (COVID-19) model-based predictions using complex noisy data greatly affects the accuracy of pandemic trend and state estimations. Quantifying the uncertainty of COVID-19 trends caused by different unobserved hidden variables is needed to evaluate the accuracy of the predictions for complex compartmental epidemiological models. A new approach for estimating the measurement noise covariance from real COVID-19 pandemic data has been presented based on the marginal likelihood (Bayesian evidence) for Bayesian model selection of the stochastic part of the Extended Kalman filter (EKF), with a sixth-order nonlinear epidemic model, known as the SEIQRD (Susceptible-Exposed-Infected-Quarantined-Recovered-Dead) compartmental model. This study presents a method for testing the noise covariance in cases of dependence or independence between the infected and death errors, to better understand their impact on the predictive accuracy and reliability of EKF statistical models. The proposed approach is able to reduce the error in the quantity of interest compared to the arbitrarily chosen values in the EKF estimation.

A semiparametric isotonic regression model for skewed distributions with application to DNA-RNA-protein analysis.

In this paper, we propose a semiparametric regression model that is built upon an isotonic regression model with the assumption that the random error follows a skewed distribution. We develop an expectation-maximization algorithm for obtaining the maximum likelihood estimates of the model parameters, examine the asymptotic properties of the estimators, conduct simulation studies to explore the performance of the proposed model, and apply the method to evaluate the DNA-RNA-protein relationship and identify genes that are key factors in tumor progression.

Information-Theoretic Aspects of Location Parameter Estimation under Skew-Normal Settings.

In several applications, the assumption of normality is often violated in data with some level of skewness, so skewness affects the mean's estimation. The class of skew-normal distributions is considered, given their flexibility for modeling data with asymmetry parameter. In this paper, we considered two location parameter (µ) estimation methods in the skew-normal setting, where the coefficient of variation and the skewness parameter are known. Specifically, the least square estimator (LSE) and the best unbiased estimator (BUE) for µ are considered. The properties for BUE (which dominates LSE) using classic theorems of information theory are explored, which provides a way to measure the uncertainty of location parameter estimations. Specifically, inequalities based on convexity property enable obtaining lower and upper bounds for differential entropy and Fisher information. Some simulations illustrate the behavior of differential entropy and Fisher information bounds.

A Bayesian multivariate mixture model for skewed longitudinal data with intermittent missing observations: An application to infant motor development.

In studies of infant growth, an important research goal is to identify latent clusters of infants with delayed motor development-a risk factor for adverse outcomes later in life. However, there are numerous statistical challenges in modeling motor development: the data are typically skewed, exhibit intermittent missingness, and are correlated across repeated measurements over time. Using data from the Nurture study, a cohort of approximately 600 mother-infant pairs, we develop a flexible Bayesian mixture model for the analysis of infant motor development. First, we model developmental trajectories using matrix skew-normal distributions with cluster-specific parameters to accommodate dependence and skewness in the data. Second, we model the cluster-membership probabilities using a Pólya-Gamma data-augmentation scheme, which improves predictions of the cluster-membership allocations. Lastly, we impute missing responses from conditional multivariate skew-normal distributions. Bayesian inference is achieved through straightforward Gibbs sampling. Through simulation studies, we show that the proposed model yields improved inferences over models that ignore skewness or adopt conventional imputation methods. We applied the model to the Nurture data and identified two distinct developmental clusters, as well as detrimental effects of food insecurity on motor development. These findings can aid investigators in targeting interventions during this critical early-life developmental window.

Scale mixture of skew-normal linear mixed models with within-subject serial dependence.

In longitudinal studies, repeated measures are collected over time and hence they tend to be serially correlated. These studies are commonly analyzed using linear mixed models (LMMs), and in this article we consider an extension of the skew-normal/independent LMM, where the error term has a dependence structure, such as damped exponential correlation or autoregressive correlation of order p. The proposed model provides flexibility in capturing the effects of skewness and heavy tails simultaneously when continuous repeated measures are serially correlated. For this robust model, we present an efficient EM-type algorithm for parameters estimation via maximum likelihood and the observed information matrix is derived analytically to account for standard errors. The methodology is illustrated through an application to schizophrenia data and some simulation studies. The proposed algorithm and methods are implemented in the new R package skewlmm.

Spatial skew-normal/independent models for nonrandomly missing clustered data.

Clinical studies on periodontal disease (PD) often lead to data collected which are clustered in nature (viz. clinical attachment level, or CAL, measured at tooth-sites and clustered within subjects) that are routinely analyzed under a linear mixed model framework, with underlying normality assumptions of the random effects and random errors. However, a careful look reveals that these data might exhibit skewness and tail behavior, and hence the usual normality assumptions might be questionable. Besides, PD progression is often hypothesized to be spatially associated, that is, a diseased tooth-site may influence the disease status of a set of neighboring sites. Also, the presence/absence of a tooth is informative, as the number and location of missing teeth informs about the periodontal health in that region. In this paper, we develop a (shared) random effects model for site-level CAL and binary presence/absence status of a tooth under a Bayesian paradigm. The random effects are modeled using a spatial skew-normal/independent (S-SNI) distribution, whose dependence structure is conditionally autoregressive (CAR). Our S-SNI density presents an attractive parametric tool to model spatially referenced asymmetric thick-tailed structures. Both simulation studies and application to a clinical dataset recording PD status reveal the advantages of our proposition in providing a significantly improved fit, over models that do not consider these features in a unified way.

Robust calibration of hierarchical population models for heterogeneous cell populations.

Cellular heterogeneity is known to have important effects on signal processing and cellular decision making. To understand these processes, multiple classes of mathematical models have been introduced. The hierarchical population model builds a novel class which allows for the mechanistic description of heterogeneity and explicitly takes into account subpopulation structures. However, this model requires a parametric distribution assumption for the cell population and, so far, only the normal distribution has been employed. Here, we incorporate alternative distribution assumptions into the model, assess their robustness against outliers and evaluate their influence on the performance of model calibration in a simulation study and a real-world application example. We found that alternative distributions provide reliable parameter estimates even in the presence of outliers, and can in fact increase the convergence of model calibration.

Skew-normal random-effects model for meta-analysis of diagnostic test accuracy (DTA) studies.

Hierarchical models are recommended for meta-analyzing diagnostic test accuracy (DTA) studies. The bivariate random-effects model is currently widely used to synthesize a pair of test sensitivity and specificity using logit transformation across studies. This model assumes a bivariate normal distribution for the random-effects. However, this assumption is restrictive and can be violated. When the assumption fails, inferences could be misleading. In this paper, we extended the current bivariate random-effects model by assuming a flexible bivariate skew-normal distribution for the random-effects in order to robustly model logit sensitivities and logit specificities. The marginal distribution of the proposed model is analytically derived so that parameter estimation can be performed using standard likelihood methods. The method of weighted-average is adopted to estimate the overall logit-transformed sensitivity and specificity. An extensive simulation study is carried out to investigate the performance of the proposed model compared to other standard models. Overall, the proposed model performs better in terms of confidence interval width of the average logit-transformed sensitivity and specificity compared to the standard bivariate linear mixed model and bivariate generalized linear mixed model. Simulations have also shown that the proposed model performed better than the well-established bivariate linear mixed model in terms of bias and comparable with regards to the root mean squared error (RMSE) of the between-study (co)variances. The proposed method is also illustrated using a published meta-analysis data.

SMARTp: A SMART design for nonsurgical treatments of chronic periodontitis with spatially referenced and nonrandomly missing skewed outcomes.

This paper proposes dynamic treatment regimes (DTRs) as effective individualized treatment strategies for managing chronic periodontitis. The proposed DTRs are studied via SMARTp-a two-stage sequential multiple assignment randomized trial (SMART) design. For this design, we propose a statistical analysis plan and a novel cluster-level sample size calculation method that factors in typical features of periodontal responses such as non-Gaussianity, spatial clustering, and nonrandom missingness. Here, each patient is viewed as a cluster, and a tooth within a patient's mouth is viewed as an individual unit inside the cluster, with the tooth-level covariance structure described by a conditionally autoregressive structure. To accommodate possible skewness and tail behavior, the tooth-level clinical attachment level (CAL) response is assumed to be skew-t, with the nonrandomly missing structure captured via a shared parameter model corresponding to the missingness indicator. The proposed method considers mean comparison for the regimes with or without sharing an initial treatment, where the expected values and corresponding variances or covariance for the sample means of a pair of DTRs are derived by the inverse probability weighting and method of moments. Simulation studies are conducted to investigate the finite-sample performance of the proposed sample size formulas under a variety of outcome-generating scenarios. An R package SMARTp implementing our sample size formula is available at the Comprehensive R Archive Network for free download.

A monotone data augmentation algorithm for multivariate nonnormal data: With applications to controlled imputations for longitudinal trials.

An efficient monotone data augmentation (MDA) algorithm is proposed for missing data imputation for incomplete multivariate nonnormal data that may contain variables of different types and are modeled by a sequence of regression models including the linear, binary logistic, multinomial logistic, proportional odds, Poisson, negative binomial, skew-normal, skew-t regressions, or a mixture of these models. The MDA algorithm is applied to the sensitivity analyses of longitudinal trials with nonignorable dropout using the controlled pattern imputations that assume the treatment effect reduces or disappears after subjects in the experimental arm discontinue the treatment. We also describe a heuristic approach to implement the controlled imputation, in which the fully conditional specification method is used to impute the intermediate missing data to create a monotone missing pattern, and the missing data after dropout are then imputed according to the assumed nonignorable mechanisms. The proposed methods are illustrated by simulation and real data analyses. Sample SAS code for the analyses is provided in the supporting information.

Influence analysis for skew-normal semiparametric joint models of multivariate longitudinal and multivariate survival data.

The normality assumption of measurement error is a widely used distribution in joint models of longitudinal and survival data, but it may lead to unreasonable or even misleading results when longitudinal data reveal skewness feature. This paper proposes a new joint model for multivariate longitudinal and multivariate survival data by incorporating a nonparametric function into the trajectory function and hazard function and assuming that measurement errors in longitudinal measurement models follow a skew-normal distribution. A Monte Carlo Expectation-Maximization (EM) algorithm together with the penalized-splines technique and the Metropolis-Hastings algorithm within the Gibbs sampler is developed to estimate parameters and nonparametric functions in the considered joint models. Case deletion diagnostic measures are proposed to identify the potential influential observations, and an extended local influence method is presented to assess local influence of minor perturbations. Simulation studies and a real example from a clinical trial are presented to illustrate the proposed methodologies. Copyright © 2017 John Wiley & Sons, Ltd.

Determinants of transient and persistent hospital efficiency: The case of Italy.

In this paper, we extend the 4-random-component closed skew-normal stochastic frontier model by including exogenous determinants of hospital persistent (long-run) and transient (short-run) inefficiency, separated from unobserved heterogeneity. We apply this new model to a dataset composed by 133 Italian hospitals during the period 2008-2013. We show that average total inefficiency is about 23%, higher than previous estimates; hence, a model where the different types of inefficiency and hospital unobserved characteristics are not confounded allows us to get less biased estimates of hospital inefficiency. Moreover, we find that transient efficiency is more important than persistent efficiency, as it accounts for 60% of the total one. Last, we find that ownership (for-profit hospitals are more transiently inefficient and less persistently inefficient than not-for-profit ones, whereas public hospitals are less transiently inefficient than not-for-profit ones), specialization (specialized hospitals are more transiently inefficient than general ones; i.e., there is evidence of scope economies in short-run efficiency), and size (large-sized hospitals are better than medium and small ones in terms of transient inefficiency) are determinants of both types of inefficiency, although we do not find any statistically significant effect of multihospital systems and teaching hospitals.

Bayesian quantile regression for nonlinear mixed-effects joint models for longitudinal data in the presence of mismeasured covariate errors.

Quantile regression (QR) models have recently received increasing attention in longitudinal studies where measurements of the same individuals are taken repeatedly over time. When continuous (longitudinal) responses follow a distribution that is quite different from a normal distribution, usual mean regression (MR)-based linear models may fail to produce efficient estimators, whereas QR-based linear models may perform satisfactorily. To the best of our knowledge, there have been very few studies on QR-based nonlinear models for longitudinal data in comparison to MR-based nonlinear models. In this article, we study QR-based nonlinear mixed-effects (NLME) joint models for longitudinal data with non-central location and outliers and/or heavy tails in response, and non-normality and measurement errors in covariate under Bayesian framework. The proposed QR-based modeling method is compared with an MR-based one by an AIDS clinical dataset and through simulation studies. The proposed QR joint modeling approach can be not only applied to AIDS clinical studies, but also may have general applications in other fields as long as relevant technical specifications are met.

Testing multiple hypotheses with skewed alternatives.

In many practical cases of multiple hypothesis problems, it can be expected that the alternatives are not symmetrically distributed. If it is known a priori that the distributions of the alternatives are skewed, we show that this information yields high power procedures as compared to the procedures based on symmetric alternatives when testing multiple hypotheses. We propose a Bayesian decision theoretic rule for multiple directional hypothesis testing, when the alternatives are distributed as skewed, under a constraint on a mixed directional false discovery rate. We compare the proposed rule with a frequentist's rule of Benjamini and Yekutieli (2005) using simulations. We apply our method to a well-studied HIV dataset.

Continuous endpoints in Bayesian two-stage designs.

Single-arm two-stage designs for phase II of clinical trials typically focus on a binary endpoint obtained by dichotomizing an underlying continuous measure of treatment efficacy. To avoid the resulting loss of information, we propose a two-stage design based on a Bayesian predictive approach that directly uses the original continuous endpoint. Numerical results are provided with reference to phase II cancer trials aimed at assessing tumor shrinking effect of an experimental treatment.

Modeling error distributions of growth curve models through Bayesian methods.

Growth curve models are widely used in social and behavioral sciences. However, typical growth curve models often assume that the errors are normally distributed although non-normal data may be even more common than normal data. In order to avoid possible statistical inference problems in blindly assuming normality, a general Bayesian framework is proposed to flexibly model normal and non-normal data through the explicit specification of the error distributions. A simulation study shows when the distribution of the error is correctly specified, one can avoid the loss in the efficiency of standard error estimates. A real example on the analysis of mathematical ability growth data from the Early Childhood Longitudinal Study, Kindergarten Class of 1998-99 is used to show the application of the proposed methods. Instructions and code on how to conduct growth curve analysis with both normal and non-normal error distributions using the the MCMC procedure of SAS are provided.

Semiparametric Bayesian inference on skew-normal joint modeling of multivariate longitudinal and survival data.

We propose a semiparametric multivariate skew-normal joint model for multivariate longitudinal and multivariate survival data. One main feature of the posited model is that we relax the commonly used normality assumption for random effects and within-subject error by using a centered Dirichlet process prior to specify the random effects distribution and using a multivariate skew-normal distribution to specify the within-subject error distribution and model trajectory functions of longitudinal responses semiparametrically. A Bayesian approach is proposed to simultaneously obtain Bayesian estimates of unknown parameters, random effects and nonparametric functions by combining the Gibbs sampler and the Metropolis-Hastings algorithm. Particularly, a Bayesian local influence approach is developed to assess the effect of minor perturbations to within-subject measurement error and random effects. Several simulation studies and an example are presented to illustrate the proposed methodologies.