Package 'summata'

Description

A convenience wrapper function that automatically detects the input type and routes to the appropriate specialized forest plot function. This eliminates the need to remember which forest function to call for different model types or analysis objects, making it ideal for exploratory analysis and rapid prototyping.

Usage

autoforest(x, data = NULL, title = NULL, ...)

Arguments

Details

This function provides a convenient wrapper around the specialized forest plot functions, automatically routing to the appropriate function based on the model class or result type. All parameters are passed through to the underlying function, so the full range of options remains available.

For model-specific advanced features, individual forest functions may be called directly.

Automatic Detection Logic:

The function uses the following priority order for detection:

1. uniscreen results: Detected by class "uniscreen_result" or presence of attributes outcome, predictors, model_type, and model_scope = "Univariable". Routes to uniforest().

2. multifit results: Detected by presence of attributes predictor, outcomes, model_type, and raw_data. Routes to multiforest().

3. Cox models: Classes coxph or clogit. Routes to coxforest().

4. GLM models: Class glm. Routes to glmforest().

5. Linear models: Class lm (but not glm). Routes to lmforest().

Value

A ggplot object containing the complete forest plot. The plot can be:

• Displayed directly: print(plot)

• Saved to file: ggsave("forest.pdf", plot, width = 12, height = 8)

• Further customized with ggplot2 functions

The returned object includes an attribute "rec_dims" accessible via attr(plot, "rec_dims"), which is a list containing:

width

Numeric. Recommended plot width in specified units

height

Numeric. Recommended plot height in specified units

These recommendations are automatically calculated based on the number of variables, text sizes, and layout parameters, and are printed to console if plot_width or plot_height are not specified.

See Also

glmforest for GLM forest plots, coxforest for Cox model forest plots, lmforest for linear model forest plots, uniforest for univariable screening forest plots, multiforest for multi-outcome forest plots, fit for single-model regression, fullfit for combined univariable/multivariable regression, uniscreen for univariable screening, multifit for multi-outcome analysis

Other visualization functions: coxforest(), glmforest(), lmforest(), multiforest(), uniforest()

Examples

data(clintrial)
data(clintrial_labels)
library(survival)

# Create example model
glm_model <- glm(surgery ~ age + sex + bmi + smoking,
                 family = binomial, data = clintrial)

# Example 1: Logistic regression model
p <- autoforest(glm_model, data = clintrial)
# Automatically detects GLM and routes to glmforest()



# Example 2: Cox proportional hazards model
cox_model <- coxph(Surv(os_months, os_status) ~ age + sex + treatment + stage,
                   data = clintrial)

plot2 <- autoforest(cox_model, data = clintrial)
# Automatically detects coxph and routes to coxforest()

# Example 3: Linear regression model
lm_model <- lm(biomarker_x ~ age + sex + bmi + treatment, data = clintrial)

plot3 <- autoforest(lm_model, data = clintrial)
# Automatically detects lm and routes to lmforest()

# Example 4: With custom labels and formatting options
plot4 <- autoforest(
    cox_model,
    data = clintrial,
    labels = clintrial_labels,
    title = "Prognostic Factors for Overall Survival",
    zebra_stripes = TRUE,
    indent_groups = TRUE
)

# Example 5: From fit() result - data and labels extracted automatically
fit_result <- fit(
    data = clintrial,
    outcome = "surgery",
    predictors = c("age", "sex", "bmi", "treatment"),
    labels = clintrial_labels
)

plot5 <- autoforest(fit_result)
# No need to pass data or labels - extracted from fit_result

# Save with recommended dimensions
dims <- attr(plot5, "rec_dims")
ggplot2::ggsave(file.path(tempdir(), "forest.pdf"),
                plot5, width = dims$width, height = dims$height)

Description

Automatically detects the output format based on file extension and exports the table using the appropriate specialized function. Provides a unified interface for table export across all supported formats.

Usage

autotable(table, file, ...)

Arguments

Details

This function provides a convenient wrapper around format-specific export functions, automatically routing to the appropriate function based on the file extension. All parameters are passed through to the underlying function, so the full range of format-specific options remains available.

For format-specific advanced features, you may prefer to use the individual export functions directly:

• PDF exports support orientation, paper size, margins, and auto-sizing

• DOCX/PPTX/RTF support font customization and flextable formatting

• HTML supports CSS styling, responsive design, and custom themes

• TeX generates standalone LaTeX source with booktabs styling

Value

Invisibly returns the file path. Called primarily for its side effect of creating the output file.

See Also

table2pdf, table2docx, table2pptx, table2html, table2rtf, table2tex

Other export functions: table2docx(), table2html(), table2pdf(), table2pptx(), table2rtf(), table2tex()

Examples

# Create example data
data(clintrial)
data(clintrial_labels)
tbl <- desctable(clintrial, by = "treatment",
    variables = c("age", "sex"), labels = clintrial_labels)

# Auto-detect format from extension
if (requireNamespace("xtable", quietly = TRUE)) {
  autotable(tbl, file.path(tempdir(), "example.html"))
}


# Load example data
data(clintrial)
data(clintrial_labels)

# Create a regression table
results <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment"),
    labels = clintrial_labels
)

# Test that LaTeX can actually compile (needed for PDF export)
has_latex <- local({
  if (!nzchar(Sys.which("pdflatex"))) return(FALSE)
  test_tex <- file.path(tempdir(), "summata_latex_test.tex")
  writeLines(c("\\documentclass{article}",
               "\\usepackage{booktabs}",
               "\\begin{document}", "test",
               "\\end{document}"), test_tex)
  result <- tryCatch(
    system2("pdflatex", c("-interaction=nonstopmode",
            paste0("-output-directory=", tempdir()), test_tex),
            stdout = FALSE, stderr = FALSE),
    error = function(e) 1L)
  result == 0L
})

# Export automatically detects format from extension
autotable(results, file.path(tempdir(), "results.html"))  # Creates HTML file
autotable(results, file.path(tempdir(), "results.docx"))  # Creates Word document
autotable(results, file.path(tempdir(), "results.pptx"))  # Creates PowerPoint slide
autotable(results, file.path(tempdir(), "results.tex"))   # Creates LaTeX source
autotable(results, file.path(tempdir(), "results.rtf"))   # Creates RTF document
if (has_latex) {
  autotable(results, file.path(tempdir(), "results.pdf")) # Creates PDF
}

# Pass format-specific parameters
if (has_latex) {
  autotable(results, file.path(tempdir(), "results.pdf"), 
             orientation = "landscape",
             paper = "a4",
             font_size = 10)
}

autotable(results, file.path(tempdir(), "results.docx"),
           caption = "Table 1: Logistic Regression Results",
           font_family = "Times New Roman",
           condense_table = TRUE)

autotable(results, file.path(tempdir(), "results.html"),
           zebra_stripes = TRUE,
           dark_header = TRUE,
           bold_significant = TRUE)

# Works with any summata table output
desc <- desctable(clintrial,
                  by = "treatment",
                  variables = c("age", "sex", "bmi"))
if (has_latex) {
  autotable(desc, file.path(tempdir(), "demographics.pdf"))
}

comparison <- compfit(
    data = clintrial,
    outcome = "os_status",
    model_list = list(
        base = c("age", "sex"),
        full = c("age", "sex", "treatment", "stage")
    )
)
autotable(comparison, file.path(tempdir(), "model_comparison.docx"))

Description

A simulated dataset from a hypothetical multi-center oncology clinical trial comparing two experimental drugs against control. Designed to demonstrate the full capabilities of descriptive and regression analysis functions.

Usage

clintrial

Format

A data frame with 850 observations and 32 variables:

patient_id

Unique patient identifier (character)

age

Age at enrollment in years (numeric: 18-90)

sex

Biological sex (factor: Female, Male)

race

Self-reported race (factor: White, Black, Asian, Other)

ethnicity

Hispanic ethnicity (factor: Non-Hispanic, Hispanic)

bmi

Body mass index in kg/m$^2$ (numeric)

smoking

Smoking history (factor: Never, Former, Current)

hypertension

Hypertension diagnosis (factor: No, Yes)

diabetes

Diabetes diagnosis (factor: No, Yes)

ecog

ECOG performance status (factor: 0, 1, 2, 3)

creatinine

Baseline creatinine in mg/dL (numeric)

hemoglobin

Baseline hemoglobin in g/dL (numeric)

biomarker_x

Serum biomarker A in ng/mL (numeric)

biomarker_y

Serum biomarker B in U/L (numeric)

site

Enrolling site (factor: Site Alpha through Site Kappa)

grade

Tumor grade (factor: Well/Moderately/Poorly differentiated)

stage

Disease stage at diagnosis (factor: I, II, III, IV)

treatment

Randomized treatment (factor: Control, Drug A, Drug B)

surgery

Surgical resection (factor: No, Yes)

any_complication

Any post-operative complication (factor: No, Yes)

wound_infection

Post-operative wound infection (factor: No, Yes)

icu_admission

ICU admission required (factor: No, Yes)

readmission_30d

