[R-sig-ME] Has anyone tried Whit Armstrong's CppMCMC library?

Whit Armstrong armstrong.whit at gmail.com
Wed Nov 3 19:55:03 CET 2010


Below is a minimal example of a simple linear model using CppBugs,
Rcpp, and inline.

The cppbugs dir needs to be on the include path of the compiler for
this example to work.  Otherwise, you will need to change the #include
statements to use the local path of your cppbugs install.

Other dependencies you need are Armadillo and Boost.

feedback welcome.

-Whit

require(inline)
require(Rcpp)

cppbugs.model <- '
#include <armadillo>
#include <cppbugs/cppbugs.hpp>

using namespace arma;
using namespace cppbugs;

class TestModel: public MCModel {
public:
  const mat& y; // given
  const mat& X; // given

  Stochastic<vec> b;
  Stochastic<double> tau_y;
  Deterministic<mat> y_hat;
  Stochastic<mat> likelihood;
  Deterministic<double> rsq;

  TestModel(const mat& y_,const mat& X_): y(y_),
X(X_),b(randn<vec>(X_.n_cols)),tau_y(1),y_hat(X*b.value),likelihood(y_,true),rsq(0)
  {
    add(b);
    add(tau_y);
    add(y_hat);
    add(likelihood);
    add(rsq);
  }

  void update() {
    y_hat.value = X*b.value;
    rsq.value = as_scalar(1 - var(y - y_hat.value) / var(y));
    b.dnorm(0.0, 0.0001);
    tau_y.dunif(0,100);
    likelihood.dnorm(y_hat.value,tau_y.value);
  }
};
'

src <- '
mat X(REAL(XR),::Rf_nrows(XR),::Rf_ncols(XR));
mat y(REAL(yr),::Rf_nrows(yr),::Rf_ncols(yr));
int iterations_ = as<int>(iterations);
int burn_ = as<int>(burn);
int thin_ = as<int>(thin);

TestModel m(y,X);
m.sample(iterations_, burn_, thin_);
return Rcpp::List::create(Rcpp::Named("b", m.b.mean()),
Rcpp::Named("tau.y", m.tau_y.mean()), Rcpp::Named("ar",
m.acceptance_ratio()), Rcpp::Named("rsq", m.rsq.mean()));
'

NR <- 100
NC <- 2
icept <- 10
y = rnorm(NR,icept)
X = matrix(rnorm(NR*NC),ncol=NC)
X[,1] <- 1.0
X[,2] <- y + - icept + rnorm(NR);

bayes.reg <- cxxfunction(signature(XR="numeric",
yr="numeric",iterations="integer",burn="integer",thin="integer"),
body=src, include=cppbugs.model, plugin="Rcpp")
ans <- bayes.reg(X,y,1e5L,1e4L,10L)
print(ans)


On Tue, Nov 2, 2010 at 1:02 PM, Whit Armstrong <armstrong.whit at gmail.com> wrote:
> Thanks for the the post, Professor Bates.
>
> CppBugs can be found here: http://github.com/armstrtw/CppBugs
>
> It's a combination of WinBUGS and PyMC.
>
> The api is still very unstable as I'm attempting to navigate between
> good design and maintaining familiar WinBUGS conventions.
>
> Rcpp is the best way to use it from R.  I'll add an example of calling
> it this way on the front page so that everyone can try it.
>
> My thinking is that since it's already very inconvenient to run a BUGS
> model from R (writing the BUGS script, providing inits, passing data
> back and forth), it's not really asking too much for users to write a
> small c++ class inline in an R script that is the CppBugs equivalent
> of a BUGS model.
>
> I'm very keen to have help on this project, so if anyone is interested
> please contact me.
>
> I'll send a follow up post highlighting a few of the features of
> CppBugs and examples of WinBUGS vs PyMC conventions that I've followed
> in the design.
>
> -Whit
>
>
>
>
> On Tue, Nov 2, 2010 at 11:08 AM, Douglas Bates <bates at stat.wisc.edu> wrote:
>> Whit Armstrong has written a C++ library called CppMCMC that provides
>> BUGS-like functionality in compiled code.  I see that one of his
>> examples uses the contagious bovine pleuropneumonia data from the cbpp
>> data set in the lme4 package, fitting an overdispersed binomial-like
>> model to it.
>>
>> Does anyone have experience with these classes?  This is not directly
>> an R question as, at present, I don't know of a bridge between CppMCMC
>> and R (although it would be an interesting project to use Rcpp for
>> this).
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>




More information about the R-sig-mixed-models mailing list