Search

Despite multivariate spatio-temporal counts often containing many zeroes, zero-inflated multinomial models for space-time data have not been considered. We are interested in comparing the transmission dynamics of several co-circulating infectious diseases across space and time where some can be absent for long periods. We first assume there is a baseline disease that is well-established and always present in the region. The other diseases switch between periods of presence and absence in each area through a series of coupled Markov chains, which account for long periods of disease absence, disease interactions and disease spread from neighboring areas. Since we are mainly interested in comparing the diseases, we assume the cases of the present diseases in an area jointly follow an autoregressive multinomial model. We use the multinomial model to investigate whether there are associations between certain factors, such as temperature, and differences in the transmission intensity of the diseases. Inference is performed using efficient Bayesian Markov chain Monte Carlo methods based on jointly sampling all presence indicators. We apply the model to spatio-temporal counts of dengue, Zika and chikungunya cases in Rio de Janeiro, during the first triple epidemic.

Spatial confounding, often regarded as a major concern in epidemiological studies, relates to the difficulty of recovering the effect of an exposure on an outcome when these variables are associated with unobserved factors. This issue is particularly challenging in spatio-temporal analyses, where it has been less explored so far. To study the effects of air pollution on mortality in Italy, we argue that a model that simultaneously accounts for spatio-temporal confounding and for the non-linear form of the effect of interest is needed. To this end, we propose a Bayesian dynamic generalized linear model, which allows for a non-linear association and for a decomposition of the exposure effect into two components. This decomposition accommodates associations with the outcome at fine and coarse temporal and spatial scales of variation. These features, when combined, allow reducing the spatio-temporal confounding bias and recovering the true shape of the association, as demonstrated through simulation studies. The results from the real-data application indicate that the exposure effect seems to have different magnitudes in different seasons, with peaks in the summer. We hypothesize that this could be due to possible interactions between the exposure variable with air temperature and unmeasured confounders., Comment: 40 pages, 11 figures

Spatio-temporal disease mapping models are commonly used to estimate the relative risk of a disease over time and across areas. For each area and time point, the disease count is modelled with a Poisson distribution whose mean is the product of an offset and the disease relative risk. This relative risk is commonly decomposed in the log scale as the sum of fixed and latent effects. The Rushworth model allows for spatio-temporal autocorrelation of the random effects. We build on the Rushworth model to accommodate and identify potentially outlying areas with respect to their disease relative risk evolution, after taking into account the fixed effects. An area may display outlying behaviour at some points in time but not all. At each time point, we assume the latent effects to be spatially structured and include scaling parameters in the precision matrix, to allow for heavy-tails. Two prior specifications are considered for the scaling parameters: one where they are independent across space and one with spatial autocorrelation. We investigate the performance of the different prior specifications of the proposed model through simulation studies and analyse the weekly evolution of the number of COVID-19 cases across the 33 boroughs of Montreal and the 96 French departments during the second wave. In Montreal, 6 boroughs are found to be potentially outlying. In France, the model with spatially structured scaling parameters identified 21 departments as potential outliers. We find that these departments tend to be close to each other and within common French regions.

This paper proposes a new approach to address the problem of unmeasured confounding in spatial designs. Spatial confounding occurs when some confounding variables are unobserved and not included in the model, leading to distorted inferential results about the effect of an exposure on an outcome. We show the relationship existing between the confounding bias of a non-spatial model and that of a semi-parametric model that includes a basis matrix to represent the unmeasured confounder conditional on the exposure. This relationship holds for any basis expansion, however it is shown that using the semi-parametric approach guarantees a reduction in the confounding bias only under certain circumstances, which are related to the spatial structures of the exposure and the unmeasured confounder, the type of basis expansion utilized, and the regularization mechanism. To adjust for spatial confounding, and therefore try to recover the effect of interest, we propose a Bayesian semi-parametric regression model, where an expansion matrix of principal spline basis functions is used to approximate the unobserved factor, and spike-and-slab priors are imposed on the respective expansion coefficients in order to select the most important bases. From the results of an extensive simulation study, we conclude that our proposal is able to reduce the confounding bias with respect to the non-spatial model, and it also seems more robust to bias amplification than competing approaches., Comment: 30 pages, 5 figures, 1 table