Hospital readmission within 30 days (factor: No, Yes)

pain_score

Pain score at discharge (numeric: 0-10)

recovery_days

Days to functional recovery (numeric)

los_days

Hospital length of stay in days (numeric)

ae_count

Adverse event count (integer). Overdispersed count suitable for negative binomial or quasipoisson regression.

fu_count

Follow-up visit count (integer). Equidispersed count suitable for standard Poisson regression.

pfs_months

Progression-Free Survival Time (months)

pfs_status

Progression or Death Event

os_months

Overall survival time in months (numeric)

os_status

Death indicator (numeric: 0=censored, 1=death)

Details

This dataset includes realistic correlations between variables: - Survival is worse with higher stage, ECOG, age, and biomarker_x - Treatment effects show Drug B > Drug A > Control - ae_count is overdispersed (variance > mean) for negative binomial demos - fu_count is equidispersed (variance $\approx$ mean) for Poisson demos - Approximately 2% of values are missing at random - Median follow-up is approximately 30 months

Source

Simulated data for demonstration purposes

See Also

Other sample data: clintrial_labels

Examples

data(clintrial)
data(clintrial_labels)

# Descriptive statistics by treatment arm
desctable(clintrial,
        by = "treatment", 
        variables = c("age", "sex", "stage", "ecog", 
                     "biomarker_x", "Surv(os_months, os_status)"),
        labels = clintrial_labels)


# Poisson regression for equidispersed counts
fit(clintrial,
    outcome = "fu_count",
    predictors = c("age", "stage", "treatment"),
    model_type = "glm",
    family = "poisson",
    labels = clintrial_labels)

# Negative binomial for overdispersed counts
fit(clintrial,
    outcome = "ae_count",
    predictors = c("age", "treatment", "diabetes"),
    model_type = "negbin",
    labels = clintrial_labels)

# Complete analysis pipeline
fullfit(clintrial,
        outcome = "Surv(os_months, os_status)",
        predictors = c("age", "sex", "stage", "grade", "ecog",
                      "smoking", "biomarker_x", "biomarker_y", "treatment"),
        method = "screen",
        p_threshold = 0.20,
        model_type = "coxph",
        labels = clintrial_labels)

Description

A named character vector providing descriptive labels for all variables in the clinical_trial dataset. Use with labels parameter in functions.

Usage

clintrial_labels

Format

Named character vector with 24 elements

See Also

clintrial

Other sample data: clintrial

Description

Fits multiple regression models and provides a comprehensive comparison table with model quality metrics, convergence diagnostics, and selection guidance. Computes a composite score combining multiple quality metrics to facilitate rapid model comparison and selection.

Usage

compfit(
  data,
  outcome,
  model_list,
  model_names = NULL,
  interactions_list = NULL,
  random = NULL,
  model_type = "auto",
  family = "binomial",
  conf_level = 0.95,
  p_digits = 3,
  include_coefficients = FALSE,
  scoring_weights = NULL,
  labels = NULL,
  number_format = NULL,
  verbose = NULL,
  ...
)

Arguments

    options(summata.number_format = "eu")
  

Details

This function fits all specified models and computes comprehensive quality metrics for comparison. It generates a Composite Model Score (CMS) that combines multiple metrics: lower AIC/BIC (information criteria), higher concordance (discrimination), and model convergence status.

For GLMs, McFadden's pseudo-R-squared is calculated as 1 - (logLik/logLik_null). For survival models, the global p-value comes from the log-rank test.

Models that fail to converge are flagged and penalized in the composite score.

Interaction Terms:

When interactions_list is provided, each element specifies the interaction terms for the corresponding model in model_list. This is particularly useful for testing whether adding interactions improves model fit:

• Use NULL for models without interactions

• Specify interactions using colon notation: c("age:treatment", "sex:stage")

• Main effects for all variables in interactions must be in the predictor list

• Common pattern: Compare main effects model vs model with interactions

Scoring weights can be customized based on model type:

• GLM: "convergence", "aic", "concordance", "pseudo_r2", "brier"

• Cox: "convergence", "aic", "concordance", "global_p"

• Linear: "convergence", "aic", "pseudo_r2", "rmse"

Default weights emphasize discrimination (concordance) and model fit (AIC).

The composite score is designed as a tool to quickly rank models by their quality metrics. It should be used alongside traditional model selection criteria rather than as a definitive model selection method.

Value

A data.table with class "compfit_result" containing:

Model

Model name/identifier

CMS

Composite Model Score for model selection (higher is better)

N

Sample size

Events

Number of events (for survival/logistic)

Predictors

Number of predictors

Converged

Whether model converged properly

AIC

Akaike Information Criterion

BIC

Bayesian Information Criterion

R$^2$ / Pseudo-R$^2$

McFadden pseudo-R-squared (GLM)

Concordance

C-statistic (logistic/survival)

Brier Score

Brier accuracy score (logistic)

Global p

Overall model p-value

Attributes include:

models

List of fitted model objects

coefficients

Coefficient comparison table (if requested)

best_model

Name of recommended model

See Also

fit for individual model fitting, fullfit for automated variable selection, table2pdf for exporting results

Other regression functions: fit(), fullfit(), multifit(), print.compfit_result(), print.fit_result(), print.fullfit_result(), print.multifit_result(), print.uniscreen_result(), uniscreen()

Examples

# Load example data
data(clintrial)
data(clintrial_labels)

# Example 1: Compare nested logistic regression models
models <- list(
    base = c("age", "sex"),
    clinical = c("age", "sex", "smoking", "diabetes"),
    full = c("age", "sex", "smoking", "diabetes", "stage", "ecog")
)

comparison <- compfit(
    data = clintrial,
    outcome = "os_status",
    model_list = models,
    model_names = c("Base", "Clinical", "Full")
)
comparison



# Example 2: Compare Cox survival models
library(survival)
surv_models <- list(
    simple = c("age", "sex"),
    clinical = c("age", "sex", "stage", "grade")
)

surv_comparison <- compfit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    model_list = surv_models,
    model_type = "coxph"
)
surv_comparison

# Example 3: Test effect of adding interaction terms
interaction_models <- list(
    main = c("age", "treatment", "sex"),
    interact = c("age", "treatment", "sex")
)

interaction_comp <- compfit(
    data = clintrial,
    outcome = "os_status",
    model_list = interaction_models,
    model_names = c("Main Effects", "With Interaction"),
    interactions_list = list(
        NULL,
        c("treatment:sex")
    )
)
interaction_comp

# Example 4: Include coefficient comparison table
detailed <- compfit(
    data = clintrial,
    outcome = "os_status",
    model_list = models,
    include_coefficients = TRUE,
    labels = clintrial_labels
)

# Access coefficient table
coef_table <- attr(detailed, "coefficients")
coef_table

# Example 5: Access fitted model objects
fitted_models <- attr(comparison, "models")
names(fitted_models)

# Example 6: Get best model recommendation
best <- attr(comparison, "best_model")
cat("Recommended model:", best, "\n")

Description

Generates a publication-ready forest plot that combines a formatted data table with a graphical representation of hazard ratios from a Cox proportional hazards survival model. The plot integrates variable names, group levels, sample sizes, event counts, hazard ratios with confidence intervals, p-values, and model diagnostics in a single comprehensive visualization designed for manuscripts and presentations.

Usage

coxforest(
  x,
  data = NULL,
  title = "Cox Proportional Hazards Model",
  effect_label = "Hazard Ratio",
  digits = 2,
  p_digits = 3,
  conf_level = 0.95,
  font_size = 1,
  annot_size = 3.88,
  header_size = 5.82,
  title_size = 23.28,
  plot_width = NULL,
  plot_height = NULL,
  table_width = 0.6,
  show_n = TRUE,
  show_events = TRUE,
  indent_groups = FALSE,
  condense_table = FALSE,
  bold_variables = FALSE,
  center_padding = 4,
  zebra_stripes = TRUE,
  ref_label = "reference",
  labels = NULL,
  color = "#8A61D8",
  qc_footer = TRUE,
  units = "in",
  number_format = NULL
)

Arguments

    options(summata.number_format = "eu")
  

Details

Survival-Specific Features:

The Cox forest plot includes several survival analysis-specific components:

• Event counts: Number of events (deaths, failures) shown for each predictor category, critical for assessing statistical power

• Hazard ratios: Always exponentiated coefficients (never raw), interpreted as the multiplicative change in hazard

• Log scale: Forest plot uses log scale for HR (reference line at 1)

• Model diagnostics: Includes concordance (C-index), global log-rank test p-value, and AIC

Plot Components:

1. Title: Centered at top

2. Data Table (left): Contains:

   • Variable and Group columns

   • n: Sample sizes by group

   • Events: Event counts by group (critical for survival)

   • aHR (95% CI); p-value: Adjusted hazard ratios with CIs and p-values

3. Forest Plot (right):

   • Point estimates (squares sized by sample size)

   • 95% confidence intervals

   • Reference line at HR = 1

   • Log scale for hazard ratios

