# Sampling from a Polya Gamma distribution?

**URL:** https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067
**Category:** Modeling
**Tags:** specification
**Created:** [March 14, 2019, 6:43pm UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067 "2019-03-14T18:43:35Z")
**Posts on this page:** 6
**Page:** 1

<div class="post-metadata">

### Author: ![Red-Portal](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/red-portal/32/5349_2.png) [@Red-Portal](https://discourse.mc-stan.org/u/Red-Portal)
#### Post date: [March 14, 2019, 6:43pm UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/1 "2019-03-14T18:43:35Z")

</div>

Hello everyone.  
Anybody tried to implement the Polya-Gamma distribution?  
Or is there any plan for supporting the Polya-Gamma distribution?

> **[Bayesian inference for logistic models using Polya-Gamma latent variables](https://arxiv.org/abs/1205.0310)**
>
> We propose a new data-augmentation strategy for fully Bayesian inference in
> models with binomial likelihoods. The approach appeals to a new class of
> Polya-Gamma distributions, which are constructed in detail. A variety of
> examples are presented to...

  
This paper claims that sampling from PG distribution is beneficial for logistic models.  
Wonder if anybody tried this on Stan.

---

<div class="post-metadata">

### Author: ![bgoodri](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bgoodri/32/4451_2.png) [@bgoodri](https://discourse.mc-stan.org/u/bgoodri)
#### Post date: [March 14, 2019, 6:49pm UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/2 "2019-03-14T18:49:28Z")

</div>

This is quite useful for using a Gibbs sampler to estimate the coefficients of logit models using data augmentation, and in fact, it can be proven that this particular Gibbs sampler is uniformally ergodic, which is rare. But the approach isn’t that useful for Stan, which uses No U-turn Sampling and joint proposals rather than full-conditional proposals and usually samples from the posterior distribution of the coefficients of a logit model (without data augmentation) quite efficiently as well.

The Polya Gamma distribution might be useful for other purposes, although I haven’t seen any and I don’t know of anyone who has tried to implement it in Stan in part because that infinite sum with alternating signs would be challenging to do in a numerically stable way.

---

<div class="post-metadata">

### Author: ![Red-Portal](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/red-portal/32/5349_2.png) [@Red-Portal](https://discourse.mc-stan.org/u/Red-Portal)
#### Post date: [March 15, 2019, 3:09am UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/3 "2019-03-15T03:09:01Z")

</div>

Thanks for the answer. I was trying to implement [1, 2] which they used PG for modeling the likelihood of stick-breaking processes. Maybe it’s not beneficial in this case too?

[1] Poddar, Lahari, Wynne Hsu, and Mong Li Lee. “Quantifying aspect bias in ordinal ratings using a bayesian approach.” _arXiv preprint arXiv:1705.05098_ (2017).

[2] Chen, Jianfei, et al. “Scalable inference for logistic-normal topic models.” _Advances in Neural Information Processing Systems_ . 2013.

---

<div class="post-metadata">

### Author: ![ergou\_chen](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/ergou_chen/32/6483_2.png) [@ergou\_chen](https://discourse.mc-stan.org/u/ergou_chen)
#### Post date: [November 19, 2019, 1:38am UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/4 "2019-11-19T01:38:58Z")

</div>

Continuing the discussion from [Sampling from a Polya Gamma distribution?](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/3):

Hello,I’m also interested in this problem. And I wonder if you have solved the problem? Did you implement the code in the articles

---

<div class="post-metadata">

### Author: ![rok\_cesnovar](https://avatars.discourse-cdn.com/v4/letter/r/7bcc69/32.png) [@rok\_cesnovar](https://discourse.mc-stan.org/u/rok_cesnovar)
#### Post date: [November 19, 2019, 7:36am UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/5 "2019-11-19T07:36:47Z")

</div>

As far as I know there arent any publicly available implementations for Stan.

Not sure if this helps and apologize if you already aware of this but the [BayesLogit](https://cran.r-project.org/web/packages/BayesLogit/index.html) R packages has the implementation of PolyaGamma distribution. Its done in C so its faster than if you would reimplement PG in R but its really not super efficient.

3+ years ago me and Erik rewrote their implementation in C and then also made a GPU variant of it that is roughly a 100x faster altogether. Its also in a CRAN R package tho that is not being actively maintained and not sure if I will ever come back to that ([https://cran.r-project.org/web/packages/bayesCL/index.html](https://cran.r-project.org/web/packages/bayesCL/index.html)).

But if anyone wants to rewrite my C code in C++ for Stan Math I am happy to provide more info or help on that.

---

<div class="post-metadata">

### Author: ![Red-Portal](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/red-portal/32/5349_2.png) [@Red-Portal](https://discourse.mc-stan.org/u/Red-Portal)
#### Post date: [November 19, 2019, 12:04pm UTC](https://discourse.mc-stan.org/t/sampling-from-a-polya-gamma-distribution/8067/6 "2019-11-19T12:04:34Z")

</div>

@ergou_chen Yes I implemented a similar model on Stan.  
Instead of using the Polya-gamma sampling scheme,  
I simply typed in the model and ran the default NUTS sampler and it just worked.
