ETC5523: Communicating with Data

Tutorial 5

Author

Michael Lydeamore

Published

1 May 2001

🎯 Objectives

  • extract model information with accessor functions and broom
  • build, critique and refine regression and descriptive tables with gtsummary and kableExtra
  • apply communication choices for precision, units, alignment and labels
  • select and interpret the results that matter for a stated audience

Install the R-packages (once per machine):

install.packages(c("broom", "carData", "gtsummary", "kableExtra", "tidyverse"))

The Salaries data are in the carData package — no separate download is needed.

Data: Academic salaries

The Salaries data set in carData contains 9-month salaries for 397 college faculty in the US (carData::Salaries). Variables:

  • salary — 9-month salary (USD)
  • rank — AsstProf, AssocProf, Prof
  • discipline — A (theoretical) or B (applied)
  • yrs.since.phd — years since PhD
  • yrs.service — years of service
  • sex — Female / Male

We use one model throughout the tutorial so you can see the full workflow from model object to communicative table.

library(tidyverse)
library(carData)
library(broom)
library(gtsummary)
library(kableExtra)

data(Salaries, package = "carData")
glimpse(Salaries)
Rows: 397
Columns: 6
$ rank          <fct> Prof, Prof, AsstProf, Prof, Prof, AssocProf, Prof, Prof,…
$ discipline    <fct> B, B, B, B, B, B, B, B, B, B, B, B, B, B, B, B, B, A, A,…
$ yrs.since.phd <int> 19, 20, 4, 45, 40, 6, 30, 45, 21, 18, 12, 7, 1, 2, 20, 1…
$ yrs.service   <int> 18, 16, 3, 39, 41, 6, 23, 45, 20, 18, 8, 2, 1, 0, 18, 3,…
$ sex           <fct> Male, Male, Male, Male, Male, Male, Male, Male, Male, Fe…
$ salary        <int> 139750, 173200, 79750, 115000, 141500, 97000, 175000, 14…

👥 Exercise 5A

From model object to communicative regression table

fit <- lm(salary ~ rank + yrs.since.phd + sex, data = Salaries)

fit contains much more than a printed summary. Your audience is a technically competent, time-poor head of school who wants to know whether rank and sex are associated with salary after adjusting for time since PhD.

A1. Extract the right pieces

Using fit, compare what you get from:

  • coef(fit) and confint(fit) / summary(fit)
  • broom::tidy(fit, conf.int = TRUE) and broom::glance(fit) and broom::augment(fit) (first 6 rows)

Which broom verb gives you:

  1. one row per term,
  2. one row per model,
  3. one row per observation?

Why is tidy(fit, conf.int = TRUE) a better starting point for a table than copying values from summary(fit)?

  (Intercept) rankAssocProf      rankProf yrs.since.phd       sexMale 
  76944.48523   14012.79755   47635.89764     -92.10413    5146.60894 
                   2.5 %     97.5 %
(Intercept)   68184.4728 85704.4977
rankAssocProf  5474.7541 22550.8410
rankProf      38960.4649 56311.3304
yrs.since.phd  -347.0786   162.8704
sexMale       -2793.9305 13087.1484

Call:
lm(formula = salary ~ rank + yrs.since.phd + sex, data = Salaries)

Residuals:
   Min     1Q Median     3Q    Max 
-67230 -15338  -1530  12163 105318 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    76944.5     4455.7  17.269  < 2e-16 ***
rankAssocProf  14012.8     4342.8   3.227  0.00136 ** 
rankProf       47635.9     4412.7  10.795  < 2e-16 ***
yrs.since.phd    -92.1      129.7  -0.710  0.47801    
sexMale         5146.6     4038.9   1.274  0.20332    
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 23630 on 392 degrees of freedom
Multiple R-squared:  0.3973,    Adjusted R-squared:  0.3912 
F-statistic: 64.61 on 4 and 392 DF,  p-value: < 2.2e-16
# A tibble: 5 × 7
  term          estimate std.error statistic  p.value conf.low conf.high
  <chr>            <dbl>     <dbl>     <dbl>    <dbl>    <dbl>     <dbl>
1 (Intercept)    76944.      4456.    17.3   4.26e-50   68184.    85704.
2 rankAssocProf  14013.      4343.     3.23  1.36e- 3    5475.    22551.
3 rankProf       47636.      4413.    10.8   5.83e-24   38960.    56311.
4 yrs.since.phd    -92.1      130.    -0.710 4.78e- 1    -347.      163.
5 sexMale         5147.      4039.     1.27  2.03e- 1   -2794.    13087.
# A tibble: 1 × 12
  r.squared adj.r.squared  sigma statistic  p.value    df logLik   AIC   BIC
      <dbl>         <dbl>  <dbl>     <dbl>    <dbl> <dbl>  <dbl> <dbl> <dbl>