4. Model Statistics (footer):

   • Events analyzed (with percentage of total)

   • Global log-rank test p-value

   • Concordance (C-index) with standard error

   • AIC

Interpreting Hazard Ratios:

• HR = 1: No effect on hazard (reference)

• HR > 1: Increased hazard (worse survival)

• HR < 1: Decreased hazard (better survival)

• Example: HR = 2.0 means twice the hazard of the event at any time

Event Counts:

The "Events" column is particularly important in survival analysis:

• Indicates the number of actual events (not censored observations) in each group

• Essential for assessing statistical power

• Categories with very few events may have unreliable HR estimates

• The footer shows total events analyzed and percentage of all events in the original data

Concordance (C-index):

The concordance statistic displayed in the footer indicates discrimination:

• Range: 0.5 to 1.0

• 0.5 = random prediction (coin flip)

• 0.7-0.8 = acceptable discrimination

• > 0.8 = excellent discrimination

• Standard error provided for confidence interval calculation

Global Log-Rank Test:

The global p-value tests the null hypothesis that all coefficients are zero:

• Significant p-value (< 0.05) indicates the model as a whole predicts survival

• Non-significant global test doesn't preclude significant individual predictors

• Based on the score (log-rank) test

Stratification and Clustering:

If the model includes stratification (strata()) or clustering (cluster()):

• Stratified variables are not shown in the forest plot (they don't have HRs)

• Clustering affects standard errors but not point estimates

• Both are handled automatically by the function

Proportional Hazards Assumption:

The forest plot assumes proportional hazards (constant HR over time). Users should verify this assumption using:

• cox.zph(model) for testing

• Stratification for variables violating the assumption

• Time-dependent coefficients if needed

Value

A ggplot object containing the complete forest plot. The plot can be:

• Displayed directly: print(plot)

• Saved to file: ggsave("forest.pdf", plot, width = 12, height = 8)

• Further customized with ggplot2 functions

The returned object includes an attribute "rec_dims" accessible via attr(plot, "rec_dims"), which is a list containing:

width

Numeric. Recommended plot width in specified units

height

Numeric. Recommended plot height in specified units

These recommendations are automatically calculated based on the number of variables, text sizes, and layout parameters, and are printed to console if plot_width or plot_height are not specified.

See Also

autoforest for automatic model detection, glmforest for logistic/GLM forest plots, lmforest for linear model forest plots, uniforest for univariable screening forest plots, multiforest for multi-outcome forest plots, coxph for fitting Cox models, fit for regression modeling

Other visualization functions: autoforest(), glmforest(), lmforest(), multiforest(), uniforest()

Examples

data(clintrial)
data(clintrial_labels)
library(survival)

# Create example model
model1 <- coxph(
    survival::Surv(os_months, os_status) ~ age + sex + treatment,
    data = clintrial)

# Example 1: Basic Cox model forest plot
p <- coxforest(model1, data = clintrial)



old_width <- options(width = 180)

# Example 2: With custom labels and title
plot2 <- coxforest(
    x = model1,
    data = clintrial,
    title = "Prognostic Factors for Overall Survival",
    labels = clintrial_labels
)

# Example 3: Comprehensive model with indented layout
model3 <- coxph(
    Surv(os_months, os_status) ~ age + sex + bmi + smoking + 
        treatment + stage + grade,
    data = clintrial
)

plot3 <- coxforest(
    x = model3,
    data = clintrial,
    labels = clintrial_labels,
    indent_groups = TRUE,
    zebra_stripes = TRUE
)

# Example 4: Condensed layout for many binary predictors
model4 <- coxph(
    Surv(os_months, os_status) ~ age + sex + smoking + 
        hypertension + diabetes + surgery,
    data = clintrial
)

plot4 <- coxforest(
    x = model4,
    data = clintrial,
    condense_table = TRUE,
    labels = clintrial_labels
)

# Example 5: Stratified Cox model
model5 <- coxph(
    Surv(os_months, os_status) ~ age + sex + treatment + strata(site),
    data = clintrial
)

plot5 <- coxforest(
    x = model5,
    data = clintrial,
    title = "Stratified by Study Site",
    labels = clintrial_labels
)

# Example 6: Save with recommended dimensions
dims <- attr(plot5, "rec_dims")
ggplot2::ggsave(file.path(tempdir(), "survival_forest.pdf"),
                plot5, width = dims$width, height = dims$height)

options(old_width)

Description

Generates comprehensive descriptive statistics tables with automatic variable type detection, group comparisons, and appropriate statistical testing. This function is designed to create "Table 1"-style summaries commonly used in clinical and epidemiological research, with full support for continuous, categorical, and time-to-event variables.

Usage

desctable(
  data,
  by = NULL,
  variables,
  stats_continuous = c("median_iqr"),
  stats_categorical = "n_percent",
  digits = 1,
  p_digits = 3,
  conf_level = 0.95,
  p_per_stat = FALSE,
  na_include = FALSE,
  na_label = "Unknown",
  na_percent = FALSE,
  test = TRUE,
  test_continuous = "auto",
  test_categorical = "auto",
  total = TRUE,
  total_label = "Total",
  labels = NULL,
  number_format = NULL,
  ...
)

Arguments

    options(summata.number_format = "eu")
  

Details

Variable Type Detection:

The function automatically detects variable types and applies appropriate summaries:

• Continuous: Numeric variables (integer or double) receive statistics specified in stats_continuous

• Categorical: Character, factor, or logical variables receive frequency counts and percentages

• Time-to-Event: Variables specified as Surv(time, event) display median survival with confidence intervals (level controlled by conf_level)

Statistical Testing:

When test = TRUE and by is specified:

• Continuous with "auto": Parametric tests (t-test, ANOVA) for mean-based statistics; non-parametric tests (Wilcoxon, Kruskal-Wallis) for median-based statistics

• Categorical with "auto": Fisher exact test when any expected cell frequency < 5; $\chi^2$ test otherwise

• Survival: Log-rank test for comparing survival curves

• Range statistics: No p-value computed (ranges are descriptive)

Missing Data Handling:

Missing values are handled differently by variable type:

• Continuous: NAs excluded from calculations; optionally shown as count when na_include = TRUE

• Categorical: NAs can be included as a category when na_include = TRUE. The na_percent parameter controls whether percentages are calculated with or without NAs in the denominator

• Survival: NAs in time or event excluded from analysis

Formatting Conventions:

All numeric output respects the number_format parameter. Separators within ranges and confidence intervals adapt automatically to avoid ambiguity:

• Mean $\pm$ SD: "45.2 \eqn{\pm} 12.3" (US) or "45,2 \eqn{\pm} 12,3" (EU)

• Median [IQR]: "38.0 [28.0-52.0]" (US) or "38,0 [28,0-52,0]" (EU, en-dash separator)

• Range: "18.0-75.0" (positive, US), "-5.0 to 10.0" (when bounds are negative)

• Survival: "24.5 (21.2-28.9)" (US) or "24,5 (21,2-28,9)" (EU)

• Counts $\ge$ 1000: "1,234" (US) or "1.234" (EU)

• p-values: "< 0.001" (US) or "< 0,001" (EU)

Value

A data.table with S3 class "desctable" containing formatted descriptive statistics. The table structure includes:

Variable

Variable name or label (from labels)

Group

For continuous variables: statistic type (e.g., "Mean $\pm$ SD", "Median [IQR]"). For categorical variables: category level. Empty for variable name rows.

Total

Statistics for the total sample (if total = TRUE)

Group columns

Statistics for each group level (when by is specified). Column names match group levels.

p-value

Formatted p-values from statistical tests (when test = TRUE and by is specified)

The first row always shows sample sizes for each column. All numeric output (counts, statistics, p-values) respects the number_format setting for locale-appropriate formatting.

The returned object includes the following attributes accessible via attr():

raw_data

A data.table containing unformatted numeric values suitable for further statistical analysis or custom formatting. Includes additional columns for standard deviations, quartiles, etc.

by_variable

The grouping variable name used (value of by)

variables

The variables analyzed (value of variables)

See Also

survtable for detailed survival summary tables, fit for regression modeling, table2pdf for PDF export, table2docx for Word export, table2html for HTML export

Other descriptive functions: print.survtable(), survtable()

Examples

# Load example clinical trial data
data(clintrial)

# Example 1: Basic descriptive table without grouping
desctable(clintrial,
        variables = c("age", "sex", "bmi"))



# Example 2: Grouped comparison with default tests
desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex", "race", "bmi"))

# Example 3: Customize continuous statistics
desctable(clintrial,
        by = "treatment",
        variables = c("age", "bmi", "creatinine"),
        stats_continuous = c("median_iqr", "range"))

# Example 4: Change categorical display format
desctable(clintrial,
        by = "treatment",
        variables = c("sex", "race", "smoking"),
        stats_categorical = "n")  # Show counts only

# Example 5: Include missing values
desctable(clintrial,
        by = "treatment",
        variables = c("age", "smoking", "hypertension"),
        na_include = TRUE,
        na_label = "Missing")

# Example 6: Disable statistical testing
desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex", "bmi"),
        test = FALSE)

