**1. Introduction**

Long-track speed skating is a skillful sport where athletes glide on a 400-m ice rink at a speed of more than 50 km/h. The athletes accelerate their body using the ground reaction force exerted by the ice through an approximately 1 mm wide blade attached to the bottom of the skate shoe. Various studies have identified kinematic features of different movement phases, such as changes in the knee and trunk angles during races [1,2], that may influence the performance of long-track speed skating athletes. Additionally, a smaller push-off angle, which is the angle of the shank with respect to the floor in the frontal plane, has been shown to be associated with a greater power output [3] and skating velocity during a 5000-m race [1,2], although such a relationship was not observed during a 1500-m race [2]. Changes in the blade tilt angle during a 4000-m long-distance skating event have also been reported [4]. Several studies have also demonstrated the benefit of a greater knee flexion angle before the push-off to generate increased kinetic energy [3,5,6].

The majority of the studies that have investigated kinematic features during long-track speed skating has primarily used video analysis [1–6]. However, conventional kinematic measurements using video analysis have several limitations [7]. First, video analysis is largely influenced by the visibility of body landmarks. However, landmark visibility is often interfered by people or objects during long-track speed skating competitions or training sessions. Therefore, researchers need to synchronize measurements with multiple cameras and/or incorporate special environmental conditions to quantify kinematics with good body landmark visibility. Consequently, most kinematic analyses of long-track speed

**Citation:** Tomita, Y.; Iizuka, T.; Irisawa, K.; Imura, S. Detection of Movement Events of Long-Track Speed Skating Using Wearable Inertial Sensors. *Sensors* **2021**, *21*, 3649. https://doi.org/10.3390/s21113649

Academic Editors: Fabien Buisseret, Liesbet Van der Perre and Frédéric Dierick

Received: 14 April 2021 Accepted: 20 May 2021 Published: 24 May 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

skating only measure an isolated segment of an entire race [1,2]. Although one study showed significant variability in the knee joint angles of athletes with similar performance levels [8], the source of the variability may partly be explained by the limited precision of the measurements.

Second, video analysis requires an enormous amount of time for data processing, which includes the identification of movement phases and digitization of body landmark positions. This limits the use of acquired data for immediate feedback to the athletes. Therefore, feedback regarding skating techniques must be provided predominantly based on visual or qualitative assessments from observers, without the benefit of a quantitative assessment.

In an effort to solve these limitations, inertial measurement units (IMUs) have been used increasingly in recent years as an alternative method to measure kinematics in various sports, such as speed skating, running and skiing. An IMU utilizes three axial accelerometers, gyroscopes and geomagnetometers [7]. The validity of IMUs for gait event detection has been well established [9,10], while evidence is limited for movement phase detection during sports performance. IMUs have several advantages over conventional video analyses. First, IMUs are not restricted by the visibility of body landmarks because the IMU system does not use positional data to compute kinematic outcomes. This allows for a kinematic measurement in a crowd and for a wide range of performance areas. Second, a kinematic analysis with IMUs does not require the digitization of body landmarks, allowing for a real-time display of kinematic features, including the angular velocity, acceleration and joint and segment angles. Therefore, the IMU is a promising tool for use in kinematic data acquisition in various situations, including sporting events and clinical rehabilitation.

The validity of joint angles derived from IMUs during gait and running has been widely examined by comparing them with gold standard measurements (e.g., an optical 3D motion capture system and a magnetic motion capture system) [11–13]. However, the validity of IMUs for the identification of movement phase classifications during long-track speed skating remains unknown. While no standardized movement phase classifications exist for long-track speed skating, foot contact and foot-off of each leg represent the start of the stance and swing, respectively, and both are important features for characterizing skating strokes [4]. Therefore, in this study, we focused on the identification of foot contact and foot-off. The objective of this study was to estimate the accuracy of IMUs for identification of foot contact and foot-off in competitive speed skaters during long-track speed skating by comparing the method with phase identification using the foot pressure sensor system. The validation of IMUs would advance the applicability of the system for use during long-track speed skating competitions and training sessions to allow for comprehensive kinematic analyses in flexible environments and instant feedback to athletes.

#### **2. Materials and Methods**

#### *2.1. Participants*

Twelve healthy competitive athletes on a university long-track speed skating team participated in the study after signing an informed consent form. All the participants had more than 10 years of speed skating experience. The demographic characteristics of the participants are shown in Table 1. This study was approved by the ethics committee of the Takasaki University of Health and Welfare (approval number: 1904) in accordance with the Declaration of Helsinki. The participants had no musculoskeletal or neurological pathologies that could affect task performance.