In the last two decades, the linear model of coregionalization (LMC) has been widely used to model multivariate spatial processes. From a computational standpoint, the LMC is a substantially easier model to work with than other multidimensional alternatives. Up to now, this fact has been largely overlooked in the literature. Starting from an analogy with matrix normal models, we propose a reformulation of the LMC likelihood that highlights the linear, rather than cubic, computational complexity as a function of the dimension of the response vector. Further, we describe in detail how those simplifications can be included in a Gaussian hierarchical model. In addition, we demonstrate in two examples how the disentangled version of the likelihood we derive can be exploited to improve Markov chain Monte Carlo (MCMC) based computations when conducting Bayesian inference. The first is an interwoven approach that combines samples from centered and whitened parametrizations of the latent LMC distributed random fields. The second is a sparsity-inducing method that introduces structural zeros in the coregionalization matrix in an attempt to reduce the number of parameters in a principled way. It also provides a new way to investigate the strength of the correlation among the components of the outcome vector. Both approaches come at virtually no additional cost and are shown to significantly improve MCMC performance and predictive performance respectively. We apply our methodology to a dataset comprised of air pollutant measurements in the state of California.

Accurately forecasting bus travel time and passenger occupancy with uncertainty is essential for both travelers and transit agencies/operators. However, existing approaches to forecasting bus travel time and passenger occupancy mainly rely on deterministic models, providing only point estimates. In this paper, we develop a Bayesian Markov regime-switching vector autoregressive model to jointly forecast both bus travel time and passenger occupancy with uncertainty. The proposed approach naturally captures the intricate interactions among adjacent buses and adapts to the multimodality and skewness of real-world bus travel time and passenger occupancy observations. We develop an efficient Markov chain Monte Carlo (MCMC) sampling algorithm to approximate the resultant joint posterior distribution of the parameter vector. With this framework, the estimation of downstream bus travel time and passenger occupancy is transformed into a multivariate time series forecasting problem conditional on partially observed outcomes. Experimental validation using real-world data demonstrates the superiority of our proposed model in terms of both predictive means and uncertainty quantification compared to the Bayesian Gaussian mixture model.

In epidemiological studies, zero-inflated and hurdle models are commonly used to handle excess zeros in reported infectious disease cases. However, they can not model the persistence (from presence to presence) and reemergence (from absence to presence) of a disease separately. Covariates can sometimes have different effects on the reemergence and persistence of a disease. Recently, a zero-inflated Markov switching negative binomial model was proposed to accommodate this issue. We present a Markov switching negative binomial hurdle model as a competitor of that approach, as hurdle models are often also used as alternatives to zero-inflated models for accommodating excess zeroes. We begin the comparison by inspecting the underlying assumptions made by both models. Hurdle models assume perfect detection of the disease cases while zero-inflated models implicitly assume the case counts can be under-reported, thus we investigate when a negative binomial distribution can approximate the true distribution of reported counts. A comparison of the fit of the two types of Markov switching models is undertaken on chikungunya cases across the neighborhoods of Rio de Janeiro. We find that, among the fitted models, the Markov switching negative binomial zero-inflated model produces the best predictions and both Markov switching models produce remarkably better predictions than more traditional negative binomial hurdle and zero-inflated models.

We consider random forests and LASSO methods for model-based small area estimation when the number of areas with sampled data is a small fraction of the total areas for which estimates are required. Abundant auxiliary information is available for the sampled areas, from the survey, and for all areas, from an exterior source, and the goal is to use auxiliary variables to predict the outcome of interest. We compare areal-level random forests and LASSO approaches to a frequentist forward variable selection approach and a Bayesian shrinkage method. Further, to measure the uncertainty of estimates obtained from random forests and the LASSO, we propose a modification of the split conformal procedure that relaxes the assumption of identically distributed data. This work is motivated by Ghanaian data available from the sixth Living Standard Survey (GLSS) and the 2010 Population and Housing Census. We estimate the areal mean household log consumption using both datasets. The outcome variable is measured only in the GLSS for 3\% of all the areas (136 out of 5019) and more than 170 potential covariates are available from both datasets. Among the four modelling methods considered, the Bayesian shrinkage performed the best in terms of bias, MSE and prediction interval coverages and scores, as assessed through a cross-validation study. We find substantial between-area variation, the log consumption areal point estimates showing a 1.3-fold variation across the GAMA region. The western areas are the poorest while the Accra Metropolitan Area district gathers the richest areas., Comment: 40 pages, 10 figures

We propose and discuss a Bayesian procedure to estimate the average treatment effect (ATE) for multilevel observations in the presence of confounding. We focus on situations where the confounders may be latent (e.g., spatial latent effects). This work is motivated by an interest in determining the causal impact of directly observed therapy (DOT) on the successful treatment of Tuberculosis (TB); the available data correspond to individual-level information observed across different cities in a state in Brazil. We focus on propensity score regression and covariate adjustment to balance the treatment (DOT) allocation. We discuss the need to include latent local-level random effects in the propensity score model to reduce bias in the estimation of the ATE. A simulation study suggests that accounting for the multilevel nature of the data with latent structures in both the outcome and propensity score models has the potential to reduce bias in the estimation of causal effects.

