<p><div><div><div><div data-show-faces="false" data-height="21" data-width="100" data-layout="button_count" data-send="true" data-href="https://www.r-bloggers.com/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/"></div></div><div></div></div><div></div></div>
<div>
(This article was first published on <strong><a href="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/"> R – Strenge Jacke!</a></strong>, and kindly contributed to <a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers)</a> </div>
<p>Another quick preview of my <a rel="nofollow" href="https://github.com/strengejacke" target="_blank">R-packages</a>, especially <strong>sjPlot</strong>, which now also support <code>brmsfit</code>-objects from the great <a rel="nofollow" href="https://cran.r-project.org/package=brms" target="_blank">brms-package</a>. To demonstrate the new features, I load all my „core“-packages at once, using the <em>strengejacke</em>-package, which is only available from GitHub. This package simply loads four packages (<em>sjlabelled, sjmisc, sjstats and sjPlot</em>).</p>
<p>First, I fit two sample models, one with <strong>brms</strong> and one with <strong>rstanarm</strong>:</p>
<pre><code><span>#> install pkg &quot;strengejacke&quot; via GitHub:
#> devtools::install_github(&quot;strengejacke/strengejacke&quot;)</span>
library(strengejacke)
library(ggplot2)
library(lme4)
library(glmmTMB)
library(brms)
library(rstanarm)

data(Owls)
m1 &lt;- brm(
  SiblingNegotiation ~ FoodTreatment + ArrivalTime + SexParent + (1 | Nest), 
  data = Owls,
  family = zero_inflated_poisson(link = &quot;log&quot;, link_zi = &quot;logit&quot;)
)
m2 &lt;- stan_glmer.nb(
  SiblingNegotiation ~ FoodTreatment + ArrivalTime + SexParent + (1 | Nest), 
  data = Owls
)
m3 &lt;- stan_lmer(Reaction ~ Days + (Days | Subject), data = sleepstudy)</code></pre>
<p>Next, some summary statistics, at first simply the <em>High Density Interval</em>, obtained from the <code>hdi()</code>-function. You can define the probabilites of the interval, however, here I simply use the defaults (which is the 90% HDI).</p>
<pre><code>hdi(m1)
<span>#> # A tibble: 6 x 3
#> term hdi.low hdi.high
#> #> 1 b_Intercept 3.798693e+00 4.584998e+00
#> 2 b_FoodTreatmentSatiated -2.795021e-01 -1.586102e-01
#> 3 b_ArrivalTime -9.768530e-02 -6.616026e-02
#> 4 b_SexParentMale -7.766154e-02 3.794376e-02
#> 5 zi 2.304234e-01 2.900676e-01
#> 6 lp__ -2.020651e+03 -2.002387e+03</span>

hdi(m2)
<span>#> # A tibble: 5 x 3
#> term hdi.low hdi.high
#> #> 1 (Intercept) 3.82688888 5.93003148
#> 2 FoodTreatmentSatiated -0.86327579 -0.51376279
#> 3 ArrivalTime -0.15414266 -0.07134095
#> 4 SexParentMale -0.09371495 0.24440447
#> 5 reciprocal_dispersion 0.77345433 1.00889671</span></code></pre>
<p>Next, I implemented an own tidyr-function, similar to <strong>broom’s</strong> <code>tidy()</code>, called <code>tidy_stan()</code>. Unlike broom’s tidy-function, <code>tidy_stan()</code> computes the HDI instead of credibility intervals, and you may return fixed or random effects only, or both. Furthermore, <code>tidy_stan()</code> also returns the ratio of effective numbers of samples, <em>n_eff</em>, and Rhat statistics. </p>
<pre><code>tidy_stan(m1)
<span>#> # A tibble: 6 x 7
#> term estimate std.error hdi.low hdi.high n_eff Rhat
#> #> 1 b_Intercept 4.206 0.241 3.819 4.581 1.000 1.000
#> 2 b_FoodTreatmentSatiated -0.220 0.037 -0.279 -0.162 1.000 1.000
#> 3 b_ArrivalTime -0.082 0.009 -0.096 -0.065 1.000 1.000
#> 4 b_SexParentMale -0.021 0.037 -0.077 0.036 1.000 1.000
#> 5 zi 0.259 0.019 0.230 0.288 1.000 0.999
#> 6 lp__ -2011.050 5.537 -2020.651 -2002.925 0.186 1.001</span>

