How to implement user-defined target and derivative function from c in stan

[Lu Zhang]

I am trying to use HMC to do MCMC sampling with a self-defined target function. I had my HMC code, but it doesn't work well. So I want to implement my code with "stan".

I have already coded the target function as well as the derivative function in C. It involves for-loop, and functions needed for Cholesky decomposition and other linear algebras from Lapack. So complicated that I don't want to recode it in "stan". Is there any way to use my functions directly? 

I browsed the user's guide but didn't see any code for defining a user-defined derivative. Can users define their own derivative?

I would also appreciate it if anyone can tell me where is the code for HMC updating in the source code. Maybe I can combine my code with C++ code. 

Thank you!

 

[Andrew Gelman]

Lu:

What happens if you just write the log-density directly in Stan and let Stan do the gradients for you?

A

--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Bob Carpenter]

Here's a pointer to the root of the MCMC code on GitHub:

https://github.com/stan-dev/stan/tree/develop/src/stan/mcmc

Stan can handle for loops and Cholesky decompositions---the
gradients are computed automatically by Stan.

And you don't want to use static Euclidean HMC---much better
off using our adaptive NUTS algorithm (see the Hoffman and Gelman
JMLR paper). The reason HMC wasn't very popular before was
that (a) derivatives are a pain to code and error prone for
complex models, and (b) HMC itself is very sensitive to its
tuning parameters and even when well tuned, usually worse
than NUTS.

Given that you seem capable of some pretty serious coding,
how about contributing to Stan? We can always use developers
with good C, MCMC, matrix, and calc skills.

- Bob

[Lu Zhang]

Thanks. I just finished the code in Stan. Now I am waiting for the result. The language is much easier than I thought. And Stan user group is very helpful! 

[Bob Carpenter]

I don't know how complicated your model is, but Stan isn't
scalable to really complicated densities in high dimensions, which
it sounds like you might be dealing with.

I'd suggest setting the number of iterations low to start (say
100 total), see how things shape up, then extend to more iterations;
you usually need at least 200 to get adaptation working.

Anyway, let us know if you want us to take a look at your Stan
code and provide performance tips.

Also, rather than doing Cholesky decompositions in the code,
we provide Cholesky-factor parameter types for both correlation
and covariance matrices (and don't recommend Wishart priors;
we have an LKJ prior on correlation matrices that can be scaled
with a prior on scales---the LKJ prior has a parameter for degree
of concentration around a unit correlation matrix).

- Bob

[Lu Zhang]

    I am still testing my code. Maybe there is a bug in my code.  I am working on a spatial model which has five parameters. I tried 2500 iterations, but 3 of the 95% CIs failed to cover the true value. Thank you so much for your warm help. I will let you know if I make any progress or if I get stuck.

Best Wishes,

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+unsubscribe@googlegroups.com.

> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--

You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/ZDiVgFMZ7KU/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[Lu Zhang]

Hi Prof. Carpenter and Prof. Gelman,

    There was a mistake in my code and I fixed it. Stan did an impressive job on my model. It turns out that the resulting samples of two parameters are highly correlated (|rho| > 0.85). I think that's the reason why my vanilla HMC fails to have good performance. 

    Actually, I spent a lot of time on coding in this project. In the beginning, R was the only programming language I was familiar with. But its' low speed and poor memory management annoyed me a lot. So I self-studied C, command lines and makefile, and how to use Lapack and Blas in c. My code works pretty well on a simpler model, but things changed a lot after including one more parameter. So I didn't expect that I can finish Stan code within one day! And it gave me a much better performance than my code in C!  

    I know the changing of step size and step length does not break the detail balance, but I never came up with any good idea of using that fact until I saw the "adaptive NUTS" paper this Monday. That is delicate and beautiful! I am also trying to understand RMHMC. Now I can't see any reason for not using Stan.

   It would be wonderful if you would provide performance tips for my code. But I want to apply Stan to a higher dimension model, and I also need to polish my code a little bit. So maybe I will contact you later. ^_^

Thank you for everything!

Sincerely,

Lu

[Bob Carpenter]

> On Aug 25, 2016, at 7:18 AM, Lu Zhang <ny...@g.ucla.edu> wrote:
>
> Hi Prof. Carpenter and Prof. Gelman,
>
> There was a mistake in my code and I fixed it. Stan did an impressive job on my model. It turns out that the resulting samples of two parameters are highly correlated (|rho| > 0.85). I think that's the reason why my vanilla HMC fails to have good performance.

That can sometimes be mitigated with a change of parameterization.
See the manual description of Radford Neal's funnel example for
the simple case.

> Actually, I spent a lot of time on coding in this project. In the beginning, R was the only programming language I was familiar with. But its' low speed and poor memory management annoyed me a lot. So I self-studied C, command lines and makefile, and how to use Lapack and Blas in c. My code works pretty well on a simpler model, but things changed a lot after including one more parameter. So I didn't expect that I can finish Stan code within one day! And it gave me a much better performance than my code in C!

So glad to hear it. That's pretty much exactly where Andrew
was at with R and why he hired me and Matt in the first place
to do the heavy coding (we both knew C). Being computer
scientists, we thought we'd try to solve the problem a little
more generally. Turns out it took a few years and a whole lot
of help from other great programmers and statisticians and
applied mathematicians (see our dev list).

> I know the changing of step size and step length does not break the detail balance,

You can mess with anything during adaptation before you start
the MCMC chains---that's how we do mass matrix and step size.
The real beauty of NUTS is how it changes the number of steps (i.e.,
integration time divided by step size).

> but I never came up with any good idea of using that fact until I saw the "adaptive NUTS" paper this Monday. That is delicate and beautiful! I am also trying to understand RMHMC. Now I can't see any reason for not using Stan.

Matt's a genius! It did take him a couple years and a lot
of input from Andrew to come up with the algorithm.

- Bob

[Lu Zhang]

    Thank you so much for your generous help! I am a student, and I hope I can become a good statistician and a computer scientist like you guys. ^_^

--

You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/ZDiVgFMZ7KU/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

[Bob Carpenter]

I'd recommend good colleagues and practice. :-)

- Bob

[Lu Zhang]

Hi Dr. Bob,

    I posted a topic about my project in Stan user list. Here is the link: https://groups.google.com/forum/#!category-topic/stan-users/general/JwmsztFcqWM. Thanks. ^_^

Sincerely,

Lu