modelbased

modelbased is a lightweight package helping with model-based estimations, used in the computation of marginal means, contrast analysis and predictions.

Installation

Run the following:

install.packages("devtools")
devtools::install_github("easystats/modelbased")

library("modelbased")

Documentation

Click on the buttons above to access the package documentation and the easystats blog, and check-out these vignettes:

• Visualisation matrix
• Marginal means
• Contrast analysis
• Use a model to make predictions
• Describe non-linear curves

Features

The package is built around 5 main functions:

• estimate_means(): Estimates the average values at each factor levels
• estimate_contrasts(): Estimates and tests contrasts between different factor levels
• estimate_slopes(): Estimates the slopes of numeric predictors at different factor levels
• estimate_response(): Predict the response variable using the model
• estimate_smooth(): Describes a non-linear term (e.g. in GAMs) by its linear parts

These functions are powered by the visualisation_matrix() function, a smart tool for guessing the appropriate reference grid.

The package currently only supports rstanarm models, but will be expanded to cover a large variety of frequentist and Bayesian models.

Examples

Create smart grids to represent complex interactions

Check-out this vignette to create this plot:

Estimate marginal means

Check-out this vignette to create this plot:

library(rstanarm)

model <- stan_glm(Sepal.Width ~ Species, data = iris)

estimate_means(model)

## Species    | Mean |       95% CI
## --------------------------------
## setosa     | 3.43 | [3.33, 3.52]
## versicolor | 2.77 | [2.68, 2.87]
## virginica  | 2.97 | [2.88, 3.07]

Contrast analysis

Check-out this vignette to create this plot:

estimate_contrasts(model)
## Level1     |     Level2 | Difference |         95% CI |     pd | % in ROPE | Difference (std.)
## ----------------------------------------------------------------------------------------------
## setosa     | versicolor |       0.66 | [ 0.52,  0.78] |   100% |        0% |              1.50
## setosa     |  virginica |       0.45 | [ 0.31,  0.58] |   100% |        0% |              1.04
## versicolor |  virginica |      -0.20 | [-0.34, -0.07] | 99.78% |     6.95% |             -0.47

Check the contrasts at different points of another linear predictor

model <- stan_glm(Sepal.Width ~ Species * Petal.Length, data = iris)

estimate_contrasts(model, modulate = "Petal.Length", length = 3)

## Level1     |     Level2 | Petal.Length | Difference |        95% CI |     pd | % in ROPE | Difference (std.)
## ------------------------------------------------------------------------------------------------------------
## setosa     | versicolor |         1.00 |       1.67 | [ 1.06, 2.31] |   100% |        0% |              3.83
## setosa     |  virginica |         1.00 |       1.34 | [ 0.59, 2.05] | 99.95% |     0.05% |              3.08
## versicolor |  virginica |         1.00 |      -0.33 | [-1.25, 0.63] | 76.10% |    13.33% |             -0.75
## setosa     | versicolor |         3.95 |       1.65 | [ 0.76, 2.64] | 99.98% |     0.02% |              3.78
## setosa     |  virginica |         3.95 |       1.71 | [ 0.71, 2.70] | 99.98% |        0% |              3.92
## versicolor |  virginica |         3.95 |       0.06 | [-0.22, 0.35] | 66.27% |    47.70% |              0.14
## setosa     | versicolor |         6.90 |       1.62 | [-0.41, 3.68] | 93.70% |     2.30% |              3.71
## setosa     |  virginica |         6.90 |       2.07 | [ 0.03, 4.04] | 97.80% |     1.05% |              4.76
## versicolor |  virginica |         6.90 |       0.43 | [-0.11, 0.97] | 94.47% |     8.08% |              0.99

Find a predictor’s slopes at each factor level

estimate_slopes(model)
## Species    | Median |       95% CI |     pd | % in ROPE | Median (std.)
## -----------------------------------------------------------------------
## setosa     |   0.36 | [0.01, 0.75] | 97.10% |     8.00% |          1.44
## versicolor |   0.36 | [0.19, 0.56] |   100% |     0.22% |          1.46
## virginica  |   0.23 | [0.08, 0.39] | 99.78% |     5.10% |          0.94

Generate predictions from your model to compare it with original data

Check-out this vignette to create this plot:

estimate_response(model)

Estimate the link between the response and a predictor

See this vignette to create this plot:

model <- stan_glm(Sepal.Width ~ poly(Petal.Length, 2), data=iris)

estimate_link(model)

Describe the smooth term by its linear parts

estimate_smooth(model)
## Part | Start |  End |   Size | Trend | Linearity
## ------------------------------------------------
## 1    |  1.00 | 4.08 | 52.50% | -0.01 |      0.94
## 2    |  4.08 | 6.90 | 47.50% |  0.01 |      0.93