Prepare Inputs Without Leakage

1. Different models need different input forms

Some algorithms cannot handle missing predictors. Some cannot directly work with text categories. Preprocessing may include imputation (filling missing predictors with estimates), encoding categories, rescaling, and removing constant columns. This is a modeling decision, not a license to invent measured facts.

Rule: A preprocessing statistic (like the training median) must be learned only on training data. During CV it must be learned within each analysis fold, not on all training rows before folding.

Example: take five typed systolic readings where one is missing. A model that cannot handle NA needs a number there. The median of the four known values is the simple choice.

sbp_small <- c(118, 125, NA, 140, 132)
median(sbp_small, na.rm = TRUE)

[1] 128.5

ifelse(is.na(sbp_small), median(sbp_small, na.rm = TRUE), sbp_small)

[1] 118.0 125.0 128.5 140.0 132.0

Reading the output: the median of 118, 125, 132 and 140 is 128.5, and the NA is replaced by that value while the others stay as they were. In a real project that median would be calculated from the training rows only.

2. A recipe is a plan of transformations

The recipes package (part of tidymodels) describes preprocessing. A workflow attaches this plan to a model so each resample follows the same leak-safe procedure.

library(tidymodels)
ehr <- readr::read_csv("data/synthetic_ehr.csv", show_col_types = FALSE) |>
  mutate(
    diabetes = factor(diabetes, levels = c("Yes", "No")),
    sex = factor(sex),
    clinic_site = factor(clinic_site),
    family_history = factor(family_history)
  )
set.seed(500)
split <- initial_split(ehr, prop = 0.8, strata = diabetes)
train <- training(split)
test <- testing(split)

Example: check what the split gave us before building the recipe.

c(train_rows = nrow(train), test_rows = nrow(test))

train_rows  test_rows 
       480        120 

colSums(is.na(train[c("age", "bmi", "sbp", "pulse")]))

  age   bmi   sbp pulse 
    0    27    12    23 

Reading the output: train has 480 rows and test has 120. In train, age has no missing values, while bmi, sbp and pulse each have some. Those are the gaps the recipe in the next section will fill.

3. Define exactly which features are allowed

rec <- recipe(
  diabetes ~ age + bmi + sbp + pulse + family_history + sex + clinic_site,
  data = train
) |>
  step_impute_median(all_numeric_predictors()) |>
  step_impute_mode(all_nominal_predictors()) |>
  step_dummy(all_nominal_predictors()) |>
  step_zv(all_predictors())

summary(rec)

# A tibble: 8 x 4
  variable       type      role      source  
  <chr>          <list>    <chr>     <chr>   
1 age            <chr [2]> predictor original
2 bmi            <chr [2]> predictor original
3 sbp            <chr [2]> predictor original
4 pulse          <chr [2]> predictor original
5 family_history <chr [3]> predictor original
6 sex            <chr [3]> predictor original
7 clinic_site    <chr [3]> predictor original
8 diabetes       <chr [3]> outcome   original

Reading the output: summary() lists the variables the recipe knows about and the role of each: diabetes is the outcome and the seven allowed columns are predictors. The steps have only been written down; nothing is calculated until prep() in the next section.

Line by line:

1. The formula uses outcome ~ allowed predictors. “~” means “model outcome using these inputs.”
2. step_impute_median() replaces missing numeric predictor values using a median learned in the relevant training subset.
3. step_impute_mode() uses the most frequent value for missing categorical predictors.
4. step_dummy() turns categorical levels into numeric indicator columns.
5. step_zv() removes variables with no variation, which cannot help distinguish outcomes.

We did not include patient_id, visit_date, or the outcome itself among the inputs. Missing outcome labels are fundamentally different from missing predictors and require separate research decisions.

4. See how training transformations are learned

trained_recipe <- prep(rec, training = train)
processed_train <- bake(trained_recipe, new_data = train)
processed_test <- bake(trained_recipe, new_data = test)
names(processed_train)

[1] "age"               "bmi"               "sbp"              
[4] "pulse"             "diabetes"          "family_history_X1"
[7] "sex_M"             "clinic_site_North" "clinic_site_South"

colSums(is.na(processed_train))

              age               bmi               sbp             pulse 
                0                 0                 0                 0 
         diabetes family_history_X1             sex_M clinic_site_North 
                0                 0                 0                 0 
