Skip to main content

Pre-processing MEDS Data

Before building models, we'll pre-process the MIMIC-IV demo with MEDS-Transforms. We'll write two short pipeline configurations and run them: the first cleans the data (it filters subjects, filters rare codes and removes outlying numeric values), and the second normalizes the cleaned data and maps its codes to integer indices for a neural network. Aggregate counts after each stage show what it changed. A short optional section shows MEDS-Reader. See the notebook below, or open it on Google Colab or on our GitHub Repository.

note

This step was added after the ML4H 2025 event, so it is not part of the event's tutorial repository. Its notebook is kept in this website's repository. Like steps 3, 5 and 6, it starts by downloading the MIMIC-IV demo in MEDS format, and the model steps that follow start from that download too, so you can skip this step.

Pre-processing MEDS Data

Open In Colab

In the previous steps, we converted the MIMIC-IV demo to MEDS and extracted labels for our task: will a patient's ICU stay last longer than 3 days? Before training a model on MEDS data, most projects run the same kinds of pre-processing: drop subjects with too little data, drop codes that are too rare to learn from, remove implausible numeric values, and normalize what remains.

In this notebook, we'll do all of this with MEDS-Transforms, the MEDS ecosystem's tool for building pre-processing pipelines out of reusable, tested stages. Instead of writing pre-processing code, we'll describe a pipeline in a short YAML file and run it with one command. We'll:

  1. Download the MIMIC-IV demo in MEDS format (the same files the other steps use).
  2. Write and run a pipeline that cleans the data, and check what each stage changed.
  3. Write and run a second pipeline that normalizes the data and maps codes to integers, ready for a neural network.
  4. See how this relates to the pre-processing in the neural network step, and take a short look at MEDS-Reader.

Setup

We'll pin the versions this notebook was run with. MEDS-Transforms needs Python 3.11 or newer. MEDS-Reader is only used in the optional section at the end.

[1]: 
1 
%pip install --quiet "MEDS-transforms==0.7.1" "meds-reader==0.2.0"
Note: you may need to restart the kernel to use updated packages.
[2]: 
1 
2
3
4
from importlib.metadata import version

for package in ["MEDS-transforms", "meds", "meds-reader", "polars"]:
    print(f"{package}: {version(package)}")
MEDS-transforms: 0.7.1
meds: 0.4.1
meds-reader: 0.2.0
polars: 2.1.0

Downloading the data

As in the other steps, we download the MIMIC-IV demo in MEDS format from PhysioNet. It is open access, so no credentials are needed.

