Clean, Filter, and Summarize

1. Cleaning is a series of decisions

Real EHR exports may contain duplicate visits, mixed units, impossible ages, missing vital signs, and conflicting records. Cleaning is not just deleting inconvenient rows. Every transformation needs a reason, and a record of what changed.

In this synthetic dataset, missing BMI, SBP, and pulse were deliberately inserted for teaching. Do not “correct” them using their underlying simulation formulas in an analysis. That would give the model information unavailable to a real analyst.

Example: suppose a visit has sbp = NA and age = 212. Three decisions are possible for each problem: keep the row and leave NA as unknown, set the impossible age to NA with a note saying why, or drop the row and report how many rows were dropped. Whatever you choose, write it down. In R, a first look at how many rows are affected needs two counts:

raw <- readr::read_csv("data/synthetic_ehr.csv", show_col_types = FALSE)
sum(is.na(raw$sbp))                    # visits with no recorded SBP

[1] 20

sum(raw$age < 0 | raw$age > 110)  # impossible ages

[1] 0

The first count is 20: twenty visits have no recorded SBP, so each needs a decision. The second is 0, so in this table no age is impossible and nothing needs to be set to NA. The rule is the same one a real analyst follows: count first, decide second, and write the decision down.

2. Load dplyr and the example data

library(dplyr)
library(readr)
library(tidyr)
ehr <- read_csv("data/synthetic_ehr.csv", show_col_types = FALSE)

A pipe (|>) passes the result from one instruction into the next. Read it as “and then.” Use one transformation per line at first.

Example: these two lines count the visits of patients aged 80 or older. They give the same answer; the pipe version reads left to right.

nrow(filter(ehr, age >= 80))

[1] 36

ehr |> filter(age >= 80) |> nrow()

[1] 36

Both lines print 36. The pipe reads “take ehr, then keep ages 80 or more, then count the rows”, and the first form has to be read from the inside out.

select() drops columns, filter() drops rows, mutate() adds a column, and summarise() collapses many rows into one.

3. Choose columns with select()

vitals_only <- ehr |>
  select(patient_id, age, bmi, sbp, pulse)
head(vitals_only, 3)

# A tibble: 3 x 5
  patient_id   age   bmi   sbp pulse
  <chr>      <dbl> <dbl> <dbl> <dbl>
1 S0001         55  20.4   120    79
2 S0002         26  22.7   110    92
3 S0003         63  24.7   133    65

Choosing columns does not edit the CSV. The new object vitals_only contains the result.

When a table has many columns, naming each one is slow. select() has helpers that choose columns by pattern or by type: starts_with(), ends_with(), contains() and where().

ehr |> select(starts_with("p")) |> names()       # Names that begin with "p"

[1] "patient_id" "pulse"     

ehr |> select(ends_with("e")) |> names()         # Names that end with "e"

[1] "visit_date"  "clinic_site" "age"         "pulse"      

ehr |> select(contains("_")) |> names()          # Names containing an underscore

[1] "patient_id"     "visit_date"     "clinic_site"    "family_history"

ehr |> select(where(is.numeric)) |> names()      # Only numeric columns

[1] "age"            "bmi"            "sbp"            "pulse"         
[5] "family_history"

The first line returns patient_id and pulse, the two names starting with “p”. where(is.numeric) tests each column and keeps those for which the answer is TRUE, so the text columns (patient_id, visit_date, clinic_site, sex, diabetes) are left out.

A minus sign drops a column, a:b takes a run of neighboring columns, and everything() stands for all the columns not yet named, which is handy for moving a column to the front.

ehr |> select(-visit_date, -sex) |> names()          # Everything except two columns

[1] "patient_id"     "clinic_site"    "age"            "bmi"           
[5] "sbp"            "pulse"          "family_history" "diabetes"      

ehr |> select(age:pulse) |> names()                  # From age to pulse, inclusive

[1] "age"   "sex"   "bmi"   "sbp"   "pulse"

ehr |> select(diabetes, everything()) |> names()     # diabetes first, then the rest

 [1] "diabetes"       "patient_id"     "visit_date"     "clinic_site"   
 [5] "age"            "sex"            "bmi"            "sbp"           
 [9] "pulse"          "family_history"

The first line lists eight names, because two of the ten columns are dropped. The second line returns age, sex, bmi, sbp and pulse, since they sit side by side in the file. The third starts with diabetes and then lists the other nine in their old order.

4. Choose rows with filter()

adults <- ehr |>
  filter(age >= 18)