1     0.397         0.391 23633.      64.6 6.18e-42     4 -4559. 9130. 9153.
# ℹ 3 more variables: deviance <dbl>, df.residual <int>, nobs <int>
# A tibble: 6 × 10
  salary rank      yrs.since.phd sex   .fitted  .resid    .hat .sigma   .cooksd
   <int> <fct>             <int> <fct>   <dbl>   <dbl>   <dbl>  <dbl>     <dbl>
1 139750 Prof                 19 Male  127977.  11773. 0.00658 23656. 0.000331 
2 173200 Prof                 20 Male  127885.  45315. 0.00604 23552. 0.00450  
3  79750 AsstProf              4 Male   81723.  -1973. 0.0158  23663. 0.0000227
4 115000 Prof                 45 Male  125582. -10582. 0.0121  23657. 0.000499 
5 141500 Prof                 40 Male  126043.  15457. 0.00791 23651. 0.000688 
6  97000 AssocProf             6 Male   95551.   1449. 0.0192  23663. 0.0000150
# ℹ 1 more variable: .std.resid <dbl>
  • tidy — one row per term (estimate, SE, CI, p-value). This is the table-ready form.
  • glance — one row per model (R², AIC, n, etc.).
  • augment — one row per observation (fitted, residual, influence).

tidy() returns a data frame you can filter, rename and pipe into a table without manual copying. It preserves units and uncertainty together and works consistently across model classes. summary() is for console inspection, not for communication.

A2. Start from a default table

Build the default regression table with gtsummary:

tbl_regression(fit, intercept = TRUE)

Render it (run it yourself). Then critique it for the head-of-school audience:

  1. Which row answers the substantive question? Is the intercept useful here?
  2. Can the labels be understood without seeing the R code?
  3. Does every column (e.g. p-values) need to be shown?
  4. Is the displayed precision appropriate for salaries in dollars?
tbl_regression(fit, intercept = TRUE)
Characteristic Beta 95% CI p-value
(Intercept) 76,944 68,184, 85,704 <0.001
rank


    AsstProf — —
    AssocProf 14,013 5,475, 22,551 0.001
    Prof 47,636 38,960, 56,311 <0.001
yrs.since.phd -92 -347, 163 0.5
sex


    Female — —
    Male 5,147 -2,794, 13,087 0.2
Abbreviation: CI = Confidence Interval

Critique: the intercept (expected salary when rank=AsstProf, yrs.since.phd=0, sex=Female) is not substantively interesting and distracts from the comparisons of interest. Labels such as rank and yrs.since.phd need plain English and units. P-values add little for this audience; the confidence interval already conveys uncertainty around the size of each difference. Precision should be to the nearest dollar (or nearest $100) with comma grouping, not many decimal places.

A3. Refine the table for its purpose

Refine the table by adapting the lecture workflow:

  • hide the intercept (intercept = FALSE),
  • give plain labels with units,
  • round to an appropriate precision,
  • hide p-values for this audience,
  • embolden variable labels.

Render the refined table. What changed and why?

Characteristic Beta 95% CI
Rank

    AsstProf — —
    AssocProf 14,013 5,475, 22,551
    Prof 47,636 38,960, 56,311
Years since PhD -92 -347, 163
Sex

    Female — —
    Male 5,147 -2,794, 13,087
Abbreviation: CI = Confidence Interval

Changes: intercept removed, labels are audience-readable with no code names, estimates shown as whole dollars with 1,000 grouping (salaries do not warrant decimals), p-value hidden because the interval is the relevant uncertainty for this purpose, and bold labels aid scanning. The table now prioritises the comparisons the head of school asked about: rank and sex adjusted for time since PhD.

You could also use estimate_fun = label_style_sigfig(digits = 3) if you prefer significant-figure control; for salaries on the order of $50k–$150k, 0 decimal places and a comma is the key choice.

A4. Write the take-away sentence

Using the refined table (or tidy(fit, conf.int = TRUE)), write one sentence for the head of school that:

  • names the comparison (e.g. Prof vs AsstProf, or Male vs Female, holding yrs.since.phd constant),
  • states the size in dollars with its 95% confidence interval,
  • avoids claiming causation (“associated with”, not “caused by”).

One possible sentence (using Prof vs AsstProf; numbers will match your run):

After adjusting for years since PhD and sex, full professors earned an estimated $47,636 more (95% CI: $38,960 to $56,311) than assistant professors, on average.

Or for sex (note the interval includes zero):

After adjusting for rank and years since PhD, male faculty earned an estimated $5,147 more than female faculty (95% CI: –$2,794 to $13,087), on average — an imprecise difference that is compatible with no association after adjusting for these factors.