**Table 1.** The demographics of the study participants (n = 12).


136

#### *2.2. Data Acquisition*

Data recording was performed during a full-speed 1000-m skating event from a static start position on a 400-m, two-lane indoor oval (Meiji Hokkaido-Tokachi Oval, Obihiro, Hokkaido, Japan).

The kinematic data were acquired using eight IMU sensors at a sampling rate of 100 Hz (myoMOTION, Noraxon, Scottsdale, AZ, USA). A data logger was embedded in each sensor, allowing data recording over a wide area. The IMU sensors were attached with double-sided tape to the skin at standardized locations on the lower thorax and pelvis and bilaterally on the thighs, shanks and feet. The specific sensor locations are shown in Table 2. Subsequently, the subjects wore compressive racing suits designed specifically for the body shape of each of the individual participants, which ensured that displacement of the sensors did not occur while recording was taking place. The foot sensors on the skating shoes were also stabilized with tape. The sensor locations were marked on the skin or skating shoes with a pen, as each sensor was attached, and we verified that there were no changes in the sensor locations before and after the data recordings.

**Table 2.** Sensor locations.


The kinetic data were acquired using a portable foot pressure measurement system at a sampling rate of 100 Hz (F-Scan System, TeckScan, South Boston, MA, USA). The system consists of two sensor sheets, two cuff units, one data logger unit and two cables connecting the cuff units and the data logger. Two sensor sheets were trimmed to the participant's foot size and a sensor sheet was inserted and attached to the sole of each of the skate shoes using double-sided tape. The cuff units were stabilized at the middle shanks and the data logger was attached to the back waist using Velcro tape. The IMU and foot pressure systems were synchronized using an electrical synch signal.

#### *2.3. Data Analysis*

We excluded the data from the first and last straights (first and last 50 m) and the first curve (100-m) from the analysis because the skating technique during these segments differs substantially from the remaining segment. Therefore, we analyzed the data for the remaining 800-m (400-m straight and 400-m curve) segment. The data from each side (left and right) were analyzed separately for both the straight and the curve. We adopted three types of analytical methods to detect foot contact and foot-off (Table 3).

**Table 3.** Overview of the three analytical methods used to detect foot contact and foot-off.


IMU: inertial measurement unit.

The first detection (the kinetic detection method) was based on the foot pressure data. Foot contact and foot-off were defined as the moments in which the foot pressure exceeded (foot contact; the vertical solid line in Figure 1A) and diminished below (foot-off; the vertical dotted line in Figure 1A) 20% of the peak foot pressure, which was calculated from all evaluated strokes (the horizontal dashed line in Figure 1A). Based on our empirical observations of the data obtained from all the participants in this study, 20% peak foot

pressure was high enough to avoid false detections due to noise, but low enough to detect both the foot contact and foot-off. *Sensors* **2021**, *21*, x FOR PEER REVIEW 5 of 10

**Figure 1.** The timing identifications using three different detection methods. The vertical solid and dotted lines in each panel show the foot contact and foot-off timing, respectively. We calculated the stance time by computing the duration of the foot contact and foot-off for each skating stroke (intervals within horizontal arrows). (**A**) **Kinetic detection using the foot pressure.** The horizontal dotted line indicates the threshold level (20% peak) for the identification of foot contact and footoff. (**B**) **Acceleration detection using the sagittal foot acceleration.** Gray line: raw sagittal acceleration. Red line: high-pass filtered sagittal acceleration. Blue line: low-pass filtered sagittal acceleration. (**C**) **Integrated detection using both the foot sagittal acceleration and the knee flexion angle**. **Figure 1.** The timing identifications using three different detection methods. The vertical solid and dotted lines in each panel show the foot contact and foot-off timing, respectively. We calculated the stance time by computing the duration of the foot contact and foot-off for each skating stroke (intervals within horizontal arrows). (**A**) **Kinetic detection using the foot pressure**. The horizontal dotted line indicates the threshold level (20% peak) for the identification of foot contact and foot-off. (**B**) **Acceleration detection using the sagittal foot acceleration**. Gray line: raw sagittal acceleration. Red line: high-pass filtered sagittal acceleration. Blue line: low-pass filtered sagittal acceleration. (**C**) **Integrated detection using both the foot sagittal acceleration and the knee flexion angle**.