north <- ehr |>
  filter(clinic_site == "North")
nrow(north)

[1] 200

Use == when testing a condition. In a real clinic, filtering to adult patients might be part of eligibility criteria, but those criteria should be chosen before looking at which subgroup gives the best results.

filter() accepts several conditions. A comma or & means and (both must be true), | means or, %in% tests membership in a set, and between(x, low, high) is a short way to write x >= low & x <= high.

ehr |> filter(age >= 60 & sbp >= 140) |> nrow()               # both

[1] 55

ehr |> filter(clinic_site == "North" | clinic_site == "South") |> nrow()   # either

[1] 400

ehr |> filter(clinic_site %in% c("North", "South")) |> nrow() # same, shorter

[1] 400

ehr |> filter(between(age, 40, 49)) |> nrow()                 # 40 up to 49

[1] 91

The second and third lines give the same count, because %in% is a compact version of repeated |. filter() keeps only rows where the condition is TRUE, so a row with a missing sbp is dropped from the first line as well.

5. Sort rows and pick rows by position

arrange() sorts rows by a column, smallest first. Wrap the column in desc() for largest first. slice() picks rows by position, slice_head() takes the first rows, and slice_max() and slice_min() take the rows with the largest or smallest value of a column.

ehr |> arrange(desc(age)) |> select(patient_id, age) |> head(3)

# A tibble: 3 x 2
  patient_id   age
  <chr>      <dbl>
1 S0041         83
2 S0107         83
3 S0173         83

ehr |> slice(1:2) |> select(patient_id, age)

# A tibble: 2 x 2
  patient_id   age
  <chr>      <dbl>
1 S0001         55
2 S0002         26

ehr |> slice_head(n = 2) |> select(patient_id, age)

# A tibble: 2 x 2
  patient_id   age
  <chr>      <dbl>
1 S0001         55
2 S0002         26

ehr |> slice_max(sbp, n = 3, with_ties = FALSE) |> select(patient_id, sbp)

# A tibble: 3 x 2
  patient_id   sbp
  <chr>      <dbl>
1 S0230        154
2 S0187        152
3 S0287        151

ehr |> slice_min(bmi, n = 3, with_ties = FALSE) |> select(patient_id, bmi)

# A tibble: 3 x 2
  patient_id   bmi
  <chr>      <dbl>
1 S0490       16.7
2 S0130       16.9
3 S0519       16.9

ehr |> slice_tail(n = 2) |> select(patient_id, age)

# A tibble: 2 x 2
  patient_id   age
  <chr>      <dbl>
1 S0599         71
2 S0600         42

set.seed(1)                                   # Makes the random choice repeatable
ehr |> slice_sample(n = 3) |> select(patient_id, age)

# A tibble: 3 x 2
  patient_id   age
  <chr>      <dbl>
1 S0129         39
2 S0509         41
3 S0471         21

slice_tail(n = 2) returns the last two rows, and slice_sample(n = 3) picks three rows at random; set.seed() fixes the random draw so the same rows come back each time you run the code.

slice(1:2) and slice_head(n = 2) return the same two rows. The slice_max() and slice_min() lines return exactly three rows each because of with_ties = FALSE; without it, rows that tie with the third value are kept as well, so you can get more than three. A missing value sorts last, so it never appears among the top rows here.

6. Add calculated columns with mutate()

ehr2 <- ehr |>
  mutate(
    high_sbp = if_else(!is.na(sbp) & sbp >= 140, TRUE, FALSE, missing = FALSE),
    age_group = if_else(age < 50, "Under 50", "50 or older")
  )
ehr2 |> select(age, sbp, high_sbp, age_group) |> head(6)

# A tibble: 6 x 4
    age   sbp high_sbp age_group  
  <dbl> <dbl> <lgl>    <chr>      
1    55   120 FALSE    50 or older
2    26   110 FALSE    Under 50   
3    63   133 FALSE    50 or older
4    34   120 FALSE    Under 50   
5    71   141 TRUE     50 or older
6    42   126 FALSE    Under 50   

Important caveat: The high_sbp example codes missing SBP as FALSE; for a real clinical flag, missing should usually remain unknown, not “not high.” A safer approach is:

ehr2 <- ehr |>
  mutate(high_sbp = case_when(
    is.na(sbp) ~ NA,
    sbp >= 140 ~ TRUE,
    TRUE ~ FALSE
  ))
table(ehr2$high_sbp, useNA = "ifany")

FALSE  TRUE  <NA> 
  522    58    20 

