# Example of an rhat of 6e+13

**URL:** <https://discourse.mc-stan.org/t/example-of-an-rhat-of-6e-13/27982>\
**Category:** General\
**Tags:** brms\
**Created:** [June 28, 2022, 10:03am UTC](https://discourse.mc-stan.org/t/example-of-an-rhat-of-6e-13/27982 "2022-06-28T10:03:29Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![scholz](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/scholz/32/14298_2.png) [@scholz](https://discourse.mc-stan.org/u/scholz)\
**Post date:** [June 28, 2022, 10:03am UTC](https://discourse.mc-stan.org/t/example-of-an-rhat-of-6e-13/27982/1 "2022-06-28T10:03:30Z")

</div>

This post is based on @avehtari asking for examples of high rhat values on twitter. Turns out my simulation study is also producing those :)  
**This is not a question for how to make the model work!**  
This is just a fun showcase of really bad rhat values for Aki (and everyone else interested.)

Here is the brms fit object for you to play around with :)  
[fit\_object.rdata](https://discourse.mc-stan.org/uploads/short-url/6NC2b9L3xuWe9LAraHne8ifDRu7.rdata) (1.6 MB)

## Dataset

The dataset of size 100 is generated from the following model:  
y \sim lognormal(\mu\_y, 2)  
\mu\_y = 0.25x + 0.5z\_1 + 0.8z\_2  
x \sim N(\mu\_x, 0.5)  
\mu\_x = 0.5z\_1 + 0.8z\_3  
z\_1, z\_2, z\_3 \sim N(0,0.5)  
z4 \sim N(\mu\_{z4}, 0.5)  
\mu\_{z4} = y + 0.5x

The DAG looks like this:

 ![DAG-1](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/f/f899c5dd05893bdcf3415fe8a03b591c2c4eaaf0.png)

## Model

The dataset is then fit with the following model:  
y \sim Gompertz(\mu, \beta)  
log(\mu) = 1 + x + z1 + z2 + z4  
With flat priors on everything, 2 chains, iter = 2500, warmup = 500

## Software

I am using R (4.2), cmdstan(r) (2.29.2 / 0.5.2) and brms (2.17.4) and my own Gompertz custom family you can find [here](https://github.com/sims1253/bayesim/blob/master/R/gompertz.R).

## Fit

```no-highlight
summary(fit)
Family: gompertz 
  Links: mu = log; beta = identity 
Formula: y ~ x + z1 + z2 + z4 
   Data: dataset (Number of observations: 100) 
  Draws: 2 chains, each with iter = 2500; warmup = 500; thin = 1;
         total post-warmup draws = 4000

Population-Level Effects: 
          Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept -1.47 0.56 -2.03 -0.91 NA 2 NA
x -1.27 0.57 -1.84 -0.70 66142364047444.67 2 NA
z1 -0.11 0.33 -0.44 0.22 NA 2 NA
z2 0.36 1.54 -1.18 1.91 NA 2 NA
z4 0.02 0.00 0.02 0.03 66142364047444.67 2 NA

Family Specific Parameters: 
     Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
beta 0.54 0.12 0.42 0.66 66142364047444.67 2 NA

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).
Warning messages:
1: Parts of the model have not converged (some Rhats are > 1.05). Be careful when analysing the results! We recommend running more iterations and/or setting stronger priors. 
2: There were 1029 divergent transitions after warmup. Increasing adapt_delta above 0.8 may help. See http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup

```

---

<div class="post-metadata">

**Author:** ![avehtari](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/avehtari/32/5935_2.png) [@avehtari](https://discourse.mc-stan.org/u/avehtari)\
**Post date:** [June 28, 2022, 12:59pm UTC](https://discourse.mc-stan.org/t/example-of-an-rhat-of-6e-13/27982/2 "2022-06-28T12:59:03Z")

</div>

rhat should be NA for all parameters. There are 2 chains and both of them are stuck, and each chain is just a constant. It’s strange that the code is ending up estimating within-chain-variance to be non-zero. Need to check the code that constant chains are detected.