dCRP active cluster groups

51 views
Skip to first unread message

Paul vdb

unread,
Nov 23, 2020, 12:04:16 PM11/23/20
to nimble-users

Hello,
The DPMM model I'm implementing I want to be able to count the number of observations in each cluster and subset for just the actively assigned clusters in order to estimate my Poisson intensity.

Essentially for:
z[1:N] ~ dCRP(alpha, size = N)

On each iteration of the MCMC I'd like to be able to get:
ni = # of observations assigned to i = 1:M
and
yi = 0/1 if cluster i = 1:M is active.

my model is very roughly
for(j 1:N)
f(tj) ~ lambda[z[j]]/Lambda[z[j]]
for(i 1:M)
ni ~ Poisson(Lambda[i]*yi)

Is this something I have to calculate within my own nimble distribution when passing z into it? Or can I do this in a simple step within the nimble code?

I've implemented my own sampler in C++ so I know everything is kosher but it would be super cool to be able to do this in standard software. Let me know if this isn't enough information to go on.
Thanks so much!
Paul

Claudia Wehrhahn Cortes

unread,
Nov 23, 2020, 1:06:56 PM11/23/20
to Paul vdb, nimble-users
Hi Paul,

nice that you are using the dCRP function!

To count the number of observations per cluster and to know which are the active clusters, you would have to do it as a post processing, by looking at the MCMC’s output. This information will be in the ‘z’ vector. ‘z’ contains the assignment of observations to active clusters.  Remember to add it in the ‘monitors’ argument. 

Let me ask you something: do you want to do this to do inference or to sample the intensity? I think is for the later, and if so, you shouldn’t need to do any of these things. They are already done in nimble when sampling ‘z’. Let me know if this is what you want to do so I can provide more details. Also, having more details about the nimble code would be helpful.

I hope this is useful.

Best,
Claudia

--
You received this message because you are subscribed to the Google Groups "nimble-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to nimble-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/nimble-users/8dd3cd99-519f-46ce-9356-ba7fa9456370n%40googlegroups.com.

Paul vdb

unread,
Nov 23, 2020, 1:32:18 PM11/23/20
to nimble-users
Thanks for the awesomely quick response Claudia. I love the popularity that the DPMM is getting in ecology and stoked you guys have included it in the sampler. Thanks for that! I did see that once I've clustered I can access all that information conveniently. I haven't written it in Nimble yet as I don't really know how to get past this part and am avoiding investing completely before I know if it's possible.

Here is hopefully a more digestible form of the idea in standard form:
I have a very specific form of what lambda_i, the observation intensity is based on latent locations x_i, lambda(x_i).

For a known cluster K, the likelihood is from a non-homogeneous Poisson process:
exp(-sum_k=1^K (int_0^T lambda_k)) * prod_1^N lambda_i(x_i) * prod_1^K f(x_k)

However I don't observe K or the latent locations that create some unknown intensity which leads me to x_i ~ G.
For K active clusters,
The likelihood under the DPMM is then exp(-sum_k=1^K (int_0^T lambda_k)) * prod_1^N int lambda_i(x_i) dG0(x_i).

So I kind of need that extra exponential term in the MCMC to make inference on the lambda and x terms.
When I define f(t_i) as lambda(x_i)/ Lambda(x_i) where Lambda(x_i) = int_0^T lambda(x_i) dt then after I've clustered I can just call that extra term as proportional to a Poisson.

I suppose a better question is how can I incorporate the active clusters directly into my likelihood so that I can estimate that additional exponential term that is outside of the direct clustering from the CRP? Maybe this can be thought of as a hierarchical problem where there are shared parameters between the CRP and the additional terms.
Thanks so much,
Paul

Claudia Wehrhahn Cortes

unread,
Dec 7, 2020, 1:38:58 PM12/7/20
to Paul vdb, nimble-users
Hi Paul,

Sorry for the late response about this and I am also sorry that I'll not give you good news.

NIMBLE’s mcmc does not handle sampling when one or more nodes in the model are functions of a node having the CRP distribution (dCRP). In your case, the functions would compute "ni" and "yi" in your previous model. So, if you need to get the unique labels and counts for the likelihood, you would not be able to sample the model.  

There are still many things that we should do to use NIMBLE in more general BNP models.
Let me know if there is something else I could help with.

Best,
Claudia


Paul vdb

unread,
Dec 7, 2020, 6:00:58 PM12/7/20
to nimble-users
Thanks for the response Claudia.
I figured as much. At least all my effort to get my version working is not time wasted since I'm a real buggy coder. I've implemented the similar protocols that Nimble uses for block Random Walk samplers and appreciate all your infrastructure and documentation to make these decisions transparent. I'll keep an eye out for when these things might change.
Cheers

Daniel Turek

unread,
Dec 7, 2020, 7:38:48 PM12/7/20
to Paul vdb, nimble-users
Paul, just one quick additional point.  Another possible approach to calculating the ni and yi derived quantities could be writing a separate MCMC sampler or two, serving only to calculate these quantities during the course of MCMC sampling.  This would be done independently of the updating of the CRP group-membership variables (z, in your model), and the rest of your model, and perhaps would avoid the other problems you're running into.

But, that said, I'm not entirely sure what you've implemented in C++, or how much of NIMBLE's base MCMC infrastructure you're making use of, so I'm not sure how applicable this suggestion is.  But, if some of my guesses are right, then it might easily sidestep your problem.  Happy to follow up more, if that's helpful.

Cheers,
Daniel


Paul vdb

unread,
Dec 8, 2020, 1:56:10 PM12/8/20
to nimble-users
Sounds like there is some real potential to get it to work. I'd be very keen to hear more about that. My email is paul dot vandambates at gmail dot com.
Thanks,
Paul

Daniel Turek

unread,
Dec 8, 2020, 5:01:53 PM12/8/20
to Paul vdb, nimble-users
Sounds good, Paul, I'll follow-up with you directly by email.

Reply all
Reply to author
Forward
0 new messages