# Example 7: Force specific tests
desctable(clintrial,
        by = "surgery",
        variables = c("age", "sex"),
        test_continuous = "t",      # t-test instead of auto
        test_categorical = "fisher") # Fisher test instead of auto

# Example 8: Adjust decimal places
desctable(clintrial,
        by = "treatment",
        variables = c("age", "bmi"),
        digits = 2,    # 2 decimals for continuous
        p_digits = 4)  # 4 decimals for p-values

# Example 9: Custom variable labels
labels <- c(
    age = "Age (years)",
    sex = "Sex",
    bmi = "Body Mass Index (kg/m\u00b2)",
    treatment = "Treatment Arm"
)

desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex", "bmi"),
        labels = labels)

# Example 10: Position total column last
desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex"),
        total = "last")

# Example 11: Exclude total column
desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex"),
        total = FALSE)

# Example 12: Survival analysis
desctable(clintrial,
        by = "treatment",
        variables = "Surv(os_months, os_status)")

# Example 13: Multiple survival endpoints
desctable(clintrial,
        by = "treatment",
        variables = c(
            "Surv(pfs_months, pfs_status)",
            "Surv(os_months, os_status)"
        ),
        labels = c(
            "Surv(pfs_months, pfs_status)" = "Progression-Free Survival",
            "Surv(os_months, os_status)" = "Overall Survival"
        ))

# Example 14: Mixed variable types
desctable(clintrial,
        by = "treatment",
        variables = c(
            "age", "sex", "race",           # Demographics
            "bmi", "creatinine",            # Labs
            "smoking", "hypertension",      # Risk factors
            "Surv(os_months, os_status)"    # Survival
        ))

# Example 15: Three or more groups
desctable(clintrial,
        by = "stage",  # Assuming stage has 3+ levels
        variables = c("age", "sex", "bmi"))
# Automatically uses ANOVA/Kruskal-Wallis and chi-squared

# Example 16: Access raw unformatted data
result <- desctable(clintrial,
                  by = "treatment",
                  variables = c("age", "bmi"))
raw_data <- attr(result, "raw_data")
print(raw_data)
# Raw data includes unformatted numbers, SDs, quartiles, etc.

# Example 17: Check which grouping variable was used
result <- desctable(clintrial,
                  by = "treatment",
                  variables = c("age", "sex"))
attr(result, "by_variable")  # "treatment"

# Example 18: NA percentage calculation options
# Include NAs in percentage denominator (all sum to 100%)
desctable(clintrial,
        by = "treatment",
        variables = "smoking",
        na_include = TRUE,
        na_percent = TRUE)

# Exclude NAs from denominator (non-missing sum to 100%)
desctable(clintrial,
        by = "treatment",
        variables = "smoking",
        na_include = TRUE,
        na_percent = FALSE)

# Example 19: Passing additional test arguments
# Equal variance t-test
desctable(clintrial,
        by = "sex",
        variables = "age",
        test_continuous = "t",
        var.equal = TRUE)

# Example 20: European number formatting
desctable(clintrial,
        by = "treatment",
        variables = c("age", "sex", "bmi"),
        number_format = "eu")

# Example 21: Complete Table 1 for publication
table1 <- desctable(
    data = clintrial,
    by = "treatment",
    variables = c(
        "age", "sex", "race", "ethnicity", "bmi",
        "smoking", "hypertension", "diabetes",
        "ecog", "creatinine", "hemoglobin",
        "site", "stage", "grade",
        "Surv(os_months, os_status)"
    ),
    labels = clintrial_labels,
    stats_continuous = c("median_iqr", "range"),
    total = TRUE,
    na_include = FALSE
)
print(table1)

Description

Provides a unified interface for fitting various types of regression models with automatic formatting of results for publication. Supports generalized linear models, linear models, survival models, and mixed-effects models with consistent syntax and output formatting. Handles both univariable and multivariable models automatically.

Usage

fit(
  data = NULL,
  outcome = NULL,
  predictors = NULL,
  model = NULL,
  model_type = "glm",
  family = "binomial",
  random = NULL,
  interactions = NULL,
  strata = NULL,
  cluster = NULL,
  weights = NULL,
  conf_level = 0.95,
  reference_rows = TRUE,
  show_n = TRUE,
  show_events = TRUE,
  digits = 2,
  p_digits = 3,
  labels = NULL,
  keep_qc_stats = TRUE,
  exponentiate = NULL,
  conf_method = NULL,
  number_format = NULL,
  verbose = NULL,
  ...
)

Arguments

    options(summata.number_format = "eu")
  

Details

Model Scope Detection:

The function automatically detects whether the model is:

• Univariable: Single predictor (e.g., predictors = "age"). Effect estimates are labeled as unadjusted ("OR", "HR", etc.), representing crude (unadjusted) association

• Multivariable: Multiple predictors (e.g., predictors = c("age", "sex", "treatment")) Effect estimates are labeled as adjusted ("aOR", "aHR", etc.), representing associations adjusted for confounding

Interaction Terms:

Interactions are specified using colon notation and added to the model:

• interactions = c("age:treatment") creates interaction between age and treatment

• Main effects for both variables are automatically included

• Multiple interactions can be specified: c("age:sex", "treatment:stage")

• For interactions between categorical variables, separate terms are created for each combination of levels

Stratification (Cox/Conditional Logistic):

The strata parameter creates separate baseline hazards:

• Allows baseline hazard to vary across strata without estimating stratum effects

• Useful when proportional hazards assumption violated across strata

• Example: strata = "center" for multicenter studies

• Stratification variable is not included as a predictor

Clustering (Cox Models):

The cluster parameter computes robust standard errors:

• Accounts for within-cluster correlation (e.g., multiple observations per patient)

• Uses sandwich variance estimator

• Does not change point estimates, only standard errors and p-values

Weighting:

The weights parameter enables weighted regression:

• For survey data with sampling weights

• Inverse probability weighting for causal inference

• Frequency weights for aggregated data

• Weights should be in a column of data

Mixed-Effects Models (lmer/glmer/coxme):

Mixed effects models handle hierarchical or clustered data:

• Use model_type = "lmer" for continuous/normal outcomes

• Use model_type = "glmer" with appropriate family for GLM outcomes

• Use model_type = "coxme" for survival outcomes with clustering

• Random effects are specified in predictors using lme4 syntax:

  • "(1|site)" - Random intercepts by site

  • "(treatment|site)" - Random slopes for treatment by site

  • "(1 + treatment|site)" - Both random intercepts and slopes

• Include random effects as part of the predictors vector

• Example: predictors = c("age", "treatment", "(1|site)")

Effect Measures by Model Type:

• Logistic (family = "binomial"/"quasibinomial"): Odds ratios (OR/aOR)

• Cox (model_type = "coxph"): Hazard ratios (HR/aHR)

• Poisson/Count (family = "poisson"/"quasipoisson"): Rate ratios (RR/aRR)

• Negative binomial (model_type = "negbin"): Rate ratios (RR/aRR)

• Gamma/Log-link: Ratios (multiplicative effects)

• Linear/Gaussian: Raw coefficient estimates (additive effects)

Confidence Intervals:

Confidence interval computation is tailored to each model class using the best available method:

• GLM and negative binomial: Profile likelihood intervals via MASS::confint.glm(), which invert the profile deviance and account for asymmetry in the likelihood surface. More accurate than the Wald approximation when subgroup sizes are small or estimates are near boundary values. Quasi-likelihood families (quasibinomial, quasipoisson) fall back to Wald intervals because they lack a true likelihood function.

• Linear models: Exact t-distribution intervals via confint.lm(), based on the known sampling distribution under normality.

• Cox proportional hazards: Wald intervals (i.e., coefficient $\pm$ z $\times$ SE), the standard approach in the survival analysis literature.

• Mixed-effects models (lmer, glmer, coxme): Wald intervals. Profile likelihood is available for lme4 models via confint(model, method = "profile") but can be prohibitively slow for complex random-effects structures and is not used by default.

If profile likelihood computation fails for any reason (e.g., non-convergence during profiling), the function falls back silently to Wald intervals.

Value

A data.table with S3 class "fit_result" containing formatted regression results. The table structure includes:

Variable

Character. Predictor name or custom label

Group

Character. For factor variables: category level. For interactions: interaction term. For continuous: typically empty

n

Integer. Total sample size (if show_n = TRUE)

n_group

Integer. Sample size for this factor level

events

Integer. Total number of events (if show_events = TRUE)

events_group

Integer. Events for this factor level

OR/HR/RR/Coefficient or aOR/aHR/aRR/Adj. Coefficient (95% CI)

Character. Formatted effect estimate with confidence interval. Column name depends on model type and scope. Univariable models use: OR, HR, RR, Coefficient. Multivariable models use adjusted notation: aOR, aHR, aRR, Adj. Coefficient

p-value

Character. Formatted p-value from Wald test

The returned object includes the following attributes accessible via attr():

model

