Embedded HandbookApps
All sections › XIV. Machine learning on a microcontroller
STM32F407 · level intermediate · 18 min read
Open in the interactive handbook

117а. Preparing data for the model: collecting from the board into CSV, exploratory data analysis (EDA) in pandas, numpy and matplotlib, cleaning, class balance and labeling

Section 117 cut the signal into windows and computed features from them, section 118 split the recordings into training and test parts and trained the model. Between them lies a step that is usually skipped, after which people spend a long time hunting for a bug in the model: looking at the data itself. How many recordings arrived intact, whether each one has a label, whether the sensor hit the end of the ADC scale, whether the zero shifted between days, whether one file got into the set twice, and whether the class is given away by something unrelated — the device number or the day of recording. This kind of examination is called exploratory data analysis.

Here the whole path is walked on one data set: the board writes frames, the PC receives them and lays them out in a table, pandas, numpy and matplotlib find the troubles in it one after another, and the output is a clean set of windows and a short data report. The data are synthetic, with seed 1171, but built like real ones: the same motor vibration as in sections 117 and 118.

Basis: the stm-eda:1 image — pandas 3.1.0, numpy 2.3.3, matplotlib 3.11.2, seaborn 0.13.2, scikit-learn 1.7.2 (built on stm-ml:1); the frame code in C is built by gcc on a PC and by armclang from Keil MDK for Cortex-M4. The port and the CRC on the Python side are section 115в, exchange formats are section 115г, a USB logger with a PC program is section 198б, the SD card and FatFS are sections 41а and 42.

Terms used in the documentation

Term Meaning
exploratory data analysis (EDA) examining the data before modelling
data logger a device that records measurements
frame, sync word a frame, the sync word at its start
DataFrame, Series a pandas table and column
dtype the type of a column
merge, join joining tables by a key
resample, rolling window recomputing at another time step, a sliding window
vectorization, boolean mask whole-array operations, a true/false mask
clipping, saturation the signal hit the end of the scale
outlier, spike an outlier, a single-sample burst
sensor drift slow change of the sensor's zero or gain
class imbalance some classes have far fewer examples
data augmentation multiplying data: shift, noise
data leakage the model sees something it will not have on the device
data card, data report a report on the data set

117а.1 Where the data come from and how the board writes them

Your own sensor on your own board is the main source: only then do the data look like what the device will see at work. They can be delivered in two ways. A stream over UART or a USB virtual COM port (CDC) — the PC stands next to it and writes a file. A logger on an SD card — the device works alone on site, and the files are collected later. Open data sets are good for trying out a method, but they were taken with someone else's sensor on someone else's machine, and a model trained only on them almost always makes mistakes on your own device.

Whatever the channel, you need to record more than the samples. You need a time stamp, the sampling rate and information about the device: which one it is, which firmware, how many bits the ADC has and what the reference voltage is. Without them, a month later you cannot tell why two recordings differ. The board in this section writes two kinds of frames: the first frame of a recording is the string dev=3A0F1C22;fw=1.2.0;fs=1024;n=64;adc=12;vref=3300, then come frames of 64 samples. Each frame is the sync word A5 5A, type, length, number, the time of the first sample in microseconds and a CRC-16:

size_t frame_seal(uint8_t *out, uint8_t type, uint16_t seq, uint32_t t_us, uint8_t len)
{
    out[0] = FR_SYNC0;
    out[1] = FR_SYNC1;
    out[2] = type;
    out[3] = len;
    put16(out + 4, seq);
    put16(out + 6, (uint16_t)t_us);
    put16(out + 8, (uint16_t)(t_us >> 16));
    put16(out + FR_HEAD + len, crc16_ccitt(out + 2, FR_HEAD - 2u + len));
    return FR_HEAD + len + 2u;
}