[3]: 
1 
2
3
4
5
6
%%bash
wget -r -N -c -np -q -nH --cut-dirs=3 -P MEDS_data https://physionet.org/files/mimic-iv-demo-meds/0.0.1/
rm MEDS_data/*.* # We just need directories
find MEDS_data -name "index.html" -type f -delete # wget insists on these.

find MEDS_data -type f | sort
MEDS_data/data/held_out/0.parquet
MEDS_data/data/train/0.parquet
MEDS_data/data/tuning/0.parquet
MEDS_data/metadata/codes.parquet
MEDS_data/metadata/dataset.json
MEDS_data/metadata/subject_splits.parquet

Recall the layout of a MEDS dataset: a data/ folder of parquet shards, here already split into train, tuning and held_out, and a metadata/ folder with the code metadata (codes.parquet), the dataset metadata (dataset.json) and the assignment of subjects to splits (subject_splits.parquet). Every row of data/ is one measurement: a subject_id, a time (null for static measurements), a code, and optionally a numeric_value, plus any extra columns the ETL kept.

The demo has 100 subjects, split as follows:

[4]: 
1 
2
3
import polars as pl

pl.read_parquet("MEDS_data/metadata/subject_splits.parquet").group_by("split").len().sort("split")
[4]: 
shape: (3, 2)
splitlen
stru32
"held_out"10
"train"80
"tuning"10

To see what each pre-processing stage does, we'll summarize a data folder with a few counts: the number of subjects, measurements (rows), distinct codes, and measurements with a numeric value.

[5]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
from pathlib import Path


def summarize(data_dir: str) -> dict:
    """Counts subjects, measurements, distinct codes and numeric values over all splits of a data folder."""
    df = pl.scan_parquet(Path(data_dir) / "*" / "*.parquet")  # One sub-folder per split
    counts = df.select(
        pl.col("subject_id").n_unique().alias("subjects"),
        pl.len().alias("measurements"),
        pl.col("code").n_unique().alias("distinct codes"),
        pl.col("numeric_value").is_not_null().sum().alias("numeric values"),
    )
    return counts.collect().row(0, named=True)


summaries = {"downloaded data": summarize("MEDS_data/data")}
pl.DataFrame([{"data": name, **counts} for name, counts in summaries.items()])
[5]: 
shape: (1, 5)
datasubjectsmeasurementsdistinct codesnumeric values
stri64i64i64i64
"downloaded data"1009161667035405236

MEDS-Transforms pipelines

A MEDS-Transforms pipeline is a YAML file with an input_dir (a MEDS dataset), an output_dir, and an ordered list of stages. Each stage reads MEDS data and writes MEDS data, and comes in one of two kinds:

  • Data stages change the rows in data/, for example by filtering measurements. They run on every split.
  • Metadata stages compute statistics over the data and write them to metadata/codes.parquet, for example how often each code occurs. In a pipeline, they are fit on the train split only, so statistics such as means and standard deviations never see the tuning or held-out subjects. This is what keeps pre-processing from leaking information from your evaluation data into your model.

A data stage uses the metadata written by the most recent metadata stage before it, so pipelines follow a "fit a statistic, then apply it" pattern. You can also run the same stage more than once under different names by setting _base_stage; we'll use this to compute different statistics with the aggregate_code_metadata stage. See the MEDS-Transforms documentation for the full list of built-in stages.

Pipeline 1: cleaning the data

Our first pipeline keeps the data in MEDS format (codes stay strings, and extra columns are kept), so its output can be used by any MEDS tool. It has five stages:

  1. filter_subjects drops subjects with fewer than 10 events. In MEDS-Transforms, an event is a unique timestamp of a subject, while a measurement is a row; a subject with 3 lab results at the same time has 3 measurements but 1 event. A subject with only a handful of events gives a model little to learn from.
  2. count_codes is an aggregate_code_metadata stage that counts, for each code, how often it occurs (code/n_occurrences) and in how many subjects (code/n_subjects) in the train split.
  3. filter_measurements drops every measurement whose code occurs in fewer than 2 training subjects. A code that only one training subject has can't teach a model anything that carries over to other patients, and each extra code is one more embedding or feature to learn. Codes that never occur in the train split have no counts at all, so their measurements are dropped from every split.
  4. fit_outlier_detection is another aggregate_code_metadata stage. It computes, for each code, the number, sum and sum of squares of its numeric values in the train split, from which the mean and standard deviation follow.
  5. occlude_outliers removes numeric values more than 4.5 standard deviations from their code's mean. Occluding keeps the measurement but sets its numeric_value to null, and a new column, numeric_value/is_inlier, records which values were kept (true) and which were removed (false).

The thresholds (10 events, 2 subjects, 4.5 standard deviations) are choices, not rules: 4.5 is the stage's default, and the other two are deliberately mild for a dataset of only 80 training subjects. On your own data, choose them for your dataset and model.

[6]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
clean_pipeline = """
input_dir: MEDS_data
output_dir: cleaned_MEDS_data

description: Clean the MIMIC-IV demo in MEDS format

stages:
  - filter_subjects:
      min_events_per_subject: 10
  - count_codes:
      _base_stage: aggregate_code_metadata
      aggregations:
        - code/n_occurrences
        - code/n_subjects
  - filter_measurements:
      min_subjects_per_code: 2
  - fit_outlier_detection:
      _base_stage: aggregate_code_metadata
      aggregations:
        - values/n_occurrences
        - values/sum
        - values/sum_sqd
  - occlude_outliers:
      stddev_cutoff: 4.5
"""

_ = Path("clean_pipeline.yaml").write_text(clean_pipeline)

Running the whole pipeline is a single command:

[7]: 
1 
2
%%bash
MEDS_transform-pipeline clean_pipeline.yaml

The runner writes the output of each stage to its own sub-folder of output_dir, and puts the output of the last data stage in data/ and of the last metadata stage in metadata/. Keeping every stage's output lets the runner resume after a failure, and lets us see exactly what each stage did. (The runner also writes logs to hidden .logs folders, which we skip here.)

[8]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
def show_tree(root: str):
    """Prints the files and folders under `root`, skipping the hidden log folders."""
    root = Path(root)
    print(f"{root}/")
    for path in sorted(root.rglob("*")):
        relative = path.relative_to(root)
        if any(part.startswith(".") for part in relative.parts):
            continue
        indent = "    " * (len(relative.parts) - 1)
        print(f"{indent}    {path.name}{'/' if path.is_dir() else ''}")


show_tree("cleaned_MEDS_data")
cleaned_MEDS_data/
    count_codes/
        codes.parquet
        train/
            0.parquet
    data/
        held_out/
            0.parquet
        train/
            0.parquet
        tuning/
            0.parquet
    filter_measurements/
        held_out/
            0.parquet
        train/
            0.parquet
        tuning/
            0.parquet
    filter_subjects/
        held_out/
            0.parquet
        train/
            0.parquet
        tuning/
            0.parquet
    fit_outlier_detection/
        train/
            0.parquet
    metadata/
        codes.parquet

What changed?

Let's summarize the data after each data stage:

[9]: 
1 
2
3
4
5
summaries["after filter_subjects"] = summarize("cleaned_MEDS_data/filter_subjects")
summaries["after filter_measurements"] = summarize("cleaned_MEDS_data/filter_measurements")
summaries["after occlude_outliers"] = summarize("cleaned_MEDS_data/data")

pl.DataFrame([{"data": name, **counts} for name, counts in summaries.items()])
[9]: 
shape: (4, 5)
datasubjectsmeasurementsdistinct codesnumeric values
stri64i64i64i64
"downloaded data"1009161667035405236
"after filter_subjects"1009161667035405236
"after filter_measurements"1009051753854402540
"after occlude_outliers"1009051753854401186

Each stage did something different:

  • filter_subjects changed nothing. Every subject in the demo has a long record (the shortest has 86 events), so none falls below 10 events. On a full EHR dataset, which includes many patients with only a few visits, the same stage usually does remove subjects.
  • filter_measurements removed nearly half of the distinct codes but only about 1% of the measurements. This is typical of EHR data: most codes are rare, and a few codes (vital signs, common labs) make up most of the rows. Dropping the rare codes shrinks the vocabulary a model has to learn a lot, at a small cost in data.
  • occlude_outliers removed a small share of the numeric values, while keeping every measurement.

Let's look more closely at the last two.

Rare codes

The counts that filter_measurements used are in the metadata written by count_codes. They only cover the train split, so codes that only occur in the tuning or held-out split are not in it:

[10]: 
1 
2
3
4
5
6
7
8
9
code_counts = pl.read_parquet("cleaned_MEDS_data/count_codes/codes.parquet")
train_codes = set(code_counts["code"])
eval_codes = set(
    pl.scan_parquet("MEDS_data/data/tuning/*.parquet").select("code").collect()["code"]
) | set(pl.scan_parquet("MEDS_data/data/held_out/*.parquet").select("code").collect()["code"])

print(f"Codes counted in the train split: {len(train_codes):,}")
print(f"  ... of which in fewer than 2 training subjects: {code_counts.filter(pl.col('code/n_subjects') < 2).height:,}")
print(f"Codes in the tuning or held-out split but not in the train split: {len(eval_codes - train_codes):,}")
Codes counted in the train split: 6,619
  ... of which in fewer than 2 training subjects: 2,765
Codes in the tuning or held-out split but not in the train split: 416

Both groups were dropped. Dropping codes that only occur outside the train split is also the right thing for evaluation: a model trained on the train split has never seen them, so it has nothing to do with them anyway.

Outliers

The statistics that occlude_outliers used are in the final metadata/codes.parquet, from which we can compute each code's mean, standard deviation and the range of values that were kept. As an example, let's look at heart rate, which this dataset records with the code LAB//220045//bpm (MIMIC-IV item 220045, in beats per minute):

[11]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
HEART_RATE = "LAB//220045//bpm"
CUTOFF = 4.5

stats = (
    pl.read_parquet("cleaned_MEDS_data/metadata/codes.parquet")
    .filter(pl.col("code") == HEART_RATE)
    .select(
        "code",
        (pl.col("values/sum") / pl.col("values/n_occurrences")).alias("mean"),
        (
            (pl.col("values/sum_sqd") / pl.col("values/n_occurrences"))
            - (pl.col("values/sum") / pl.col("values/n_occurrences")) ** 2
        )
        .sqrt()
        .alias("std"),
    )
    .with_columns(
        (pl.col("mean") - CUTOFF * pl.col("std")).alias("lowest kept"),
        (pl.col("mean") + CUTOFF * pl.col("std")).alias("highest kept"),
    )
)
stats
[11]: 
shape: (1, 5)
codemeanstdlowest kepthighest kept
strf32f32f32f32
"LAB//220045//bpm"91.60488119.0802175.743904177.465851
[12]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
def value_summary(data_dir: str, code: str) -> dict:
    """Summarizes the numeric values of one code over all splits of a data folder."""
    values = pl.scan_parquet(Path(data_dir) / "*" / "*.parquet").filter(pl.col("code") == code)
    return (
        values.select(
            pl.len().alias("measurements"),
            pl.col("numeric_value").is_not_null().sum().alias("numeric values"),
            pl.col("numeric_value").min().alias("min"),
            pl.col("numeric_value").max().alias("max"),
        )
        .collect()
        .row(0, named=True)
    )


pl.DataFrame(
    [
        {"data": "downloaded data", **value_summary("MEDS_data/data", HEART_RATE)},
        {"data": "after occlude_outliers", **value_summary("cleaned_MEDS_data/data", HEART_RATE)},
    ]
)
[12]: 
shape: (2, 5)
datameasurementsnumeric valuesminmax
stri64i64f64f64
"downloaded data"13913139130.0200.0
"after occlude_outliers"139131390925.0174.0

The downloaded data has heart rates of 0 bpm, which are almost certainly recording artifacts rather than real measurements. They and a few very high values fall outside the kept range, so their values were removed, while the measurements themselves (and so the fact that a heart rate was charted at that time) remain. Over the whole dataset, the new numeric_value/is_inlier column shows how many values were kept and removed; it is null for measurements that have no numeric value:

[13]: 
1 
2
3
4
5
6
7
(
    pl.scan_parquet("cleaned_MEDS_data/data/*/*.parquet")
    .group_by("numeric_value/is_inlier")
    .len()
    .sort("numeric_value/is_inlier", nulls_last=True)
    .collect()
)
[13]: 
shape: (3, 2)
numeric_value/is_inlierlen
boolu32
false1354
true401186
null502635

