I’ve a likelihood function over some continuous parameter, x. I can calculate the likelihood of different x by: discretizing x-space, letting each x be a state in a markov chain process and finding the stationary frequencies of that (linear algebra involved, specifically: matrix inversion + matrix dot product). In principle I could interpolate between the different x’es to get a continuous distribution. Now, my problem is that this likelihood function is not implemented in Stan, and as far as I can see not implementable directly in Stan?
I currently have the code in a python function, but I could in principle implement that in C++ if necessary. In either case (if not implementable directly in Stan), could you provide me with a pointer how I would go about attaching a custom python/c++ likelihood function to my Stan code?
Is Stan unusually slow in the regular continuous space? Was it evaluating of the likelihood or slow parameter exploration that was bugging you? How many parameters are you talking about interpolating over?
Custom probabilities in Stan are most easily implemented by just manually incrementing the log likelihood (Chapter 21 in the manual).
Apparently, both matrix inverse and product is implemented in Stan, so perhaps I don’t need an external function after all.
Perhaps it was not entirely clear from my first explanation, so let me add a bit more detail on what I’m dealing with. With my current model, my likelihood function depends on three parameters: LikelihoodF(x | param1, param2, param3). I have no analytical expression for the likelihood function as you do for the Triangle in chapter 21. The relevant range of x values is perhaps [20.0 ; 40.0]. Using some linear algebra, I can calculate the probability density function over a discretized space of x’s (for instance 20.0, 20.1, 20.2 … 40.0) for a given set of (param1, param2, param3). If the observed x’s can be rounded (and they can) to the nearest 0.1 decimal, I would perhaps not have to interpolate.
Oh woops, I think I got things backwards, my bad :D. You’re interpolating across x, not your parameters? Would this mean you’d need to rebuild your approximation at every HMC step (every time the parameters changed)?
Purely out of curiosity, what’re you modeling (short desc. is fine)?
Yes. Every time I change my parameters I have to recalculate my pdf.
My project: The link between between protein stability and protein fitness is unknown. I’ve a model of evolution (some unpleasant markov chain stuff + evolutionary dynamics), which defines the pdf over protein stabilities given a protein fitness function. The protein fitness function has three parameters and I want to know the posterior of them given a bunch of protein stability measurements (x). It is pretty basic science and part of my PhD at Lund University :)
Ah, biology stuff is hard. Just yesterday a biochem guy came by the lab asking for advice on a modelling problem. Two hours later and a very informative presentation for him and the only suggestion I could offer was check out ggplot it makes pretty lines :/.
Thanks for your time. I understand that ABC can be used when you can’t calculate the likelihood function. However I my case, I’m actually able to calculate the likelihood function (for every 0.1 increments in x, which is also the precision the data comes in), so wouldn’t it be better to use STAN?
You can do the linear algebra directly in Stan and define your own probability function say
mylikelihood_pmf
since your data is discrete you only need to calculate at the discrete points. I think you can vectorize this. The thing you need to be careful of is whether mylikelihood_pmf is a continuous deterministic differentiable function of the parameters
I assumed that “some unpleasant markov chain stuff + evolutionary dynamics” would mean that you can only approximate the likelihood calculation, in which case there is the choice between 1) exact simulation and ABC and 2) approximated likelihood and Stan.
I didn’t understand yet the restriction to discretized x. Maybe you could tell a bit more about your model?
You are right. I don’t have an analytical likelihood function, but I can calculate the likelihood for discrete values with very high accuracy. In principle I could interpolate between different grid points in x-space to make it continuous. I suppose, that I’m then interested in doing approximated likelihood with Stan?