A frame with samples takes 140 bytes and carries 62 500 µs of signal. At 1024 Hz this is 2240 bytes per second — 19.4 % of the 11 520 bytes per second that a UART passes at 115 200 baud. The same recording as text, a "time,code" line per sample, took 17 bytes per sample — 17 408 bytes per second: it no longer fits into a UART at 115 200 baud. CSV is convenient on a PC, but from the board a binary frame pays off.

armclang for Cortex-M4 with -Os built the four frame functions into 278 bytes of code; with -O2 they take 956 bytes: the compiler unrolls the eight CRC steps. The CRC here is bitwise, without a table: a table of 256 values would take another 512 bytes of flash. How many cycles a frame takes was not measured.


117а.2 Receiving on the PC: frames, CRC, gaps

The set is 40 files of 4 s each, every one with 4096 samples: 64 frames with samples and a frame with information about the device. The frames are built by the same C as on the board — frame117b.c, built by gcc on a PC — and an encoder in Python repeats the format:

def frame(typ, seq, t_us, pay):
    body = struct.pack("<BBHI", typ, len(pay), seq, t_us) + pay
    return b"\xa5\x5a" + body + struct.pack("<H", crc16(body))

All 40 files matched byte for byte: Python understands the frame the same way C writes it. This is the first check of any collection rig. If the struct format string has the byte order or a field size wrong, the error shows up here, and not a month later in the model weights. crc16 is one line, binascii.crc_hqx(b, 0xFFFF), and it is the same CRC-16/CCITT-FALSE that the board computes.

Between the board and the file, in the experiment three frames were lost — the board did not manage to send them — and in two, interference flipped a bit. The reader looks for the sync word, checks the CRC and on an error moves on by one byte:

while (i := buf.find(b"\xa5\x5a", i)) >= 0 and i + 12 <= len(buf):
    typ, n, seq, t_us = struct.unpack_from("<BBHI", buf, i + 2)
    end = i + 12 + n
    if end > len(buf) or crc16(buf[i + 2:end - 2]) != struct.unpack_from("<H", buf, end - 2)[0]:

2555 frames out of 2560 were received, CRC errors 2. All 163 520 received samples matched what the board sent: a broken frame is dropped whole, and not a single wrong number got into the table. A gap is seen from the frame time: the step is 62 500 µs, anything else is a loss. But the board's time is 32-bit, the microsecond counter overflows once every 71.6 min, and on one device this happened in the middle of a recording. A plain difference found 6 gaps, the difference modulo 2³² — the real 5:

gaps += int((np.diff(t) % 2**32 != DT_US).sum())   # время 32-битное: разность — по модулю 2^32

What is received is written to CSV: columns rec, t_s, adc, 3 588 343 bytes for 163 520 rows. The same signal in frames took 360 920 bytes. The columnar Parquet format would be both smaller and faster, but pandas writes it through the pyarrow or fastparquet library, and the image has neither: to_parquet answers ImportError. They were not installed here — CSV is enough for a set of a few megabytes.


117а.3 pandas: the table, types, metadata and the log

read_csv guesses the types itself: rec is strings (str), t_s is float64, adc is int64. Eight bytes for a 12-bit sample is wasteful; the table took 11 283 012 bytes in memory, and after changing the types — 1 146 892:

df = df.astype({"rec": "category", "t_s": "float32", "adc": "uint16"})

The class label is not in the board's data: the board does not know whether the motor is healthy. The labels come from the operator's log: recording, motor, class, day, who took it. The recordings table is the information from the first frames plus the log. The join is best done with checks: validate stops the work if a recording got into the log twice, and indicator shows recordings without a line in the log:

rec_tab = rec_tab.merge(labels, on="rec", how="left", indicator=True, validate="one_to_one")

This is how recording r017 was found: the operator did not log it, and it has no class. Guessing the class from the signal is not allowed — that is labeling "by eye", see section 118. The recording is set aside. A duplicate is found not by name but by content: a hash of the samples of each recording and duplicated:

