Original link:tecdat.cn/?p=11974
Original source:Tuo End number according to the tribe public number
The R2WinBUGS package provides a convenient way to call WinBUGS from R. It automatically writes data and scripts in a Winbug-readable format for batch processing (since version 1.4). Once the WinBUGS process is complete, the resulting data can be read into R through the package itself (which provides a compact graphical summary of inference and convergence diagnostics), or the output can be further analyzed using the capabilities of the Coda package.
The WinBUGS software is available free at http://www.mrc-bsu.cam.ac.uk/bugs/.
R, a “language for Data Analysis and graphics processing,” is an open source and freely available statistical package that implements the language. See http://www.R-project.org/. R and R2WinBUGS are available from CRAN, http://CRAN.R-Project.org, or one of its mirrors. If an Internet connection is available, you can install R2WinBUGS by typing install.Packages (” R2WinBUGS “) at the R command prompt. Don’t forget to use library(R2WinBUGS)
example
School data
The Scholastic Aptitude Test (SAT) measures high school students’ abilities to help colleges make admissions decisions. Our data came from an experiment conducted in the late 1970s with eight high school SAT V (Scholastic Aptitude Test Language) tests. Sat-v is a standard multiple choice test administered by the Educational Testing Service. The service is interested in the effectiveness of the coaching program at each school in the selected school.
implementation
The implementation of the R2WinBUGS package is very simple. The main “function bugs () is called by the user. In principle, it is a package for several other functions that are progressively called init, as follows: 1.bugs.data.inits () writes data files’ data.txt ‘and’ inits1.txt ‘, ‘inits2.txt’… Go to the working directory. 2. bugs.script () writes to the file “script.txt” that WinBUGS use for batch processing. 3. bugs.run () updates the WinBUGS registry, calls WinBUGS, and runs it in batch mode using ‘script.txt’. 4. Bugs.sims () only called if the codaPkg argument is set to false (the default). Otherwise, bugs () returns the file name in which the data was stored. For example, these can be imported via packaged Coda, which provides convergence diagnostics, monte Carlo estimation calculations, trace plots, and more. The bug.sims () function reads simulations from WinBUGS into R, formats them, monitors convergence, performs convergence checks, and calculates medians and quantiles. It also prepares output for bugs () itself. These functions are not directly invoked by the user. Arguments are passed from bugs () to other functions.
example
We applied the functionality provided by R2WinBUGS to the sample data and analyzed the output.
School data
Sample data:
> schools
Copy the code
To model these data, we use the hierarchical model proposed by Gelman et al. We assume that the observation estimates for each school have a normal distribution with the mean theta and variance tau.y, and the inverse variance is 1 =σ.y2, whose prior distribution is uniform on (0,1000). For mean theta, we use another normal distribution with mean mu. Theta and inverse variance tau. Theta. For its prior distribution, see the following WinBUGS code:
model { for (j in 1:J) { y[j] ~ dnorm (theta[j], tau.y[j]) theta[j] ~ dnorm (mu.theta, Tau. Theta) tau. Y [j] < -pow (sigma. Y [j], -2)} mu. -2) sigma.theta ~ dunif (0, 1000) }Copy the code
This model must be stored in a separate file, such as ‘schools.bug’2, in the appropriate directory, such as C: / schools /. In R, the user must prepare the data input required by the Bugs () function. This can be a list containing the names of each data vector, for example
> J <- nrow(schools)
Copy the code
Using these data and model files, we can run the MCMC simulation to get estimates of theta, mu. Theta, and Sigma. Theta. Before running, the user must determine how many chains to run (N.chain = 3) and the number of iterations (N.iter = 1000). In addition, the user must specify the initial value of the chain, for example by writing a function:
> inits <- function(){ + list(theta = rnorm(J, 0, 100), mu.theta = rnorm(1, 0, 100), + sigma.theta = runif(1, 0, 100)) +}Copy the code
To begin MCMC emulation, the argument bug.directory in R must point to the installation directory of WinBUGS. You can easily print the results in the school.sim object by printing (schools.sim). For this example, you get similar results
Inference for Bugs model at "c:/schools/schools.bug" 3 chains, Each with 1000 iterations (first 500 discarded) n.ims = 1500 iterations saved mean SD 2.5% 25% 50% 75% 97.5% Rhat n.fest Theta [1] 11.1 9.1-3.0 5.0 10.0 16.0 31.8 1.1 39 Theta [2] 7.6 6.6-4.7 3.3 7.8 11.6 21.1 1.1 42 Theta [3] 5.7 8.4-12.5 0.6 6.1 10.8 21.8 1.0 150 Theta [4] 7.1 7.0-6.6 2.7 7.2 11.5 21.0 1.1 42 Theta [5] 5.1 6.8-9.5 0.7 5.2 9.7 18.1 1.0 83 Theta [6] 5.7 7.3-9.7 1.0 6.2 10.2 20.0 1.0 56 Theta [7] 10.4 7.3-2.1 5.3 9.8 15.3 25.5 1.1 27 Theta [8] 8.3 8.4-6.6 2.8 8.1 12.7 26.2 1.0 64 Mu. Theta 7.6 5.9 -3.0 3.7 8.0 11.0 19.5 1.1 35 sigma Deviance 60.8 57.0 59.1 60.2 62.1 66.6 1.0 170 pD = 3 and DIC = 63.8 (using the rule, pD = var(deviance)/2) For each parameter, n.eff is a crude measure of effective sample size, and Rhat is the potential scale reduction factor (at convergence, Rhat=1). DIC is an estimate of expected predictive error (lower deviance is better).Copy the code
In addition, users can generate a resulting graph by typing plot (schools.sim). The resulting graph is shown in the figure below. In the figure, the left column shows a quick summary of the following:
Inference and convergence (Rb of all parameters is close to 1.0, indicating a good mixture of the three chains and therefore approximate convergence); The right column shows the inferences for each set of parameters. As you can see from the right column, R2WinBUGS uses the parameter names in WinBUGS to construct the output as scalar, vector, and parameter arrays.
Most welcome insight
1. Matlab uses Bayesian optimization for deep learning
2. Matlab Bayesian hidden Markov HMM model implementation
3. Simple Bayesian linear regression simulation of Gibbs sampling in R language
4. Block Gibbs Sampling Bayesian multiple linear regression in R language
5. Bayesian model of MCMC sampling for Stan probabilistic programming in R language
6.Python implements Bayesian linear regression model with PyMC3
7.R language uses Bayesian hierarchical model for spatial data analysis
8.R language random search variable selection SSVS estimation Bayesian vector autoregression (BVAR) model
9. Matlab Bayesian hidden Markov HMM model implementation