The fitted model object (glm, lm, coxph, etc.). Access for diagnostics, predictions, or further analysis

raw_data

data.table. Unformatted numeric results with columns for coefficients, standard errors, confidence bounds, quality statistics, etc.

outcome

Character. The outcome variable name

predictors

Character vector. The predictor variable names

formula_str

Character. The complete model formula as a string

model_scope

Character. "Univariable" (one predictor) or "Multivariable" (multiple predictors)

model_type

Character. The regression model type used

interactions

Character vector (if interactions specified). The interaction terms included

strata

Character (if stratification used). The stratification variable

cluster

Character (if clustering used). The cluster variable

weights

Character (if weighting used). The weights variable

significant

Character vector. Names of predictors with p-value below 0.05, suitable for downstream variable selection workflows

See Also

uniscreen for univariable screening of multiple predictors, fullfit for complete univariable-to-multivariable workflow, compfit for comparing multiple models, m2dt for model-to-table conversion

Other regression functions: compfit(), fullfit(), multifit(), print.compfit_result(), print.fit_result(), print.fullfit_result(), print.multifit_result(), print.uniscreen_result(), uniscreen()

Examples

# Load example data
data(clintrial)
data(clintrial_labels)
library(survival)

# Example 1: Univariable logistic regression
uni_model <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = "age"
)
print(uni_model)
# Labeled as "Univariable OR"



# Example 2: Multivariable logistic regression
multi_model <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "bmi", "treatment"),
    labels = clintrial_labels
)
print(multi_model)

# Example 3: Cox proportional hazards model
cox_model <- fit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    predictors = c("age", "sex", "treatment", "stage"),
    model_type = "coxph",
    labels = clintrial_labels
)
print(cox_model)

# Example 4: Model with interaction terms
interact_model <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "treatment", "sex"),
    interactions = c("age:treatment"),
    labels = clintrial_labels
)
print(interact_model)

# Example 5: Cox model with stratification
strat_model <- fit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    predictors = c("age", "sex", "treatment"),
    model_type = "coxph",
    strata = "site",  # Separate baseline hazards by site
    labels = clintrial_labels
)
print(strat_model)

# Example 6: Cox model with clustering
cluster_model <- fit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    predictors = c("age", "treatment"),
    model_type = "coxph",
    cluster = "site",  # Robust SEs accounting for site clustering
    labels = clintrial_labels
)
print(cluster_model)

# Example 7: Linear regression
linear_model <- fit(
    data = clintrial,
    outcome = "bmi",
    predictors = c("age", "sex", "smoking"),
    model_type = "lm",
    labels = clintrial_labels
)
print(linear_model)

# Example 8: Poisson regression for equidispersed count data
# fu_count has variance ~= mean, appropriate for standard Poisson
poisson_model <- fit(
    data = clintrial,
    outcome = "fu_count",
    predictors = c("age", "stage", "treatment", "surgery"),
    model_type = "glm",
    family = "poisson",
    labels = clintrial_labels
)
print(poisson_model)
# Returns rate ratios (RR/aRR)

# Example 9: Negative binomial regression for overdispersed counts
# ae_count has variance > mean (overdispersed), use negbin or quasipoisson
if (requireNamespace("MASS", quietly = TRUE)) {
    nb_result <- fit(
        data = clintrial,
        outcome = "ae_count",
        predictors = c("age", "treatment", "diabetes", "surgery"),
        model_type = "negbin",
        labels = clintrial_labels
    )
    print(nb_result)
}

# Example 10: Gamma regression for positive continuous outcomes
gamma_model <- fit(
    data = clintrial,
    outcome = "los_days",
    predictors = c("age", "treatment", "surgery"),
    model_type = "glm",
    family = Gamma(link = "log"),
    labels = clintrial_labels
)
print(gamma_model)

# Example 11: Access the underlying fitted model
result <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "bmi")
)

# Get the model object
model_obj <- attr(result, "model")
summary(model_obj)

# Model diagnostics
plot(model_obj)

# Predictions
preds <- predict(model_obj, type = "response")

# Example 12: Access raw numeric data
raw_data <- attr(result, "raw_data")
print(raw_data)
# Contains unformatted coefficients, SEs, CIs, AIC, BIC, etc.

# Example 13: Multiple interactions
complex_model <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment", "bmi"),
    interactions = c("age:treatment", "sex:bmi"),
    labels = clintrial_labels
)
print(complex_model)

# Example 14: Customize output columns
minimal <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment"),
    show_n = FALSE,
    show_events = FALSE,
    reference_rows = FALSE
)
print(minimal)

# Example 15: Different confidence levels
ci90 <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "treatment"),
    conf_level = 0.90  # 90% confidence intervals
)
print(ci90)

# Example 16: Force coefficient display instead of OR
coef_model <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "bmi"),
    exponentiate = FALSE  # Show log odds instead of OR
)
print(coef_model)

# Example 17: Confidence interval method
# Default: profile likelihood CIs for GLM (more accurate)
profile_result <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "treatment"),
    p_digits = 4,
    conf_method = "profile"
)
print(profile_result)

# Wald CIs (faster, suitable for simulation or exploratory work)
wald_result <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "treatment"),
    p_digits = 4,
    conf_method = "wald"
)
print(wald_result)

# Example 18: Check model quality statistics
result <- fit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment", "stage"),
    keep_qc_stats = TRUE
)

raw <- attr(result, "raw_data")
cat("AIC:", raw$AIC[1], "\n")
cat("BIC:", raw$BIC[1], "\n")
cat("C-statistic:", raw$c_statistic[1], "\n")

# Example 19: Interaction effects - treatment effect modified by stage
interaction_model <- fit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    predictors = c("age", "treatment", "stage"),
    interactions = c("treatment:stage"),
    model_type = "coxph",
    labels = clintrial_labels
)
print(interaction_model)
# Shows main effects plus all treatment×stage interaction terms

# Example 20: Multiple interactions in logistic regression
multi_interaction <- fit(
    data = clintrial,
    outcome = "readmission_30d",
    predictors = c("age", "sex", "surgery", "diabetes"),
    interactions = c("surgery:diabetes", "age:sex"),
    labels = clintrial_labels
)
print(multi_interaction)

# Example 21: Quasipoisson for overdispersed count data
# Alternative to negative binomial when MASS not available
quasi_model <- fit(
    data = clintrial,
    outcome = "ae_count",
    predictors = c("age", "treatment", "diabetes", "surgery"),
    model_type = "glm",
    family = "quasipoisson",
    labels = clintrial_labels
)
print(quasi_model)
# Adjusts standard errors for overdispersion

# Example 22: Quasibinomial for overdispersed binary data
quasi_logistic <- fit(
    data = clintrial,
    outcome = "any_complication",
    predictors = c("age", "bmi", "diabetes", "surgery"),
    model_type = "glm",
    family = "quasibinomial",
    labels = clintrial_labels
)
print(quasi_logistic)

# Example 23: Gamma regression with identity link for additive effects
gamma_identity <- fit(
    data = clintrial,
    outcome = "los_days",
    predictors = c("age", "treatment", "surgery", "any_complication"),
    model_type = "glm",
    family = Gamma(link = "identity"),
    labels = clintrial_labels
)
print(gamma_identity)
# Shows additive effects (coefficients) instead of multiplicative (ratios)

# Example 24: Inverse Gaussian regression for highly skewed data
inverse_gaussian <- fit(
    data = clintrial,
    outcome = "recovery_days",
    predictors = c("age", "surgery", "pain_score"),
    model_type = "glm",
    family = inverse.gaussian(link = "log"),
    labels = clintrial_labels
)
print(inverse_gaussian)

# Example 25: Linear mixed effects with random intercepts
# Accounts for clustering of patients within sites
if (requireNamespace("lme4", quietly = TRUE)) {
    lmer_model <- fit(
        data = clintrial,
        outcome = "los_days",
        predictors = c("age", "treatment", "stage", "(1|site)"),
        model_type = "lmer",
        labels = clintrial_labels
    )
    print(lmer_model)
}

# Example 26: Generalized linear mixed effects (logistic with random effects)
if (requireNamespace("lme4", quietly = TRUE)) {
    glmer_model <- fit(
        data = clintrial,
        outcome = "readmission_30d",
        predictors = c("age", "surgery", "los_days", "(1|site)"),
        model_type = "glmer",
        family = "binomial",
        labels = clintrial_labels
    )
    print(glmer_model)
}

# Example 27: Cox mixed effects for clustered survival data
if (requireNamespace("coxme", quietly = TRUE)) {
    coxme_model <- fit(
        data = clintrial,
        outcome = "Surv(os_months, os_status)",
        predictors = c("age", "treatment", "stage", "(1|site)"),
        model_type = "coxme",
        labels = clintrial_labels
    )
    print(coxme_model)
}

# Example 28: Random slopes - treatment effect varies by site
if (requireNamespace("lme4", quietly = TRUE)) {
    random_slopes <- fit(
        data = clintrial,
        outcome = "los_days",
        predictors = c("age", "treatment", "stage", "(treatment|site)"),
        model_type = "lmer",
        labels = clintrial_labels
    )
    print(random_slopes)
}

