tidy() confuses response and group values for brms-mv-models
- Dominant language
- R
- Stars
- 231
- Forks
- 28
- PR merge metrics
- No merged PRs in 30d
Description
For the random-effect intercepts, the values in the columns `response` and `group` of `tidy()` are inverted, and the values in the `term`-column contain `NA` in its description.
I think these models are very special, and hence this issue probably can't be fixed quickly, so I would not add this to the wish-list for the next update...
``` r
library(sjstats)
library(brms)
library(broom.mixed)
data(efc)
f3 <- bf(neg_c_7 ~ e42dep + c12hour + c172code + (1 |ID| e15relat))
f4 <- bf(c12hour ~ c172code + (1 |ID| e15relat))
b14 <- brm(f3 + f4 + set_rescor(FALSE), data = efc, iter = 500, chains = 1)
b14
#> Family: MV(gaussian, gaussian)
#> Links: mu = identity; sigma = identity
#> mu = identity; sigma = identity
#> Formula: neg_c_7 ~ e42dep + c12hour + c172code + (1 | ID | e15relat)
#> c12hour ~ c172code + (1 | ID | e15relat)
#> Data: efc (Number of observations: 834)
#> Samples: 1 chains, each with iter = 500; warmup = 250; thin = 1;
#> total post-warmup samples = 250
#>
#> Group-Level Effects:
#> ~e15relat (Number of levels: 8)
#> Estimate Est.Error l-95% CI
#> sd(negc7_Intercept) 0.59 0.33 0.13
#> sd(c12hour_Intercept) 23.67 6.84 14.91
#> cor(negc7_Intercept,c12hour_Intercept) 0.52 0.36 -0.36
#> u-95% CI Eff.Sample Rhat
#> sd(negc7_Intercept) 1.32 96 1.02
#> sd(c12hour_Intercept) 40.34 154 1.05
#> cor(negc7_Intercept,c12hour_Intercept) 0.97 39 1.04
#>
#> Population-Level Effects:
#> Estimate Est.Error l-95% CI u-95% CI Eff.Sample Rhat
#> negc7_Intercept 6.78 0.58 5.68 7.86 250 1.00
#> c12hour_Intercept 41.44 9.73 22.92 60.36 173 1.01
#> negc7_e42dep 1.34 0.15 1.05 1.60 250 1.00
#> negc7_c12hour 0.01 0.00 0.00 0.01 250 1.00
#> negc7_c172code 0.33 0.20 -0.06 0.70 250 1.00
#> c12hour_c172code -3.67 2.63 -8.77 0.97 150 1.00
#>
#> Family Specific Parameters:
#> Estimate Est.Error l-95% CI u-95% CI Eff.Sample Rhat
#> sigma_negc7 3.58 0.09 3.43 3.76 250 1.00
#> sigma_c12hour 46.33 1.15 43.96 48.79 250 1.00
#>
#> Samples were drawn using sampling(NUTS). For each parameter, Eff.Sample
#> is a crude measure of effective sample size, and Rhat is the potential
#> scale reduction factor on split chains (at convergence, Rhat = 1).
tidy_stan(b14, type = "all")
#>
#> # Summary Statistics of Stan-Model
#>
#> ## Fixed effects for response: negc7
#>
#> estimate std.error HDI(89%) ratio rhat mcse
#> Intercept 6.75 0.57 [ 5.91 7.69] 1 1 0.04
#> e42dep 1.34 0.16 [ 1.12 1.57] 1 1 0.01
#> c12hour 0.01 0.00 [ 0.00 0.01] 1 1 0.00
#> c172code 0.32 0.20 [ 0.05 0.69] 1 1 0.01
#>
#> ## Random effect (Intercept: e15relat) for response negc7
#>
#> estimate std.error HDI(89%) ratio rhat mcse
#> e15relat.1 0.63 0.45 [-0.01 1.46] 0.23 1.03 0.06
#> e15relat.2 0.33 0.35 [-0.18 0.88] 0.33 1.00 0.04
#> e15relat.3 -0.08 0.31 [-0.84 0.54] 0.68 1.00 0.03
#> e15relat.4 -0.02 0.29 [-0.57 0.61] 0.37 1.01 0.04
#> e15relat.5 -0.18 0.37 [-1.04 0.46] 0.48 1.00 0.04
#> e15relat.6 -0.55 0.51 [-1.47 0.16] 0.49 1.00 0.05
#> e15relat.7 0.15 0.49 [-0.65 1.10] 1.00 1.02 0.04
#> e15relat.8 -0.21 0.31 [-0.74 0.28] 0.31 1.00 0.04
#>
#> ## Fixed effects for response: c12hour
#>
#> estimate std.error HDI(89%) ratio rhat mcse
#> Intercept 40.99 9.93 [26.15 57.00] 0.69 1.01 0.74
#> c172code -3.59 2.74 [-7.60 0.25] 0.60 1.00 0.21
#>
#> ## Random effect (Intercept: e15relat) for response c12hour
#>
#> estimate std.error HDI(89%) ratio rhat mcse
#> e15relat.1 47.18 9.36 [ 35.02 62.55] 0.50 1.02 0.78
#> e15relat.2 2.58 8.06 [-11.09 14.77] 0.52 1.03 0.72
#> e15relat.3 0.12 10.27 [-15.99 15.00] 0.50 1.02 0.87
#> e15relat.4 -0.50 9.26 [-13.31 15.10] 0.65 1.02 0.70
#> e15relat.5 -15.90 12.68 [-36.28 1.29] 1.00 1.00 0.77
#> e15relat.6 -19.46 12.56 [-37.71 -0.25] 1.00 1.00 0.77
#> e15relat.7 20.65 16.87 [ -4.98 49.87] 1.00 1.01 1.07
#> e15relat.8 -13.11 9.56 [-27.17 2.29] 0.67 1.02 0.69
broom.mixed::tidy(b14, conf.int = FALSE)
#> # A tibble: 11 x 7
#> response effect component group term estimate std.error
#>
#> 1 negc7 fixed cond (Intercept) 6.78e+0 0.578
#> 2 c12hour fixed cond (Intercept) 4.14e+1 9.73
#> 3 negc7 fixed cond e42dep 1.34e+0 0.146
#> 4 negc7 fixed cond c12hour 6.26e-3 0.00278
#> 5 negc7 fixed cond c172code 3.29e-1 0.204
#> 6 c12hour fixed cond c172code -3.67e+0 2.63
#> 7 e15relat ran_pa~ cond negc7_Int~ sd__NA.negc7_~ 5.90e-1 0.329
#> 8 e15relat ran_pa~ cond c12hour_I~ sd__NA.c12hou~ 2.37e+1 6.84
#> 9 e15relat ran_pa~ cond negc7_Int~ cor__c12hour_~ 5.20e-1 0.356
#> 10 negc7 ran_pa~ cond Residual sd__Observati~ 3.58e+0 0.0870
#> 11 c12hour ran_pa~ cond Residual sd__Observati~ 4.63e+1 1.15
```
Created on 2018-10-14 by the [reprex package](https://reprex.tidyverse.org) (v0.2.1)
Contributor guide
No contributing guide indexed for this repository
Research direction
Start with broom.mixed::tidy() and reproduce the brms multivariate model shown in the issue. Trace the random-effect rows for response and group, including the term descriptions containing NA; done means these fields correctly identify the response and grouping values without the malformed term text.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- r
- Domain
- data
- Issue type
- Bug
- Difficulty
- 4/5
- Estimated time
- 3-5 days
- Activity status
- Stale
- Clarity
- Mostly clear
- Newbie friendliness
- 35/100