rec_tab["sha1"] = samples.groupby("rec").adc.agg(lambda a: hashlib.sha1(a.to_numpy().tobytes()).hexdigest()).to_numpy()

Recordings r008 and r033 are identical down to the sample: one file was uploaded twice and logged twice. A duplicate is more dangerous than it seems: if the copies land in training and in the test part, the model "recognises" a recording it has already seen.

value_counts over recordings — 24 normal, 10 imbalance, 5 bearing. crosstab of device by class shows what cannot be seen in the log by eye:

Device Firmware normal imbalance bearing
3A0F1C22 1.2.0 24 10 0
7C41D905 1.3.1 0 0 5

All the bearings were taken with another device — on a rig with a known bad bearing. We will come back to this in 117а.6.

The samples are unfolded into a grid "sample number × recording" — pivot and reindex over 0…4095. Where a frame was lost, NaN appears, and the type of the columns becomes float64: by default an integer in pandas cannot be empty. There are 320 missing values, exactly five frames of 64 samples. Recomputing at a step of 250 ms shows which windows lack them:

cnt = wide.set_axis(pd.to_timedelta(wide.index / FS, unit="s")).resample("250ms").count()

There are 16 windows of 250 ms per recording, incomplete ones — 5, one for each loss.


117а.4 Plots: look first

Six exploratory plots from the run in the stm-eda:1 image: recording r021 of strong imbalance hits codes 0 and 4095 with flat tops; the rolling RMS of recording r005 grows during the first half second of run-up and breaks at the lost frame; histograms of window RMS by class — normal and bearing lie together, imbalance further right; on the scatter plot of RMS and the energy share in 100–248 Hz the bearing stands apart above 0.25; in the correlation matrix the band energy share is linked with the crest factor by 0.89 and with the zero offset by minus 0.56; the mean window spectrum of all classes has a 25 Hz peak, the bearing also has a hump near 160 Hz

The figure is assembled from six plots of one run, matplotlib and seaborn, each one or two lines:

ax[1, 0].scatter(w.rms, w.band, s=7, color=COL[c], alpha=0.6, label=c, lw=0)

1. Time series. A piece of recording r021: the tops of the sine are cut at codes 0 and 4095. This is saturation — the sensor or amplifier went beyond the ADC scale. The RMS and spectrum of such a window lie: part of the signal is simply lost.

2. Rolling RMS. Recording r005: for the first half second the RMS grows — the motor is running up — then it stays flat. The gap near 1.19 s is lost frame number 20.

3. Histograms by class. By RMS alone normal cannot be separated from bearing: their windows lie in the same range. Imbalance is further right.

4. Scatter plot. RMS against the energy share in the 100–248 Hz band (section 117): the bearing stands apart. Something else is visible too: the bearing points lie in five tight clusters, one per recording. Five clusters are five motors, not five hundred.

5. Correlation matrix. The band energy share and the crest factor are linked by 0.89: both features see the bearing impacts, and the second adds little to the first. The alarming cell is another one: zero offset and band — −0.56. The zero offset has nothing to do with the defect, and a strong link means that the class is tied to the device.

6. Mean spectrum. All classes have a 25 Hz peak — the rotation frequency, the bearing also has a hump near 160 Hz. This is the band on which the feature of section 117 is built.


117а.5 numpy: masks instead of loops

Checks for saturation, outliers and run-up are comparisons over all samples at once. In numpy they are written as boolean masks: one line per rule, and no Python loop over 163 840 samples.

Saturation — three samples in a row at the edge of the scale:

rail = (X <= 0) | (X >= 4095)
flat = rail[:-2] & rail[1:-1] & rail[2:]                # три отсчёта подряд на краю шкалы — плоская вершина

Flat tops were found in r021 and r024 — two recordings of strong imbalance, 3053 samples at the edge of the scale. There is nothing to repair such windows with, they are set aside. But the conclusion is wider: exactly the strongest imbalances are set aside, and the model will not see them. The right cure is on the device: a sensor with a larger range or less gain in front of the ADC (section 14).

