Skip to content

Commit 683259c

Browse files
committed
Start reviewing section 3.
1 parent 19cb67c commit 683259c

File tree

2 files changed

+251
-213
lines changed

2 files changed

+251
-213
lines changed

_02_Material.qmd

Lines changed: 35 additions & 213 deletions
Original file line numberDiff line numberDiff line change
@@ -235,7 +235,7 @@ An additional preprocessing step is required as part of the synchronization betw
235235
236236
```{r}
237237
#| label: tbl-raw-gaitrite-data
238-
#| tbl-cap: "Example (33rd session) of an extraction of the occurence of the key gait events from the GAITRite data *before* suppressing uncorrectly labelled events."
238+
#| tbl-cap: "Extraction of the occurence of the key gait events from the GAITRite data *before* suppressing uncorrectly labelled events (33rd session)."
239239
#| tbl-subcap: ["First half", "Second half"]
240240
#| tbl-pos: "H"
241241
#| layout-ncol: 2
@@ -333,7 +333,7 @@ gaitrite_data[[33]] |>
333333
334334
```{r}
335335
#| label: tbl-preprocessed-gaitrite-data
336-
#| tbl-cap: "Example (33rd session) of an extraction of the occurence of the key gait events from the GAITRite data *after* suppressing uncorrectly labelled events."
336+
#| tbl-cap: "Extraction of the occurence of the key gait events from the GAITRite data *after* suppressing uncorrectly labelled events (33rd session)."
337337
#| tbl-subcap: ["First half", "Second half"]
338338
#| tbl-pos: "H"
339339
#| layout-ncol: 2
@@ -397,11 +397,12 @@ gaitrite_data[[33]] |>
397397
398398
## Labelled IMU data {#sec-labelled-imu-data}
399399
400-
@sec-sec-imu-data presented the preprocessed IMU data from which we aim at identifying the four key gait events (RHS, LTO, LHS, RTO) and @sec-gaitrite-data presented the preprocessed GAITRite data which provides the ground truth labels for each timepoint. We can now perform a table join use the timepoints as key variable to create the labelled IMU data that will serve as the basis for building the feature space for training supervised classification models.
401-
402-
Before doing this join, we observe that as the obervations are ordered timepoints, we can think of including data from the past as features in addition to data from a given timepoint. This means that there is no need to keep IMU data collected *after* the last event detected in the GAITRite data. Similarly, we only keep the IMU data collected starting $10$ timepoints *before* the first event detected in the GAITRite data, which provides a lower bound to the lag interval that we will explore when training and tuning the classification models later on. The following piece of code
400+
@sec-imu-data presented the preprocessed IMU data from which we aim at identifying the four key gait events (RHS, LTO, LHS, RTO) and @sec-gaitrite-data presented the preprocessed GAITRite data which provides the ground truth labels for each timepoint. We can now perform a table join use the timepoints as key variable to create the labelled IMU data that will serve as the basis for building the feature space for training supervised classification models. First, we observe that, as the obervations are ordered timepoints, it makes sense to include data *anterior* to the considered timepoint as features in addition to data from the considered timepoint but we cannot include data recorded *posterior* to the considered timepoint. This means that there is no need to keep IMU data collected *after* the last event detected in the GAITRite data. Similarly, we only keep the IMU data collected starting $10$ timepoints *before* the first event detected in the GAITRite data. This choice prescribes a lower bound to the lag interval that we will explore when training and tuning the classification models later on. Once this filtering step is achieved, we can perform the table join to add the event labels to the timepoints in the IMU data. @fig-labelled-imu-data produced by the code below shows the filtered IMU data with events of interest clearly identified by joining the GAITRite data.
403401
404402
```{r}
403+
#| label: fig-labelled-imu-data
404+
#| fig-cap: "Filtered IMU data (33rd session) with key gait events superimposed with colored points."
405+
#| fig-pos: "H"
405406
min_timepoints <- gaitrite_data |>
406407
purrr::map("event_time") |>
407408
purrr::map_dbl(min)
@@ -420,215 +421,36 @@ imu_data <- purrr::pmap(
420421
}
421422
)
422423
423-
source("scripts/utils-viz.R")
424-
425-
gaitrite33 <- gaitrite_data[gaitrite_data$session == 33, ]
426-
rhs33 <- gaitrite33 |>
427-
dplyr::filter(event_type == "RHS") |>
428-
dplyr::pull(event_time)
429-
lto33 <- gaitrite33 |>
430-
dplyr::filter(event_type == "LTO") |>
431-
dplyr::pull(event_time)
432-
lhs33 <- gaitrite33 |>
433-
dplyr::filter(event_type == "LHS") |>
434-
dplyr::pull(event_time)
435-
rto33 <- gaitrite33 |>
436-
dplyr::filter(event_type == "RTO") |>
437-
dplyr::pull(event_time)
438-
plot_ts(imu_data[[33]], rhs33, rto33, lhs33, lto33)
439-
```
440-
441-
## Feature space {#sec-feature-space}
442-
443-
In this work, the objective is to identify which timepoints of unit QTS correspond to the RHS, LTO, LHS and RTO events. For this purpose, we consider the timepoints as statistical units (observations) and we aim at labelling them by means of supervised classification models. We therefore need to define a so-called *feature space* which consists of a data table listing the timepoints by row and collecting a number of features for each of them. A first important feature is the actual label that we want to predict with the trained model. @sec-labelled-gait-data details the elaboration of what
444-
445-
### Labelled gait data {#sec-labelled-gait-data}
446-
447-
In this view, we can first create the data set that we will use for training. The following code achieves this task by binding together all timepoints from all walking sessions while attaching to each timepoint:
448-
449-
- an `event_type` which affects it to one of the five gait events defined in @sec-gaitrite-data;
450-
- a `phase_type` which affects it one of the four gait pahses defined in @sec-gaitrite-data.
451-
452-
```{r}
453-
events_to_phases <- function(events) {
454-
events_of_interest <- events != "None"
455-
first_event <- events[events_of_interest][1]
456-
phase_durations <- diff(c(0, sort(which(events_of_interest))))
457-
n_phases <- length(phase_durations)
458-
phase_names <- switch(
459-
first_event,
460-
"RHS" = rep(
461-
c("Swing", "Pre-Stance", "Stance", "Pre-Swing"),
462-
times = n_phases
463-
)[1:n_phases],
464-
"LTO" = rep(
465-
c("Pre-Stance", "Stance", "Pre-Swing", "Swing"),
466-
times = n_phases
467-
)[1:n_phases],
468-
"LHS" = rep(
469-
c("Stance", "Pre-Swing", "Swing", "Pre-Stance"),
470-
times = n_phases
471-
)[1:n_phases],
472-
"RTO" = rep(
473-
c("Pre-Swing", "Swing", "Pre-Stance", "Stance"),
474-
times = n_phases
475-
)[1:n_phases]
476-
)
477-
purrr::map2(
478-
phase_names,
479-
phase_durations,
480-
\(phase_name, phase_duration) rep(phase_name, times = phase_duration)
481-
) |>
482-
purrr::list_c()
483-
}
484-
485-
labelled_gait_data <- purrr::map(1:nrow(bhg), \(session_index) {
486-
gaitrite_data <- gaitrite_data |>
487-
dplyr::filter(session == session_index) |>
488-
dplyr::select(-session)
489-
imu_data[[session_index]] |>
490-
dplyr::left_join(gaitrite_data, by = c("time" = "event_time")) |>
491-
dplyr::mutate(
492-
event_type = dplyr::if_else(is.na(event_type), "None", event_type),
493-
phase_type = events_to_phases(event_type)
494-
)
495-
}) |>
496-
dplyr::bind_rows(.id = "session") |>
497-
dplyr::mutate(
498-
session = as.numeric(session),
499-
event_type = factor(
500-
event_type,
501-
levels = c("RHS", "LTO", "LHS", "RTO", "None")
502-
),
503-
phase_type = factor(
504-
phase_type,
505-
levels = c("Pre-Stance", "Stance", "Pre-Swing", "Swing")
506-
)
507-
)
508-
class(labelled_gait_data) <- class(labelled_gait_data)[-1]
509-
head(labelled_gait_data)
510-
```
511-
512-
The feature space is an important piece of machine learning models as it defines the data that will be used to train them. In our application, we constructed a feature space from the raw QTS data recorded by the IMU sensor. First, we need to define what is an observation in our context. An observation corresponds to the data recorded at a given time point $t_j$. Since QTS are ordered sets of unit quaternions, for the $j$-*th* observation (row) of the feature space, we can then use features computed from the data observed at time $t_j$ or any other time points preceding $t_j$. One feature is of course the label of the observation that we get from the gold standard. This variable is available for the train/test process but will not be available at prediction time since it corresponds to the events that we want to predict. The other features are predictors that we compute from the sensor data. The ultimate goal is to design a feature space to label each observation $t_j$ with the gait event happening at that time (if any) using only the predictors computed from the QTS at time points $t_k \le t_j$, with $1 \le k \le j$. In the remainder of this section, we describe the predictors that we computed to build our feature space.
513-
514-
Angular velocity and acceleration vectors
515-
516-
: If we can have access to the first and second time derivatives of a QTS, we can compute the *angular velocity vector* $\pmb{\Omega}$ and *angular acceleration vector* $\dot{\pmb{\Omega}}$ [@narayan2017]. These vectors are aligned with the axis of rotation and carry in their norm the angular velocity and acceleration respectively. The angular velocity vector is computed as follows:
517-
518-
$$
519-
\begin{bmatrix}
520-
0 \\
521-
\pmb{\Omega}
522-
\end{bmatrix}
523-
= 2 \mathbf{q}^{-1} \dot{\mathbf{q}} \hspace{3mm} \text{with} \hspace{3mm} \dot{\mathbf{q}} = \frac{d \mathbf{q}}{dt} = \frac{1}{2} \mathbf{q} \begin{bmatrix}
524-
0 \\
525-
\pmb{\Omega}
526-
\end{bmatrix}.
527-
$$ {#eq-angular-vel}
528-
529-
Recall that $\mathbf{q} = \exp \mathbf{v}$ where $\mathbf{v}$ is the logarithm of $\mathbf{q}$ introduced in @eq-smoothed-qts. Then, we have after simple mathematical computations that $\dot{\mathbf{q}} = \mathbf{q} \dot{\mathbf{v}}$. Therefore, the angular velocity vector is nothing but twice the vector part of the temporal derivative $\dot{\mathbf{v}}$ of the log-QTS. This is implemented in the `squat::qts2avvts()` function, which returns a tibble with columns `time`, `x`, `y` and `z` storing the coordinates of the angular velocity vector at each time point as illustrated by @fig-avvts:
530-
531-
```{r}
532-
#| label: fig-avvts
533-
#| fig-cap: "An illustration (33rd session) of the angular velocity vector time series representation."
534-
#| fig-pos: "H"
535-
avvts <- squat::qts2avvts(imu_data[[33]], spar = spar)
536-
avvts |>
537-
dplyr::rename(`v[x]` = x, `v[y]` = y, `v[z]` = z) |>
538-
tidyr::pivot_longer(
539-
cols = c(`v[x]`, `v[y]`, `v[z]`),
540-
names_to = "component",
541-
values_to = "angular_velocity"
542-
) |>
543-
ggplot(aes(x = time, y = angular_velocity)) +
544-
geom_line() +
545-
facet_wrap(~component, ncol = 1, scales = "free_y", labeller = label_parsed) +
546-
theme_bw() +
547-
labs(title = "", x = "Time (seconds)", y = "Angular velocity (rad/s)")
548-
```
549-
550-
The angular acceleration vector $\dot{\pmb{\Omega}}$ is then computed as the derivative of the angular velocity vector:
551-
552-
$$
553-
\begin{bmatrix}
554-
0 \\
555-
\dot{\pmb{\Omega}}
556-
\end{bmatrix}
557-
= 2 \left( \mathbf{q}^{-1} \ddot{\mathbf{q}} - (\mathbf{q}^{-1}\dot{\mathbf{q}})^2 \right) \hspace{3mm} \text{with} \hspace{3mm} \ddot{\mathbf{q}} = \frac{d^2 \mathbf{q}}{dt^2} = \frac{1}{2} \left( \dot{\mathbf{q}} \begin{bmatrix}
558-
0 \\
559-
\pmb{\Omega}
560-
\end{bmatrix} + \mathbf{q} \begin{bmatrix}
561-
0 \\
562-
\dot{\pmb{\Omega}}
563-
\end{bmatrix} \right)
564-
$$ {#eq-angular-acc}
565-
566-
From $\dot{\mathbf{q}} = \mathbf{q} \dot{\mathbf{v}}$, we can write the second temporal derivative of $\mathbf{q}$ as $\ddot{\mathbf{q}} = \mathbf{q} \left( \dot{\mathbf{v}} \dot{\mathbf{v}} + \ddot{\mathbf{v}} \right)$. Therefore, the angular acceleration vector is nothing but twice the vector part of the second temporal derivative $\ddot{\mathbf{v}}$ of the log-QTS. This is implemented in the `squat::qts2aavts()` function, which returns a tibble with columns `time`, `x`, `y` and `z` storing the coordinates of the angular acceleration vector at each time point as illustrated by @fig-aavts:
567-
568-
```{r}
569-
#| label: fig-aavts
570-
#| fig-cap: "An illustration (33rd session) of the angular acceleration vector time series representation."
571-
#| fig-pos: "H"
572-
aavts <- squat::qts2aavts(imu_data[[33]], spar = spar)
573-
aavts |>
574-
dplyr::rename(`a[x]` = x, `a[y]` = y, `a[z]` = z) |>
575-
tidyr::pivot_longer(
576-
cols = c(`a[x]`, `a[y]`, `a[z]`),
577-
names_to = "component",
578-
values_to = "angular_acceleration"
579-
) |>
580-
ggplot(aes(x = time, y = angular_acceleration)) +
581-
geom_line() +
582-
facet_wrap(~component, ncol = 1, scales = "free_y", labeller = label_parsed) +
583-
theme_bw() +
584-
labs(title = "", x = "Time (seconds)", y = "Angular acceleration (rad/s²)")
585-
```
586-
587-
Euler angles
588-
589-
: The angles named *Roll*, *Pitch* and *Yaw* represent rotations around the three principal axes. They are computed from a unit quaternion $\mathbf{q} = (q_w, q_x, q_y, q_z)$ using the following formulas:
424+
labelled_imu_data <- purrr::map2(
425+
imu_data,
426+
gaitrite_data,
427+
\(.imu, .gaitrite) {
428+
out <- .imu |>
429+
dplyr::left_join(.gaitrite, by = c("time" = "event_time")) |>
430+
dplyr::mutate(
431+
event_type = dplyr::if_else(is.na(event_type), "None", event_type),
432+
event_type = factor(
433+
event_type,
434+
levels = c("RHS", "LTO", "LHS", "RTO", "None")
435+
)
436+
)
437+
class(out) <- class(out)[-1]
438+
out
439+
}
440+
)
590441
591-
$$
592-
\begin{bmatrix}
593-
\mathrm{Roll} \\
594-
\mathrm{Pitch} \\
595-
\mathrm{Yaw}
596-
\end{bmatrix}
597-
=
598-
\begin{bmatrix}
599-
\arctan\!2 \left(2(q_w q_x + q_y q_z), 1-2(q_x^2 + q_y^2) \right) \\
600-
\arcsin \left(2(q_w q_y - q_x q_z) \right) \\
601-
\arctan\!2 \left(2(q_w q_z + q_x q_y), 1-2(q_y^2 + q_z^2) \right)
602-
\end{bmatrix}
603-
$$ {#eq-rpy}
604-
605-
where $\arctan\!2(y, x)$ computes the angle $\theta$ (in radians) between the positive $x$-axis and the ray from the origin to the point $(x, y)$ in the Cartesian plane. This is implemented in the `squat::qts2rpyts()` function, which returns a tibble with columns `time`, `roll`, `pitch` and `yaw` storing the roll, pitch and yaw angles at each time point as illustrated by @fig-rpyts:
442+
df <- tidyr::pivot_longer(labelled_imu_data[[33]], cols = -c(time, event_type))
606443
607-
```{r}
608-
#| label: fig-rpyts
609-
#| fig-cap: "An illustration (33rd session) of the roll-pitch-yaw time series representation."
610-
#| fig-pos: "H"
611-
rpyts <- squat::qts2rpyts(imu_data[[33]])
612-
rpyts |>
613-
dplyr::rename(Roll = roll, Pitch = pitch, Yaw = yaw) |>
614-
tidyr::pivot_longer(
615-
cols = c(Roll, Pitch, Yaw),
616-
names_to = "component",
617-
values_to = "angle_values"
618-
) |>
619-
ggplot(aes(x = time, y = angle_values)) +
444+
df |>
445+
ggplot(aes(x = time, y = value)) +
446+
geom_point(
447+
data = dplyr::filter(df, event_type != "None"),
448+
mapping = aes(color = event_type),
449+
size = 2
450+
) +
620451
geom_line() +
621-
facet_wrap(~component, ncol = 1, scales = "free_y", labeller = label_parsed) +
452+
facet_wrap(vars(name), ncol = 1, scales = "free") +
622453
theme_bw() +
623-
labs(title = "", x = "Time (seconds)", y = "Angle (rad)")
624-
```
625-
626-
Walking speed
627-
628-
: It is likely that the shape of the QTS can be quite different according to the walking speed. We therefore included this information in the feature space from the GAITRite® mat output. Since this predictor comes from the gold standard, it is not available when predicting on new time series where patients only used the wearable sensor. To counter this issue, we can estimate the walking speed from the mean angular velocity with a simple linear regression.
629-
630-
Hyper-parameters
631-
632-
: The feature space depends on a number of hyper-parameters. The first one is a smoothness parameter called *spar* that represents the amount of smoothing that we seek to achieve when fitting cubic splines to the original QTS. Since time derivatives are used to compute some predictors, it is indeed useful to smooth the QTS sufficiently so the derivative values remain stable and not dominated by noise. The *spar* parameter takes its values in the interval $]0,1]$. The second hyper-parameter is the *lag* parameter. To keep information from the past as we predict at a given time point $t_j$, we include features computed at $t_j$ as well as the same features computed at time points $t_{j-1}, \dots, t_{j-\mathrm{lag}}$. The size of the feature space therefore depends on the *lag* parameter as it contains $10 + 9 \times \mathrm{lag}$ predictors. The ten first predictors are the three-component angular speed vector, the three-component angular acceleration vector, the three angles Roll, Pitch and Yaw and the walking speed. We then add nine predictors per *lag*: the angular speed and acceleration and the Roll-Pitch-Yaw angles from the previous time.
633-
634-
Finally, a last hyper-parameter is used when we label the observations with the gait event times from the gold standard. It is called *k* and controls the number of points around a gait event that we label as part of the event. In @sec-classification-strategies, we elicit two strategies for building segmentation models, one of which does not use this parameter. Therefore, we differ other details about this parameter in the dedicated section.
454+
labs(x = "Time (s)", y = "", color = "Event type") +
455+
theme(legend.position = "top")
456+
```

0 commit comments

Comments
 (0)