Bayesian Inference for Diffusion-Driven Mixed-Effects Models

Following its introduction over 30 years ago, the Markov cohort state-transition model has been used extensively to model population trajectories over time in health decision modeling and cost-effectiveness analysis studies. We recently showed that a cohort model represents the average of a continuous-time stochastic process on a multidimensional integer lattice governed by a master equation, which represents the time-evolution of the probability function of an integer-valued random vector. By leveraging this theoretical connection, this study introduces an alternative modeling method using a stochastic differential equation (SDE) approach, which captures not only the mean behavior but also the variance of the population process. We show the derivation of an SDE model from first principles, describe an algorithm to construct an SDE and solve the SDE via simulation for use in practice, and demonstrate the two applications of an SDE in detail. The first example demonstrates that the population trajectories, and their mean and variance, from the SDE and other commonly used methods in decision modeling match. The second example shows that users can readily apply the SDE method in their existing works without the need for additional inputs beyond those required for constructing a conventional cohort model. In addition, the second example demonstrates that the SDE model is superior to a microsimulation model in terms of computational speed. In summary, an SDE model provides an alternative modeling framework which includes information on variance, can accommodate for time-varying parameters, and is computationally less expensive than a microsimulation for a typical cohort modeling problem.

Although the governing equations of many systems, when derived from first principles, may be viewed as known, it is often too expensive to numerically simulate all the interactions they describe. Therefore, researchers often seek simpler descriptions that describe complex phenomena without numerically resolving all the interacting components. Stochastic differential equations (SDEs) arise naturally as models in this context. The growth in data acquisition, both through experiment and through simulations, provides an opportunity for the systematic derivation of SDE models in many disciplines. However, inconsistencies between SDEs and real data at short time scales often cause problems, when standard statistical methodology is applied to parameter estimation. The incompatibility between SDEs and real data can be addressed by deriving sufficient statistics from the time-series data and learning parameters of SDEs based on these. Here, we study sufficient statistics computed from time averages, an approach that we demonstrate to lead to sufficient statistics on a variety of problems and that has the secondary benefit of obviating the need to match trajectories. Following this approach, we formulate the fitting of SDEs to sufficient statistics from real data as an inverse problem and demonstrate that this inverse problem can be solved by using ensemble Kalman inversion. Furthermore, we create a framework for non-parametric learning of drift and diffusion terms by introducing hierarchical, refinable parameterizations of unknown functions, using Gaussian process regression. We demonstrate the proposed methodology for the fitting of SDE models, first in a simulation study with a noisy Lorenz ’63 model, and then in other applications, including dimension reduction in deterministic chaotic systems arising in the atmospheric sciences, large-scale pattern modeling in climate dynamics and simplified models for key observables arising in molecular dynamics. The results confirm that the proposed methodology provides a robust and systematic approach to fitting SDE models to real data.

The tenth special issue of the Japanese Journal of Statistics and Data Science (JJSD) focuses on statistics for stochastic processes. With the advancement of information and communication technology and measurement technology, a vast amount of diverse time series data have been readily obtained. Consequently, the significance of stochastic processes as mathematical models for random phenomena evolving in time or space has increased dramatically, and statistics for various types of stochastic processes have been rapidly developing. In fact, a multitude of studies on statistics for stochastic processes such as ARMA, CARMA, COGARCH, diffusion processes, fractional Brownian motions, GARCH, jump-diffusion processes, Lèvy processes, Markov processes, point processes, random fields, stochastic differential equations and stochastic partial differential equations have been conducted. Nonetheless, because of the mathematical complexity of statistics for stochastic processes, numerous theoretical and computational challenges remain to be addressed. This special feature contains ten original papers contributed by leading mathematical statisticians and data scientists on current topics in statistics for stochastic processes. The first paper examines the efficient particle Markov chain Monte Carlo method. The second and third papers respectively address the estimation problems of Hawkes processes and Cox proportional excess hazard models. The fourth and fifth papers investigate parametric estimation for probability partial differential equation models. The latter five papers focus on statistical inference for stochastic differential equation models driven by various stochastic processes, including the Wiener process, fractional Brownian motion, β-stable process and Lèvy process. The ten papers in the special feature are summarized as follows: Awaya and Omori propose an efficient particle Markov chain Monte Carlo methodology for the rolling-window estimation of state space models. They provide theoretical justifications for their proposed method and demonstrate the efficiency of B Masayuki Uchida

