Search

Growing anthropogenic pressures have increased the need for robust predictive models. Meeting this demand requires approaches that can handle bigger data to yield forecasts that capture the variability and underlying uncertainty of ecological systems. Bayesian models are especially adept at this and are growing in use in ecology. Yet many ecologists today are not trained to take advantage of the bigger ecological data needed to generate more flexible robust models. Here we describe a broadly generalizable workflow for statistical analyses and show how it can enhance training in ecology. Building on the increasingly computational toolkit of many ecologists, this approach leverages simulation to integrate model building and testing for empirical data more fully with ecological theory. In turn this workflow can fit models that are more robust and well-suited to provide new ecological insights -- allowing us to refine where to put resources for better estimates, better models, and better forecasts.

Survival models capture the relationship between an accumulating hazard and the occurrence of a singular event stimulated by that accumulation. When the model for the hazard is sufficiently flexible survival models can accommodate a wide range of behaviors. If the hazard model is less flexible, for example when it is constrained by an external physical process, then the resulting survival model can be much too rigid. In this paper I introduce a modified survival model that generalizes the relationship between accumulating hazard and event occurrence with particular emphasis on capturing thresholding behavior. Finally I demonstrate the utility of this approach on a physiological application., Comment: 28 pages, 9 figures

Derivative-based algorithms are ubiquitous in statistics, machine learning, and applied mathematics. Automatic differentiation offers an algorithmic way to efficiently evaluate these derivatives from computer programs that execute relevant functions. Implementing automatic differentiation for programs that incorporate implicit functions, such as the solution to an algebraic or differential equation, however, requires particular care. Contemporary applications typically appeal to either the application of the implicit function theorem or, in certain circumstances, specialized adjoint methods. In this paper we show that both of these approaches can be generalized to any implicit function, although the generalized adjoint method is typically more effective for automatic differentiation. To showcase the relative advantages and limitations of the two methods we demonstrate their application on a suite of common implicit functions.

This short note reviews the basic theory for quantifying both the asymptotic and preasymptotic convergence of Markov chain Monte Carlo estimators., Comment: 19 pages, 10 figures

A complete recipe of measure-preserving diffusions in Euclidean space was recently derived unifying several MCMC algorithms into a single framework. In this paper, we develop a geometric theory that improves and generalises this construction to any manifold. We thereby demonstrate that the completeness result is a direct consequence of the topology of the underlying manifold and the geometry induced by the target measure $P$; there is no need to introduce other structures such as a Riemannian metric, local coordinates, or a reference measure. Instead, our framework relies on the intrinsic geometry of $P$ and in particular its canonical derivative, the deRham rotationnel, which allows us to parametrise the Fokker--Planck currents of measure-preserving diffusions using potentials. The geometric formalism can easily incorporate constraints and symmetries, and deliver new important insights, for example, a new complete recipe of Langevin-like diffusions that are suited to the construction of samplers. We also analyse the reversibility and dissipative properties of the diffusions, the associated deterministic flow on the space of measures, and the geometry of Langevin processes. Our article connects ideas from various literature and frames the theory of measure-preserving diffusions in its appropriate mathematical context.

Inferences about hypotheses are ubiquitous in the cognitive sciences. Bayes factors provide one general way to compare different hypotheses by their compatibility with the observed data. Those quantifications can then also be used to choose between hypotheses. While Bayes factors provide an immediate approach to hypothesis testing, they are highly sensitive to details of the data/model assumptions. Moreover it's not clear how straightforwardly this approach can be implemented in practice, and in particular how sensitive it is to the details of the computational implementation. Here, we investigate these questions for Bayes factor analyses in the cognitive sciences. We explain the statistics underlying Bayes factors as a tool for Bayesian inferences and discuss that utility functions are needed for principled decisions on hypotheses. Next, we study how Bayes factors misbehave under different conditions. This includes a study of errors in the estimation of Bayes factors. Importantly, it is unknown whether Bayes factor estimates based on bridge sampling are unbiased for complex analyses. We are the first to use simulation-based calibration as a tool to test the accuracy of Bayes factor estimates. Moreover, we study how stable Bayes factors are against different MCMC draws. We moreover study how Bayes factors depend on variation in the data. We also look at variability of decisions based on Bayes factors and how to optimize decisions using a utility function. We outline a Bayes factor workflow that researchers can use to study whether Bayes factors are robust for their individual analysis, and we illustrate this workflow using an example from the cognitive sciences. We hope that this study will provide a workflow to test the strengths and limitations of Bayes factors as a way to quantify evidence in support of scientific hypotheses. Reproducible code is available from https://osf.io/y354c/.