tidy_stan(m2)
<span>#> # A tibble: 5 x 7
#> term estimate std.error hdi.low hdi.high n_eff Rhat
#> #> 1 (Intercept) 4.873 0.652 3.847 5.880 1 1.000
#> 2 FoodTreatmentSatiated -0.693 0.112 -0.860 -0.522 1 1.000
#> 3 ArrivalTime -0.115 0.025 -0.154 -0.074 1 1.000
#> 4 SexParentMale 0.071 0.102 -0.088 0.241 1 0.999
#> 5 reciprocal_dispersion 0.883 0.073 0.776 1.005 1 1.000</span></code></pre>
<p>Finally, the new and generic <code>plot_model()</code> function in <strong>sjPlot</strong> also supports brmsfit or stanreg objects, and plots estimates, random effects and marginal effects.</p>
<p>The coefficients-plot was already shown <a rel="nofollow" href="https://strengejacke.wordpress.com/2017/09/15/quick-sjplot-status-update-rstats-rstanarm-ggplot2/" target="_blank">in a previous post</a>. It shows the 50% and 89% HDI and the posterior median. The style can be changed using the <code>bpe</code> and <code>bpe.style</code> arguments.</p>
<pre><code>theme_set(theme_sjplot())
plot_model(m1, axis.lim = c(.6, 1.4))
plot_model(m2, bpe = &quot;mean&quot;, bpe.style = &quot;dot&quot;, colors = &quot;grey30&quot;)</code></pre>
<p><strong>Coefficients, brms-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot.png" data-orig-size="700,400" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot.png?w=300" data-attachment-id="1413" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot.png?w=450&h=375&fit=656%2C375" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot-6/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot.png?w=450&amp;h=375&amp;fit=656%2C375" alt="" data-image-title="Rplot" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1413" data-recalc-dims="1"></a></p>
<p><strong>Coefficients, rstanarm-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot02.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot02.png" data-orig-size="700,400" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot02.png?w=300" data-attachment-id="1415" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot02.png?w=450&h=375&fit=656%2C375" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot02-4/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot02.png?w=450&amp;h=375&amp;fit=656%2C375" alt="" data-image-title="Rplot02" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1415" data-recalc-dims="1"></a></p>
<p>Random effects are displayed in a similar way. For models with random slopes and intercepts, random effects are displayed in a grid layout (use <code>grid = FALSE</code> to create a separate plot for each random effect).</p>
<pre><code>plot_model(m1, type = &quot;re&quot;)
plot_model(m3, type = &quot;re&quot;)</code></pre>
<p><strong>Random effects, brms-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot03.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot03.png" data-orig-size="700,600" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot03.png?w=300" data-attachment-id="1417" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot03.png?w=450&h=562&fit=656%2C562" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot03/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot03.png?w=450&amp;h=562&amp;fit=656%2C562" alt="" data-image-title="Rplot03" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1417" data-recalc-dims="1"></a></p>
<p><strong>Random effects, rstanarm-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot04.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot04.png" data-orig-size="700,600" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot04.png?w=300" data-attachment-id="1418" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot04.png?w=450&h=562&fit=656%2C562" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot04/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot04.png?w=450&amp;h=562&amp;fit=656%2C562" alt="" data-image-title="Rplot04" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1418" data-recalc-dims="1"></a></p>
<p>Marginal effects in <strong>sjPlot</strong> are based on the <a rel="nofollow" href="https://cran.r-project.org/package=ggeffects" target="_blank">ggeffects-package</a>. By default, predictions are computed with <code>rstantools::posterior_linpred()</code>.</p>
<pre><code>plot_model(m1, type = &quot;pred&quot;, terms = c(&quot;ArrivalTime&quot;, &quot;FoodTreatment&quot;))
plot_model(m2, type = &quot;pred&quot;, terms = c(&quot;ArrivalTime&quot;, &quot;FoodTreatment&quot;))</code></pre>
<p><strong>Marginal effects, brms-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot05.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot05.png" data-orig-size="700,500" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot05.png?w=300" data-attachment-id="1419" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot05.png?w=450&h=469&fit=656%2C469" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot05/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot05.png?w=450&amp;h=469&amp;fit=656%2C469" alt="" data-image-title="Rplot05" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1419" data-recalc-dims="1"></a></p>
<p><strong>Marginal effects, rstanarm-model</strong><br>
<a rel="nofollow" href="https://strengejacke.files.wordpress.com/2017/10/rplot1.png" target="_blank"><img data-orig-file="https://strengejacke.files.wordpress.com/2017/10/rplot1.png" data-orig-size="700,500" data-medium-file="https://strengejacke.files.wordpress.com/2017/10/rplot1.png?w=300" data-attachment-id="1423" data-image-meta="{" data-image-description="" src="https://strengejacke.files.wordpress.com/2017/10/rplot1.png?w=450&h=469&fit=656%2C469" data-permalink="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/rplot-7/" data-large-file="https://strengejacke.files.wordpress.com/2017/10/rplot1.png?w=450&amp;h=469&amp;fit=656%2C469" alt="" data-image-title="Rplot" sizes="(max-width: 656px) 100vw, 656px" data-comments-opened="1" class="wp-image-1423" data-recalc-dims="1"></a></p>
<p>The work on sjPlot for the next CRAN release is almost done… I hope to submit the sjPlot-update in the course of the next week! In the meantime, you may try out the new features using the <a rel="nofollow" href="https://github.com/strengejacke" target="_blank">GitHub-versions of my R-packages</a> – at least <em>sjstats</em> and <em>ggeffects</em> are required.</p>
<p> Tagged: <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/brms/" target="_blank">brms</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/data-visualization/" target="_blank">data visualization</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/ggplot/" target="_blank">ggplot</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/r/" target="_blank">R</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/rstan/" target="_blank">rstan</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/rstanarm/" target="_blank">rstanarm</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/rstats/" target="_blank">rstats</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/sjplot/" target="_blank">sjPlot</a>, <a rel="nofollow" href="https://strengejacke.wordpress.com/tag/stan/" target="_blank">Stan</a> <a rel="nofollow" href="http://feeds.wordpress.com/1.0/gocomments/strengejacke.wordpress.com/1409/" target="_blank"><img border="0" alt="" src="http://feeds.wordpress.com/1.0/comments/strengejacke.wordpress.com/1409/"></a> <img border="0" alt="" src="https://i1.wp.com/pixel.wp.com/b.gif?resize=1%2C1&ssl=1" data-recalc-dims="1"></p>
	<p>