A cut-off in standard deviations is a blunt tool: it assumes values are roughly symmetric around their mean, and it can't know clinically plausible ranges. For a real project, look at the removed values for the codes you care about, or replace this stage with clinical range checks.

Pipeline 2: preparing the data for a neural network

The cleaned data is still a general-purpose MEDS dataset. A neural network needs two more things: codes as integer indices (to look up embeddings) and numeric values on a common scale. Our second pipeline does this in three stages:

  1. fit_normalization computes code counts and the mean and standard deviation of each code's numeric values on the train split, again. We fit these on the cleaned data, because removing outliers changes the mean and standard deviation.
  2. fit_vocabulary_indices assigns each code an integer index (code/vocab_index in codes.parquet), in alphabetical order of the codes.
  3. normalization replaces each code with its index and each numeric value with its z-score, (value - mean) / std.

The output keeps only the subject_id, time, code and numeric_value columns, with an integer code. It is no longer a standard MEDS dataset, which requires string codes, so it is meant as input to a model, not to other MEDS tools.

[14]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
model_ready_pipeline = """
input_dir: cleaned_MEDS_data
output_dir: model_ready_MEDS_data

description: Normalize the cleaned MIMIC-IV demo and map its codes to vocabulary indices

stages:
  - fit_normalization:
      _base_stage: aggregate_code_metadata
      aggregations:
        - code/n_occurrences
        - code/n_subjects
        - values/n_occurrences
        - values/sum
        - values/sum_sqd
  - fit_vocabulary_indices
  - normalization
"""