clinic_site_South 
                0 

The trained_recipe object holds preprocessing rules learned from train; bake applies those same rules elsewhere. The test-set values never chose the imputation medians. Note that this explicit prep() example is explanatory; in later CV, workflows re-prep the recipe inside each fold.

5. Should we always impute?

No. Median imputation is a simple baseline with assumptions. Missingness may be informative or systematic; a single median can distort variation and subgroup behavior. Compare strategies on training folds, document missingness by site/group, and consider missingness indicator features when justified. Methods like MICE require careful outcome/time handling and resampling integration; do not run a global imputation on the entire dataset before splitting.

Example: look at where BMI is missing before deciding that a single median is fine.

train |>
  group_by(clinic_site) |>
  summarise(n = n(), missing_bmi = sum(is.na(bmi)),
            share_missing = round(mean(is.na(bmi)), 3))

# A tibble: 3 x 4
  clinic_site     n missing_bmi share_missing
  <fct>       <int>       <int>         <dbl>
1 Mobile        166           9         0.054
2 North         163           9         0.055
3 South         151           9         0.06 

Reading the output: the table gives, for each site, how many training rows it has and how many have no BMI. Here the shares are 5.4%, 5.5% and 6.0%, so missingness is about the same at every site and one overall median is a reasonable simple choice. If one site had a much higher share missing than the others, the missingness would be tied to the site, and one overall median could distort that site. Roughly equal shares would be less worrying, though still worth documenting.

6. Scaling and outliers

Some algorithms depend strongly on variable scales (K-means, PCA, regularized regression); many tree models do not. We have not added step_normalize() because it is not necessary for our basic logistic regression or tree examples. If needed, add it inside the recipe, with its parameters learned on training folds.

Example: add scaling to a copy of the recipe and look at the result for age.

rec_scaled <- rec |> step_normalize(all_numeric_predictors())
scaled_train <- bake(prep(rec_scaled, training = train), new_data = train)
summary(train$age)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  18.00   34.00   50.00   50.17   66.25   83.00 

summary(scaled_train$age)

     Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
-1.697169 -0.853036 -0.008903  0.000000  0.848419  1.732121 

Reading the output: the first summary is in years (roughly 18 to 83). The second is in standard deviations, centered at 0, which is the scale that K-means or regularized regression would see. The ordering of patients has not changed; only the units have.

Never remove clinically unusual measurements simply because they do not match an average. Check units, plausible ranges, and source records with authorized staff.

7. Skewed variables and transformations

A variable is right-skewed when most values sit at the low end and a long tail stretches to the high end. Lab values such as triglycerides, costs, and counts of visits often look like this. Methods that work with means and variances (linear regression, PCA, K-means) can be pulled around by a few very large values, so a transformation that pulls the tail in is sometimes worth trying.

The synthetic ehr table has no strongly skewed column, so here we make up a small lab table of triglyceride values (mg/dL) with a long right tail.

set.seed(12)
lab <- data.frame(triglycerides = round(rlnorm(400, meanlog = log(120), sdlog = 0.6)))
par(mfrow = c(1, 2))
hist(lab$triglycerides, main = "Raw", xlab = "Triglycerides")
hist(log(lab$triglycerides), main = "After log()", xlab = "log(triglycerides)")

↔ Swipe the plot sideways to see all of it.

par(mfrow = c(1, 1))

The left histogram leans to the left with a long right tail; the right one, after log(), is close to symmetric.

Skewness is a single number for this lean. Standardize the variable (subtract the mean, divide by the standard deviation), cube each value, and take the mean. A value near 0 means symmetric, a clearly positive value means a right tail, and a negative value means a left tail.

skew <- function(x) mean(((x - mean(x)) / sd(x))^3)
c(raw = skew(lab$triglycerides),
  log = skew(log(lab$triglycerides)))

       raw        log 
 1.6334226 -0.0134994 

The raw values give a clearly positive number (about 1.6 on this seed) and the logged values give a number near 0. The e1071 package has e1071::skewness(). Its default (type = 3) gives the same number as this formula; type = 1 divides by n in the standard deviation and gives a slightly different number.

