Tuesday, June 4, 2013

New R Package for BEST (Bayesian ESTimation supersedes the t test)

A completely re-packaged version of the BEST software (from the article, "Bayesian estimation supersedes the t test") has been prepared by Michael E. Meredith. Mike is a key member of the Wildlife Conservation Society in Malaysia. For his new R package, Mike included additional MCMC diagnostic information, combined the two-group and one-group cases into a single function, made additional plot options and numerical-summary output, made the whole thing match conventional R style for plot commands, wrote new documentation, etc. Mike did all this completely independently of any prompt or help from me. I am very grateful for his efforts. The package looks great, and here's how to get it:

Open R. At the R command prompt, simply type:
install.packages("BEST")
Or, in RStudio, go through menu Tools,  Install Packages. You only need to install the BEST package once.

[Updated 05 June: When I initially posted this announcement yesterday, the binaries were not yet available on CRAN, and so I provided a set of instructions for how to install the package from its source code. Since then, the binaries have been posted, and installation takes only the single step above. Binaries for MacOS might be delayed another day or so.]

At the command line in R, type
library(BEST)
The functions for BEST are now in R's working memory. You need to do this command every time you invoke R.

Type
?BEST
to display a help file with a complete example of using the functions. Further documentation is available in this PDF that Mike prepared. Please note that the syntax of Mike's R package is a bit different than the original BEST programs. Perhaps the best way to learn the new syntax is to run the examples provided in the help file and PDF.

If you like the package, please let Mike know; his e-mail address is shown in the help files.