Outliers — a single sample far from the median of five neighbours. The threshold is 12 times the typical deviation of the same recording, but not less than 400 codes, otherwise the strong 160 Hz bearing tone would itself pass for an outlier:

med = wide.rolling(5, center=True).median().to_numpy()
thr = np.maximum(12 * np.nanmedian(np.abs(r_), axis=0), 400)
spike = big & ~near                                     # выброс — одиночный: соседи в норме, до шкалы далеко

8 samples in 4 recordings were found — exactly the interference bursts that were put into the experiment; among them the copy r033, which has the same two bursts as r008. An outlier is replaced by the median of its neighbours. The same calculation written as a loop over samples gave the same 163 840 numbers. A loop is fine for checking a mask, but the working code is masks.

Run-up — the rolling RMS over 256 samples and the first moment it reached 90 % of the steady value:

runup = (rms_roll >= 0.9 * steady).idxmax() / FS        # первая секунда, где СКЗ дошло до 90 % установившегося

The median is 0.502 s, the latest is 0.562 s: the first three windows of 256 samples in each recording are run-up, not operation, and they are discarded.

Drift between days — groupby by day and device:

Day, device Mean sensor voltage
21.10, 3A0F1C22 1.650 V
22.10, 3A0F1C22 1.698 V
22.10, 7C41D905 1.610 V

The zero of the first device moved by 48 mV in a day — about 60 ADC codes. RMS and spectrum are not affected: the window is cleared of the constant component anyway (section 117). But any feature that looks at the raw level will tell days apart, not motors. Sensitivity drift is more dangerous — it is dealt with by the mount calibration in section 118.


117а.6 Class imbalance, leakage and the data report

After cleaning, 466 of 640 windows from 36 recordings remain. Discarded are 32 windows of the two set-aside recordings, 114 run-up windows, 3 windows with a gap and 25 windows with saturation. By class: normal — 298 windows, imbalance — 103, bearing — 65. describe of the RMS of the clean windows: 466 windows, mean 0.395 V, spread 0.169 V, from 0.172 to 0.869 V.

Class imbalance. There are more than four times fewer bearings than normal. What to do about it — in order of usefulness. First — collect more: five bearing recordings are five motors, and no processing will add a sixth. Second — class weights in training (section 118). For this set compute_class_weight("balanced") gave 0.52 for normal, 1.51 for imbalance and 2.39 for bearing. Third — augmentation: the window is shifted by 64 samples and weak noise is added to a copy. From 65 bearing windows this gave 455. These are the same five recordings, and only the training part may be augmented, already after the split by recordings. Otherwise a shifted copy of a test window ends up in training — the same leakage as in section 118.

Leakage through metadata. A one-question tree trained only on the "device" column finds all 65 bearing windows without a single false alarm. By the zero offset alone — 51 of 65, false alarms 9:

p = DecisionTreeClassifier(max_depth=1, random_state=117).fit(col, yb).predict(col)

A model given such a column learns "the rig with the bad bearing", and a line device will have no such rig. Device number, firmware, day and raw signal level are not taken as features. But even without them a set with this bias is dangerous: the device and the day of the bearing leave traces in the signal itself. The only cure is collection — bearing recordings from line devices.

The data report is what is attached to the set together with the window file clean117b.csv. It is written by the same script, and it shows what is wrong with the set before any training:

ОТЧЁТ О ДАННЫХ — вибрация моторов, 1024 Гц, окно 256
файлов 40, кадров принято 2555 из 2560, ошибок CRC 2, разрывов времени 5
без строки в журнале: r017 — отложена
дубликат: r008 = r033 — оставлена первая
упор в шкалу АЦП: r021, r024 — окна с плоскими вершинами отложены
выбросы: 8 отсчётов в 4 записях заменены медианой соседей
разгон: до 0,56 с — первые 3 окна записи отброшены
смещение нуля: 10-21 3A0F1C22 — 1,650 В, 10-22 3A0F1C22 — 1,698 В, 10-22 7C41D905 — 1,610 В
окон всего 640, в чистом наборе 466 из 36 записей
по классам: норма 298, дисбаланс 103, подшипник 65
ВНИМАНИЕ: все записи подшипника — прибор 7C41D905, прошивка 1.3.1, день 10-22;
  dev, fw, день и смещение нуля в признаки не брать; нужны записи подшипника с приборов линии

From here section 118 works with this set: labels on whole recordings, split by recordings, cross-validation.


117а.7 Pitfalls

The model is good on the PC, nonsense on the board, although the feature code is the same. The frame reader parses something other than what the board writes: byte order, field size, sign. Your own Python encoder must produce the same bytes as C — check before collecting data.

A "gap" in a recording out of nowhere. The 32-bit microsecond counter overflowed, once every 71.6 min. The time difference is taken modulo 2³².

The ADC code column suddenly became float64. NaN from lost frames appeared in it. This is not an error but missing values; count them with isna and do not fill them with zeros — zero is a legal ADC value.

win.clip & … fails with an obscure error. The column is named the same as a table method: clip, count, mean. Through the dot, pandas returns the method. Write win["clip"] or name the column differently.

"Outliers" every other sample on the bearing tone. The outlier threshold is set in codes without regard to the signal itself. The threshold comes from the spread of the same recording, and an outlier is a single sample.

The class is guessed suspiciously easily. It is given away by the device, day, firmware or zero offset. crosstab of the class with every metadata column — before training.

Many windows but few motors. Shift and noise multiply windows, not motors. Count recordings, not windows, and augment only the training part.

The strongest defects fell out of the set. Their signal hit the ADC scale, and cleaning set them aside. Look at what exactly was discarded, and fix the range on the device.


117а.8 Summary

Full sources

frame117b.c — the recording frame on the board and on a PC. The rest — sim117b.py (recordings), eda117b.py (reception, cleaning, report) and fig117b.py (figure) — is in _verify/ml117b/.

/* Раздел 117а: кадр записи вибрации — плата пишет, ПК читает (eda117b.py).
 * Кадр: A5 5A | тип | длина | номер (2) | время первого отсчёта, мкс (4) | данные | CRC-16 (2).
 * Числа — младшим байтом вперёд; CRC-16/CCITT-FALSE от типа до конца данных.
 * На плате кадр собирается из половины буфера DMA АЦП и уходит в UART или USB CDC.
 * С -DFRAME_HOST main() читает записи из stdin и пишет по файлу кадров на запись —
 * так проверка сверяет C с читателем на Python. Проверка — _verify/ml117b_check.py. */
#include <stdint.h>
#include <stddef.h>

#define FR_SYNC0 0xA5u
#define FR_SYNC1 0x5Au
#define FR_META  1u                     /* данные — строка key=value: прибор, прошивка, частота */
#define FR_DATA  2u                     /* данные — FR_N отсчётов АЦП по 2 байта */
#define FR_N     64u                    /* отсчётов в кадре: 62,5 мс при 1024 Гц */
#define FR_HEAD  10u                    /* синхро 2 + тип 1 + длина 1 + номер 2 + время 4 */
#define FR_MAX   (FR_HEAD + 2u * FR_N + 2u)

uint16_t crc16_ccitt(const uint8_t *p, size_t n)
{
    uint16_t c = 0xFFFFu;
    while (n--) {
        c ^= (uint16_t)(*p++ << 8);
        for (int i = 0; i < 8; i++)
            c = (c & 0x8000u) ? (uint16_t)((c << 1) ^ 0x1021u) : (uint16_t)(c << 1);
    }
    return c;
}

static void put16(uint8_t *p, uint16_t v)
{
    p[0] = (uint8_t)v;
    p[1] = (uint8_t)(v >> 8);
}