log() works only for values above 0. If zeros are possible (counts, for example), log1p(x), which is log(1 + x), is the usual fix. If the variable can be zero or negative, or you do not want to choose the transformation yourself, step_YeoJohnson() in a recipe picks the best power from the training data. Like every step in this chapter, it learns its setting (here a power called lambda) from the training rows only and then applies the same setting to new rows.

idx <- sample(nrow(lab), 300)
lab_train <- lab[idx, , drop = FALSE]
lab_test  <- lab[-idx, , drop = FALSE]

rec_yj <- recipe(~ triglycerides, data = lab_train) |>
  step_YeoJohnson(triglycerides)
prepped_yj <- prep(rec_yj, training = lab_train)
tidy(prepped_yj, number = 1)            # the lambda learned from training rows

# A tibble: 1 x 3
  terms          value id              
  <chr>          <dbl> <chr>           
1 triglycerides 0.0373 YeoJohnson_GTxr1

test_yj <- bake(prepped_yj, new_data = lab_test)
c(before = skew(lab_test$triglycerides),
  after  = skew(test_yj$triglycerides))

   before     after 
1.7880443 0.2977128 

The value column is the lambda learned from lab_train. Applied to lab_test, the same lambda reduces the skewness from about 1.8 to about 0.3, although the test rows did not choose it.

When not to transform. Leave a variable alone when it is already roughly symmetric (like age, bmi, sbp and pulse here), when the model does not care (trees split on order, so any monotone transform gives the same tree), or when you need results in the original units, such as a coefficient per mg/dL. A transformed variable is also harder to explain to clinicians: say clearly when a result is on the log scale. And do not log a variable to hide a few extreme values that are really data-entry errors; check those first.

c(train_missing_bmi = sum(is.na(train$bmi)),
  test_missing_bmi  = sum(is.na(test$bmi)),
  train_missing_sbp = sum(is.na(train$sbp)),
  test_missing_sbp  = sum(is.na(test$sbp)))

train_missing_bmi  test_missing_bmi train_missing_sbp  test_missing_sbp 
               27                 8                12                 8 

rec_norm <- rec |>
  step_normalize(all_numeric_predictors())
prepped_norm <- prep(rec_norm, training = train)
norm_train <- bake(prepped_norm, new_data = train)
norm_test  <- bake(prepped_norm, new_data = test)

c(train_mean = mean(norm_train$age), train_sd = sd(norm_train$age),
  test_mean = mean(norm_test$age), test_sd = sd(norm_test$age))

   train_mean      train_sd     test_mean       test_sd 
-1.502302e-16  1.000000e+00  8.210510e-02  1.023569e+00 

c(missing_bmi = sum(is.na(ehr$bmi)),
  missing_outcome = sum(is.na(ehr$diabetes)))

    missing_bmi missing_outcome 
             35               0 

names(processed_train)

[1] "age"               "bmi"               "sbp"              
[4] "pulse"             "diabetes"          "family_history_X1"
[7] "sex_M"             "clinic_site_North" "clinic_site_South"

grep("^(sex|clinic_site|family_history)_", names(processed_train), value = TRUE)

[1] "family_history_X1" "sex_M"             "clinic_site_North"
[4] "clinic_site_South"

c(full_data_median = median(ehr$bmi, na.rm = TRUE),
  train_only_median = median(train$bmi, na.rm = TRUE))

 full_data_median train_only_median 
             30.2              30.5 

skew <- function(x) mean(((x - mean(x)) / sd(x))^3)
c(age = skew(na.omit(train$age)), bmi = skew(na.omit(train$bmi)))

        age         bmi 
0.021465863 0.002659681 

set.seed(7)
v <- data.frame(v = rexp(300))
v_train <- v[1:200, , drop = FALSE]
v_test  <- v[201:300, , drop = FALSE]

yj <- prep(recipe(~ v, data = v_train) |> step_YeoJohnson(v), training = v_train)
tidy(yj, number = 1)

# A tibble: 1 x 3
  terms  value id              
  <chr>  <dbl> <chr>           
1 v     -0.820 YeoJohnson_06BjG

c(raw = skew(v_test$v),
  log1p = skew(log1p(v_test$v)),
  yeo_johnson = skew(bake(yj, new_data = v_test)$v))

        raw       log1p yeo_johnson 
 1.29300663  0.50322467 -0.09373028