Both make the unit (USD), the reference category, and the uncertainty explicit, and do not claim the model has established a causal effect. See week5/index.qmd “Interpret the result before formatting it” for the lecture wording template.

🛠️ Exercise 5B

Descriptive tables and communication polish

B1. A descriptive summary that makes comparison easy

Produce a descriptive summary by sex using the same workflow as the lecture’s by = cyl example. Adapt it for Salaries:

Salaries |>
  select(salary, yrs.since.phd, yrs.service, rank, discipline, sex) |>
  tbl_summary(
    by = sex,
    label = list(
      salary ~ "Salary (USD)",
      yrs.since.phd ~ "Years since PhD",
      yrs.service ~ "Years of service"
    ),
    statistic = all_continuous() ~ "{mean} ({sd})",
    digits = all_continuous() ~ 0,
    missing = "no"
  ) |>
  bold_labels()

Questions:

  1. Why might "{mean} ({sd})" be reasonable for salary here, and when would "{median} ({p25}, {p75})" be preferred?
  2. What does missing = "no" do, and when would you set it to "ifany"?
  3. How do the column labels and digits choices reflect the communication principles for precision and units?
Characteristic Female
N = 391		 Male
N = 3581
Salary (USD) 101,002 (25,952) 115,090 (30,437)
Years since PhD 17 (10) 23 (13)
Years of service 12 (9) 18 (13)
rank

    AsstProf 11 (28%) 56 (16%)
    AssocProf 10 (26%) 54 (15%)
    Prof 18 (46%) 248 (69%)
discipline

    A 18 (46%) 163 (46%)
    B 21 (54%) 195 (54%)
1 Mean (SD); n (%)
  1. Mean (SD) is reasonable if distributions are roughly symmetric; salary is somewhat right-skewed but not extremely so (see summary(Salaries$salary)). If skew or outliers were strong, median (IQR) would better represent centre/spread without being pulled by high salaries — the same point made in week5/index.qmd for mtcars (“Mean (SD) is a choice — not an automatic default”).
  2. missing = "no" suppresses a “Unknown” row when there are no missings. Use "ifany" when missings exist and the audience needs to see how much data are missing — otherwise you hide a data-quality signal.
  3. Labels put units once in the header (“Salary (USD)”), not repeated in every cell. digits = 0 matches the precision of the raw measure (whole dollars) and avoids implying sub-dollar accuracy. This mirrors the lecture’s “Numerical precision”, “Labels within tables”, and “Column alignment” sections.

B2. Polish a table for publication

The code below produces a final-product table badly. It repeats units in every cell, shows 5 decimal places, lacks comma grouping, and misaligns numbers.

bad <- Salaries |>
  slice_head(n = 4) |>
  transmute(
    Rank = rank,
    Salary = paste0("$", salary, " USD"),
    `Years since PhD` = as.character(round(yrs.since.phd + rnorm(4)/1000, 5))
  )

bad |>
  kbl(align = "lcc", caption = "First four faculty (poorly formatted)") |>
  kable_classic(full_width = FALSE)
First four faculty (poorly formatted)
Rank Salary Years since PhD
Prof $139750 USD 19.00077
Prof $173200 USD 20.00106
AsstProf $79750 USD 3.99922
Prof $115000 USD 44.99897

Task: Fix the table using kableExtra (as in the lecture). Your polished table should:

  • put units once in the column header (e.g. Salary (USD)), not in each cell,
  • show salaries as whole dollars with comma grouping (e.g. 139,750),
  • display years with at most 1 decimal and trailing zeroes where needed,
  • right-align numeric columns and left-align text, centre any spanner if you add one,
  • apply kable_classic() styling and add a clear caption.

Hint: scales::comma(), formatC() or kbl(digits=, format.args=) + column_spec() help. The lecture’s “Numerical precision”, “Column alignment”, and “Labels within tables” slides show the exact patterns.

First four faculty — polished for publication
Rank Salary (USD) Years since PhD
Prof 139,750 19
Prof 173,200 20
AsstProf 79,750 4
Prof 115,000 45

Alternatives that also earn full marks:

Key checks: units appear once in headers, numbers are right-aligned, commas aid reading for a 5–6 digit salary (lecture: “1,000,000 is easier than 1000000”), and no “USD” is repeated 397 times.

Take-away checklist (from lecture)
  • Start with the audience’s question, not everything the model produced.
  • Show estimates with uncertainty (CI), units, reference categories, and sample/context — not just p-values.
  • Choose precision, alignment, and labels deliberately; polish with gtsummary + kableExtra.

Workshop 5 continues this by turning a statistical output into a full results paragraph and a self-contained caption.