Derived variables. A column computed from other columns is a derived variable. It can carry more meaning than its parts. The classic example is body mass index: height alone and weight alone say little about body size, but weight divided by height squared compares people of different heights. This dataset already contains bmi. A second example uses two vital signs. The shock index is pulse divided by systolic blood pressure; clinicians use it as a quick warning sign because a fast pulse with a low pressure raises it.

ehr3 <- ehr |>
  mutate(shock_index = round(pulse / sbp, 2))
ehr3 |> select(pulse, sbp, shock_index) |> head(4)

# A tibble: 4 x 3
  pulse   sbp shock_index
  <dbl> <dbl>       <dbl>
1    79   120        0.66
2    92   110        0.84
3    65   133        0.49
4    77   120        0.64

summary(ehr3$shock_index)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.     NAs 
  0.400   0.540   0.615   0.620   0.700   0.990      50 

Each shock_index is one pulse divided by the matching sbp and rounded to two decimals. Where either input is missing the result is NA, and summary() counts those under its NAs heading (50 of them).

Two relatives of mutate():

• transmute() keeps only the new columns. The modern way to write the same thing is mutate(.keep = "none"); transmute() is still available but is no longer the recommended form.
• min_rank() turns a column into ranks, where 1 is the smallest and tied values share the lowest rank. Wrap the column in desc() to rank from largest.

ehr |> transmute(patient_id, shock_index = round(pulse / sbp, 2)) |> head(3)

# A tibble: 3 x 2
  patient_id shock_index
  <chr>            <dbl>
1 S0001             0.66
2 S0002             0.84
3 S0003             0.49

ehr |> mutate(shock_index = round(pulse / sbp, 2), .keep = "none") |> head(3)

# A tibble: 3 x 1
  shock_index
        <dbl>
1        0.66
2        0.84
3        0.49

min_rank(c(10, 20, 20, 5))

[1] 2 3 3 1

The two ehr lines print the same shock-index column; the transmute() one also keeps patient_id because we named it. In min_rank(c(10, 20, 20, 5)), 5 is rank 1, 10 is rank 2, and both 20s share rank 3.

Ranks help inside a table: mutate(sbp_rank = min_rank(desc(sbp))) gives each patient a position from highest pressure down. The filter below asks for ranks 3 or better, yet four rows print: two patients tie at 151 and both get rank 3.

ehr |>
  mutate(sbp_rank = min_rank(desc(sbp))) |>
  filter(sbp_rank <= 3) |>
  select(patient_id, sbp, sbp_rank) |>
  arrange(sbp_rank)

# A tibble: 4 x 3
  patient_id   sbp sbp_rank
  <chr>      <dbl>    <int>
1 S0230        154        1
2 S0187        152        2
3 S0287        151        3
4 S0371        151        3

7. Group and summarize

site_summary <- ehr |>
  group_by(clinic_site) |>
  summarise(
    n_visits = n(),
    mean_age = mean(age),
    missing_bmi = sum(is.na(bmi)),
    mean_sbp = mean(sbp, na.rm = TRUE),
    .groups = "drop"
  )
site_summary

# A tibble: 3 x 5
  clinic_site n_visits mean_age missing_bmi mean_sbp
  <chr>          <int>    <dbl>       <int>    <dbl>
1 Mobile           200     51.5          12     123.
2 North            200     49.5          11     122.
3 South            200     50.4          12     122.

Each group becomes one row. Always say whether your denominator is visits, unique patients, or observed measurements.

summarise() can use many summary functions. Besides mean() and sum(), these come up constantly: n() (rows in the group), n_distinct() (how many different values), median(), sd() (standard deviation), var() (variance, the square of the standard deviation), min() and max(). Most need na.rm = TRUE when the column has missing values.

ehr |>
  group_by(clinic_site) |>
  summarise(
    n = n(),
    n_ids = n_distinct(patient_id),
    median_age = median(age),
    sd_sbp = round(sd(sbp, na.rm = TRUE), 1),
    var_sbp = round(var(sbp, na.rm = TRUE), 1),
    min_bmi = min(bmi, na.rm = TRUE),
    max_bmi = max(bmi, na.rm = TRUE)
  )

# A tibble: 3 x 8
  clinic_site     n n_ids median_age sd_sbp var_sbp min_bmi max_bmi
  <chr>       <int> <int>      <dbl>  <dbl>   <dbl>   <dbl>   <dbl>