_ = Path("model_ready_pipeline.yaml").write_text(model_ready_pipeline)
[15]: 
1 
2
%%bash
MEDS_transform-pipeline model_ready_pipeline.yaml

Let's check the result. The data now has integer codes, and the vocabulary has one index per code that survived cleaning:

[16]: 
1 
2
3
4
5
model_ready = pl.scan_parquet("model_ready_MEDS_data/data/*/*.parquet")
print(model_ready.collect_schema())

vocabulary = pl.read_parquet("model_ready_MEDS_data/metadata/codes.parquet", columns=["code", "code/vocab_index"])
print(f"Vocabulary size: {vocabulary.height:,} codes, with indices {vocabulary['code/vocab_index'].min()} to {vocabulary['code/vocab_index'].max()}")
Schema({'subject_id': Int64, 'time': Datetime(time_unit='us', time_zone=None), 'code': UInt16, 'numeric_value': Float32})
Vocabulary size: 3,854 codes, with indices 1 to 3854

And the heart rates are now z-scores. In the train split, where the statistics were fit, they have a mean of (almost exactly) 0 and a standard deviation of 1:

[17]: 
1 
2
3
4
5
6
7
8
9
10
11
12
heart_rate_index = vocabulary.filter(pl.col("code") == HEART_RATE)["code/vocab_index"].item()

