## 2 Materials and Methods
### Principal Motion Analysis to Estimate MoS Values from Kinematic Data  
We employ PMA [28] to estimate the MoS values from the triaxial translational velocities of a body feature. PMA encompasses a suite of linear time-series data analysis techniques. It involves decomposing a multidimensional time-series data sample into a set of independent bases, which are also multidimensional time series. The unsupervised variant, known as PMA0, extends the principal component analysis to time-series data and has been applied in the analysis of human motion [35; 36; 37; 38; 39]. For this study, we used a supervised variant, PMA1 [28], which is a time-series extension of partial least-squares regression analysis. In PMA1, a time-series data sample is represented as a linear combination of multiple independent time-series bases. These bases are selected to maximize the covariance between the objective variable and principal motion scores, which quantify the extent to which each sample incorporates the bases. Herein, we introduce PMA1 by first explaining partial least-squares regression analysis, followed by the adaptation of its explanatory variables to time-series data. Let \(x\in\mathbb{R}^{p\times 1}\) and \(y\) be the explanatory vector and objective variable, respectively. \(x_{i}\) and \(y_{i}\) are the values from the \(i\)th sample (\(i=1,\ldots,n\)). These variables are centered such that their mean values are zero: \(\overline{x}=\mathbf{0}\) and \(\overline{y}=0\). The explanatory vectors of all the samples are displayed in a matrix (\(\mathbf{X}\in\mathbb{R}^{n\times p}\)):  
\[\mathbf{X}= \begin{bmatrix}\mathbf{x}_{1}^{\mathrm{T}}\\ \vdots\\ \mathbf{x}_{n}^{\mathrm{T}}\end{bmatrix} . \tag{5}\]  
Similarly, the objective values of all the samples are stored in a vector \(\mathbf{y}\in\mathbb{R}^{n\times 1}\):  
\[\mathbf{y}=[y_{1}\ \ldots y_{n}]^{\mathrm{T}}. \tag{6}\]  
\(\mathbf{X}\) and \(\mathbf{y}\) are approximated as a linear summation of \(a^{\prime}\) (\(\in\mathbb{N}\)) principal components:  
\[\mathbf{X} =\sum_{a=1}^{a^{\prime}}\mathbf{t}_{a}\mathbf{I}_{a}^{\mathrm{T}}+\mathbf{E} \tag{7}\] \[\mathbf{y} =\sum_{a=1}^{a^{\prime}}b_{a}\mathbf{t}_{a}+\mathbf{e}, \tag{8}\]  
where \(\mathbf{t}_{a}\in\mathbb{R}^{n\times 1}\) and \(\mathbf{I}_{a}\in\mathbb{R}^{p\times 1}\) are the scores and loadings of the \(a\)th principal component (\(a=1,\ldots,a^{\prime}\)). \(\mathbf{E}\) and \(\mathbf{e}\) are the error matrix and vectors. Let \(\mathbf{X}_{a}\) and \(\mathbf{y}_{a}\) be  
\[\mathbf{X}_{a} =\mathbf{X}_{a-1}-\mathbf{t}_{a-1}\mathbf{I}_{a-1}^{\mathrm{T}} \tag{9}\] \[\mathbf{y}_{a} =\mathbf{y}_{a-1}-\mathbf{b}_{a-1}\mathbf{t}_{a-1}. \tag{10}\]  
Here, \(\mathbf{X}_{1}=\mathbf{X}\) and \(\mathbf{y}_{1}=\mathbf{y}\). \(\mathbf{t}_{a}\) and \(\mathbf{I}_{a}\) are computed as  
\[\mathbf{t}_{a} =\mathbf{X}_{a}\frac{\mathbf{X}_{a}^{\mathrm{T}}\mathbf{y}_{a}}{|\mathbf{X}_{a}^{ \mathrm{T}}\mathbf{y}_{a}|} \tag{11}\] \[\mathbf{I}_{a} =\mathbf{X}_{a}^{\mathrm{T}}\frac{\mathbf{t}_{a}}{\mathbf{t}_{a}^{\mathrm{T}} \mathbf{t}_{a}}\] (12) \[b_{a} =\frac{\mathbf{t}_{a}^{\mathrm{T}}}{\mathbf{t}_{a}^{\mathrm{T}}\mathbf{t}_{a} }\mathbf{y}_{a}, \tag{13}\]  
where \(|\circ|\) is the L2 norm. \(a^{\prime}\) is the number of principal components to approximate the original data and is determined by cross-validation. In this study, \(a^{\prime}=3\), which is close to the typical value of two or three for approximating the gait motions of healthy people [26; 38; 39; 40]. To expand the above-mentioned partial least-squares regression to time-series data, \(x\) is replaced with a time-series vector. In the case of triaxial kinematic data for gait analysis, each sample includes the translational velocities \(v\) along the x, y, and z-directions for one gait cycle, that is, two succeeding steps. The time-series data are discretized by each percentage of the gait cycle.