We propose a multivariate GARCH model for non-stationary health time series by modifying the variance of the observations of the standard state space model. The proposed model provides an intuitive way of dealing with heteroskedastic data using the conditional nature of state space models. We follow the Bayesian paradigm to perform the inference procedure. In particular, we use Markov chain Monte Carlo methods to obtain samples from the resultant posterior distribution. Due to the natural temporal correlation structure induced on model parameters, we use the forward filtering backward sampling algorithm to efficiently obtain samples from the posterior distribution. The proposed model also handles missing data in a fully Bayesian fashion. We validate our model on synthetic data, and then use it to analyze a data set obtained from an intensive care unit in a Montreal hospital. We further show that our proposed models offer better performance, in terms of WAIC, than standard state space models. The proposed model provides a new way to model multivariate heteroskedastic non-stationary time series data and the simplicity in applying the WAIC allows us to compare competing models.

Recurrent COVID-19 outbreaks have placed immense strain on the hospital system in Quebec. We develop a Bayesian three-state coupled Markov switching model to analyze COVID-19 outbreaks across Quebec based on admissions in the 30 largest hospitals. Within each catchment area, we assume the existence of three states for the disease: absence, a new state meant to account for many zeroes in some of the smaller areas, endemic and outbreak. Then we assume the disease switches between the three states in each area through a series of coupled nonhomogeneous hidden Markov chains. Unlike previous approaches, the transition probabilities may depend on covariates and the occurrence of outbreaks in neighboring areas, to account for geographical outbreak spread. Additionally, to prevent rapid switching between endemic and outbreak periods we introduce clone states into the model which enforce minimum endemic and outbreak durations. We make some interesting findings, such as that mobility in retail and recreation venues had a positive association with the development and persistence of new COVID-19 outbreaks in Quebec. Based on model comparison our contributions show promise in improving state estimation retrospectively and in real-time, especially when there are smaller areas and highly spatially synchronized outbreaks. Furthermore, our approach offers new and interesting epidemiological interpretations, such as being able to estimate the effect of covariates on disease extinction., Comment: The accepted version in The Annals of Applied Statistics

For many infectious disease outbreaks, the at-risk population changes their behavior in response to the outbreak severity, causing the transmission dynamics to change in real-time. Behavioral change is often ignored in epidemic modeling efforts, making these models less useful than they could be. We address this by introducing a novel class of data-driven epidemic models which characterize and accurately estimate behavioral change. Our proposed model allows time-varying transmission to be captured by the level of "alarm" in the population, with alarm specified as a function of the past epidemic trajectory. We investigate the estimability of the population alarm across a wide range of scenarios, applying both parametric functions and non-parametric functions using splines and Gaussian processes. The model is set in the data-augmented Bayesian framework to allow estimation on partially observed epidemic data. The benefit and utility of the proposed approach is illustrated through applications to data from real epidemics., Comment: 20 pages, 10 figures

Many models for point process data are defined through a thinning procedure where locations of a base process (often Poisson) are either kept (observed) or discarded (thinned). In this paper, we go back to the fundamentals of the distribution theory for point processes and provide a colouring theorem that characterizes the joint density of thinned and observed locations in any of such models. In practice, the marginal model of observed points is often intractable, but thinned locations can be instantiated from their conditional distribution and typical data augmentation schemes can be employed to circumvent this problem. Such approaches have been employed in recent publications, but conceptual flaws have been introduced in this literature. We concentrate on an example: the so-called sigmoidal Gaussian Cox process. We apply our general theory to resolve what are contradicting viewpoints in the data augmentation step of the inference procedures therein. Finally, we provide a multitype extension to this process and conduct Bayesian inference on data consisting of positions of 2 different species of trees in Lansing Woods, Illinois.

