Go from fitted models to publication-ready tables without
hand-formatting effect estimates. gtregression supports
logistic, log-binomial, Poisson, robust Poisson, negative binomial, Cox
survival, parametric survival, and linear regression.
library(gtregression)
library(dplyr)
data("data_birthwt", package = "gtregression")
birthwt_data <- data_birthwt |>
mutate(
race = factor(race, levels = c(1, 2, 3),
labels = c("White", "Black", "Other")),
smoke = factor(smoke, levels = c(0, 1), labels = c("No", "Yes")),
ht = factor(ht, levels = c(0, 1), labels = c("No", "Yes")),
ui = factor(ui, levels = c(0, 1), labels = c("No", "Yes")),
low = factor(low, levels = c(0, 1), labels = c("Normal BW", "Low BW")),
ptl_cat = factor(ifelse(ptl > 0, "Yes", "No"), levels = c("No", "Yes")),
ftv_cat = factor(case_when(
ftv == 0 ~ "None",
ftv == 1 ~ "One",
ftv >= 2 ~ "Two or more"
), levels = c("None", "One", "Two or more"))
)
birthwt_exposures <- c(
"age", "lwt", "race", "smoke", "ht", "ui", "ptl_cat", "ftv_cat"
)
attr(birthwt_data$age, "label") <- "Maternal age"
attr(birthwt_data$lwt, "label") <- "Maternal weight"
attr(birthwt_data$race, "label") <- "Maternal race"
attr(birthwt_data$smoke, "label") <- "Smoking during pregnancy"
attr(birthwt_data$ht, "label") <- "Hypertension"
attr(birthwt_data$ui, "label") <- "Uterine irritability"
attr(birthwt_data$ptl_cat, "label") <- "Previous preterm labour"
attr(birthwt_data$ftv_cat, "label") <- "First trimester visits"uni_reg() fits one model per exposure and returns a
table ready for reports. Variable labels set with
attr(x, "label") or labelled::var_label() are
used automatically in the displayed table, while raw column names remain
in $table_body.
birthwt_uni <- uni_reg(
data = birthwt_data,
outcome = "low",
exposures = birthwt_exposures,
approach = "logit",
theme = clinical
)
birthwt_uni## a flextable object.
## col_keys: `Characteristic`, `N`, `OR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 22 row(s)
## original dataset sample:
## 'data.frame': 22 obs. of 4 variables:
## $ Characteristic: chr "Maternal age" "Maternal weight" "Maternal race" "White" ...
## $ N : chr "189" "189" "189" "" ...
## $ OR (95% CI) : chr "0.95 (0.89-1.01)" "0.99 (0.97-1.00)" "" "Ref." ...
## $ p-value : chr "0.105" "0.023" "" "" ...With adjust_for = NULL, multi_reg() fits
all supplied exposures in one multivariable model. This is the usual
fully adjusted model when every exposure listed should appear in the
same formula.
birthwt_full <- multi_reg(
data = birthwt_data,
outcome = low,
exposures = c("age", "lwt", "race", "smoke", "ht", "ui"),
approach = logit,
theme = clinical
)
birthwt_full## a flextable object.
## col_keys: `Characteristic`, `Adjusted OR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 15 row(s)
## original dataset sample:
## 'data.frame': 15 obs. of 3 variables:
## $ Characteristic : chr "Maternal age" "Maternal weight" "Maternal race" "White" ...
## $ Adjusted OR (95% CI): chr "0.98 (0.92–1.05)" "0.98 (0.97–1.00)" "" "Ref." ...
## $ p-value : chr "0.606" "0.018" "" "" ...Use adjust_for when you want one adjusted model per
exposure, each adjusted for the same core covariate set. This is useful
for screening several clinically important exposures while keeping the
adjustment strategy explicit.
birthwt_multi <- multi_reg(
data = birthwt_data,
outcome = low,
exposures = c("smoke", "ht", "ui", "ptl_cat", "ftv_cat"),
adjust_for = c("age", "lwt", "race"),
approach = logit,
theme = striped
)
birthwt_multi## a flextable object.
## col_keys: `Characteristic`, `Adjusted OR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 16 row(s)
## original dataset sample:
## 'data.frame': 16 obs. of 3 variables:
## $ Characteristic : chr "Smoking during pregnancy" "No" " Yes" "Hypertension" ...
## $ Adjusted OR (95% CI): chr "" "Ref." "2.87 (1.36–6.04)" "" ...
## $ p-value : chr "" "" "0.006" "" ...The adjustment variables are recorded in a compact table footnote, so the result is ready for manuscript-style reporting without making the table unnecessarily tall.
Publication tables should stay readable. When you need model-fit
information, set model_stats = TRUE and inspect the
returned object’s $model_stats element. This keeps AIC,
BIC, log-likelihood, deviance, pseudo R-squared, and linear-model
R-squared values available without adding clutter to the main table.
birthwt_uni_stats <- uni_reg(
data = birthwt_data,
outcome = low,
exposures = birthwt_exposures,
approach = logit,
model_stats = TRUE
)
birthwt_uni_stats$model_stats## model AIC BIC logLik deviance null_deviance pseudo_r2
## 1 age 235.9120 242.3955 -115.9560 231.9120 234.672 0.01176126
## 2 lwt 232.6907 239.1742 -114.3453 228.6907 234.672 0.02548803
## 3 race 235.6616 245.3869 -114.8308 229.6616 234.672 0.02135051
## 4 smoke 233.8046 240.2881 -114.9023 229.8046 234.672 0.02074128
## 5 ht 234.6499 241.1334 -115.3249 230.6499 234.672 0.01713938
## 6 ui 233.5959 240.0794 -114.7979 229.5959 234.672 0.02163060
## 7 ptl_cat 225.8978 232.3812 -110.9489 221.8978 234.672 0.05443445
## 8 ftv_cat 238.0851 247.8103 -116.0425 232.0851 234.672 0.01102346
## r_squared adj_r_squared n
## 1 NA NA 189
## 2 NA NA 189
## 3 NA NA 189
## 4 NA NA 189
## 5 NA NA 189
## 6 NA NA 189
## 7 NA NA 189
## 8 NA NA 189For adjusted mode, multi_reg() returns one row of
statistics per adjusted exposure-specific model.
birthwt_multi_stats <- multi_reg(
data = birthwt_data,
outcome = low,
exposures = c("smoke", "ht", "ui", "ptl_cat"),
adjust_for = c("age", "lwt", "race"),
approach = logit,
model_stats = TRUE
)
birthwt_multi_stats$model_stats## model AIC BIC logLik deviance null_deviance pseudo_r2
## 1 smoke 226.5772 246.0277 -107.2886 214.5772 234.672 0.08562914
## 2 ht 227.7487 247.1991 -107.8743 215.7487 234.672 0.08063739
## 3 ui 231.0194 250.4699 -109.5097 219.0194 234.672 0.06669985
## 4 ptl_cat 222.4421 241.8926 -105.2211 210.4421 234.672 0.10324998
## r_squared adj_r_squared n
## 1 NA NA 189
## 2 NA NA 189
## 3 NA NA 189
## 4 NA NA 189Switch the approach to change the estimand.
uni_reg(
data = birthwt_data,
outcome = low,
exposures = c("smoke", "ht", "ui", "ptl_cat"),
approach = logbinomial
)## a flextable object.
## col_keys: `Characteristic`, `N`, `RR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 12 row(s)
## original dataset sample:
## 'data.frame': 12 obs. of 4 variables:
## $ Characteristic: chr "Smoking during pregnancy" "No" " Yes" "Hypertension" ...
## $ N : chr "189" "" "" "189" ...
## $ RR (95% CI) : chr "" "Ref." "1.61 (1.06-2.44)" "" ...
## $ p-value : chr "" "" "0.026" "" ...Use approach = firth when a binary-outcome logistic
model has sparse cells, very wide intervals, or separation concerns. The
output remains an odds-ratio table, but the model is fitted with Firth
penalized logistic regression.
The built-in data_endometrial dataset is a useful
teaching example because neovascularization is completely absent among
low-grade cases, a pattern that can make ordinary logistic regression
unstable.
data("data_endometrial", package = "gtregression")
endometrial_data <- data_endometrial |>
mutate(
HG = factor(HG, levels = c(0, 1),
labels = c("Low grade", "High grade")),
NV = factor(NV, levels = c(0, 1), labels = c("Absent", "Present"))
)
multi_reg(
data = endometrial_data,
outcome = HG,
exposures = c(NV, PI, EH),
approach = firth
)## a flextable object.
## col_keys: `Characteristic`, `Adjusted OR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 5 row(s)
## original dataset sample:
## 'data.frame': 5 obs. of 3 variables:
## $ Characteristic : chr "NV" "Absent" " Present" "PI" ...
## $ Adjusted OR (95% CI): chr "" "Ref." "18.71 (1.84–2,577.65)" "0.97 (0.88–1.04)" ...
## $ p-value : chr "" "" "0.009" "0.387" ...cox_reg() uses direct time and
event arguments and returns hazard ratios. Without
adjust_for, the table shows crude HRs. With
adjust_for, the table shows adjusted HRs.
data("data_lungcancer", package = "gtregression")
lung_data <- data_lungcancer |>
mutate(
trt = factor(trt, levels = c(1, 2),
labels = c("Standard treatment", "Test treatment")),
prior = factor(prior, levels = c(0, 10), labels = c("No", "Yes"))
)
attr(lung_data$trt, "label") <- "Treatment group"
attr(lung_data$celltype, "label") <- "Cancer cell type"
attr(lung_data$karno, "label") <- "Karnofsky performance score"
attr(lung_data$age, "label") <- "Age"
attr(lung_data$prior, "label") <- "Prior therapy"lung_hr <- cox_reg(
data = lung_data,
time = time,
event = status,
exposures = c("trt", "celltype", "karno", "age"),
theme = clinical
)
lung_hr## a flextable object.
## col_keys: `Characteristic`, `N`, `HR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 10 row(s)
## original dataset sample:
## 'data.frame': 10 obs. of 4 variables:
## $ Characteristic: chr "Treatment group" "Standard treatment" " Test treatment" "Cancer cell type" ...
## $ N : chr "137" "" "" "137" ...
## $ HR (95% CI) : chr "" "Ref." "1.02 (0.71–1.45)" "" ...
## $ p-value : chr "" "" "0.922" "" ...lung_adj_hr <- cox_reg(
data = lung_data,
time = time,
event = status,
exposures = c(trt, celltype, prior),
adjust_for = c(age, karno),
model_stats = TRUE,
theme = striped
)
lung_adj_hr## a flextable object.
## col_keys: `Characteristic`, `Adjusted HR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 11 row(s)
## original dataset sample:
## 'data.frame': 11 obs. of 3 variables:
## $ Characteristic : chr "Treatment group" "Standard treatment" " Test treatment" "Cancer cell type" ...
## $ Adjusted HR (95% CI): chr "" "Ref." "1.21 (0.84–1.74)" "" ...
## $ p-value : chr "" "" "0.307" "" ...## model AIC BIC logLik concordance events n
## 1 trt 973.7560 982.3121 -483.8780 0.7119491 128 137
## 2 celltype 961.0882 975.3484 -475.5441 0.7349500 128 137
## 3 prior 974.7459 983.3019 -484.3729 0.7134257 128 137Use interaction = exposure*modifier for planned effect
modification. In the default exposure-by-exposure workflow, provide one
exposure.
cox_interaction <- cox_reg(
data = lung_data,
time = time,
event = status,
exposures = trt,
adjust_for = c(age, karno),
interaction = trt*prior
)
cox_interaction## a flextable object.
## col_keys: `Characteristic`, `Adjusted HR (95% CI)`, `p-value`
## header has 1 row(s)
## body has 4 row(s)
## original dataset sample:
## 'data.frame': 4 obs. of 3 variables:
## $ Characteristic : chr "Treatment group" "Standard treatment" " Test treatment" " Treatment group: Test treatment x Prior therapy: Yes"
## $ Adjusted HR (95% CI): chr "" "Ref." "1.52 (0.98–2.36)" "0.47 (0.21–1.06)"
## $ p-value : chr "" "" "0.061" "0.070"surv_reg() uses the same time,
event, exposures, and adjust_for
grammar as cox_reg(), but fits parametric survival models
with survival::survreg(). The table reports time ratios
rather than hazard ratios. A time ratio above 1 suggests longer survival
time; below 1 suggests shorter survival time, conditional on the chosen
distribution.
lung_time_ratio <- surv_reg(
data = lung_data,
time = time,
event = status,
exposures = c("trt", "celltype", "karno", "age"),
distribution = weibull,
theme = clinical
)
lung_time_ratio## a flextable object.
## col_keys: `Characteristic`, `N`, `Time Ratio (95% CI)`, `p-value`
## header has 1 row(s)
## body has 10 row(s)
## original dataset sample:
## 'data.frame': 10 obs. of 4 variables:
## $ Characteristic : chr "Treatment group" "Standard treatment" " Test treatment" "Cancer cell type" ...
## $ N : chr "137" "" "" "137" ...
## $ Time Ratio (95% CI): chr "" "Ref." "1.05 (0.70–1.58)" "" ...
## $ p-value : chr "" "" "0.818" "" ...lung_adj_time_ratio <- surv_reg(
data = lung_data,
time = time,
event = status,
exposures = c(trt, celltype, prior),
adjust_for = c(age, karno),
distribution = lognormal,
model_stats = TRUE,
theme = striped
)
lung_adj_time_ratio## a flextable object.
## col_keys: `Characteristic`, `Adjusted Time Ratio (95% CI)`, `p-value`
## header has 1 row(s)
## body has 11 row(s)
## original dataset sample:
## 'data.frame': 11 obs. of 3 variables:
## $ Characteristic : chr "Treatment group" "Standard treatment" " Test treatment" "Cancer cell type" ...
## $ Adjusted Time Ratio (95% CI): chr "" "Ref." "0.87 (0.60–1.27)" "" ...
## $ p-value : chr "" "" "0.468" "" ...## model distribution AIC BIC logLik scale events n
## 1 trt lognormal 1451.264 1465.864 -720.6318 1.110644 128 137
## 2 celltype lognormal 1444.379 1464.819 -715.1894 1.064351 128 137
## 3 prior lognormal 1451.768 1466.368 -720.8840 1.112741 128 137surv_interaction <- surv_reg(
data = lung_data,
time = time,
event = status,
exposures = trt,
adjust_for = c(age, karno),
interaction = trt*prior,
distribution = weibull
)
surv_interaction## a flextable object.
## col_keys: `Characteristic`, `Adjusted Time Ratio (95% CI)`, `p-value`
## header has 1 row(s)
## body has 4 row(s)
## original dataset sample:
## 'data.frame': 4 obs. of 3 variables:
## $ Characteristic : chr "Treatment group" "Standard treatment" " Test treatment" " Treatment group: Test treatment x Prior therapy: Yes"
## $ Adjusted Time Ratio (95% CI): chr "" "Ref." "0.68 (0.44–1.04)" "2.29 (1.03–5.09)"
## $ p-value : chr "" "" "0.073" "0.042"Linear regression outputs beta coefficients and keeps diagnostics
under $reg_check.
birthwt_linear <- multi_reg(
data = birthwt_data,
outcome = bwt,
exposures = c("age", "lwt", "race", "smoke", "ht", "ui"),
approach = linear
)
birthwt_linear## a flextable object.
## col_keys: `Characteristic`, `Adjusted Beta (95% CI)`, `p-value`
## header has 1 row(s)
## body has 15 row(s)
## original dataset sample:
## 'data.frame': 15 obs. of 3 variables:
## $ Characteristic : chr "Maternal age" "Maternal weight" "Maternal race" "White" ...
## $ Adjusted Beta (95% CI): chr "-4.67 (-22.97–13.63)" "4.40 (1.05–7.74)" "" "Ref." ...
## $ p-value : chr "0.617" "0.011" "" "" ...$table: publication-ready table.$table_body: numeric estimates behind the display.$models: fitted model objects.$model_summaries: model-level summaries.$model_stats: optional model-fit statistics when
model_stats = TRUE.$variable_labels: labels used in publication
output.$reg_check: diagnostics for linear models.