Approximate Bayesian computation (ABC) has become one of the major tools of likelihood-free statistical inference in complex mathematical models. Simultaneously, stochastic differential equations (SDEs) have developed to an established tool for modelling time-dependent, real-world phenomena with underlying random effects. When applying ABC to stochastic models, two major difficulties arise: First, the derivation of effective summary statistics and proper distances is particularly challenging, since simulations from the stochastic process under the same parameter configuration result in different trajectories. Second, exact simulation schemes to generate trajectories from the stochastic model are rarely available, requiring the derivation of suitable numerical methods for the synthetic data generation. To obtain summaries that are less sensitive to the intrinsic stochasticity of the model, we propose to build up the statistical method (e.g. the choice of the summary statistics) on the underlying structural properties of the model. Here, we focus on the existence of an invariant measure and we map the data to their estimated invariant density and invariant spectral density. Then, to ensure that these model properties are kept in the synthetic data generation, we adopt measure-preserving numerical splitting schemes. The derived property-based and measure-preserving ABC method is illustrated on the broad class of partially observed Hamiltonian type SDEs, both with simulated data and with real electroencephalography data. The derived summaries are particularly robust to the model simulation, and this fact, combined with the proposed reliable numerical scheme, yields accurate ABC inference. In contrast, the inference returned using standard numerical methods (Euler–Maruyama discretisation) fails. The proposed ingredients can be incorporated into any type of ABC algorithm and directly applied to all SDEs that are characterised by an invariant distribution and for which a measure-preserving numerical method can be derived.

Stochastic differential equations and Copula theories are important topics that have many advantages for applications in almost every discipline. Many studies in forestry collect longitudinal, multi-dimensional, and discrete data for which the amount of measurement of individual variables does not match. For example, during sampling experiments, the diameters of all trees, the heights of approximately 10% of the trees, and the tree crown base height and crown width for a significantly smaller number of trees are measured. In this study, for estimating five-dimensional dependencies, we used a normal copula approach, where the dynamics of individual tree variables (diameter, potentially available area, height, crown base height, and crown width) are described by a stochastic differential equation with mixed-effect parameters. The approximate maximum likelihood method was used to obtain parameter estimates of the presented stochastic differential equations, and the normal copula dependence parameters were estimated using the pseudo-maximum likelihood method. This study introduced the normalized multi-dimensional interaction information index based on differential entropy to capture dependencies between state variables. Using conditional copula-type probability density functions, the exact form equations defining the links among the diameter, potentially available area, height, crown base height, and crown width were derived. All results were implemented in the symbolic algebra system MAPLE.

Data-Driven Discovery of Stochastic Differential Equations

Stochastic differential equations (SDEs) or diffusions are continuous-valued continuous-time stochastic processes widely used in the applied and mathematical sciences. Simulating paths from these processes is usually an intractable problem, and typically involves time-discretization approximations. We propose an exact Markov chain Monte Carlo sampling algorithm that involves no such time-discretization error. Our sampler is applicable to the problem of prior simulation from an SDE, posterior simulation conditioned on noisy observations, as well as parameter inference given noisy observations. Our work recasts an existing rejection sampling algorithm for a class of diffusions as a latent variable model, and then derives an auxiliary variable Gibbs sampling algorithm that targets the associated joint distribution. At a high level, the resulting algorithm involves two steps: simulating a random grid of times from an inhomogeneous Poisson process, and updating the SDE trajectory conditioned on this grid. Our work allows the vast literature of Monte Carlo sampling algorithms from the Gaussian process literature to be brought to bear to applications involving diffusions. We study our method on synthetic and real datasets, where we demonstrate superior performance over competing methods. Supplementary materials for this article are available online.