We study Bayesian approaches to causal inference via propensity score regression. Much of the Bayesian literature on propensity score methods have relied on approaches that cannot be viewed as fully Bayesian in the context of conventional `likelihood times prior' posterior inference; in addition, most methods rely on parametric and distributional assumptions, and presumed correct specification. We emphasize that causal inference is typically carried out in settings of mis-specification, and develop strategies for fully Bayesian inference that reflect this. We focus on methods based on decision-theoretic arguments, and show how inference based on loss-minimization can give valid and fully Bayesian inference. We propose a computational approach to inference based on the Bayesian bootstrap which has good Bayesian and frequentist properties.

Spatio-temporal processes in environmental applications are often assumed to follow a Gaussian model, possibly after some transformation. However, heterogeneity in space and time might have a pattern that will not be accommodated by transforming the data. In this scenario, modelling the variance laws is an appealing alternative. This work adds flexibility to the usual Multivariate Dynamic Gaussian model by defining the process as a scale mixture between a Gaussian and log-Gaussian processes. The scale is represented by a process varying smoothly over space and time which is allowed to depend on covariates. State-space equations define the dynamics over time for both mean and variance processes resulting infeasible inference and prediction. Analysis of artificial datasets show that the parameters are identifiable and simpler models are well recovered by the general proposed model. The analyses of two important environmental processes, maximum temperature and maximum ozone, illustrate the effectiveness of our proposal in improving the uncertainty quantification in the prediction of spatio-temporal processes.

In disease mapping, the relative risk of a disease is commonly estimated across different areas within a region of interest. The number of cases in an area is often assumed to follow a Poisson distribution whose mean is decomposed as the product between an offset and the logarithm of the disease's relative risk. The log risk may be written as the sum of fixed effects and latent random effects. The BYM2 model decomposes each latent effect into a weighted sum of independent and spatial effects. We build on the BYM2 model to allow for heavy-tailed latent effects and accommodate potentially outlying risks, after accounting for the fixed effects. We assume a scale mixture structure wherein the variance of the latent process changes across areas and allows for outlier identification. We propose two prior specifications for this scale mixture parameter. These are compared through simulation studies and in the analysis of Zika cases from the first (2015-2016) epidemic in Rio de Janeiro city, Brazil. The simulation studies show that, in terms of the model assessment criterion WAIC and outlier detection, the two proposed parametrisations perform better than the model proposed by Congdon (2017) to capture outliers. In particular, the proposed parametrisations are more efficient, in terms of outlier detection, than Congdon's when outliers are neighbours. Our analysis of Zika cases finds 19 out of 160 districts of Rio as potential outliers, after accounting for the socio-development index. Our proposed model may help prioritise interventions and identify potential issues in the recording of cases., Comment: 50 pages, 33 figures, 7 tables

Spatio-temporal counts of infectious disease cases often contain an excess of zeros. With existing zero inflated count models applied to such data it is difficult to quantify space-time heterogeneity in the effects of disease spread between areas. Also, existing methods do not allow for separate dynamics to affect the reemergence and persistence of the disease. As an alternative, we develop a new zero-state coupled Markov switching negative binomial model, under which the disease switches between periods of presence and absence in each area through a series of partially hidden nonhomogeneous Markov chains coupled between neighboring locations. When the disease is present, an autoregressive negative binomial model generates the cases with a possible 0 representing the disease being undetected. Bayesian inference and prediction is illustrated using spatio-temporal counts of dengue fever cases in Rio de Janeiro, Brazil.

The class of $\alpha$-stable distributions enjoys multiple practical applications in signal processing, finance, biology and other areas because it allows to describe interesting and complex data patterns, such as asymmetry or heavy tails, in contrast with the simpler and widely used Gaussian distribution. The density associated with a general $\alpha$-stable distribution cannot be obtained in closed form, which hinders the process of estimating its parameters. A nonlinear population Monte Carlo (NPMC) scheme is applied in order to approximate the posterior probability distribution of the parameters of an $\alpha$-stable random variable given a set of random realizations of the latter. The approximate posterior distribution is computed by way of an iterative algorithm and it consists of a collection of samples in the parameter space with associated nonlinearly-transformed importance weights. A numerical comparison of the main existing methods to estimate the $\alpha$-stable parameters is provided, including the traditional frequentist techniques as well as a Markov chain Monte Carlo (MCMC) and a likelihood-free Bayesian approach. It is shown by means of computer simulations that the NPMC method outperforms the existing techniques in terms of parameter estimation error and failure rate for the whole range of values of $\alpha$, including the smaller values for which most existing methods fail to work properly. Furthermore, it is shown that accurate parameter estimates can often be computed based on a low number of observations. Additionally, numerical results based on a set of real fish displacement data are provided., Comment: To appear in Computational Statistics and Data Analysis

The Brazilian Mathematical Olympiads for Public Schools (OBMEP) is held every year since 2005. In the 2013 edition there were over 47,000 schools registered involving nearly 19.2 million students. The Brazilian public educational system is structured into three administrative levels: federal, state and municipal. Students participating in the OBMEP come from three educational levels, two in primary and one in secondary school. We aim at studying the performance of Brazilian public schools which have been taking part of the OBMEP from 2006 until 2013. We propose a standardization of the mean scores of schools per year and educational level which is modeled through a hierarchical dynamic beta regression model. Both the mean and precision of the beta distribution are modeled as a function of covariates whose effects evolve smoothly with time. Results show that, regardless of the educational level, federal schools have better performance than municipal or state schools. The mean performance of schools increases with the human development index (HDI) of the municipality the school is located in. Moreover, the difference in mean performance between federal and state or municipal schools tends to increase with the HDI. Schools with higher proportion of boys tend to have better mean performance in the second and third educational levels of OBMEP.