Gradient-based techniques are becoming increasingly critical in quantitative fields, notably in statistics and computer science. The utility of these techniques, however, ultimately depends on how efficiently we can evaluate the derivatives of the complex mathematical functions that arise in applications. In this paper we introduce a discrete adjoint method that efficiently evaluates derivatives for functions of discrete sequences., Comment: 16 pages, 1 figure

Recent research on accelerated gradient methods of use in optimization has demonstrated that these methods can be derived as discretizations of dynamical systems. This, in turn, has provided a basis for more systematic investigations, especially into the geometric structure of those dynamical systems and their structure--preserving discretizations. In this work, we introduce dynamical systems defined through a contact geometry which are not only naturally suited to the optimization goal but also subsume all previous methods based on geometric dynamical systems. As a consequence, all the deterministic flows used in optimization share an extremely interesting geometric property: they are invariant under contact transformations. In our main result, we exploit this observation to show that the celebrated Bregman Hamiltonian system can always be transformed into an equivalent but separable Hamiltonian by means of a contact transformation. This in turn enables the development of fast and robust discretizations through geometric contact splitting integrators. As an illustration, we propose the Relativistic Bregman algorithm, and show in some paradigmatic examples that it compares favorably with respect to standard optimization algorithms such as classical momentum and Nesterov's accelerated gradient., Comment: 15 pages

Reparameterizing a probabilisitic system is common advice for improving the performance of a statistical algorithm like Markov chain Monte Carlo, even though in theory such reparameterizations should leave the system, and the performance of any algorithm, invariant. In this paper I show how the reparameterizations common in practice are only incomplete reparameterizations which result in different interactions between a target probabilistic system and a given algorithm. I then consider how these changing interactions manifest in the context of Markov chain Monte Carlo algorithms defined on Riemannian manifolds. In particular I show how any incomplete reparameterization is equivalent to modifying the metric geometry directly., Comment: 34 pages, 15 figures

Experiments in research on memory, language, and in other areas of cognitive science are increasingly being analyzed using Bayesian methods. This has been facilitated by the development of probabilistic programming languages such as Stan, and easily accessible front-end packages such as brms. The utility of Bayesian methods, however, ultimately depends on the relevance of the Bayesian model, in particular whether or not it accurately captures the structure of the data and the data analyst's domain expertise. Even with powerful software, the analyst is responsible for verifying the utility of their model. To demonstrate this point, we introduce a principled Bayesian workflow (Betancourt, 2018) to cognitive science. Using a concrete working example, we describe basic questions one should ask about the model: prior predictive checks, computational faithfulness, model sensitivity, and posterior predictive checks. The running example for demonstrating the workflow is data on reading times with a linguistic manipulation of object versus subject relative clause sentences. This principled Bayesian workflow also demonstrates how to use domain knowledge to inform prior distributions. It provides guidelines and checks for valid data analysis, avoiding overfitting complex models to noise, and capturing relevant data structure in a probabilistic model. Given the increasing use of Bayesian methods, we aim to discuss how these methods can be properly employed to obtain robust answers to scientific questions. All data and code accompanying this paper are available from https://osf.io/b2vx9/., Comment: 75 pages, 19 figures

First-order automatic differentiation is a ubiquitous tool across statistics, machine learning, and computer science. Higher-order implementations of automatic differentiation, however, have yet to realize the same utility. In this paper I derive a comprehensive, differential geometric treatment of automatic differentiation that naturally identifies the higher-order differential operators amenable to automatic differentiation as well as explicit procedures that provide a scaffolding for high-performance implementations., Comment: 55 pages, 10 figures