Research Highlights: Today’s approaches to modeling of forest stands are in most cases based on that the regression models and they are constructed as static sub-models describing individual stands variables. The disadvantages of this method; it is laborious because too many different equations need to be assessed and empirical choices of candidate equations make the results subjective; it does not relate to the stand variables dynamics against the age dimension (time); and does not consider the underlying covariance structure driving changes in the stand variables. In this study, the dynamical model defined by a fixed-and mixed effect parameters trivariate stochastic differential equation (SDE) is introduced and described how such a model can be used to model quadratic mean diameter, mean height, number of trees per hectare, self-thinning line, stand basal area, stand volume per hectare and much more. Background and Objectives: New developed marginal and conditional trivariate probability density functions, combining information generated from an age-dependent variance-covariance matrix of quadratic mean diameter, mean height and number of trees per hectare, improve stand growth prediction, and forecast (in forecast the future is completely unavailable and must only be estimated from historical patterns) accuracies. Materials and Methods: Fixed-and mixed effect parameters SDE models were harmonized to predict and forecast the dynamics of quadratic mean diameter, mean height, number of trees per hectare, basal area, stand volume per hectare, and their current and mean increments. The results and experience from applying the SDE concepts and techniques in an extensive whole stand growth and yield analysis are described using a Scots pine (Pinus sylvestris L.) experimental dataset in Lithuania. Results: The mixed effects scenario SDE model showed high accuracy, the percentage root mean square error values for quadratic mean diameter, mean height, number of trees per hectare, stand basal area and stand volume per hectare predictions (forecasts) were 3.37% (10.44%), 1.82% (2.07%), 1.76% (2.93%), 6.65% (10.41%) and 6.50% (8.93%), respectively. In the same way, the quadratic mean diameter, mean height, number of trees per hectare, stand basal area and stand volume per hectare prediction (forecast) relationships had high values of the coefficient of determination, R2, 0.998 (0.987), 0.997 (0.992), 0.997 (0.988), 0.968 (0.984) and 0.966 (0.980), respectively. Conclusions: The approach presented in this paper can be used for developing a new generation stand growth and yield models.

We studied isometric stochastic flows of a Stratonovich stochastic differential equation on spheres, i.e., on the standard sphere and Gromoll-Meyer exotic sphere . In this case, and are homeomorphic but not diffeomorphic. The standard sphere can be constructed as the quotient manifold with the so-called -action of S3, whereas the Gromoll-Meyer exotic sphere as the quotient manifold with respect to the so-called -action of S3. The corresponding continuous-time stochastic process and its properties on the Gromoll-Meyer exotic sphere can be obtained by constructing a homeomorphism . The stochastic flow can be regarded as the same stochastic flow on S7, but viewed in Gromoll-Meyer differential structure. The flow on and the corresponding flow on constructed in this paper have the same regularities. There is no difference between the stochastic flow's appearance on S7 viewed in standard differential structure and the appearance of the same stochastic flow viewed in the Gromoll-Meyer differential structure. Furthermore, since the inverse mapping h-1 is differentiable on , the Riemannian metric tensor on , i.e., the pull-back of the Riemannian metric tensor G on the standard sphere , is also differentiable. This fact implies, for instance, the fact that the Fokker-Planck equation associated with the stochastic flow and the Fokker-Planck equation associated with the stochastic differential equation have the same regularities provided that the function β is C1-differentiable. Therefore both differential structures on S7 give the same description of the dynamics of the distribution function of the stochastic process understudy on seven spheres.

There is a history of simple error growth models designed to capture the key properties of error growth in operational numerical weather prediction models. We propose here such a scalar model that relies on the previous ones, but captures the effect of small scales on the error growth via additive noise in a nonlinear stochastic differential equation (SDE). We nondimensionalize the equation and study its behavior with respect to the error saturation value, the growth rate of small errors, and the magnitude of noise. We show that the addition of noise can change the curvature of the error growth curve. The SDE model seems to improve substantially the fit to operational error growth curves, compared to the deterministic counterparts.