/* дописать заголовок и CRC к данным, уже лежащим в out + FR_HEAD; вернуть длину кадра */
size_t frame_seal(uint8_t *out, uint8_t type, uint16_t seq, uint32_t t_us, uint8_t len)
{
    out[0] = FR_SYNC0;
    out[1] = FR_SYNC1;
    out[2] = type;
    out[3] = len;
    put16(out + 4, seq);
    put16(out + 6, (uint16_t)t_us);
    put16(out + 8, (uint16_t)(t_us >> 16));
    put16(out + FR_HEAD + len, crc16_ccitt(out + 2, FR_HEAD - 2u + len));
    return FR_HEAD + len + 2u;
}

size_t frame_data(uint8_t *out, uint16_t seq, uint32_t t_us, const uint16_t *adc)
{
    for (unsigned i = 0; i < FR_N; i++)
        put16(out + FR_HEAD + 2u * i, adc[i]);
    return frame_seal(out, FR_DATA, seq, t_us, 2u * FR_N);
}

size_t frame_meta(uint8_t *out, uint16_t seq, uint32_t t_us, const char *kv)
{
    uint8_t n = 0;
    while (kv[n] && n < 2u * FR_N) {
        out[FR_HEAD + n] = (uint8_t)kv[n];
        n++;
    }
    return frame_seal(out, FR_META, seq, t_us, n);
}

#ifdef FRAME_HOST
#include <stdio.h>
#include <stdlib.h>
static int rd(void *p, size_t n) { return fread(p, 1, n, stdin) == n; }
static uint32_t le32(const uint8_t *b) { return b[0] | b[1] << 8 | b[2] << 16 | (uint32_t)b[3] << 24; }

int main(int argc, char **argv)
{
    uint8_t len, b8[8], fr[FR_MAX];
    char kv[256], name[512];
    int rec = 0, frames = 0;
    while (rd(&len, 1)) {
        if (!rd(kv, len) || !rd(b8, 8))
            return 1;
        kv[len] = 0;
        uint32_t t0 = le32(b8), ns = le32(b8 + 4);
        uint16_t *adc = malloc(ns * 2u);
        for (uint32_t i = 0; i < ns; i++) {
            uint8_t s[2];
            if (!rd(s, 2))
                return 1;
            adc[i] = (uint16_t)(s[0] | s[1] << 8);
        }
        snprintf(name, sizeof name, "%s/r%03d.bin", argc > 1 ? argv[1] : ".", ++rec);
        FILE *f = fopen(name, "wb");
        if (!f)
            return 1;
        fwrite(fr, 1, frame_meta(fr, 0, t0, kv), f);
        for (uint32_t k = 0; k < ns / FR_N; k++)      /* время первого отсчёта: 64 отсчёта = 62 500 мкс */
            fwrite(fr, 1, frame_data(fr, (uint16_t)(k + 1), t0 + k * 62500u, adc + k * FR_N), f);
        frames += 1 + (int)(ns / FR_N);
        fclose(f);
        free(adc);
    }
    printf("recs=%d\nframes=%d\nframe_bytes=%u\n", rec, frames, (unsigned)FR_MAX);
    return 0;
}
#else
/* на плате: половина буфера DMA готова — собрать кадр и отдать передатчику */
volatile uint16_t g_adc[FR_N];
uint8_t g_tx[FR_MAX];
volatile size_t g_len;
static uint16_t seq;
void SystemInit(void) {}
int main(void)
{
    g_len = frame_meta(g_tx, seq++, 0, "dev=3A0F1C22;fw=1.2.0;fs=1024;n=64;adc=12;vref=3300");
    for (uint32_t t_us = 0;; t_us += 62500u)
        g_len = frame_data(g_tx, seq++, t_us, (const uint16_t *)g_adc);
}
#endif

Next: section 118 "Training a small model in a container" — from windows with features to weights in the firmware.

Open in the interactive handbook