# Example 29: Format a pre-fitted model (model-based workflow)
# Useful for models fitted outside of fit()
pre_fitted <- glm(os_status ~ age + sex + treatment,
                  family = binomial, data = clintrial)
result <- fit(model = pre_fitted, 
              data = clintrial,
              labels = clintrial_labels)
print(result)

Description

Executes a comprehensive regression analysis pipeline that combines univariable screening, automatic/manual variable selection, and multivariable modeling in a single function call. This function is designed to streamline the complete analytical workflow from initial exploration to final adjusted models, with publication-ready formatted output showing both univariable and multivariable results side-by-side if desired.

Usage

fullfit(
  data,
  outcome,
  predictors,
  method = "screen",
  multi_predictors = NULL,
  p_threshold = 0.05,
  columns = "both",
  model_type = "glm",
  family = "binomial",
  random = NULL,
  conf_level = 0.95,
  reference_rows = TRUE,
  show_n = TRUE,
  show_events = TRUE,
  digits = 2,
  p_digits = 3,
  labels = NULL,
  metrics = "both",
  return_type = "table",
  keep_models = FALSE,
  exponentiate = NULL,
  conf_method = NULL,
  parallel = TRUE,
  n_cores = NULL,
  number_format = NULL,
  verbose = NULL,
  ...
)

Arguments

    options(summata.number_format = "eu")
  

Details

Analysis Workflow:

The function implements a complete regression analysis pipeline:

1. Univariable screening: Fits separate models for each predictor (outcome ~ predictor). Each predictor is tested independently to understand crude associations.

2. Variable selection: Based on the method parameter:

   • "screen": Automatically selects predictors with univariable p $\le$ p_threshold

   • "all": Includes all predictors (no selection)

   • "custom": Uses predictors specified in multi_predictors

3. Multivariable modeling: Fits a single model with selected predictors (outcome ~ predictor1 + predictor2 + ...). Estimates are adjusted for all other variables in the model.

4. Output formatting: Combines results into publication-ready table with appropriate effect measures and formatting.

Variable Selection Strategies:

"Screen" Method (method = "screen"):

• Uses p-value threshold for automatic selection

• Liberal thresholds (e.g., 0.20) cast a wide net to avoid missing important predictors

• Stricter thresholds (e.g., 0.05) focus on strongly associated predictors

• Helps reduce overfitting and multicollinearity

• Common in exploratory analyses and when sample size is limited

"All" Method (method = "all"):

• No variable selection - includes all predictors

• Appropriate when all variables are theoretically important

• Risk of overfitting with many predictors relative to sample size

• Useful for confirmatory analyses with pre-specified models

"Custom" Method (method = "custom"):

• Manual selection based on subject matter knowledge

• Runs univariable analysis for all predictors (for comparison)

• Includes only specified predictors in multivariable model

• Ideal for theory-driven model building

• Allows comparison of unadjusted vs adjusted effects for all variables

Interpreting Results:

When columns = "both" (default), tables show:

• Univariable columns: Crude associations, unadjusted for other variables. Labeled as "OR/HR/RR/Coefficient (95% CI)" and "Uni p"

• Multivariable columns: Adjusted associations, accounting for all other predictors in the model. Labeled as "aOR/aHR/aRR/Adj. Coefficient (95% CI)" and "Multi p" ("a" = adjusted)

• Variables not meeting selection criteria show "-" in multivariable columns

Comparing univariable and multivariable results helps identify:

• Confounding: Large changes in effect estimates

• Independent effects: Similar univariable and multivariable estimates

• Mediation: Attenuated effects in multivariable model

• Suppression: Effects that emerge only after adjustment

Sample Size Considerations:

Rule of thumb for multivariable models:

• Logistic regression: $\ge$ 10 events per predictor variable

• Cox regression: $\ge$ 10 events per predictor variable

• Linear regression: $\ge$ 10-20 observations per predictor

Use screening methods to reduce predictor count when these ratios are not met.

Value

Depends on return_type parameter:

When return_type = "table" (default): A data.table with S3 class "fullfit_result" containing:

Variable

Character. Predictor name or custom label

Group

Character. Category level for factors, empty for continuous

n/n_group

Integer. Sample sizes (if show_n = TRUE). For variables included in the multivariable model, reflects the complete-case sample size from the fitted model (listwise deletion across all included predictors). For variables not selected into the multivariable model, reflects the per-variable sample size from the univariable analysis. This follows STROBE guideline item 12, which recommends reporting the number of participants included at each stage of analysis.

events/events_group

Integer. Event counts (if show_events = TRUE). Same complete-case convention as n: multivariable rows show events from the fitted model, univariable-only rows show per-variable counts.

OR/HR/RR/Coefficient (95% CI)

Character. Unadjusted effect (if columns includes "uni" and metrics includes "effect")

Uni p

Character. Univariable p-value (if columns includes "uni" and metrics includes "p")

aOR/aHR/aRR/Adj. Coefficient (95% CI)

Character. Adjusted effect (if columns includes "multi" and metrics includes "effect")

Multi p

Character. Multivariable p-value (if columns includes "multi" and metrics includes "p")

When return_type = "model": The fitted multivariable model object (glm, lm, coxph, etc.).

When return_type = "both": A list with two elements:

table

The formatted results data.table

model

The fitted multivariable model object

The table includes the following attributes:

outcome

Character. The outcome variable name

model_type

Character. The regression model type

method

Character. The variable selection method used

columns

Character. Which columns were displayed

model

The multivariable model object (if fitted)

uni_results

The complete univariable screening results

n_multi

Integer. Number of predictors in multivariable model

screened

Character vector. Names of predictors that passed univariable screening at the specified p-value threshold

significant

Character vector. Names of variables with p < 0.05 in the multivariable model (or univariable if multivariable was not fitted)

See Also

uniscreen for univariable screening only, fit for fitting a single multivariable model, compfit for comparing multiple models, desctable for descriptive statistics

Other regression functions: compfit(), fit(), multifit(), print.compfit_result(), print.fit_result(), print.fullfit_result(), print.multifit_result(), print.uniscreen_result(), uniscreen()

Examples

# Load example data
data(clintrial)
data(clintrial_labels)

# Example 1: Basic screening with p < 0.05 threshold
result1 <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "bmi", "smoking",
                   "hypertension", "diabetes",
                   "treatment", "stage"),
    method = "screen",
    p_threshold = 0.05,
    labels = clintrial_labels
)
print(result1)
# Shows both univariable and multivariable results
# Only significant univariable predictors in multivariable model



# Example 2: Include all predictors (no selection)
result2 <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment", "stage"),
    method = "all",
    labels = clintrial_labels
)
print(result2)

# Example 3: Custom variable selection
result3 <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "bmi", "smoking", "treatment", "stage"),
    method = "custom",
    multi_predictors = c("age", "treatment", "stage"),
    labels = clintrial_labels
)
print(result3)
# Univariable for all, multivariable for selected only

# Example 4: Cox regression with screening
library(survival)
cox_result <- fullfit(
    data = clintrial,
    outcome = "Surv(os_months, os_status)",
    predictors = c("age", "sex", "treatment", "stage"),
    model_type = "coxph",
    method = "screen",
    p_threshold = 0.10,
    labels = clintrial_labels
)
print(cox_result)

# Example 5: Linear regression without screening
linear_result <- fullfit(
    data = clintrial,
    outcome = "bmi",
    predictors = c("age", "sex", "smoking", "creatinine"),
    model_type = "lm",
    method = "all",
    labels = clintrial_labels
)
print(linear_result)

# Example 6: Poisson regression for count outcomes
poisson_result <- fullfit(
    data = clintrial,
    outcome = "fu_count",
    predictors = c("age", "stage", "treatment", "surgery"),
    model_type = "glm",
    family = "poisson",
    method = "all",
    labels = clintrial_labels
)
print(poisson_result)

# Example 7: Show only multivariable results
multi_only <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment", "stage"),
    method = "all",
    columns = "multi",
    labels = clintrial_labels
)
print(multi_only)

# Example 8: Return both table and model object
both <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "sex", "treatment", "stage"),
    method = "all",
    return_type = "both"
)
print(both$table)
summary(both$model)

# Example 9: Keep univariable models for diagnostics
with_models <- fullfit(
    data = clintrial,
    outcome = "os_status",
    predictors = c("age", "bmi", "creatinine"),
    keep_models = TRUE
)
uni_results <- attr(with_models, "uni_results")
uni_models <- attr(uni_results, "models")
summary(uni_models[["age"]])

# Example 10: Linear mixed effects with site clustering
if (requireNamespace("lme4", quietly = TRUE)) {
    lmer_result <- fullfit(
        data = clintrial,
        outcome = "los_days",
        predictors = c("age", "treatment", "surgery", "stage"),
        random = "(1|site)",
        model_type = "lmer",
        method = "all",
        labels = clintrial_labels
    )
    print(lmer_result)
}