In this article we test the ability of the stochastic differential model proposed by Fatone et al. [Maximum likelihood estimation of the parameters of a system of stochastic differential equations that models the returns of the index of some classes of hedge funds, J. Inv. Ill-Posed Probl. 15 (2007), pp. 329–362] of forecasting the returns of a long-short equity hedge fund index and of a market index, that is, of the Hedge Fund Research performance Index (HFRI)-Equity index and of the S&P 500 (Standard & Poor 500 New York Stock Exchange) index, respectively. The model is based on the assumptions that the value of the variation of the log-return of the hedge fund index (HFRI-Equity) is proportional up to an additive stochastic error to the value of the variation of the log-return of a market index (S&P 500) and that the log-return of the market index can be satisfactorily modelled using the Heston stochastic volatility model. The model consists of a system of three stochastic differential equations, two of them are the Heston stochastic volatility model and the third one is the equation that models the behaviour of the hedge fund index and its relation with the market index. The model is calibrated on observed data using a method based on filtering and maximum likelihood proposed by Mariani et al. [Maximum likelihood estimation of the Heston stochastic volatility model using asset and option prices: An application of nonlinear filtering theory, Opt. Lett., 2 (2008), pp. 177–222] and further developed in Fatone et al. [Maximum likelihood estimation of the parameters of a system of stochastic differential equations that models the returns of the index of some classes of hedge funds, J. Inv. Ill-Posed Probl. 15 (2007), pp. 329–362; The calibration of the Heston stochastic volatility model using filtering and maximum likelihood methodsin Proceedings of Dynamic Systems and Applications, Vol. 5, G.S. Ladde, N.G. Medhin, C. Peng, and M. Sambandham, eds., Dynamic Publishers, Atlanta, USA, 2008, pp. 170–181]. That is, an inverse problem for the stochastic dynamical system representing the model is solved using the calibration procedure. The data analysed is from January 1990 to June 2007, and are monthly data. For each observation time, they consist of the value at the observation time of the log-returns of the HFRI-Equity and of the S&P 500 indices. The calibration procedure uses appropriate subsets of data, that is the data observed in a 6 months time period. The 6 months data time period used in the calibration is rolled through the time series generating a sequence of calibration problems. The values of the HFRI-Equity and S&P 500 indices forecasted using the calibrated models are compared to the values of the indices observed. The result of the comparison is very satisfactory. The website http://www.econ.univpm.it/recchioni/finance/w8 contains some auxiliary material including some animations that helps the understanding of this article. A more general reference to the work of some of the authors and of their coauthors in mathematical finance is the website: http://www.econ.univpm.it/recchioni/finance.

We considered the issues associated with methods for solving the problems of representation of continuous Markov processes in discrete-time by the example of radar signals, which are subject to the simultaneous effect of additive and multiplicative noise described by stochastic differential equations. We considered the algorithm, which allows us to carry out computer simulation of fading amplitude, phase and information parameter of the signal in discrete observation time using the example of a continuous frequency-modulated signal, the amplitude of which adheres to the Nakagami distribution that is typical for radar signals under simultaneous effect of additive and multiplicative noise. It is shown that with the help of computer simulation it is possible to create not only optimal algorithms, but also their corresponding optimal structural diagrams, which makes it possible to elevate the processing of radar signals to higher standards.

Mixing soil samples across experimental units ignores uncertainty and generates incorrect estimates of soil biota effects on plants: Response to Cahill etal. (2017) 'No silver bullet: different soil handling techniques are useful for different research questions, exhibit differential type I and II error rates, and are sensitive to sampling intensity'.

Many fluctuation phenomena, in physics and other fields, can be modeled byFokker–Planck or stochastic differential equations whose coefficients, associated with driftand diffusion components, may be estimated directly from the observed time series. Itscorrect characterization is crucial to determine the system quantifiers. However, due to thefinite sampling rates of real data, the empirical estimates may significantly differ from theirtrue functional forms. In the literature, low-order corrections, or even no corrections, havebeen applied to the finite-time estimates. A frequent outcome consists of lineardrift and quadratic diffusion coefficients. For this case, exact corrections havebeen recently found, from Itô–Taylor expansions. Nevertheless, model validationconstitutes a necessary step before determining and applying the appropriatecorrections. Here, we exploit the consequences of the exact theoretical results obtainedfor the linear–quadratic model. In particular, we discuss whether the observedfinite-time estimates are actually a manifestation of that model. The relevanceof this analysis is put into evidence by its application to two contrasting realdata examples in which finite-time linear drift and quadratic diffusion coefficientsare observed. In one case the linear–quadratic model is readily rejected while inthe other, although the model constitutes a very good approximation, low-ordercorrections are inappropriate. These examples give warning signs about the properinterpretation of finite-time analysis even in more general diffusion processes.

A minimum Hellinger distance estimator for stochastic differential equations: An application to statistical inference for continuous time interest rate models