*2.4. Statistical Analysis*  The stance time as calculated based on the kinetic detection method was considered as the gold standard measure in our study. The difference in the stance times among the three detection methods was examined by a repeated measures analysis of variance (ANOVA). The Tukey honestly significant difference test was performed for the post-hoc pairwise comparisons. We used the intra-class coefficient (ICC) to examine the similarity between the kinetic detection and acceleration/integrated detection methods by computing the ICC (2,1) and their 95% confidence intervals (95% CIs). To assess the validity of the proposed detection method, Bland-Altman plots and limits of agreement were calculated for both the acceleration and integrated detection methods, where we estimated the level of agreement between the proposed methods and the gold standard measure (i.e., the kinetic detection method). The bias between the proposed methods and the gold standard measure was calculated as the mean difference between the measurements from each method. The upper and lower limits of agreement, which The second detection (the acceleration detection method) was based on the foot sagittal (anterior-posterior direction of the foot) acceleration data obtained from the IMU sensors on each foot. The sagittal acceleration signal was used because it showed consistent changes at both the foot contact and foot-off throughout the entire 1000-m of skating. The measured sagittal acceleration signals were filtered and decomposed to their highfrequency (Butterworth high-pass filter at a cut-off frequency of 20 Hz; the red line in Figure 1B) and low-frequency (Butterworth low-pass filter at a cut-off frequency of 10 Hz; the blue line in Figure 1B) components. The high-frequency component represents instant acceleration changes and clearly shows foot contact (the vertical solid line in Figure 1B) and foot-off (the vertical dotted line in Figure 1B). The low-frequency component represents slower acceleration changes and shows the swing movement of the leg. We divided the data such that each segment comprised the start and end of the swing movement based on the low-frequency component of the acceleration. We then looked for the moment at which the high-frequency component of the foot sagittal acceleration reached its peak. The first and second peaks were set as the foot contact and foot-off, respectively.

defined the margin in which 95% of the differences between the methods were expected to lie, were calculated as a bias of ± 1.96 SD. The precision of the limits of agreement is reported as the 95% confidence interval. SPSS ver. 21 (IBM, Armonk, NY, USA) was used The third detection (the integrated detection method) was based on the combination of the foot acceleration and the knee flexion angle. Raw data were automatically filtered using

for the statistical analyses. A statistical significance level of *p* < 0.05 was used for all the

tests.

a robust fusion algorithm (Kalman filter) optimized for IMU data by the IMU software (myoRESEARCH 3.10, Noraxon, Scottsdale, AZ, USA). Four element quaternion values were derived by combining the elemental sensor component axes to estimate the angular offset of each sensor from the calibrated position in the global coordinate [14,15]. The knee flexion angle was automatically calculated using the biomechanical model adopted by the IMU system software. The bias (normalized root mean square) of the knee flexion angle, derived by the present IMU software from the angle based on the model recommended by the International Society of Biomechanics [16], has been reported to be 16.9 ± 5.1% during gait. The knee flexion angle showed a phasic pattern within each stroke (Figure 1C), allowing for clear segmentation of the strokes. We divided the data such that each segment comprised the start and end of the phasic pattern, which constituted one stroke and swing of each leg. We then looked for the moment at which the high-frequency component of the foot sagittal acceleration reached its peak. The first and second peaks were set as the foot contact (the vertical solid line in Figure 1C) and foot-off (the vertical dotted line in Figure 1C), respectively. All the timing detections were performed by combining the automated and visual identifications.

Based on previous studies that reported the validity of IMU systems to detect gait events [17,18], we calculated the stance time for each stroke (calculated as the time from foot contact to foot-off for each stroke) separately for each leg (right and left), section (straight and curve) and detection method (kinetic, acceleration and integrated detection). The data analysis was performed using custom-made programs (MATLAB 2014a, MathWorks, Natick, MA, USA).

## *2.4. Statistical Analysis*

The stance time as calculated based on the kinetic detection method was considered as the gold standard measure in our study. The difference in the stance times among the three detection methods was examined by a repeated measures analysis of variance (ANOVA). The Tukey honestly significant difference test was performed for the post-hoc pairwise comparisons. We used the intra-class coefficient (ICC) to examine the similarity between the kinetic detection and acceleration/integrated detection methods by computing the ICC (2,1) and their 95% confidence intervals (95% CIs). To assess the validity of the proposed detection method, Bland-Altman plots and limits of agreement were calculated for both the acceleration and integrated detection methods, where we estimated the level of agreement between the proposed methods and the gold standard measure (i.e., the kinetic detection method). The bias between the proposed methods and the gold standard measure was calculated as the mean difference between the measurements from each method. The upper and lower limits of agreement, which defined the margin in which 95% of the differences between the methods were expected to lie, were calculated as a bias of ±1.96 SD. The precision of the limits of agreement is reported as the 95% confidence interval. SPSS ver. 21 (IBM, Armonk, NY, USA) was used for the statistical analyses. A statistical significance level of *p* < 0.05 was used for all the tests.