17 comments:

  1. Re installing the package: what's wrong with simply typing

    install.packages("BEST")

    ? As long as you have a working internet connection, this should work on Windows, Linux, Apple OS etc.

    ReplyDelete
  2. Dear Anonymous: Thanks for the comment. When I originally created the post, the binaries were not yet on CRAN. But now they are, and I have updated the blog post.

    ReplyDelete
  3. I've bought your book and have begun the process of converting myself over to Bayesian. I was wondering if you had any opinions on BLOG (https://sites.google.com/site/bloginference/).

    Or should I just learn JAGS instead?

    ReplyDelete
  4. Dear Anonymous June 10:

    BLOG looks interesting, but it is still in version < 1.0. On the other hand, JAGS is stable and I recommend it. (The likely successor to JAGS might be STAN, but I'd say stick to JAGS for learning and for moderately sized data applications. STAN may be the way to go if you have large data.)

    And thanks for getting the book. I hope it serves you well.

    ReplyDelete
  5. I tried to download and install BEST (for MacOS 10.6.8), but it failed:

    > install.packages("BEST")
    trying URL 'http://cran.uk.r-project.org/bin/macosx/contrib/3.0/BEST_0.1.0.tgz'
    Content type 'application/x-gzip' length 587233 bytes (573 Kb)
    opened URL
    ==================================================
    downloaded 573 Kb


    The downloaded binary packages are in
    /var/folders/r6/r6IOB+RkFwqO0yZdeBE1c++++TI/-Tmp-//RtmpLoCxqL/downloaded_packages
    > library(BEST)
    Loading required package: rjags
    Error : .onLoad failed in loadNamespace() for 'rjags', details:
    call: dyn.load(file, DLLpath = DLLpath, ...)
    error: unable to load shared object '/Library/Frameworks/R.framework/Versions/3.0/Resources/library/rjags/libs/rjags.so':
    dlopen(/Library/Frameworks/R.framework/Versions/3.0/Resources/library/rjags/libs/rjags.so, 10): Library not loaded: /usr/local/lib/libjags.3.dylib
    Referenced from: /Library/Frameworks/R.framework/Versions/3.0/Resources/library/rjags/libs/rjags.so
    Reason: image not found
    Error: package ‘rjags’ could not be loaded
    trying URL 'http://cran.uk.r-project.org/bin/macosx/contrib/3.0/foreign_0.8-54.tgz'
    Content type 'application/x-gzip' length 254326 bytes (248 Kb)
    opened URL
    ==================================================
    downloaded 248 Kb

    trying URL 'http://cran.uk.r-project.org/bin/macosx/contrib/3.0/mgcv_1.7-24.tgz'
    Content type 'application/x-gzip' length 1448309 bytes (1.4 Mb)
    opened URL
    ==================================================
    downloaded 1.4 Mb


    The downloaded binary packages are in
    /var/folders/r6/r6IOB+RkFwqO0yZdeBE1c++++TI/-Tmp-//RtmpLoCxqL/downloaded_packages

    any suggestions for how I can get it to run would be gratefully received.
    Thanks

    ReplyDelete
  6. Dear Anonymous June 17:

    From the error message you included, it looks like what is actually failing is rjags, not BEST. Have you installed rjags? Please follow all the steps at http://doingbayesiandataanalysis.blogspot.com/search/label/installation
    Hope that solves it for you.

    ReplyDelete
  7. Has the BEST package been discontinued or something? I tried to download it today and could not locate it:

    > install.packages("BEST")
    Installing package(s) into ‘C:\Users\aletheist\Documents/R/win-library/2.12’
    (as ‘lib’ is unspecified)
    Warning message:
    In getDependencies(pkgs, dependencies, available, lib) :
    package ‘BEST’ is not available

    ReplyDelete
  8. It works for me... Perhaps try a different CRAN mirror site. Thanks.

    ReplyDelete
  9. For what it's worth, it seems the problem was that I hadn't updated R on my machine in way too long. After installing the current version of R, I had no problems getting the package.

    ReplyDelete
  10. Hi,

    First thanks for making all this available as packages and the theory. It is much a appreciated!

    My question is if the BEST can be used where one would used a paired t-test? e.g. repeated measures.

    Sorry if this is in the book, I have just started reading it.

    ReplyDelete
  11. Yes, BEST can be used for paired data. The original programs and the re-packaged programs by Mike Meredith do it. See http://doingbayesiandataanalysis.blogspot.com/2012/09/one-group-version-of-best-bayesian.html
    Thanks for your interest.

    ReplyDelete
  12. Sir, I can make the plots, sorry but how to analyse the plots and write a paper. Will you please put up a model paper with such analyses. In India we have been brought up on a strict diet of ANOVA-RBD/CRD, and t tests with CV values (Analysis of Variance - Randomized Block Design / Completely Randomized design). How can I do these analyses / similar in Bayesian methods, with regards, Dr.D.K.Samuel, Ind Inst of Hort Research, Bangalore - 89

    ReplyDelete
  13. What follows is a reply to Deleep from Mike Meredith:

    Dear Duleep,

    I can't point you to papers describing a Bayesian analysis of ANOVA-type models, but I will indicate sources with advice and pass this on to John Kruscke who may have more references to hand.

    John has a blog post on the topic at http://doingbayesiandataanalysis.blogspot.com/2012/05/how-to-report-bayesian-analysis.html. This has a reproduction of the first part of section (23.1) on "Reporting a Bayesian analysis" from Kruschke, J.K. (2011) Doing Bayesian data analysis: a tutorial with R and BUGS Elsevier, Amsterdam etc.

    The main points are summarised in Kruschke, J.K., Aguinis, H., & Joo, H. (2012) The time has come: Bayesian methods for data analysis in the organizational sciences. Organizational Research Methods, 15, 722-752 in the section "Recommendations and Illustration of How to Report Bayesian Analysis Results", where they also give an example of writing up a regression analysis.

    Also Kruschke, J.K. (2013) Bayesian estimation supersedes the t test. *Journal of Experimental Psychology: General*, 142, 573-603 has a section on "Reporting the results of a Bayesian analysis". You would need to cite that paper as the description of the methods for robust Bayesian analysis.

    There are links to the two Kruschke papers at http://www.indiana.edu/~kruschke/publications.html

    In my own field (wildlife ecology) we can rarely do experiments, so hypothesis testing is scarcely relevant anyway. The focus is on estimation (of density, abundance, occupancy, survival, etc) or on modelling. In the latter context, information theoretic approaches (using Akaike's Information Criterion, AIC) have been the norm, and moving from there to a Bayesian approach is relatively painless. So I can't point you to papers reporting a Bayesian analysis of ANOVA-type data; Marc Kery's book, An Introduction to WinBUGS for Ecologists (Academic Press, 2010), shows how to do this, but not how to write it up.

    Regards, Mike.

    ReplyDelete
  14. Mike,

    I've recently started using BEST and have found it easy to set up and play with.

    However, I have a question about the 'parallel' option. I would like to run BEST with more than the 3 cores BESTmcmc appears to be defaulted too. Is there a way to run BEST on the nCores-1? Or to choose the number of cores that I can run BEST on?

    Many thanks,

    Rich

    ReplyDelete
  15. Please ask Mike Meredith directly about the innards of the BEST package. In the DBDA2E software (not the BEST package), the number of cores used is set when DBDA2E-utilities.R is sourced. The relevant code is as follows:

    library(parallel) # for detectCores().
    nCores = detectCores()
    if ( !is.finite(nCores) ) { nCores = 1 }
    if ( nCores > 4 ) {
    nChainsDefault = 4 # because JAGS has only 4 rng's.
    runjagsMethodDefault = "parallel"
    }
    if ( nCores == 4 ) {
    nChainsDefault = 3 # save 1 core for other processes.
    runjagsMethodDefault = "parallel"
    }
    if ( nCores < 4 ) {
    nChainsDefault = 3
    runjagsMethodDefault = "rjags" # NOT parallel
    }

    Notice that it maxes out at 4 because JAGS would recycle its PRNGs and produce identical chains. You can manually manipulate the initial values and PRNGs -- see the runjags manual at https://cran.r-project.org/web/packages/runjags/index.html

    ReplyDelete
  16. Ahhh! Thanks John - I didn't realise that would happen in jags, and it makes sense that the default for BEST would match that.

    I'll have a look at the runjags manual to see if there is a simple workaround.

    Thanks again!

    Rich

    ReplyDelete
  17. I'm just following up with a new question - i'm new to Bayesian Analysis, so apologies if i've missed something obvious, but is there a methodology for testing the difference of non-symmetric distributions?

    For example, how would I go about comparing the skewness of two distributions to see if they are different?

    Thanks,

    Rich

    ReplyDelete