Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

A while back I was running a bunch of JAGS models through R, using the rjags (written by Martyn Plummer) and R2jags (by Yu-Sung Su) packages. These packages provide a great interface to the JAGS software, which allows analysis of Bayesian models (written in the BUGS language) through Markov chain Monte Carlo simulation.

Running a JAGS model using these tools returns an rjags object, which when printed to the screen, summarises the posterior distribution of each monitored node, giving its mean and standard deviation, a range of quantiles, and its Gelman-Rubin convergence diagnostic statistic (Rhat), which indicates the ratio of variance within chains to that among chains. The summary is great, but when monitoring a large number of nodes, printing these to the screen can cause R to hang, and can exceed the screen buffer (not to mention making it painful to find the nodes you’re immediately interested in).

To help deal with this I wrote a couple of simple R functions:

• `jagsresults:` return a matrix containing summary results for just the nodes you are interested in (using regular expression pattern-matching, if desired).
• `rhats`: sort the output by the nodes’ Rhat values, making it easy to show the n least converged nodes.

## Examples

Consider an rjags object that has the following summary (this is the result of the example model given in ?jagsresults). Note that in the examples below I’ve rounded the output to 3 decimal places to prevent line-wrapping/scrollbars on my website (although I’ve excluded the call to `round()`).

```> jagsfit
Inference for Bugs model at "C:/Users/John/AppData/Local/Temp/RtmpIDmQOb/model112833e66869.txt", fit using jags,
3 chains, each with 10000 iterations (first 5000 discarded), n.thin = 5
n.sims = 3000 iterations saved
mu.vect sd.vect    2.5%     25%     50%     75%   97.5%  Rhat n.eff
a            0.041   0.052  -0.057   0.005   0.041   0.076   0.146 1.001  3000
beta.rain    0.997   0.055   0.888   0.960   0.999   1.034   1.105 1.001  3000
beta.temp    1.414   0.054   1.310   1.377   1.415   1.451   1.519 1.001  2400
beta.wind   -0.458   0.057  -0.570  -0.494  -0.458  -0.420  -0.342 1.001  3000
resid[1]    -0.090   0.028  -0.144  -0.109  -0.090  -0.071  -0.035 1.001  3000
resid[2]    -0.048   0.024  -0.094  -0.064  -0.048  -0.032   0.001 1.001  3000
resid[3]     0.076   0.031   0.015   0.055   0.075   0.097   0.138 1.001  3000
resid[4]     0.170   0.050   0.069   0.137   0.171   0.204   0.267 1.001  3000
resid[5]     0.110   0.028   0.054   0.091   0.109   0.129   0.165 1.001  3000
resid[6]     0.086   0.045   0.001   0.055   0.086   0.115   0.175 1.001  3000
resid[7]    -0.016   0.025  -0.066  -0.033  -0.017   0.001   0.033 1.001  3000
resid[8]    -0.226   0.027  -0.280  -0.243  -0.226  -0.207  -0.174 1.001  3000
resid[9]     0.250   0.023   0.204   0.235   0.250   0.266   0.294 1.001  3000
resid[10]   -0.017   0.022  -0.060  -0.032  -0.017  -0.002   0.027 1.001  3000
resid[11]   -0.069   0.025  -0.117  -0.085  -0.069  -0.052  -0.020 1.001  3000
resid[12]    0.145   0.022   0.102   0.131   0.145   0.159   0.187 1.001  3000
resid[13]   -0.184   0.030  -0.246  -0.204  -0.184  -0.163  -0.128 1.001  3000
resid[14]   -0.027   0.029  -0.085  -0.046  -0.027  -0.007   0.030 1.001  3000
resid[15]   -0.015   0.027  -0.067  -0.032  -0.015   0.003   0.038 1.001  3000
resid[16]   -0.176   0.028  -0.231  -0.195  -0.176  -0.158  -0.122 1.001  3000
resid[17]   -0.083   0.032  -0.146  -0.105  -0.083  -0.062  -0.021 1.001  3000
resid[18]   -0.010   0.030  -0.069  -0.030  -0.010   0.010   0.046 1.001  3000
resid[19]   -0.122   0.030  -0.180  -0.142  -0.122  -0.101  -0.063 1.001  3000
resid[20]    0.111   0.024   0.064   0.096   0.111   0.127   0.158 1.001  3000
resid[21]   -0.009   0.039  -0.085  -0.036  -0.010   0.018   0.070 1.001  3000
resid[22]   -0.171   0.026  -0.223  -0.189  -0.171  -0.154  -0.118 1.001  3000
resid[23]   -0.063   0.043  -0.147  -0.092  -0.064  -0.035   0.024 1.001  2800
resid[24]    0.069   0.037  -0.002   0.044   0.070   0.094   0.140 1.001  3000
resid[25]    0.061   0.029   0.003   0.042   0.062   0.081   0.119 1.001  3000
resid[26]    0.046   0.031  -0.017   0.025   0.045   0.066   0.106 1.001  3000
resid[27]    0.138   0.031   0.074   0.117   0.138   0.160   0.197 1.001  3000
resid[28]    0.116   0.020   0.076   0.102   0.116   0.129   0.155 1.001  3000
resid[29]   -0.069   0.035  -0.141  -0.093  -0.069  -0.045  -0.002 1.001  3000
resid[30]   -0.238   0.029  -0.294  -0.257  -0.238  -0.218  -0.183 1.001  3000
resid[31]    0.131   0.036   0.063   0.106   0.131   0.156   0.203 1.001  3000
resid[32]    0.024   0.036  -0.046   0.000   0.024   0.049   0.094 1.001  3000
resid[33]    0.073   0.036   0.002   0.049   0.073   0.097   0.140 1.001  3000
resid[34]   -0.053   0.036  -0.123  -0.077  -0.053  -0.029   0.019 1.001  2700
resid[35]    0.106   0.030   0.047   0.086   0.106   0.126   0.164 1.001  3000
resid[36]    0.154   0.028   0.102   0.136   0.154   0.173   0.211 1.001  3000
resid[37]    0.166   0.029   0.111   0.147   0.166   0.186   0.223 1.001  2100
resid[38]   -0.045   0.043  -0.129  -0.074  -0.045  -0.016   0.040 1.002  2700
resid[39]   -0.076   0.036  -0.149  -0.101  -0.076  -0.052  -0.005 1.001  3000
resid[40]    0.181   0.032   0.117   0.159   0.181   0.203   0.243 1.001  3000
resid[41]    0.135   0.051   0.033   0.101   0.136   0.169   0.233 1.001  3000
resid[42]    0.118   0.026   0.068   0.101   0.119   0.136   0.170 1.001  3000
resid[43]    0.026   0.031  -0.034   0.005   0.025   0.046   0.085 1.001  3000
resid[44]   -0.059   0.032  -0.120  -0.082  -0.059  -0.038   0.005 1.001  3000
resid[45]   -0.035   0.041  -0.117  -0.062  -0.035  -0.007   0.043 1.001  3000
resid[46]    0.058   0.035  -0.013   0.035   0.058   0.082   0.126 1.001  3000
resid[47]    0.035   0.033  -0.030   0.013   0.036   0.057   0.100 1.001  3000
resid[48]    0.087   0.024   0.041   0.070   0.086   0.103   0.135 1.001  3000
resid[49]    0.070   0.026   0.017   0.053   0.070   0.087   0.119 1.001  3000
resid[50]    0.085   0.042   0.001   0.058   0.086   0.113   0.168 1.001  3000
resid[51]    0.018   0.034  -0.047  -0.005   0.018   0.040   0.085 1.002  2000
resid[52]   -0.125   0.044  -0.210  -0.154  -0.125  -0.097  -0.040 1.001  2900
resid[53]    0.443   0.030   0.384   0.423   0.443   0.463   0.502 1.001  3000
resid[54]   -0.351   0.036  -0.421  -0.375  -0.351  -0.326  -0.281 1.001  3000
resid[55]   -0.204   0.021  -0.244  -0.218  -0.204  -0.190  -0.162 1.001  3000
resid[56]   -0.079   0.027  -0.131  -0.098  -0.079  -0.062  -0.027 1.001  3000
resid[57]    0.025   0.029  -0.034   0.005   0.025   0.044   0.081 1.000  3000
resid[58]   -0.149   0.043  -0.233  -0.177  -0.148  -0.120  -0.064 1.001  3000
resid[59]    0.058   0.041  -0.020   0.031   0.057   0.085   0.136 1.001  2200
resid[60]    0.050   0.019   0.012   0.037   0.050   0.063   0.088 1.001  3000
resid[61]    0.237   0.019   0.198   0.223   0.237   0.249   0.275 1.001  3000
resid[62]   -0.062   0.031  -0.121  -0.083  -0.062  -0.041  -0.004 1.001  3000
resid[63]   -0.196   0.030  -0.257  -0.216  -0.197  -0.175  -0.139 1.001  3000
resid[64]    0.009   0.032  -0.055  -0.012   0.010   0.031   0.073 1.001  3000
resid[65]    0.138   0.041   0.057   0.109   0.137   0.165   0.220 1.001  3000
resid[66]    0.124   0.030   0.066   0.103   0.124   0.144   0.183 1.002  1800
resid[67]    0.217   0.022   0.175   0.202   0.217   0.232   0.262 1.001  3000
resid[68]    0.043   0.036  -0.028   0.018   0.043   0.067   0.114 1.001  3000
resid[69]   -0.224   0.034  -0.290  -0.246  -0.224  -0.200  -0.158 1.001  3000
resid[70]    0.072   0.028   0.019   0.053   0.072   0.090   0.128 1.001  3000
resid[71]   -0.184   0.028  -0.238  -0.203  -0.184  -0.165  -0.129 1.001  3000
resid[72]    0.183   0.033   0.117   0.161   0.183   0.205   0.250 1.003  1600
resid[73]    0.041   0.028  -0.015   0.023   0.041   0.061   0.094 1.001  3000
resid[74]    0.131   0.025   0.079   0.114   0.131   0.147   0.180 1.001  3000
resid[75]   -0.138   0.022  -0.182  -0.153  -0.138  -0.123  -0.095 1.001  3000
resid[76]    0.341   0.023   0.295   0.326   0.341   0.357   0.387 1.001  3000
resid[77]   -0.394   0.045  -0.481  -0.424  -0.394  -0.365  -0.306 1.001  2700
resid[78]   -0.306   0.027  -0.359  -0.324  -0.306  -0.287  -0.252 1.001  3000
resid[79]    0.041   0.035  -0.029   0.017   0.041   0.065   0.110 1.001  3000
resid[80]   -0.295   0.030  -0.356  -0.315  -0.294  -0.274  -0.239 1.001  3000
resid[81]    0.159   0.036   0.087   0.135   0.160   0.183   0.226 1.002  3000
resid[82]    0.111   0.033   0.044   0.089   0.110   0.133   0.179 1.002  2000
resid[83]   -0.075   0.025  -0.125  -0.092  -0.075  -0.058  -0.028 1.001  3000
resid[84]   -0.336   0.031  -0.397  -0.356  -0.335  -0.313  -0.275 1.001  3000
resid[85]    0.203   0.019   0.167   0.191   0.203   0.216   0.241 1.001  3000
resid[86]    0.171   0.020   0.132   0.157   0.171   0.185   0.211 1.001  3000
resid[87]   -0.021   0.032  -0.083  -0.042  -0.021  -0.001   0.041 1.001  3000
resid[88]   -0.112   0.031  -0.173  -0.133  -0.112  -0.091  -0.052 1.001  3000
resid[89]    0.012   0.019  -0.026  -0.001   0.012   0.025   0.048 1.001  3000
resid[90]   -0.025   0.031  -0.086  -0.046  -0.025  -0.004   0.034 1.001  3000
resid[91]    0.194   0.037   0.119   0.170   0.194   0.218   0.264 1.001  3000
resid[92]   -0.060   0.036  -0.133  -0.084  -0.060  -0.035   0.009 1.001  3000
resid[93]   -0.073   0.022  -0.117  -0.088  -0.072  -0.058  -0.028 1.001  3000
resid[94]    0.180   0.020   0.140   0.167   0.180   0.193   0.220 1.001  2900
resid[95]    0.050   0.031  -0.012   0.029   0.050   0.071   0.110 1.001  3000
resid[96]   -0.246   0.023  -0.292  -0.261  -0.246  -0.230  -0.202 1.001  3000
resid[97]    0.020   0.037  -0.053  -0.004   0.021   0.045   0.093 1.001  3000
resid[98]   -0.298   0.035  -0.368  -0.322  -0.298  -0.274  -0.229 1.001  3000
resid[99]   -0.038   0.025  -0.089  -0.055  -0.038  -0.021   0.011 1.001  3000
resid[100]  -0.209   0.033  -0.273  -0.231  -0.209  -0.187  -0.145 1.001  3000
sd           0.159   0.012   0.138   0.151   0.158   0.166   0.183 1.001  2900
deviance   -85.830   3.245 -90.105 -88.214 -86.519 -84.102 -78.114 1.001  2500

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 info (using the rule, pD = var(deviance)/2)
pD = 5.3 and DIC = -80.6
DIC is an estimate of expected predictive error (lower deviance is better).
```