(
    pl.scan_parquet("model_ready_MEDS_data/data/train/*.parquet")
    .filter(pl.col("code") == heart_rate_index)
    .select(
        pl.lit(HEART_RATE).alias("code"),
        pl.col("numeric_value").mean().alias("mean"),
        pl.col("numeric_value").std().alias("std"),
    )
    .collect()
)
[17]: 
shape: (1, 3)
codemeanstd
strf32f32
"LAB//220045//bpm"8.5494e-81.000041

How this relates to the neural network step

In the neural network step, MTD_preprocess from MEDS TorchData prepares the data. It is itself a MEDS-Transforms pipeline: it runs fit_normalization, fit_vocabulary_indices and normalization, just like our second pipeline, followed by two stages of its own, tokenization and tensorization, which reshape each subject's measurements into the sequences of tensors that its PyTorch dataset loads.

What MTD_preprocess does not do is the cleaning from our first pipeline. Because the first pipeline's output is still a MEDS dataset, you can hand it to MEDS TorchData in place of the downloaded data:

MTD_preprocess MEDS_dataset_dir="cleaned_MEDS_data" output_dir="tensorized_MEDS_data"

One caveat: at the time of writing, the latest MEDS TorchData release (0.9.0) requires a MEDS-Transforms version older than 0.7, so install it in a separate Python environment from the MEDS-Transforms 0.7.1 we used here.

Optional: MEDS-Reader

MEDS-Reader serves a different need. Rather than running a pipeline of stages, it converts a MEDS dataset into a database that you can read subject by subject, quickly, from Python. It's useful when you want to write your own pre-processing or feature extraction as plain Python code.

First, we convert the dataset. We convert the downloaded data rather than our cleaned data, because MEDS-Reader 0.2.0 can't convert boolean columns such as numeric_value/is_inlier.

[18]: 
1 
2
%%bash
meds_reader_convert MEDS_data MEDS_reader_db

Then we can open the database and loop over subjects. Each subject has a list of events, in time order, each with the subject's MEDS fields (such as time, code and numeric_value). For example, let's count how many ICU stays each subject has:

[19]: 
1 
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
from collections import Counter

import meds_reader

database = meds_reader.SubjectDatabase("MEDS_reader_db")
print(f"{len(database)} subjects")

subjects_by_n_icu_stays = Counter()
for subject_id in database:
    subject = database[subject_id]
    n_icu_stays = sum(event.code.startswith("ICU_ADMISSION") for event in subject.events)
    subjects_by_n_icu_stays[n_icu_stays] += 1

pl.DataFrame(
    {"ICU stays": list(subjects_by_n_icu_stays), "subjects": list(subjects_by_n_icu_stays.values())}
).sort("ICU stays")
100 subjects
[19]: 
shape: (5, 2)
ICU stayssubjects
i64i64
174
216
37
42
51

On larger datasets, SubjectDatabase.map runs a function over the subjects in parallel worker processes; see the MEDS-Reader documentation.

Key takeaways

  • Common pre-processing steps (filtering subjects, filtering rare codes, removing outliers, normalizing values, building a vocabulary) are ready-made stages in MEDS-Transforms. A pipeline is a short YAML file, run with MEDS_transform-pipeline.
  • Metadata stages are fit on the train split only, so the statistics behind filtering and normalization don't leak information from the tuning and held-out splits.
  • Every stage writes its output to its own folder, so you can check what each stage changed.
  • Cleaning stages keep the data in MEDS format, so any MEDS tool can use the result; normalization produces model-ready data with integer codes.
  • The neural network step's MTD_preprocess is a MEDS-Transforms pipeline too, and you can run it on the cleaned data.

Next, we'll build models: a tabular baseline and a neural network. Both start again from the downloaded data, so you can run them without this step.