1 Mobile        200   200       51.5   13      168     17.5    43.2
2 North         200   200       49.5   13.1    173.    16.9    43.5
3 South         200   200       50.5   12.9    166.    16.7    43.7

In every site n equals n_ids, so no patient ID repeats within a site. That does not rule out the same ID appearing at two different sites; the whole-table check is anyDuplicated() in section 10. var_sbp is sd_sbp squared, up to rounding.

ungroup(). group_by() stays in effect until you remove it. After a mutate() on grouped data, later steps would still work group by group, which is a quiet source of wrong answers. Call ungroup() once the grouped step is done.

centered <- ehr |>
  group_by(clinic_site) |>
  mutate(age_minus_site_mean = age - mean(age)) |>
  ungroup()
is_grouped_df(centered)

[1] FALSE

mean(centered$age_minus_site_mean)

[1] 7.578458e-16

is_grouped_df() returns FALSE: the grouping is gone, so any step after this treats the table as a whole. The mean of the centered ages prints as about 8e-16, which is zero apart from rounding error, because each site was centered on its own mean.

8. Put the verbs together in one pipe

Real questions need several steps. Reading from top to bottom, this pipe filters the rows, adds a derived column, groups, summarizes, and sorts the result.

ehr |>
  filter(!is.na(sbp), !is.na(pulse), age >= 18) |>
  mutate(shock_index = pulse / sbp) |>
  group_by(clinic_site, sex) |>
  summarise(
    n = n(),
    mean_shock_index = round(mean(shock_index), 3),
    .groups = "drop"
  ) |>
  filter(n >= 100) |>
  arrange(desc(mean_shock_index))

# A tibble: 3 x 4
  clinic_site sex       n mean_shock_index
  <chr>       <chr> <int>            <dbl>
1 Mobile      F       110            0.723
2 South       F       109            0.621
3 North       F       111            0.516

Without the filter(n >= 100) line the result would have one row for every combination of site and sex. A filter() after summarise() works on the summary table, so here it keeps only the groups with at least 100 visits (n counts rows, and each row is one visit). The rows are ordered from the largest mean shock index to the smallest. Only rows with both sbp and pulse present are counted in n, which is a decision you should state when you report the numbers.

9. Tidy data: the shape a table should have

Many functions, including those that make plots, expect tidy data:

1. each variable is a column,
2. each observation is a row,
3. each value is a cell.

Spreadsheets are often laid out for people instead, with one column per measurement time. This small table is wide: the variable “visit” is hidden in the column names, and the variable “sbp” is spread over two columns.

messy <- data.frame(
  patient = c("A", "B", "C"),
  sbp_visit1 = c(130, 142, 118),
  sbp_visit2 = c(128, 150, 121)
)
messy

  patient sbp_visit1 sbp_visit2
1       A        130        128
2       B        142        150
3       C        118        121

pivot_longer() from tidyr turns the two sbp_ columns into one pair of columns: one holds the old column name (visit), the other holds the value (sbp).

tidy <- messy |>
  pivot_longer(
    cols = starts_with("sbp_"),
    names_to = "visit",
    names_prefix = "sbp_",
    values_to = "sbp"
  )
tidy

# A tibble: 6 x 3
  patient visit    sbp
  <chr>   <chr>  <dbl>
1 A       visit1   130
2 A       visit2   128
3 B       visit1   142
4 B       visit2   150
5 C       visit1   118
6 C       visit2   121

The table now has six rows, one per patient per visit, and each row holds exactly one blood pressure. pivot_wider() does the reverse, which is useful for a report that needs one column per visit:

tidy |>
  pivot_wider(names_from = visit, values_from = sbp, names_prefix = "sbp_")

# A tibble: 3 x 3
  patient sbp_visit1 sbp_visit2
  <chr>        <dbl>      <dbl>
1 A              130        128
2 B              142        150
3 C              118        121

This gives back the three-row table we started with.

10. Check obvious data problems

stopifnot(nrow(ehr) == 600)
stopifnot(!anyDuplicated(ehr$patient_id))
table(is.na(ehr$bmi))

FALSE  TRUE 
  565    35 

summary(ehr$sbp)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.     NAs 
   87.0   113.0   122.0   122.2   131.0   154.0      20 

table(is.na(ehr$bmi)) shows 565 FALSE and 35 TRUE, so 35 visits have no BMI, and summary(ehr$sbp) reports 20 NAs.

stopifnot() stops if a condition is not TRUE. It’s useful when you expect something essential (for example, one row per patient) and want your program to fail loudly when that assumption is broken.