The `jagsresults()` function can be used to return a matrix containing just `a` and `beta.wind`:

```> jagsresults(x=jagsfit, params=c('a', 'beta.wind'))
mean    sd   2.5%    25%    50%    75%  97.5%  Rhat n.eff
a          0.041 0.052 -0.057  0.005  0.041  0.076  0.146 1.001  3000
beta.wind -0.458 0.057 -0.570 -0.494 -0.458 -0.420 -0.342 1.001  3000
```

Or, by specifying `invert=TRUE`, we can return a matrix containing all parameters except for `resid`:

```> jagsresults(x=jagsfit, params='resid', invert=TRUE)
mean    sd    2.5%     25%     50%     75%   97.5%  Rhat n.eff
a           0.041 0.052  -0.057   0.005   0.041   0.076   0.146 1.001  3000
beta.rain   0.997 0.055   0.888   0.960   0.999   1.034   1.105 1.001  3000
beta.temp   1.414 0.054   1.310   1.377   1.415   1.451   1.519 1.001  2400
beta.wind  -0.458 0.057  -0.570  -0.494  -0.458  -0.420  -0.342 1.001  3000
deviance  -85.830 3.245 -90.105 -88.214 -86.519 -84.102 -78.114 1.001  2500
sd          0.159 0.012   0.138   0.151   0.158   0.166   0.183 1.001  2900
```