<div>
<div>To <strong>leave a comment</strong> for the author, please follow the link and comment on their blog: <strong><a href="https://strengejacke.wordpress.com/2017/10/11/more-support-for-bayesian-analysis-in-the-sj-packages-rstats-rstan-brms/"> R – Strenge Jacke!</a></strong>.</div>
<hr>
<a rel="nofollow" href="https://www.r-bloggers.com/">R-bloggers.com</a> offers <strong><a rel="nofollow" href="https://feedburner.google.com/fb/a/mailverify?uri=RBloggers">daily e-mail updates</a></strong> about <a rel="nofollow" href="https://www.r-project.org/" title="The R Project for Statistical Computing">R</a> news and <a rel="nofollow" href="https://www.r-bloggers.com/search/tutorial" title="R tutorials">tutorials</a> on topics such as: <a rel="nofollow" href="https://www.r-bloggers.com/search/data%20science" title="Data science">Data science</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Big%20Data" title="Big Data">Big Data, <a rel="nofollow" href="https://www.r-users.com/" title="R jobs">R jobs</a>, visualization (<a rel="nofollow" href="https://www.r-bloggers.com/search/ggplot2" title="ggplot and ggplot2 tutorials">ggplot2</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/boxplot" title="Boxplots using lattice and ggplot2 tutorials">Boxplots</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/map" title="Maps and gis">maps</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/animation" title="Animation in R">animation</a>), programming (<a rel="nofollow" href="https://www.r-bloggers.com/search/RStudio" title="RStudio IDE for R">RStudio</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/sweave" title="Sweave and literate programming">Sweave</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/LaTeX" title="LaTeX in R">LaTeX</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/SQL" title="SQL and databases">SQL</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/eclipse" title="Eclipse IDE for R">Eclipse</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/git" title="git and github, Version Control System">git</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/hadoop" title="Large data in R using Hadoop">hadoop</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/Web+Scraping" title="Web Scraping of google, facebook, yahoo, twitter and more using R">Web Scraping</a>) statistics (<a rel="nofollow" href="https://www.r-bloggers.com/search/regression" title="Regressions and ANOVA analysis tutorials">regression</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/PCA" title="principal component analysis tutorial">PCA</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/time+series" title="Time series">time series</a>, <a rel="nofollow" href="https://www.r-bloggers.com/search/trading" title="finance trading">trading</a>) and more...
<img alt="" width="1" src="http://feeds.feedburner.com/~r/RBloggers/~4/iuS72btPsa8" height="1"></a></div>