The absence of an exhaustive bibliography of novels published in the British Isles and Ireland during the 19th century blocks several lines of research in sociologically-inclined literary history and book history. Without a detailed account of novelistic production, it is difficult to characterize, for example, the population of individuals who pursued careers as novelists. This paper contributes to efforts to develop such an account by estimating yearly rates of new novel publication in the British Isles and Ireland between 1789 and 1919. This period witnessed, in aggregate, the publication of between 40,000 and 63,000 previously unpublished novels. The number of new novels published each year counts as essential information for researchers interested in understanding the development of the text industry between 1789 and 1919., Comment: Cultural Analytics (2021)

Verifying the correctness of Bayesian computation is challenging. This is especially true for complex models that are common in practice, as these require sophisticated model implementations and algorithms. In this paper we introduce \emph{simulation-based calibration} (SBC), a general procedure for validating inferences from Bayesian algorithms capable of generating posterior samples. This procedure not only identifies inaccurate computation and inconsistencies in model implementations but also provides graphical summaries that can indicate the nature of the problems that arise. We argue that SBC is a critical part of a robust Bayesian workflow, as well as being a useful tool for those developing computational algorithms and statistical software., Comment: 19 pages, 13 figures

As the frontiers of applied statistics progress through increasingly complex experiments we must exploit increasingly sophisticated inferential models to analyze the observations we make. In order to avoid misleading or outright erroneous inferences we then have to be increasingly diligent in scrutinizing the consequences of those modeling assumptions. Fortunately model-based methods of statistical inference naturally define procedures for quantifying the scope of inferential outcomes and calibrating corresponding decision making processes. In this paper I review the construction and implementation of the particular procedures that arise within frequentist and Bayesian methodologies., Comment: 35 pages, 10 figures

Accelerated gradient methods have had significant impact in machine learning -- in particular the theoretical side of machine learning -- due to their ability to achieve oracle lower bounds. But their heuristic construction has hindered their full integration into the practical machine-learning algorithmic toolbox, and has limited their scope. In this paper we build on recent work which casts acceleration as a phenomenon best explained in continuous time, and we augment that picture by providing a systematic methodology for converting continuous-time dynamics into discrete-time algorithms while retaining oracle rates. Our framework is based on ideas from Hamiltonian dynamical systems and symplectic integration. These ideas have had major impact in many areas in applied mathematics, but have not yet been seen to have a relationship with optimization., Comment: 20 pages, 5 figures

Bayesian data analysis is about more than just computing a posterior distribution, and Bayesian visualization is about more than trace plots of Markov chains. Practical Bayesian data analysis, like all data analysis, is an iterative process of model building, inference, model checking and evaluation, and model expansion. Visualization is helpful in each of these stages of the Bayesian workflow and it is indispensable when drawing inferences from the types of modern, high-dimensional models that are used by applied researchers., Comment: 17 pages, 11 Figures. Includes supplementary material

A key sticking point of Bayesian analysis is the choice of prior distribution, and there is a vast literature on potential defaults including uniform priors, Jeffreys' priors, reference priors, maximum entropy priors, and weakly informative priors. These methods, however, often manifest a key conceptual tension in prior modeling: a model encoding true prior information should be chosen without reference to the model of the measurement process, but almost all common prior modeling techniques are implicitly motivated by a reference likelihood. In this paper we resolve this apparent paradox by placing the choice of prior into the context of the entire Bayesian analysis, from inference to prediction to model evaluation., Comment: 13 pages

From its inception in the 1950s to the modern frontiers of applied statistics, Markov chain Monte Carlo has been one of the most ubiquitous and successful methods in statistical computing. In that time its development has been fueled by increasingly difficult problems and novel techniques from physics. In this article I will review the history of Markov chain Monte Carlo from its inception with the Metropolis method to today's state-of-the-art in Hamiltonian Monte Carlo. Along the way I will focus on the evolving interplay between the statistical and physical perspectives of the method., Comment: 13 pages

