We derive an exact and efficient Bayesian regression algorithm for piecewise constant functions of unknown segment number, boundary locations, and levels. The derivation works for any noise and segment level prior, e.g. Cauchy which can handle outliers. We derive simple but good estimates for the in-segment variance. We also propose a Bayesian regression curve as a better way of smoothing data without blurring boundaries. The Bayesian approach also allows straightforward determination of the evidence, break probabilities and error estimates, useful for model selection and significance and robustness studies. We discuss the performance on synthetic and real-world examples. Many possible extensions are discussed.

Our primary goal is to obtain a smoothed summary estimate of the magnetic field generated in and near to the Milky Way by using Faraday rotation measures (RM's). Each RM in our data set provides an integrated measure of the effect of the magnetic field along the entire line of sight to an extragalactic radio source. The ability to estimate the magnetic field generated locally by our galaxy and its environs will help astronomers distinguish local versus distant properties of the universe. RM's can be considered analogous to %prototypical of geostatistical data on a sphere. In order to model such data, we employ a Bayesian process convolution approach which uses Markov chain Monte Carlo (MCMC) for estimation and prediction. Complications arise due to contamination in the RM measurements, and we resolve these by means of a mixture prior on the errors.

MAP estimators and HPD credible sets are often criticized in the literature because of paradoxical behaviour due to a lack of invariance under reparametrization. In this paper, we propose a new version of MAP estimators and HPD credible sets that avoid this undesirable feature. Moreover, in the special case of non-informative prior, the new MAP estimators coincide with the invariant frequentist ML estimators. We also propose several adaptations in the case of nuisance parameters.

In the context of statistical analysis, elicitation is the process of translating someone's beliefs about some uncertain quantities into a probability distribution. The person's judgements about the quantities are usually fitted to some member of a convenient parametric family. This approach does not allow for the possibility that any number of distributions could fit the same judgements.

In this paper, elicitation of an expert's beliefs is treated as any other inference problem: the facilitator of the elicitation exercise has prior beliefs about the form of the expert's density function, the facilitator elicits judgements about the density function, and the facilitator's beliefs about the expert's density function are updated in the light of these judgements. This paper investigates prior beliefs about an expert's density function and shows how many different types of judgement can be handled by this method.

This elicitation method begins with the belief that the expert's density will roughly have the shape of a $t$~density. This belief is then updated through a Gaussian process model using judgements from the expert. The method gives a framework for quantifying the facilitator's uncertainty about a density given judgements about the mean and percentiles of the expert's distribution. A property of Gaussian processes can be manipulated to include judgements about the derivatives of the density, which allows the facilitator to incorporate mode judgements and judgements on the sign of the density at any given point. The benefit of including the second type of judgement is that substantial computational time can be saved.

Suppose that $S({\bf Y},\Theta)$ is a function of data ${\bf Y}$ and a model parameter $\Theta$, and suppose that the sampling distribution of $S({\bf Y},\Theta)$ is invariant when evaluated at $\Theta_0$, the "true" (i.e., data-generating) value of $\Theta$. Then $S({\bf Y},\Theta)$ is a pivotal quantity, and it follows from simple probability calculus that the distribution of $S({\bf Y},\Theta_0)$ is identical to the distribution of $S({\bf Y},\Theta_{\bf Y})$, where $\Theta_{\bf Y}$ is a value of $\Theta$ drawn from the posterior distribution given ${\bf Y}$. This fact makes it possible to define a large number of Bayesian model diagnostics having a known sampling distribution. It also facilitates the calibration of the joint sampling of model diagnostics based on pivotal quantities.

Motivated by an increasing number of Bayesian hierarchical model applications, the objective of this paper is to evaluate properties of several diagnostic techniques when the fitted model includes some hierarchical structure, but the data are from a model with additional, unknown hierarchical structure. Because there has been no apparent evaluation of Bayesian diagnostics used for this purpose, we start by studying the simple situation where the data come from a normal model with two-stage hierarchical structure while the fitted model does not have any hierarchical structure, and then extend this to the case where the fitted model has two-stage normal hierarchical structure while the data come from a model with three-stage normal structure. We use exact derivations, large sample approximations and numerical examples to evaluate the quality of the diagnostic techniques. Our investigation suggests two promising techniques: distribution of individual posterior predictive $p$ values and the conventional posterior predictive $p$ value with the $F$ statistic as a checking function. We show that (at least) for large sample sizes these $p$ values are uniformly distributed under the null model and are effective in detecting hierarchical structure not included in the null model. Finally, we apply these two techniques to examine the fit of a model for data from the Patterns of Care Study, a two-stage cluster sample of cancer patients undergoing radiation therapy.

We show how the mean of a monotone function (defined on a state space equipped with a partial ordering) can be estimated, using ergodic averages calculated from upper and lower dominating processes of a stationary irreducible Markov chain. In particular, we do not need to simulate the stationary Markov chain and we eliminate the problem of whether an appropriate burn-in is determined or not. Moreover, when a central limit theorem applies, we show how confidence intervals for the mean can be estimated by bounding the asymptotic variance of the ergodic average based on the equilibrium chain. Our methods are studied in detail for three models using Markov chain Monte Carlo methods and we also discuss various types of other models for which our methods apply.

In recent times, complex computer models have received wide attention in scientific research. However, in order to make conventional statistical statements regarding the scientific research, many expensive runs of the computer model are usually needed. New statistical theories, making their appearances, hold promise to alleviate the technical challenges. However, in cases where the underlying complex system is evolving with time, an effective theory for statistical analyses is lacking. In this paper, we propose a novel Bayesian methodology that extends the existing methodologies to the case of complex dynamic systems. The approach described in the paper exploits the recursive nature of dynamic simulation models to give a more efficient and accurate emulator.

The motivating example, although not a real model for any physical process, may be thought of as a proxy for a model representing climate change, where it is of interest to predict, over time $t$, the four-dimensional proxy time series $\mathbb{y}_t=\mbox{(temperature},$ $\mbox{ice melting rate}$, $\mbox{barren land, $\mbox{CO}_2$ emission)}$. Also available are proxy observations on deforestation, recorded over time; hence treated as known. The latter is known as forcing input, denoted by $z_t$. The computer model is treated as a black box.

Typically, Gaussian processes are used to model unknown computer models, which we adopt in our article. In order to exploit the recursive nature of dynamic computer models, we introduce a grid within the range of the unknown function where the entire dynamic sequence is expected to lie. This grid essentially defines a look-up table. %Instead of the full high-dimensional posterior %predictive distribution for the future time series, Our proposed method then assumes that conditional on the response surface on the grid, and the available training data, the future responses are approximately independent. Exploiting the properties of Gaussian process, we justify our proposal theoretically and with ample simulation studies. We also apply our proposed methodology to the motivating example.

Hierarchical Bayes models provide a natural way of incorporating covariate information into the inferential process through the elaboration of regression equations for one or more of the model parameters, with errors that are often assumed to be i.i.d. Gaussian. Unfortunately, building adequate regression models is a complicated art form that requires the practitioner to make numerous decisions along the way. Assessing the validity of the modeling decisions is often difficult.

In this article I develop a simple and effective device for ascertaining the quality of the modeling choices and detecting lack-of-fit. I specify an artificial autoregressive structure (AAR) in the probability model for the errors that incorporates the i.i.d. model as a special case. Lack-of-fit can be detected by examining the posterior distribution of AAR parameters. In general, posterior distributions that assign considerable mass to a region of the AAR parameter space away from zero provide evidence that apparent dependencies in the errors are compensating for misspecifications of some other aspects (typically conditional means) of the model. I illustrate the methodology through several examples including its application to the analysis of data on brain and body weights of mammalian species and response time data.