I’m working on the next update of my sjPlot-package, which will get a generic
plot_model() method, which plots any kind of regression model, with different plot types being supported (forest plots for estimates, marginal effects and predictions, including displaying interaction terms, …).
The package also supports
rstanarm models. Since these are typically presented in a slightly different way (e.g., „outer“ and „inner“ probability of credible intervals), I implemented a special handling for these models, for which I wanted to show a quick preview here:
Weiterlesen „Quick #sjPlot status update… #rstats #rstanarm #ggplot2“
Here’s a small preview of forthcoming features in the ggeffects-package, which are already available in the GitHub-version: For marginal effects from models fitted with
glmer.nb(), confidence intervals are now also computed.
If you want to test these features, simply install the package from GitHub:
Here are three examples:
Weiterlesen „Marginal effects for negative binomial mixed effects models (glmer.nb and glmmTMB) #rstats“
Some time ago I started working with Bayesian methods, using the great rstanarm-package. Beside the fantastic package-vignettes, and books like Statistical Rethinking or Doing Bayesion Data Analysis, I also found the ressources from Tristan Mahr helpful to both better understand Bayesian analysis and rstanarm. This motivated me to implement tools for Bayesian analysis into my packages, as well.
Due to the latest tidyr-update, I had to update some of my packages, in order to make them work again, so – beside some other features – some Bayes-stuff is now avaible in my packages on CRAN.
Weiterlesen „Going Bayes #rstats“
My sjstats-package has been updated on CRAN. The past updates introduced new functions for various purposes, e.g. predictive accuracy of regression models or improved support for the marvelous glmmTMB-package. The current update, however, added some ANOVA tools to the package.
Weiterlesen „Effect Size Statistics for Anova Tables #rstats“
I started writing my first package as collection of various functions that I needed for (almost) daily work. Meanwhile, packages were growing and bit by bit I sourced out functions to put them into new packages. Although this means more work for CRAN members when they have more packages to manage on their network, from a user-perspective it is much better if packages have a clear focus and a well defined set of functions. That’s why I now released a new package on CRAN, sjlabelled, which contains all functions that deal with labelled data. These functions use to live in the sjmisc-package, where they now are deprecated and will be removed in a future update.
My aim is not only to provide packages with a clear focus, but also with a consistent design and philosophy, making it easier and more intuitive to use (see also here) – I prefer to follow the so called „tidyverse“-approach here. It is still work in progress, but so far I think I’m on a good way…
So, what are the packages I use for (almost) daily work?
- sjlabelled – reading, writing and working with labelled data (especially since I collaborate a lot with people who use Stata or SPSS)
- sjmisc – data and variable transformation utilities (the complement to dplyr and tidyr, when it comes down from data frames to variables within the data wrangling process)
- sjstats – out-of-the-box statistics that is not already provided by base R or packages
- sjPlot – to quickly generate tables and plot
- ggeffects – to visualize regression models
The next step is creating cheat sheets for my packages. I think if you can map the scope and idea of your package (functions) on a cheat sheet, its focus is probably well defined.
Btw, if you also use some of the above packages more or less regularly, you can install the „strengejacke“-package to load them all in one step. This package is not on CRAN, because its only purpose is to load other packages.
Disclaimer: Of course I use other packages everyday as well – this posting is just focussing on my packages that I created because I frequently needed these kind of functions.
Aim of the ggeffects-package
The aim of the ggeffects-package is similar to the broom-package: transforming “untidy” input into a tidy data frame, especially for further use with ggplot. However, ggeffects does not return model-summaries; rather, this package computes marginal effects at the mean or average marginal effects from statistical models and returns the result as tidy data frame (as tibbles, to be more precisely).
Since the focus lies on plotting the data (the marginal effects), at least one model term needs to be specified for which the effects are computed. It is also possible to compute marginal effects for model terms, grouped by the levels of another model’s predictor. The package also allows plotting marginal effects for two- or three-way-interactions, or for specific values of a model term only. Examples are shown below.
Weiterlesen „ggeffects: Create Tidy Data Frames of Marginal Effects for ‚ggplot‘ from Model Outputs #rstats“
The survey-package from Thomas Lumley is a great toolkit when analyzing complex samples. It provides
svyglm(), to fit generalised linear models to data from a complex survey design.
svyglm() covers all families that are also provided by R’s
glm() – however, the survey-package has no function to fit negative binomial models, which might be useful for overdispersed count models. Yet, the package provides a generic
svymle() to fit user-specified likelihood estimations. In his book, Appendix E, Thomas Lumley describes how to write your own likelihood-function, passed to
svymle(), to fit negative binomial models for complex samples. So I wrote a small „wrapper“ and implemented a function
svyglm.nb() in my sjstats-package.
# This example reproduces the results from
# Lumley 2010, figure E.7 (Appendix E, p256)
# create survey design
des <- svydesign(
id = ~ SDMVPSU,
strat = ~ SDMVSTRA,
weights = ~ WTINT2YR,
nest = TRUE,
data = nhanes_sample
# fit negative binomial regression
fit <- svyglm.nb(total ~ factor(RIAGENDR) * (log(age) + factor(RIDRETH1)), des)
# print coefficients and standard errors
round(cbind(coef(fit), survey::SE(fit)), 2)
#> [,1] [,2]
#> theta.(Intercept) 0.81 0.05
#> eta.(Intercept) 2.29 0.16
#> eta.factor(RIAGENDR)2 -0.80 0.18
#> eta.log(age) 1.07 0.23
#> eta.factor(RIDRETH1)2 0.08 0.15
#> eta.factor(RIDRETH1)3 0.09 0.18
#> eta.factor(RIDRETH1)4 0.82 0.30
#> eta.factor(RIDRETH1)5 0.06 0.38
#> eta.factor(RIAGENDR)2:log(age) -1.22 0.27
#> eta.factor(RIAGENDR)2:factor(RIDRETH1)2 -0.18 0.26
#> eta.factor(RIAGENDR)2:factor(RIDRETH1)3 0.60 0.19
#> eta.factor(RIAGENDR)2:factor(RIDRETH1)4 0.06 0.37
#> eta.factor(RIAGENDR)2:factor(RIDRETH1)5 0.38 0.44
The functions returns an object of class
svymle, so all methods provided by the survey-package for this class work – it’s just that there are only a few, and common methods like
predict() are currently not implemented. Maybe, hopefully, future updates of the survey-package will include such features.
When describing a sample, researchers in my field often show proportions of specific characteristics as description. For instance, proportion of female persons, proportion of persons with higher or lower income etc. Since it happens often that I like to know these characteristics when exploring data, I decided to write a function,
prop(), which is part of my sjstats-package – a package dedicated to summary-functions, mostly for fit- or association-measures of regression models or descriptive statistics.
Weiterlesen „Descriptive summary: Proportions of values in a vector #rstats“
A new update of my sjPlot-package was just released on CRAN. Thanks to @c_schwemmer, it’s now possible to easily integrate the HTML-ouput of all table-functions into knitr-rmarkdown-documents.
Weiterlesen „Direct integration of sjPlot-tables in knitr-rmarkdown-documents #rstats“
Yesterday, I had the pleasure to give a talk at the 8th Hamburg R User-Group meeting. I talked about data wrangling and data transformation, and how the philosophy behind the tidyverse makes these tasks easier. If you like, you can download the slides here (or at the GitHub-Repo from the Hamburg R-User-Group). Feel free to add your comments to the slide here.