The argument `exact` can be set to `FALSE` to return a matrix containing results for all parameters whose names contain (rather than match exactly) the strings given in params:

```> jagsresults(x=jagsfit, params='beta', exact=FALSE)
mean    sd   2.5%    25%    50%    75%  97.5%  Rhat n.eff
beta.rain  0.997 0.055  0.888  0.960  0.999  1.034  1.105 1.001  3000
beta.temp  1.414 0.054  1.310  1.377  1.415  1.451  1.519 1.001  2400
beta.wind -0.458 0.057 -0.570 -0.494 -0.458 -0.420 -0.342 1.001  3000
```

Finally, `regex=TRUE` can be specified to use regular expression pattern-matching, where `param` is now the pattern to be matched. In this case, additional arguments accepted by the `grep` function can be supplied to `jagsresults`, such as `perl=TRUE`, which is required for the syntax of some patterns:

```> jagsresults(x=jagsfit, param='^(?!res)', regex=TRUE, perl=TRUE)
mean    sd    2.5%     25%     50%     75%   97.5%  Rhat n.eff
a           0.041 0.052  -0.057   0.005   0.041   0.076   0.146 1.001  3000
beta.rain   0.997 0.055   0.888   0.960   0.999   1.034   1.105 1.001  3000
beta.temp   1.414 0.054   1.310   1.377   1.415   1.451   1.519 1.001  2400
beta.wind  -0.458 0.057  -0.570  -0.494  -0.458  -0.420  -0.342 1.001  3000
deviance  -85.830 3.245 -90.105 -88.214 -86.519 -84.102 -78.114 1.001  2500
sd          0.159 0.012   0.138   0.151   0.158   0.166   0.183 1.001  2900
```

Note, though, that the above (which returns all parameters except for `resid`) can be more simply produced by the earlier example demonstrating the use of the `invert` argument.

Standard matrix subsetting is possible for the resulting object, e.g.:

```> jags.resids <- jagsresults(x=jagsfit, param='resid')
> jags.resids[c(1, 3, 10:15), c('mean', '2.5%', '97.5%')]
mean        2.5%       97.5%
resid[1]  -0.08970919 -0.14412395 -0.03469326
resid[3]   0.07606610  0.01462996  0.13840041
resid[10] -0.01695619 -0.06031680  0.02725152
resid[11] -0.06851323 -0.11697671 -0.02024998
resid[12]  0.14516502  0.10249481  0.18689894
resid[13] -0.18398696 -0.24553446 -0.12795644
resid[14] -0.02688139 -0.08512743  0.02994298
resid[15] -0.01459422 -0.06701636  0.03812518
```

## Installation

The simplest way to install the `jagstools` package is with `install_github` in the `devtools` package, which will grab the source directly from Github:

```install.packages('devtools')
library(devtools)
```