Description

Generates a publication-ready forest plot that combines a formatted data table with a graphical representation of effect estimates (odds ratios, risk ratios, or coefficients) from a generalized linear model. The plot integrates variable names, group levels, sample sizes, effect estimates with confidence intervals, p-values, and model diagnostics in a single comprehensive visualization designed for manuscripts and presentations.

Usage

glmforest(
  x,
  data = NULL,
  title = "Generalized Linear Model",
  effect_label = NULL,
  digits = 2,
  p_digits = 3,
  conf_level = 0.95,
  font_size = 1,
  annot_size = 3.88,
  header_size = 5.82,
  title_size = 23.28,
  plot_width = NULL,
  plot_height = NULL,
  table_width = 0.6,
  show_n = TRUE,
  show_events = TRUE,
  indent_groups = FALSE,
  condense_table = FALSE,
  bold_variables = FALSE,
  center_padding = 4,
  zebra_stripes = TRUE,
  ref_label = "reference",
  labels = NULL,
  color = NULL,
  exponentiate = NULL,
  qc_footer = TRUE,
  units = "in",
  number_format = NULL
)

Arguments

    options(summata.number_format = "eu")
  

Details

Plot Components:

The forest plot consists of several integrated components:

1. Title: Centered at top, describes the analysis

2. Data Table (left side): Contains columns for:

   • Variable: Predictor names (or custom labels)

   • Group: Factor levels (optional, hidden when indenting)

   • n: Sample sizes by group (optional)

   • Events: Event counts by group (optional)

   • Effect (95% CI); p-value: Formatted estimates with p-values

3. Forest Plot (right side): Graphical display with:

   • Point estimates (squares sized by sample size)

   • 95% confidence intervals (error bars)

   • Reference line (at OR/RR = 1 or coefficient = 0)

   • Log scale for odds/risk ratios

   • Labeled axis

4. Model Statistics (footer): Summary of:

   • Observations analyzed (with percentage of total data)

   • Model family (Binomial, Poisson, etc.)

   • Deviance statistics

   • Pseudo-R$^2$ (McFadden)

   • AIC

Automatic Effect Measure Selection:

When effect_label = NULL and exponentiate = NULL, the function intelligently selects the appropriate effect measure:

• Logistic regression (family = binomial(link = "logit")): Odds Ratios (OR)

• Log-link models (link = "log"): Risk Ratios (RR) or Rate Ratios

• Other exponential families: exp(coefficient)

• Identity link: Raw coefficients

Reference Categories:

For factor variables, the first level (determined by factor ordering or alphabetically for character variables) serves as the reference category:

• Displayed with the ref_label instead of an estimate

• No confidence interval or p-value shown

• Visually aligned with other categories

• When condense_table = TRUE, reference-only variables may be omitted entirely

Layout Optimization:

The function automatically optimizes layout based on content:

• Calculates appropriate axis ranges to accommodate all confidence intervals

• Selects meaningful tick marks on log or linear scales

• Sizes point markers proportional to sample size (larger = more data)

• Adjusts table width based on variable name lengths when table_width = NULL

• Recommends overall dimensions based on number of rows

Visual Grouping Options:

Three display modes are available:

1. Standard (indent_groups = FALSE, condense_table = FALSE): Separate "Variable" and "Group" columns, all categories shown

2. Indented (indent_groups = TRUE, condense_table = FALSE): Hierarchical display with groups indented under variables

3. Condensed (condense_table = TRUE): Binary variables shown in single rows, automatically indented

Zebra Striping:

When zebra_stripes = TRUE, alternating variables (not individual rows) receive light gray backgrounds. This helps visually group all levels of a factor variable together, making the plot easier to read especially with many multi-level factors.

Model Statistics Display:

The footer shows key diagnostic information:

• Observations analyzed: Total N and percentage of original data (accounting for missing values)

• Null/Residual Deviance: Model fit improvement

• Pseudo-R$^2$: McFadden R$^2$ = 1 - (log L_1 / log L_2)

• AIC: For model comparison (lower is better)

For logistic regression, concordance (C-statistic/AUC) may also be displayed if available.

Saving Plots:

Use ggplot2::ggsave() with recommended dimensions:

  p <- glmforest(model, data)
  dims <- attr(p, "rec_dims")
  ggplot2::ggsave("forest.pdf", p, width = dims$width, height = dims$height)

Or specify custom dimensions:

ggplot2::ggsave("forest.png", p, width = 12, height = 8, dpi = 300)

Value

A ggplot object containing the complete forest plot. The plot can be:

• Displayed directly: print(plot)

• Saved to file: ggsave("forest.pdf", plot, width = 12, height = 8)

• Further customized with ggplot2 functions

The returned object includes an attribute "rec_dims" accessible via attr(plot, "rec_dims"), which is a list containing:

width

Numeric. Recommended plot width in specified units

height

Numeric. Recommended plot height in specified units

These recommendations are automatically calculated based on the number of variables, text sizes, and layout parameters, and are printed to console if plot_width or plot_height are not specified.

See Also

autoforest for automatic model detection, coxforest for Cox proportional hazards forest plots, lmforest for linear model forest plots, uniforest for univariable screening forest plots, multiforest for multi-outcome forest plots, glm for fitting GLMs, fit for regression modeling

Other visualization functions: autoforest(), coxforest(), lmforest(), multiforest(), uniforest()

Examples

data(clintrial)
data(clintrial_labels)

# Create example model
model1 <- glm(os_status ~ age + sex + bmi + treatment,
              data = clintrial, family = binomial)

# Example 1: Basic logistic regression forest plot
p <- glmforest(model1, data = clintrial)



old_width <- options(width = 180)

# Example 2: With custom variable labels
plot2 <- glmforest(
    x = model1,
    data = clintrial,
    title = "Risk Factors for Mortality",
    labels = clintrial_labels
)

# Example 3: Indented layout with formatting options
plot3 <- glmforest(
    x = model1,
    data = clintrial,
    indent_groups = TRUE,
    zebra_stripes = TRUE,
    color = "#D62728",
    labels = clintrial_labels
)

# Example 4: Condensed layout for many binary variables
model4 <- glm(os_status ~ age + sex + smoking + hypertension + 
                  diabetes + surgery,
              data = clintrial,
              family = binomial)

plot4 <- glmforest(
    x = model4,
    data = clintrial,
    condense_table = TRUE,
    labels = clintrial_labels
)
# Binary variables shown in single rows

# Example 5: Poisson regression for count data
model5 <- glm(ae_count ~ age + treatment + diabetes + surgery,
               data = clintrial,
               family = poisson)

plot5 <- glmforest(
    x = model5,
    data = clintrial,
    title = "Rate Ratios for Adverse Events",
    labels = clintrial_labels
)

# Example 6: Save with recommended dimensions
dims <- attr(plot5, "rec_dims")
ggplot2::ggsave(file.path(tempdir(), "forest.pdf"),
                plot5, width = dims$width, height = dims$height)

options(old_width)

Description

Generates a publication-ready forest plot that combines a formatted data table with a graphical representation of regression coefficients from a linear model. The plot integrates variable names, group levels, sample sizes, coefficients with confidence intervals, p-values, and model diagnostics (R$^2$, F-statistic, AIC) in a single comprehensive visualization designed for manuscripts and presentations.

Usage

lmforest(
  x,
  data = NULL,
  title = "Linear Model",
  effect_label = "Coefficient",
  digits = 2,
  p_digits = 3,
  conf_level = 0.95,
  font_size = 1,
  annot_size = 3.88,
  header_size = 5.82,
  title_size = 23.28,
  plot_width = NULL,
  plot_height = NULL,
  table_width = 0.6,
  show_n = TRUE,
  indent_groups = FALSE,
  condense_table = FALSE,
  bold_variables = FALSE,
  center_padding = 4,
  zebra_stripes = TRUE,
  ref_label = "reference",
  labels = NULL,
  units = "in",
  color = "#5A8F5A",
  qc_footer = TRUE,
  number_format = NULL
)

Arguments

    options(summata.number_format = "eu")
  

Details

Linear Model-Specific Features:

The linear model forest plot differs from logistic and Cox plots in several ways:

• Coefficients: Raw regression coefficients shown (not exponentiated)

• Reference line: At coefficient = 0 (not at 1)

• Linear scale: Forest plot uses linear scale (not log scale)

• No events column: Only sample sizes shown (no event counts)

• R$^2$ statistics: Model fit assessed by R$^2$ and adjusted R$^2$

• F-test: Overall model significance from F-statistic

Plot Components:

1. Title: Centered at top

2. Data Table (left): Contains:

   • Variable: Predictor names

   • Group: Factor levels (if applicable)

   • n: Sample sizes by group

   • Coefficient (95% CI); p-value: Raw coefficients with CIs and p-values

3. Forest Plot (right):

   • Point estimates (squares sized by sample size)

   • 95% confidence intervals (error bars)

   • Reference line at coefficient = 0

   • Linear scale