You may encounter a literal string "NA" in a bad export; that is text, not necessarily R’s actual missing value. Inspect types and source documentation before deciding how to convert it.

ehr |> filter(!is.na(sbp), sbp > 140) |> nrow()

[1] 51

ehr |>
  group_by(clinic_site) |>
  summarise(median_pulse = median(pulse, na.rm = TRUE), .groups = "drop")

# A tibble: 3 x 2
  clinic_site median_pulse
  <chr>              <dbl>
1 Mobile                88
2 North                 62
3 South                 75

ehr |>
  group_by(clinic_site) |>
  summarise(pct_missing_bmi = round(100 * mean(is.na(bmi)), 1), .groups = "drop")

# A tibble: 3 x 2
  clinic_site pct_missing_bmi
  <chr>                 <dbl>
1 Mobile                  6  
2 North                   5.5
3 South                   6  

ehr |>
  mutate(flag_false = if_else(is.na(sbp), FALSE, sbp >= 140),
         flag_na = if_else(is.na(sbp), NA, sbp >= 140)) |>
  summarise(false_coding = sum(flag_false), na_coding = sum(flag_na, na.rm = TRUE),
            unknown = sum(is.na(flag_na)))

# A tibble: 1 x 3
  false_coding na_coding unknown
         <int>     <int>   <int>
1           58        58      20

analysis_table <- ehr |>
  select(patient_id, age, bmi, sbp, pulse, family_history, diabetes)
dim(analysis_table)

[1] 600   7

names(analysis_table)

[1] "patient_id"     "age"            "bmi"            "sbp"           
[5] "pulse"          "family_history" "diabetes"      

ehr |> select(where(is.numeric)) |> names()

[1] "age"            "bmi"            "sbp"            "pulse"         
[5] "family_history"

ehr |> select(starts_with("clinic"), ends_with("_id")) |> names()

[1] "clinic_site" "patient_id" 

ehr |> filter(clinic_site %in% c("Mobile", "North"), age >= 65) |> nrow()

[1] 118

ehr |> filter(between(age, 30, 39)) |> nrow()

[1] 91

ehr |> arrange(desc(age)) |> select(patient_id, age) |> head(3)

# A tibble: 3 x 2
  patient_id   age
  <chr>      <dbl>
1 S0041         83
2 S0107         83
3 S0173         83

ehr |> slice_min(bmi, n = 3, with_ties = FALSE) |> select(patient_id, bmi)

# A tibble: 3 x 2
  patient_id   bmi
  <chr>      <dbl>
1 S0490       16.7
2 S0130       16.9
3 S0519       16.9

ranked <- ehr |>
  mutate(shock_index = pulse / sbp,
         si_rank = min_rank(desc(shock_index))) |>
  filter(si_rank == 1) |>
  select(patient_id, pulse, sbp, shock_index, si_rank)
ranked

# A tibble: 1 x 5
  patient_id pulse   sbp shock_index si_rank
  <chr>      <dbl> <dbl>       <dbl>   <int>
1 S0521         91    92       0.989       1

ehr |>
  group_by(sex) |>
  summarise(n = n(),
            n_sites = n_distinct(clinic_site),
            median_age = median(age),
            sd_bmi = round(sd(bmi, na.rm = TRUE), 2),
            .groups = "drop")

# A tibble: 2 x 5
  sex       n n_sites median_age sd_bmi
  <chr> <int>   <int>      <dbl>  <dbl>
1 F       360       3       50.5   7.1 
2 M       240       3       50.5   7.11

ehr |>
  filter(diabetes == "Yes", !is.na(pulse), !is.na(sbp)) |>
  mutate(shock_index = pulse / sbp) |>
  group_by(clinic_site) |>
  summarise(mean_shock_index = round(mean(shock_index), 3), .groups = "drop")

# A tibble: 3 x 2
  clinic_site mean_shock_index
  <chr>                  <dbl>
1 Mobile                 0.649
2 North                  0.468
3 South                  0.556

wide <- data.frame(id = 1:2, t1 = c(5, 6), t2 = c(7, 8))
long <- wide |> pivot_longer(cols = c(t1, t2), names_to = "time", values_to = "value")
long

# A tibble: 4 x 3
     id time  value
  <int> <chr> <dbl>
1     1 t1        5
2     1 t2        7
3     2 t1        6
4     2 t2        8

long |> pivot_wider(names_from = time, values_from = value)

# A tibble: 2 x 3
     id    t1    t2
  <int> <dbl> <dbl>
1     1     5     7
2     2     6     8