#### **3. Results**

In total, 1036 strokes (86.3 ± 10.4 strokes per participant) were analyzed in this study. The stance times detected by the three methods are summarized in Table 4. The results of the repeated measures ANOVA showed significant differences of stance time among detection methods on the right side (straight: *F* = 15.236, *p* < 0.001; curve: *F* = 92.298, *p* < 0.001). The post-hoc analysis showed that acceleration and integrated detection methods on the right side significantly overestimated the stance time by 2.4–3.6%, compared to the kinetic detection method (*p* < 0.05; Table 4). No significant difference was found between acceleration and integrated detections.


**Table 4.** Stance time detected by kinetic, acceleration and integrated detection methods.

The differences between the acceleration and integrated detection methods and the kinetic detection methods are shown as ∆%. The proportion of cases within the limits of agreement is shown as LOA%. The *F* value and *p* value were obtained by the repeated measures analysis of variance. \* Significantly different from the kinetic detection method in the post-hoc analysis. **3. Results** 

In total, 1036 strokes (86.3 ± 10.4 strokes per participant) were analyzed in this study.

The gray-shaded areas in the figures show the limits of agreement of the two

The ICC (2,1) results for each detection method are summarized in Table 5. The ICC (2,1) was ≥0.700 for the integrated detection method for all the sections on both sides, while the ICC (2,1) was 0.657 for the acceleration detection method in the curve for the left side. The stance times detected by the three methods are summarized in Table 4. The results of the repeated measures ANOVA showed significant differences of stance time among detection methods on the right side (straight: *F* = 15.236, *p* < 0.001; curve: *F* = 92.298, *p* < 0.001). The post-hoc analysis showed that acceleration and integrated detection methods

**Table 5.** The intra-class coefficient as computed by the acceleration and integrated detection methods. on the right side significantly overestimated the stance time by 2.4–3.6%, compared to the


ICC: intra-class coefficient; 95% CI: 95% confidence interval. (right and left), and detection method (acceleration and integrated) (Figures 2 and 3).

The Bland-Altman plot is shown separately for each section (straight and curve), side (right and left), and detection method (acceleration and integrated) (Figures 2 and 3). detection methods. The proportion of cases within the limits of agreement was greater than 90% for all the measurements (LOA% in Table 4).

**Figure 2.** The Bland-Altman plot depicts the differences between the different detection methods in the straight, with 95% limits of agreement. The mean difference is shown by the dotted line. The 95% confidence intervals of the limits of agreement are also depicted (gray-shaded area). (**A**) **Figure 2.** The Bland-Altman plot depicts the differences between the different detection methods in the straight, with 95% limits of agreement. The mean difference is shown by the dotted line. The 95% confidence intervals of the limits of agreement are also depicted (gray-shaded area). (**A**) Acceleration detection on the right side. (**B**) Acceleration detection on the left side. (**C**) Integrated detection on the right side. (**D**) Integrated detection on the left side.

Acceleration detection on the right side. (**B**) Acceleration detection on the left side. (**C**) Integrated

*Sensors* **2021**, *21*, x FOR PEER REVIEW 7 of 10

**Figure 3.** The Bland-Altman plot depicts the differences between the different detection methods in the curve, with 95% limits of agreement. The mean difference is shown by the dotted line. The 95% confidence intervals of the limits of agreement are also depicted (gray-shaded area). (**A**) Acceleration detection on the right side. (**B**) Acceleration detection on the left side. (**C**) Integrated **Figure 3.** The Bland-Altman plot depicts the differences between the different detection methods in the curve, with 95% limits of agreement. The mean difference is shown by the dotted line. The 95% confidence intervals of the limits of agreement are also depicted (gray-shaded area). (**A**) Acceleration detection on the right side. (**B**) Acceleration detection on the left side. (**C**) Integrated detection on the right side. (**D**) Integrated detection on the left side.

detection on the right side. (**D**) Integrated detection on the left side. **Table 4.** Stance time detected by kinetic, acceleration and integrated detection methods. The gray-shaded areas in the figures show the limits of agreement of the two detection methods. The proportion of cases within the limits of agreement was greater than 90% for all the measurements (LOA% in Table 4).