4. Model Statistics (footer):

   • Observations analyzed (with percentage of total data)

   • R$^2$ and adjusted R$^2$

   • F-statistic with degrees of freedom and p-value

   • AIC

Interpreting Coefficients:

Linear regression coefficients represent the change in the outcome variable for a one-unit change in the predictor:

• Continuous predictors: Coefficient = change in Y per unit of X

• Binary predictors: Coefficient = difference in Y between groups

• Factor predictors: Coefficients = differences from reference category

• Sign matters: Positive = increase in Y, Negative = decrease in Y

• Zero crossing: CI crossing zero suggests no significant effect

Example: If the coefficient for "age" is 0.50 when predicting BMI, BMI increases by 0.50 kg/m$^2$ for each additional year of age.

Model Fit Statistics:

The footer displays key diagnostics:

• R$^2$: Proportion of variance explained (0 to 1)

  • 0.0-0.3: Weak explanatory power

  • 0.3-0.5: Moderate

  • 0.5-0.7: Good

  • > 0.7: Strong (rare in social/biological sciences)

• Adjusted R$^2$: R$^2$ penalized for number of predictors

  • Always $\le$ R$^2$

  • Preferred for model comparison

  • Accounts for model complexity

• F-statistic: Tests null hypothesis that all coefficients = 0

  • Degrees of freedom: df1 = # predictors, df2 = # observations - # predictors - 1

  • Significant p-value indicates model explains variance better than intercept-only

• AIC: For model comparison (lower is better)

Assumptions:

Linear regression assumes:

1. Linearity of relationships

2. Independence of observations

3. Homoscedasticity (constant variance)

4. Normality of residuals

5. No multicollinearity

Check assumptions using:

• plot(model) for diagnostic plots

• car::vif(model) for multicollinearity

• lmtest::bptest(model) for heteroscedasticity

• shapiro.test(residuals(model)) for normality

Reference Categories:

For factor variables:

• First level is the reference (coefficient = 0)

• Other levels show difference from reference

• Reference displayed with ref_label

• Relevel factors before modeling if needed: factor(x, levels = c("desired_ref", ...))

Sample Size Reporting:

The "n" column shows:

• For continuous variables: Total observations with non-missing data

• For factor variables: Number of observations in each category

• Footer shows total observations analyzed and percentage of original data (accounting for missing values)

Value

A ggplot object containing the complete forest plot. The plot can be:

• Displayed directly: print(plot)

• Saved to file: ggsave("forest.pdf", plot, width = 12, height = 8)

• Further customized with ggplot2 functions

The returned object includes an attribute "rec_dims" accessible via attr(plot, "rec_dims"), which is a list containing:

width

Numeric. Recommended plot width in specified units

height

Numeric. Recommended plot height in specified units

These recommendations are automatically calculated based on the number of variables, text sizes, and layout parameters, and are printed to console if plot_width or plot_height are not specified.

See Also

autoforest for automatic model detection, glmforest for logistic/GLM forest plots, coxforest for Cox model forest plots, uniforest for univariable screening forest plots, multiforest for multi-outcome forest plots, lm for fitting linear models, fit for regression modeling

Other visualization functions: autoforest(), coxforest(), glmforest(), multiforest(), uniforest()

Examples

data(clintrial)
data(clintrial_labels)

# Create example model
model1 <- lm(bmi ~ age + sex + smoking, data = clintrial)

# Example 1: Basic linear model forest plot
p <- lmforest(model1, data = clintrial)



old_width <- options(width = 180)

# Example 2: With custom labels and title
plot2 <- lmforest(
    x = model1,
    data = clintrial,
    title = "Predictors of Body Mass Index",
    effect_label = "Change in BMI (kg/m^2)",
    labels = clintrial_labels
)

# Example 3: Comprehensive model with indented layout
model3 <- lm(
    bmi ~ age + sex + smoking + hypertension + diabetes + creatinine,
    data = clintrial
)

plot3 <- lmforest(
    x = model3,
    data = clintrial,
    labels = clintrial_labels,
    indent_groups = TRUE,
    zebra_stripes = TRUE
)

# Example 4: Condensed layout
plot4 <- lmforest(
    x = model3,
    data = clintrial,
    condense_table = TRUE,
    labels = clintrial_labels
)

# Example 5: Different outcome (hemoglobin)
model5 <- lm(
    hemoglobin ~ age + sex + bmi + smoking + creatinine,
    data = clintrial
)

plot5 <- lmforest(
    x = model5,
    data = clintrial,
    title = "Predictors of Baseline Hemoglobin",
    effect_label = "Change in Hemoglobin (g/dL)",
    labels = clintrial_labels
)

# Example 6: Save with recommended dimensions
dims <- attr(plot5, "rec_dims")
ggplot2::ggsave(file.path(tempdir(), "linear_forest.pdf"),
                plot5, width = dims$width, height = dims$height)

options(old_width)

Description

Extracts coefficients, confidence intervals, and comprehensive model statistics from fitted regression models and converts them to a standardized data.table format suitable for further analysis or publication. This is a core utility function frequently used internally by other summata regression functions, although it can be used as a standalone function as well.

Usage

m2dt(
  data,
  model,
  conf_level = 0.95,
  keep_qc_stats = TRUE,
  include_intercept = TRUE,
  terms_to_exclude = NULL,
  reference_rows = TRUE,
  reference_label = "reference",
  skip_counts = FALSE,
  conf_method = NULL
)

Arguments

Details

This function is the core extraction utility used by fit() and other regression functions. It handles the complexities of different model classes and provides a consistent output format suitable for tables and forest plots.

Model Type Detection: The function automatically detects model type and applies appropriate:

• Effect measure naming (OR, HR, RR, Coefficient)

• Confidence interval calculation (see below)

• Event counting for binary/survival outcomes

Confidence Interval Methods: The CI method is selected per model class using stats::confint() dispatch:

• GLM/negative binomial: Profile likelihood via MASS::confint.glm(), except quasi-families which use Wald

• Linear models: Exact t-distribution via confint.lm()

• Cox PH: Wald intervals (coefficient $\pm$ z $\times$ SE)

• Mixed-effects models: Wald intervals

Falls back to Wald on profiling failure.

Mixed Effects Models: For lme4 models (glmer, lmer), the function extracts fixed effects only. Random effects variance components are not included in the output table, as they represent clustering structure rather than predictor effects.

Value

A data.table containing extracted model information with the following standard columns:

model_scope

Character. Either "Univariable" (unadjusted model with single predictor) or "Multivariable" (adjusted model with multiple predictors)

model_type

Character. Type of regression (e.g., "Logistic", "Linear", "Cox PH", "Poisson", etc.)

variable

Character. Variable name (for factor variables, the base variable name without the level)

group

Character. Group/level name for factor variables; empty string for continuous variables

n

Integer. Total sample size used in the model

n_group

Integer. Sample size for this specific variable level (factor variables only)

events

Integer. Total number of events in the model (for survival and logistic models)

events_group

Integer. Number of events for this specific variable level (for survival and logistic models with factor variables)

coefficient

Numeric. Raw regression coefficient (log odds, log hazard, etc.)

se

Numeric. Standard error of the coefficient

OR/HR/RR/Coefficient

Numeric. Effect estimate - column name depends on model type:

• OR for logistic regression (odds ratio)

• HR for Cox models (hazard ratio)

• RR for Poisson regression (rate/risk ratio)

• Coefficient for linear models or other GLMs

ci_lower

Numeric. Lower bound of confidence interval for effect estimate

ci_upper

Numeric. Upper bound of confidence interval for effect estimate

statistic

Numeric. Test statistic (z-value for GLM/Cox, t-value for LM)

p_value

Numeric. p-value for coefficient test

sig

Character. Significance markers: *** (p < 0.001), ** (p < 0.01), * (p < 0.05), . (p < 0.10).

sig_binary

Logical. Binary indicator: TRUE if p < 0.05, FALSE otherwise

reference

Character. Contains reference_label for reference category rows when reference_rows = TRUE, empty string otherwise

See Also

fit for the main regression interface, glmforest, coxforest, lmforest for forest plot visualization

Examples

# Load example data
data(clintrial)

# Example 1: Extract from logistic regression
glm_model <- glm(os_status ~ age + sex + treatment, 
                 data = clintrial, family = binomial)

glm_result <- m2dt(clintrial, glm_model)
glm_result



# Example 2: Extract from linear model
lm_model <- lm(los_days ~ age + sex + surgery, data = clintrial)

lm_result <- m2dt(clintrial, lm_model)
lm_result

# Example 3: Cox proportional hazards model
library(survival)
cox_model <- coxph(Surv(os_months, os_status) ~ age + sex + stage,
                   data = clintrial)

cox_result <- m2dt(clintrial, cox_model)
cox_result

# Example 4: Exclude intercept for cleaner tables
clean_result <- m2dt(clintrial, glm_model, include_intercept = FALSE)
clean_result

# Example 5: Change confidence level
result_90ci <- m2dt(clintrial, glm_model, conf_level = 0.90)
result_90ci