The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.
The most common function used to view an overall summary of a model
is the summary() function:
library(logitr)
model <- logitr(
data = yogurt,
outcome = "choice",
obsID = "obsID",
pars = c("price", "feat", "brand")
)
summary(model)
#> =================================================
#>
#> Model estimated on: Wed Jul 24 05:47:01 2024
#>
#> Using logitr version: 1.1.2
#>
#> Call:
#> logitr(data = yogurt, outcome = "choice", obsID = "obsID", pars = c("price",
#> "feat", "brand"))
#>
#> Frequencies of alternatives:
#> 1 2 3 4
#> 0.402156 0.029436 0.229270 0.339138
#>
#> Exit Status: 3, Optimization stopped because ftol_rel or ftol_abs was reached.
#>
#> Model Type: Multinomial Logit
#> Model Space: Preference
#> Model Run: 1 of 1
#> Iterations: 21
#> Elapsed Time: 0h:0m:0.01s
#> Algorithm: NLOPT_LD_LBFGS
#> Weights Used?: FALSE
#> Robust? FALSE
#>
#> Model Coefficients:
#> Estimate Std. Error z-value Pr(>|z|)
#> price -0.366555 0.024365 -15.0441 < 2.2e-16 ***
#> feat 0.491439 0.120062 4.0932 4.254e-05 ***
#> brandhiland -3.715477 0.145417 -25.5506 < 2.2e-16 ***
#> brandweight -0.641138 0.054498 -11.7645 < 2.2e-16 ***
#> brandyoplait 0.734519 0.080642 9.1084 < 2.2e-16 ***
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#>
#> Log-Likelihood: -2656.8878790
#> Null Log-Likelihood: -3343.7419990
#> AIC: 5323.7757580
#> BIC: 5352.7168000
#> McFadden R2: 0.2054148
#> Adj McFadden R2: 0.2039195
#> Number of Observations: 2412.0000000The summary prints out a table of the model coefficients as well as other information about the model, such as the log-likelihood, the number of observations, etc.
You can extract the coefficient from the summary table as a
data.frame using coef(summary(model)):
coefs <- coef(summary(model))
coefs
#> Estimate Std. Error z-value Pr(>|z|)
#> price -0.3665546 0.02436526 -15.044150 0.000000e+00
#> feat 0.4914392 0.12006175 4.093221 4.254226e-05
#> brandhiland -3.7154773 0.14541671 -25.550553 0.000000e+00
#> brandweight -0.6411384 0.05449794 -11.764450 0.000000e+00
#> brandyoplait 0.7345195 0.08064229 9.108366 0.000000e+00Another approach for extracting the model coefficients as a data
frame is to use the tidy() function from the {broom} package:
library(broom)
coefs <- tidy(model)
coefs
#> # A tibble: 5 × 5
#> term estimate std.error statistic p.value
#> <chr> <dbl> <dbl> <dbl> <dbl>
#> 1 price -0.367 0.0244 -15.0 0
#> 2 feat 0.491 0.120 4.09 0.0000425
#> 3 brandhiland -3.72 0.145 -25.6 0
#> 4 brandweight -0.641 0.0545 -11.8 0
#> 5 brandyoplait 0.735 0.0806 9.11 0The tidy() function returns a tibble and
provides a more standardized output and interfaces well with other
packages. You can also append a confidence interval to the data
frame:
coefs <- tidy(model, conf.int = TRUE, conf.level = 0.95)
coefs
#> # A tibble: 5 × 7
#> term estimate std.error statistic p.value conf.low conf.high
#> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 brandhiland -3.72 0.145 -25.6 0 -4.01 -3.43
#> 2 brandweight -0.641 0.0545 -11.8 0 -0.747 -0.535
#> 3 brandyoplait 0.735 0.0806 9.11 0 0.576 0.890
#> 4 feat 0.491 0.120 4.09 0.0000425 0.256 0.728
#> 5 price -0.367 0.0244 -15.0 0 -0.413 -0.318You can also extract other values of interest at the solution, such as:
The estimated coefficients
coef(model)
#> price feat brandhiland brandweight brandyoplait
#> -0.3665546 0.4914392 -3.7154773 -0.6411384 0.7345195The estimated standard errors
se(model)
#> price feat brandhiland brandweight brandyoplait
#> 0.02436526 0.12006175 0.14541671 0.05449794 0.08064229The log-likelihood
The variance-covariance matrix
vcov(model)
#> price feat brandhiland brandweight
#> price 0.0005936657 5.729619e-04 0.001851795 1.249988e-04
#> feat 0.0005729619 1.441482e-02 0.000855011 5.092011e-06
#> brandhiland 0.0018517954 8.550110e-04 0.021146019 1.490080e-03
#> brandweight 0.0001249988 5.092012e-06 0.001490080 2.970026e-03
#> brandyoplait -0.0015377721 -1.821331e-03 -0.003681036 7.779428e-04
#> brandyoplait
#> price -0.0015377721
#> feat -0.0018213311
#> brandhiland -0.0036810363
#> brandweight 0.0007779427
#> brandyoplait 0.0065031782You can also view a summary of statistics about the model using the
glance() function from the {broom} package:
Often times you will need to create summary tables that are formatted
for publication. The {gtsummary} package
offers a convenient solution that works well with
logitrmodels. For example, a formatted summary table can be
obtained using the tbl_regression() function:
| Characteristic | Beta | 95% CI1 | p-value |
|---|---|---|---|
| brand | |||
| dannon | — | — | |
| hiland | -3.7 | -4.0, -3.4 | <0.001 |
| weight | -0.64 | -0.75, -0.53 | <0.001 |
| yoplait | 0.73 | 0.58, 0.90 | <0.001 |
| feat | 0.49 | 0.26, 0.72 | <0.001 |
| price | -0.37 | -0.42, -0.32 | <0.001 |
| 1 CI = Confidence Interval | |||
The tbl_regression() function has many options
for customizing the output table. For example, you can change the
coefficient names with the label argument:
| Characteristic | Beta | 95% CI1 | p-value |
|---|---|---|---|
| Yogurt's brand | |||
| dannon | — | — | |
| hiland | -3.7 | -4.0, -3.4 | <0.001 |
| weight | -0.64 | -0.75, -0.53 | <0.001 |
| yoplait | 0.73 | 0.58, 0.89 | <0.001 |
| Newspaper ad shown? | 0.49 | 0.25, 0.73 | <0.001 |
| price | -0.37 | -0.41, -0.32 | <0.001 |
| 1 CI = Confidence Interval | |||
The {gtsummary} package supports a wide variety of output types,
including support for LaTeX. One you create the
table with tbl_regression(), you can print it a variety of
ways. For example, once you’ve created a table x,
you can print it to LaTeX with any of the following ways:
as_gt(x) |> gt::as_latex()as_kable_extra(x, format = "latex")as_hux_table(x) |> to_latex()as_kable(x, format = "latex")Multiple models can also be printed in the same table:
model1 <- model
model2 <- logitr(
data = yogurt,
outcome = "choice",
obsID = "obsID",
pars = c("price*feat", "brand")
)
# Make individual tables
t1 <- tbl_regression(model1)
t2 <- tbl_regression(model2)
# Merge tables
tbl_merge(
tbls = list(t1, t2),
tab_spanner = c("**Baseline**", "**Interaction**")
)| Characteristic | Baseline | Interaction | ||||
|---|---|---|---|---|---|---|
| Beta | 95% CI1 | p-value | Beta | 95% CI1 | p-value | |
| brand | ||||||
| dannon | — | — | — | — | ||
| hiland | -3.7 | -4.0, -3.4 | <0.001 | -3.7 | -4.0, -3.4 | <0.001 |
| weight | -0.64 | -0.75, -0.53 | <0.001 | -0.64 | -0.75, -0.53 | <0.001 |
| yoplait | 0.73 | 0.58, 0.89 | <0.001 | 0.72 | 0.57, 0.88 | <0.001 |
| feat | 0.49 | 0.26, 0.72 | <0.001 | 1.2 | 0.39, 1.9 | 0.002 |
| price | -0.37 | -0.41, -0.32 | <0.001 | -0.36 | -0.41, -0.31 | <0.001 |
| price * feat | -0.09 | -0.18, 0.01 | 0.068 | |||
| 1 CI = Confidence Interval | ||||||
Another option for obtaining a formatted table is to use the {texreg} package. This is particularly useful for obtaining tables formatted for use in LaTeX.
For example, you can print a summary to the screen using
screenreg():
library(texreg)
screenreg(model, stars = c(0.01, 0.05, 0.1))
#>
#> ============================
#> Model 1
#> ----------------------------
#> price -0.37 ***
#> (0.02)
#> feat 0.49 ***
#> (0.12)
#> brandhiland -3.72 ***
#> (0.15)
#> brandweight -0.64 ***
#> (0.05)
#> brandyoplait 0.73 ***
#> (0.08)
#> ----------------------------
#> Log Likelihood -2656.89
#> null.logLik -3343.74
#> AIC 5323.78
#> BIC 5352.72
#> R^2 0.21
#> Adj. R^2 0.20
#> nobs 2412.00
#> ============================
#> *** p < 0.01; ** p < 0.05; * p < 0.1Likewise, you can print the LaTeX code for a summary table using
texreg()
library(texreg)
texreg(model, stars = c(0.01, 0.05, 0.1))
#>
#> \begin{table}
#> \begin{center}
#> \begin{tabular}{l c}
#> \hline
#> & Model 1 \\
#> \hline
#> price & $-0.37^{***}$ \\
#> & $(0.02)$ \\
#> feat & $0.49^{***}$ \\
#> & $(0.12)$ \\
#> brandhiland & $-3.72^{***}$ \\
#> & $(0.15)$ \\
#> brandweight & $-0.64^{***}$ \\
#> & $(0.05)$ \\
#> brandyoplait & $0.73^{***}$ \\
#> & $(0.08)$ \\
#> \hline
#> Log Likelihood & $-2656.89$ \\
#> null.logLik & $-3343.74$ \\
#> AIC & $5323.78$ \\
#> BIC & $5352.72$ \\
#> R$^2$ & $0.21$ \\
#> Adj. R$^2$ & $0.20$ \\
#> nobs & $2412.00$ \\
#> \hline
#> \multicolumn{2}{l}{\scriptsize{$^{***}p<0.01$; $^{**}p<0.05$; $^{*}p<0.1$}}
#> \end{tabular}
#> \caption{Statistical models}
#> \label{table:coefficients}
#> \end{center}
#> \end{table}Similar to {gtsummary}, multiple models can be printed using
screenreg() or texreg():
model1 <- model
model2 <- logitr(
data = yogurt,
outcome = "choice",
obsID = "obsID",
pars = c("price*feat", "brand")
)
screenreg(
list(
model1,
model2
),
stars = c(0.01, 0.05, 0.1),
custom.model.names = c("Baseline", "Interaction")
)
#>
#> ==========================================
#> Baseline Interaction
#> ------------------------------------------
#> price -0.37 *** -0.36 ***
#> (0.02) (0.02)
#> feat 0.49 *** 1.16 ***
#> (0.12) (0.38)
#> brandhiland -3.72 *** -3.72 ***
#> (0.15) (0.15)
#> brandweight -0.64 *** -0.64 ***
#> (0.05) (0.05)
#> brandyoplait 0.73 *** 0.72 ***
#> (0.08) (0.08)
#> price:feat -0.09 *
#> (0.05)
#> ------------------------------------------
#> Log Likelihood -2656.89 -2655.54
#> null.logLik -3343.74 -3343.74
#> AIC 5323.78 5323.08
#> BIC 5352.72 5357.81
#> R^2 0.21 0.21
#> Adj. R^2 0.20 0.20
#> nobs 2412.00 2412.00
#> ==========================================
#> *** p < 0.01; ** p < 0.05; * p < 0.1These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.