Hamiltonian Monte Carlo has proven a remarkable empirical success, but only recently have we begun to develop a rigorous understanding of why it performs so well on difficult problems and how it is best applied in practice. Unfortunately, that understanding is confined within the mathematics of differential geometry which has limited its dissemination, especially to the applied communities for which it is particularly important. In this review I provide a comprehensive conceptual account of these theoretical foundations, focusing on developing a principled intuition behind the method and its optimal implementations rather of any exhaustive rigor. Whether a practitioner or a statistician, the dedicated reader will acquire a solid grasp of how Hamiltonian Monte Carlo works, when it succeeds, and, perhaps most importantly, when it fails., Comment: 60 pages, 42 figures

Mass is a fundamental property of galaxy groups and clusters. In theory weak gravitational lensing will enable an approximately unbiased measurement of mass, but parametric methods for extracting cluster masses from data require the additional knowledge of concentration. Measurements of both mass and concentration are limited by the degeneracy between the two parameters, particularly in low mass, high redshift systems where the signal-to-noise is low. In this paper we develop a hierarchical model of mass and concentration for mass inference we test our method on toy data and then apply it to a sample of galaxy groups and poor clusters down to masses of $\sim$ 1e13 M$_\odot$. Our fit and model gives a relationship among masses, concentrations and redshift that allow prediction of these parameters from incomplete and noisy future measurements. Additionally the underlying population can be used to infer an observationally based concentration-mass relation. Our method is equivalent to a quasi- stacking approach with the degree of stacking set by the data. We also demonstrate that mass and concentration derived from pure stacking can be offset from the population mean with differing values depending on the method of stacking., Comment: 16 pages, 10 figures, submitted to MNRAS

When properly tuned, Hamiltonian Monte Carlo scales to some of the most challenging high-dimensional problems at the frontiers of applied statistics, but when that tuning is suboptimal the performance leaves much to be desired. In this paper I show how suboptimal choices of one critical degree of freedom, the cotangent disintegration, manifest in readily observed diagnostics that facilitate the robust application of the algorithm., Comment: 17 pages, 9 figures

Throughout the different phases of a drug development program, randomized trials are used to establish the tolerability, safety, and efficacy of a candidate drug. At each stage one aims to optimize the design of future studies by extrapolation from the available evidence at the time. This includes collected trial data and relevant external data. However, relevant external data are typically available as averages only, for example from trials on alternative treatments reported in the literature. Here we report on such an example from a drug development for wet age-related macular degeneration. This disease is the leading cause of severe vision loss in the elderly. While current treatment options are efficacious, they are also a substantial burden for the patient. Hence, new treatments are under development which need to be compared against existing treatments. The general statistical problem this leads to is meta-analysis, which addresses the question of how we can combine datasets collected under different conditions. Bayesian methods have long been used to achieve partial pooling. Here we consider the challenge when the model of interest is complex (hierarchical and nonlinear) and one dataset is given as raw data while the second dataset is given as averages only. In such a situation, common meta-analytic methods can only be applied when the model is sufficiently simple for analytic approaches. When the model is too complex, for example nonlinear, an analytic approach is not possible. We provide a Bayesian solution by using simulation to approximately reconstruct the likelihood of the external summary and allowing the parameters in the model to vary under the different conditions. We first evaluate our approach using fake-data simulations and then report results for the drug development program that motivated this research., Comment: Code is available at https://github.com/wds15/baad and as supplementary material at the publisher site

We establish general conditions under which Markov chains produced by the Hamiltonian Monte Carlo method will and will not be geometrically ergodic. We consider implementations with both position-independent and position-dependent integration times. In the former case we find that the conditions for geometric ergodicity are essentially a gradient of the log-density which asymptotically points towards the centre of the space and grows no faster than linearly. In an idealised scenario in which the integration time is allowed to change in different regions of the space, we show that geometric ergodicity can be recovered for a much broader class of tail behaviours, leading to some guidelines for the choice of this free parameter in practice., Comment: 29 pages + supplement (included in arXival as Appendix), 1 figure

By leveraging the natural geometry of a smooth probabilistic system, Hamiltonian Monte Carlo yields computationally efficient Markov Chain Monte Carlo estimation. At least provided that the algorithm is sufficiently well-tuned. In this paper I show how the geometric foundations of Hamiltonian Monte Carlo implicitly identify the optimal choice of these parameters, especially the integration time. I then consider the practical consequences of these principles in both existing algorithms and a new implementation called \emph{Exhaustive Hamiltonian Monte Carlo} before demonstrating the utility of these ideas in some illustrative examples., Comment: 31 pages, 21 figures

As computational challenges in optimization and statistical inference grow ever harder, algorithms that utilize derivatives are becoming increasingly more important. The implementation of the derivatives that make these algorithms so powerful, however, is a substantial user burden and the practicality of these algorithms depends critically on tools like automatic differentiation that remove the implementation burden entirely. The Stan Math Library is a C++, reverse-mode automatic differentiation library designed to be usable, extensive and extensible, efficient, scalable, stable, portable, and redistributable in order to facilitate the construction and utilization of such algorithms. Usability is achieved through a simple direct interface and a cleanly abstracted functional interface. The extensive built-in library includes functions for matrix operations, linear algebra, differential equation solving, and most common probability functions. Extensibility derives from a straightforward object-oriented framework for expressions, allowing users to easily create custom functions. Efficiency is achieved through a combination of custom memory management, subexpression caching, traits-based metaprogramming, and expression templates. Partial derivatives for compound functions are evaluated lazily for improved scalability. Stability is achieved by taking care with arithmetic precision in algebraic expressions and providing stable, compound functions where possible. For portability, the library is standards-compliant C++ (03) and has been tested for all major compilers for Windows, Mac OS X, and Linux., Comment: 96 pages, 9 figures

The predictive performance of any inferential model is critical to its practical success, but quantifying predictive performance is a subtle statistical problem. In this paper I show how the natural structure of any inferential problem defines a canonical measure of relative predictive performance and then demonstrate how approximations of this measure yield many of the model comparison techniques popular in statistics and machine learning., Comment: 20 pages, 11 figures

Historically, light curve studies of supernovae (SNe) and other transient classes have focused on individual objects with copious and high signal-to-noise observations. In the nascent era of wide field transient searches, objects with detailed observations are decreasing as a fraction of the overall known SN population, and this strategy sacrifices the majority of the information contained in the data about the underlying population of transients. A population level modeling approach, simultaneously fitting all available observations of objects in a transient sub-class of interest, fully mines the data to infer the properties of the population and avoids certain systematic biases. We present a novel hierarchical Bayesian statistical model for population level modeling of transient light curves, and discuss its implementation using an efficient Hamiltonian Monte Carlo technique. As a test case, we apply this model to the Type IIP SN sample from the Pan-STARRS1 Medium Deep Survey, consisting of 18,837 photometric observations of 76 SNe, corresponding to a joint posterior distribution with 9,176 parameters under our model. Our hierarchical model fits provide improved constraints on light curve parameters relevant to the physical properties of their progenitor stars relative to modeling individual light curves alone. Moreover, we directly evaluate the probability for occurrence rates of unseen light curve characteristics from the model hyperparameters, addressing observational biases in survey methodology. We view this modeling framework as an unsupervised machine learning technique with the ability to maximize scientific returns from data to be collected by future wide field transient searches like LSST. \smallskip, Comment: Submitted; 10 pages, 11 figures, plus appendix with code

With its systematic exploration of probability distributions, Hamiltonian Monte Carlo is a potent Markov Chain Monte Carlo technique; it is an approach, however, ultimately contingent on the choice of a suitable Hamiltonian function. By examining both the symplectic geometry underlying Hamiltonian dynamics and the requirements of Markov Chain Monte Carlo, we construct the general form of admissible Hamiltonians and propose a particular choice with potential application in Bayesian inference., Comment: 9 pages, 1 table