Prediction of Lower Limb Kinetics and Kinematics during Walking by a Single IMU on the Lower Back Using Machine Learning

Recent studies have reported the application of artificial neural network (ANN) techniques on data of inertial measurement units (IMUs) to predict ground reaction forces (GRFs), which could serve as quantitative indicators of sports performance or rehabilitation. The number of IMUs and their measurement locations are often determined heuristically, and the rationale underlying the selection of these parameter values is not discussed. Using the dynamic relationship between the center of mass (CoM), the GRFs and joint kinetics, we propose the CoM as a single measurement location with which to predict the dynamic data of the lower limbs, using an ANN. Data from seven subjects walking on a treadmill at various speeds were collected from a single IMU worn near the sacrum. The data was segmented by step and numerically processed for integration. Six segment angles of the stance and swing leg, three joint torques, and two GRFs were estimated from the kinematics of the CoM measured from a single IMU sensor, with fair accuracy. These results indicate the importance of the CoM as a dynamic determinant of multi-segment kinetics during walking. The tradeoff between data quantity and wearable convenience can be solved by utilizing a machine learning algorithm based on the dynamic characteristics of human walking.


Introduction
Kinetic data of human motion, such as ground reaction forces (GRFs) and joint torques, serve as a quantitative indicator of sports performance or the effect of rehabilitation. GRFs and joint torques have served as indicators of injury risks and pain during running [1][2][3][4][5]. The progress of gait rehabilitation of hemiparetic patients has been monitored by the magnitude and the degree of asymmetry of the GRF [6][7][8]. The motion analysis used for these studies is performed in a laboratory with accurate motion trackers and force transducers, which often impose spatial and temporal constraints on subjects.
With the development of sensor technology, the market for wearable motion monitoring systems has rapidly grown. Wearable motion monitoring products such as Galaxy, Garmin, Fitbit, and Apple watches monitor the overall rough motion information, such as the steps, cadence, distance, etc. That rough motion information, however, makes it difficult to estimate the kinetic information for the analysis of injury risks and rehabilitation. One of the common issues of current wearables is how to increase the quantity and quality of the measured data with minimal increase in the system complexity. In particular, obtaining force information using wearables is very challenging due to the size and weight of the force transducers. Instead of direct measurement, researchers have attempted to predict unmeasured GRF information from joint kinematics information. Oh et al. predicted 3D GRFs Figure 1. Schematics of the procedures of lower limb kinetics prediction from a single inertial measurement unit (IMU) using machine learning. Processed IMU signals, such as velocity, position and time, were fed into input nodes of the network. The network has seven input nodes followed by one hidden layer of 20 nodes and 11 output nodes. The output predictions are the angle of the thigh, shank, and foot of the stance and swing leg and the joint torques of the hip, knee, and ankle and the vertical and horizontal ground reaction forces.

Finding the Location of a Single Measurement by the IMU
To examine which part of body motion should be measured to best estimate the lower limb joint kinetics, we used the dynamic characteristics of human walking. Since the outputs of the ANN are obtained from the weighted sum of the inputs passed through a nonlinear activation function, the equations of motion of the CoM during walking obtained from the recently published SLIP model [29] were approximated, to estimate the GRF and joint dynamics in terms of the weighted sum of the CoM. An approximation of the relationship between the CoM position and the state variables of the mass-spring model of human walking was formulated (Equation (1)). The GRFs and angles of multiple joints were formulated in terms of state variables (Equations (2) and (3)). Finally, the GRFs, ankle joint torque, and angles of the foot, shank, and thigh were formulated as the weighted sum of the CoM positions. These calculations led us to propose the CoM as the best measurement position for prediction of the unmeasured GRF and joint dynamics. A brief summary follows, the details of which can be found in the Appendix A.
The position of the CoM, x = (x, y), is expressed as model parameters of the SLIP, such as the leg angle , the leg length l, the radius of the curvy foot R, the offset of the ankle joint with respect to the foot d, and the ground contact position chs, with the subscript hs denoting values at the onset of the stance phase (Figure 2), as follows: Figure 1. Schematics of the procedures of lower limb kinetics prediction from a single inertial measurement unit (IMU) using machine learning. Processed IMU signals, such as velocity, position and time, were fed into input nodes of the network. The network has seven input nodes followed by one hidden layer of 20 nodes and 11 output nodes. The output predictions are the angle of the thigh, shank, and foot of the stance and swing leg and the joint torques of the hip, knee, and ankle and the vertical and horizontal ground reaction forces.

Finding the Location of a Single Measurement by the IMU
To examine which part of body motion should be measured to best estimate the lower limb joint kinetics, we used the dynamic characteristics of human walking. Since the outputs of the ANN are obtained from the weighted sum of the inputs passed through a nonlinear activation function, the equations of motion of the CoM during walking obtained from the recently published SLIP model [29] were approximated, to estimate the GRF and joint dynamics in terms of the weighted sum of the CoM. An approximation of the relationship between the CoM position and the state variables of the mass-spring model of human walking was formulated (Equation (1)). The GRFs and angles of multiple joints were formulated in terms of state variables (Equations (2) and (3)). Finally, the GRFs, ankle joint torque, and angles of the foot, shank, and thigh were formulated as the weighted sum of the CoM positions. These calculations led us to propose the CoM as the best measurement position for prediction of the unmeasured GRF and joint dynamics. A brief summary follows, the details of which can be found in the Appendix A.
The position of the CoM, x = (x, y), is expressed as model parameters of the SLIP, such as the leg angle θ, the leg length l, the radius of the curvy foot R, the offset of the ankle joint with respect to the foot d, and the ground contact position c hs , with the subscript hs denoting values at the onset of the stance phase (Figure 2), as follows: x and y are the positions of the horizontal and vertical axes, respectively, and the model state variables l and θ are the total spring length and angle, respectively, with a positive sign in the clockwise direction with respect to the vertical. The length and stiffness of the CoM-ankle are represented by la and ka, respectively. The offcentered curvy foot parameters d and R are the offset from the ankle from the center of the foot and the radius of the curvy foot, respectively. The forces at the ground contact point are presented by fx and fy in the horizontal and vertical directions, respectively. The constraint force, fa, acts on the ankle in the vertical direction of the spring force, and the constraint ankle torque, Ta, is generated to prevent foot rotation relative to the leg rotation, with a positive sign in the extension direction. With the kinematic constraint and the positions of the CoM and the ankle, (C) multibody segment angles are determined by inverse kinematics. Segment lengths and angles are represented by l and θ, respectively, with subscripts 1-3, and angles are positive in the clockwise direction.
By linearization, we can obtain the mapping equation from the CoM position x to the augmented state z =l,,l  , which includes the state parameters, l and  of the SLIP model, and its nonlinear term l as follows: Wb (1) where α = Rhs − chs + d.  (2) The relationship between the augmented state z and the augmented segment angles  = [θ1, θ2, θ3] T can be rewritten as follows: x and y are the positions of the horizontal and vertical axes, respectively, and the model state variables l and θ are the total spring length and angle, respectively, with a positive sign in the clockwise direction with respect to the vertical. The length and stiffness of the CoM-ankle are represented by l a and k a , respectively. The off-centered curvy foot parameters d and R are the offset from the ankle from the center of the foot and the radius of the curvy foot, respectively. The forces at the ground contact point are presented by f x and f y in the horizontal and vertical directions, respectively. The constraint force, f a , acts on the ankle in the vertical direction of the spring force, and the constraint ankle torque, T a , is generated to prevent foot rotation relative to the leg rotation, with a positive sign in the extension direction. With the kinematic constraint and the positions of the CoM and the ankle, (C) multibody segment angles are determined by inverse kinematics. Segment lengths and angles are represented by l and θ, respectively, with subscripts 1-3, and angles are positive in the clockwise direction.
By linearization, we can obtain the mapping equation from the CoM position x to the augmented state z = [l, θ, lθ] T , which includes the state parameters, l and θ of the SLIP model, and its nonlinear term lθ as follows: The relationship between the augmented state z and the augmented kinetics f = [f x , f y , T a ] T such as GRF and ankle torque are obtained as follows: Sensors 2020, 20, 130 5 of 29 The relationship between the augmented state z and the augmented segment angles θ = [θ 1 , θ 2 , θ 3 ] T can be rewritten as follows: In summary, the kinematics and kinetics of the lower limb can be approximated as a weighted sum of augmented states z of the spring mechanics and the bias as follows: From the approximated relationship between the CoM position x and the lower limb kinetics and kinematics [f T , θ T ] as a form of weighted sum of x and bias, we tried to predict the unmeasured lower limb kinematics and kinetics using the sacrum motion profiles, the measurable approximation of the CoM. We may use the terms 'CoM' and 'sacrum' interchangeably throughout the manuscript (see also the Appendix A for the formulation of the lower limb kinematics and kinetics as an approximation of the weighted sum of the CoM states).

Experimental Protocols and Data Collection
Seven young (25.0 ± 2.9 years) healthy male subjects with average heights and weights of 168.8 ± 7.5 cm and 72.1 ± 7.7 kg, respectively, participated in the data collection, after signing fully informed consent forms approved by KAIST Institutional Review Board on 22 January 2018 (IRB-18-023). The subjects participated in one-minute treadmill walking trials at three customized gait speeds, ranging from 1.0 m/s to 2.3 m/s. The test speeds were determined by each subject's preferred gait speed, the midpoint between his preferred and maximum speeds, and a speed slower than the preferred speed by an amount equal to approximately 20% of the difference between the preferred and maximum speeds. The label 'slow', 'moderate', and 'fast' indicates slow speed, preferred speed, and fast speed, respectively. The speed of each subject is in Table 1. Three one-minute gait trials per test speed were randomly ordered. The kinematics of the lower back, which is assumed to be close to the CoM, were measured by an IMU (EBIMU-9DOFV4) with an acceleration sensitivity of ±16 g, and resolution of 0.001 g. The width of the IMU was 15.6 mm, and the length was 18.6 mm. The IMU data were collected at a sampling frequency of 100 Hz, and used the universal asynchronous receiver/transmitter (UART) protocol for communication with a baud rate of 921,600 bps. The IMU was aligned with the subject's heading direction to set the local frame of the IMU to the body frame. The 2D acceleration data from the IMU in the sagittal plane were used for the analysis. The lower limb kinematics and ground reaction forces under each foot were measured by an optical motion capture system (Motion Analysis Hawk ® ) and force plates (Bertec, FP 6012 ® ), respectively. For the segment angle calculation, 12 markers were located at the toe, heel, ankle, knee, and anterior and posterior superior iliac spine (ASI and PSI) of both legs. The midpoint of the PSIs was used as the reference position for the IMU. Motion data and force data were collected at the sampling frequency of 100 Hz and 200 Hz, respectively, and filtered with a fifth-order low-pass Butterworth filter with a cutoff frequency of 10 Hz. The joint torque was calculated by the inverse dynamics of the seven-segmental rigid body model, which includes the hip, knee, and ankle. The mass distribution of the rigid body were based on the anatomical reports [31].

Preprocessing of IMU Data for ANN Input
The IMU measurement data that contain the double support phase must be segmented into one stance phase per leg before feeding into the ANN for GRF prediction. The stance phase is defined based on the GRF as the duration from the onset of heel strike (HS) to toe off (TO) ( Figure 3A,C). To identify the characteristic features of the IMU data which indicate the gait events of HS and TO, the IMU acceleration and GRF profiles were aligned over time, and three characteristic features were chosen: the CoM apex, HS, and TO ( Figure 3). The filter frequency of IMU acceleration was selected heuristically to best extract the characteristic features, and was chosen to be 40 Hz for identifying HS and 10 Hz for identifying TO. The CoM apex is the point at which the height of CoM is maximal in sagittal plane around the mid-stance phase, and the corresponding vertical GRF has a local minimum ( Figure 3B). The CoM apex was identified as the local minimum of the vertical IMU acceleration data filtered with a cutoff frequency of 10 Hz ( Figure 3B). The HS was identified as the first local minimum of the vertical acceleration filtered by a 40 Hz cutoff frequency, following the positive local maximum after CoM apex ( Figure 3B). The TO was identified as the first local minimum of the horizontal acceleration data filtered by a 10 Hz cutoff frequency after CoM apex ( Figure 3D). The detailed process of IMU data segmentation is shown in Figure 3E. The precedence of HS over TO was checked for validation. When HS estimation follows the TO, the HS was re-identified as the preceding local minimum of vertical acceleration compared with the one originally found. After completing segmentation of the full stance, the consecutive segmentation process was performed by identifying the subsequent CoM apex based on the acceleration profiles from the TO event, the last event of the previous step.
Another step of preprocessing the IMU data is to calculate the velocity and position of the CoM, i.e., the sacrum, from the acceleration, which includes drift removal and identification of the integral constant. To eliminate the integration drift, steady walking and random error were assumed. Then, the drift was removed linearly with time over the duration from one apex to the following apex.
The drift-removed velocity of the CoMv(t) could be obtained as follows: where v (t) is the time integration of the acceleration data for time t over the duration of the stance phase T. To satisfy the steady walking assumption, the mean velocity change was set to zero so that position drift due to the random error was simultaneously removed. The velocity and the position offsets were estimated in a heuristic manner to obtain the velocity and the position from the acceleration integration. From the observed proportionality of the gait speed v with the average magnitude of the acceleration, A, and the gait frequency, f, (0.79 and 0.73 Pearson correlation coefficients, respectively), the gait speed offset v 0 was estimated as follows: where a, b, and c are the regression coefficients obtained in a least square manner. To minimize the size effect of subjects, we normalized the speed with h · g, where h and g are the height of the subject and gravitational acceleration, respectively. The position of the sacrum was obtained by integrating the compensated velocity with the position offset, defined as the subject's leg length. All data processing was conducted by MATLAB R2018a (Mathworks, Inc., Natick, MA, USA).

Structure of the ANN and Its Training and Test Procedures
The GRFs and segment angles of the lower limbs were approximated using the weighted sum of CoM position (Equation (4)), so we hypothesized that a simple ANN that generates its output from the weighted sum of inputs could be a good candidate for a prediction network. We used a fully

Structure of the ANN and Its Training and Test Procedures
The GRFs and segment angles of the lower limbs were approximated using the weighted sum of CoM position (Equation (4)), so we hypothesized that a simple ANN that generates its output from the weighted sum of inputs could be a good candidate for a prediction network. We used a fully connected feed-forward ANN with one input layer, one hidden layer and one output layer. The seven input nodes consist of a 7 × 1 column vector consisting of a time t, and the corresponding CoM kinematics (t, x, y, v x , v y , a x , a y ), where x, y are horizontal and vertical position, respectively, and v and a are velocity and acceleration measured at a specific time, such as t = t 0 . The position of the CoM was reset to be zero at every HS, and we defined 'displacement' as the reset position. The eleven output nodes constitute an 11 × 1 column vector consisting of the GRFs and the joint kinetics and kinematics corresponding to t = t 0 , such as (GRF x , GRF y , T ank , T knee , T hip , θ foot,stance , θ shank,stance , θ thigh,stance , θ foot,swing , θ shank,swing , θ thigh,swing ), where T, and θ are joint torque and segment angle, respectively, with subscripts showing the names of the lower limb segment. All kinematic data of the sacrum, such as the acceleration, velocity and displacement, were normalized by its maximum norm such that the values ranged from 0 to 1. The hidden layer has 20 nodes, and node weights are iteratively updated to best match the data. The sigmoid and linear activation transfer functions were used in the hidden and output layers. The loss function of the network was set as the mean square error of the predicted angles of both the stance and swing legs, the joint torque of the stance leg and the GRF with respect to the observed data. The network employed an Adam optimizer for the back-propagation process. To match the size of the training dataset collected at various gait speeds of many subjects, gait data of one stance phase were interpolated into 200 points. As a quantitative measure of prediction performance, we used the normalized root mean square error (NRMSE), which was normalized by the difference between the maximum and minimum values of the outputs of each subject at each speed. The NRMSE value of each subject represents the error of a total of 30 trials at one speed for one output. The neural network programming was conducted using PyTorch 0.4.1 [32].
To test the proposed prediction methods, we used LOO cross validation. A total of 540 stance phase trials obtained from six subjects was used as the training dataset, and the results were validated using data from the remaining subject. Validation was performed for each of seven subjects. To evaluate which kinematic information contributed most to the prediction accuracy, various combinations of kinematic inputs, such as CoM kinematics with or without position, velocity, and acceleration, were used as the inputs to the ANN. To examine whether the walking speed at which the training datasets were collected affected the prediction accuracy of the GRFs at other speeds, the data collected at various walking speeds were used as training datasets of the ANN.

Results
From the single IMU measurement at the lower back, the practical approximate of the CoM, eleven joint kinetics data points were predicted by the ANN over walking speeds ranging from 1.0 to 2.3 m/s. The predicted joint data are the thigh, shank and foot angles of the swing and stance leg, the hip, knee, and ankle joint torques of the stance leg, and the vertical and A-P ground reaction forces.
Gait event detection from the IMU data showed a reasonable match with those defined from the GRF, with average detection errors of the HS and TO of 0.025 s and 0.014 s, respectively ( Table 2). The detection algorithm showed a more accurate performance for TO than that for HS, with the accuracy increasing with the gait speed. With drift removal and the estimation of the integral constant of the IMU measured acceleration, the velocity and displacement trajectories of the CoM were obtained by integrating the IMU data based on the detected gait events (Figure 4). The calculated sacrum speed and the displacement segmented by the stance phase matched well with the data obtained from optical markers, with NRMSEs of approximately 25% and 12%, respectively (Table 3). No significant tendencies were observed for the estimation errors of the sacrum kinematics as a function of the gait speed.  The prediction results (NRMSE and R 2 ) of each subjects are presented in Tables A1 and A2 in the Appendix A. The lower-limb joint kinematics and kinetics at various gait speeds for all seven subjects were fairly well predicted from the kinematics of the sacrum, such as the displacement, velocity, and acceleration, with an average NRMSE of approximately 7% (Figures 5 and 6, Table 4). The segment angles of the stance and swing leg showed a reasonably good match between the prediction and experimental data, even for the worst case test set ( Figure 6A,B, Table 5), whereas the joint torques and GRFs showed larger error, especially for the estimation of the high-frequency component of the hip joint torques ( Figure 6C,D, Table 5). The average prediction errors also showed a similar trend in the prediction accuracy ( Figure 5, Table 4). As long as the displacement information of the sacrum is provided, the velocity and acceleration input to the ANN do not significantly improve the prediction accuracy, whereas the omission of the displacement input significantly (p < 0.05) increases the prediction errors of the segment angles, torques and GRFs up to 3% (Figure 7, Table 6). These results imply that the position of the sacrum contributes most to the predictions of the lower-limb joint kinetics.
The ANN predicted the kinematics and kinetics of the lower limb at various gait speeds that were not used for the training datasets. The joint kinematics and kinetics are known to vary as a function of the walking speed. To examine the effect of the walking speed at which the training data sets are collected to the prediction accuracy, the ANN was trained by the data collected at three different speed conditions: at moderate speed only, at slow and moderate speeds, and at slow, moderate and fast speeds (Figure 8). The results showed that the inclusion of the data of a specific gait speed, which are to be predicted through the ANN, would slightly but not significantly increase the prediction performance. More importantly, the sacrum kinematics measured at the moderate  The prediction results (NRMSE and R 2 ) of each subjects are presented in Tables A1 and A2 in the Appendix A. The lower-limb joint kinematics and kinetics at various gait speeds for all seven subjects were fairly well predicted from the kinematics of the sacrum, such as the displacement, velocity, and acceleration, with an average NRMSE of approximately 7% (Figures 5 and 6, Table 4). The segment angles of the stance and swing leg showed a reasonably good match between the prediction and experimental data, even for the worst case test set ( Figure 6A,B, Table 5), whereas the joint torques and GRFs showed larger error, especially for the estimation of the high-frequency component of the hip joint torques ( Figure 6C,D, Table 5). The average prediction errors also showed a similar trend in the prediction accuracy ( Figure 5, Table 4). As long as the displacement information of the sacrum is provided, the velocity and acceleration input to the ANN do not significantly improve the prediction accuracy, whereas the omission of the displacement input significantly (p < 0.05) increases the prediction errors of the segment angles, torques and GRFs up to 3% (Figure 7, Table 6). These results imply that the position of the sacrum contributes most to the predictions of the lower-limb joint kinetics.
The ANN predicted the kinematics and kinetics of the lower limb at various gait speeds that were not used for the training datasets. The joint kinematics and kinetics are known to vary as a function of the walking speed. To examine the effect of the walking speed at which the training data sets are collected to the prediction accuracy, the ANN was trained by the data collected at three different speed conditions: at moderate speed only, at slow and moderate speeds, and at slow, moderate and fast speeds (Figure 8). The results showed that the inclusion of the data of a specific gait speed, which are to be predicted through the ANN, would slightly but not significantly increase the prediction performance. More importantly, the sacrum kinematics measured at the moderate speed could predict the eleven unmeasured joint kinetics at various walking speeds, which were not included in the training of the ANN (Figure 8). The prediction results of the joint torques and GRFs at the various gait speeds obtained from the ANN trained by trials of moderate speed showed a reasonable match with the data (Figure 9, Table 7).

Discussion
Biomechanical knowledge about walking dynamics was used to design an ANN to predict unmeasured motion data from a single IMU measurement. Equations related to CoM dynamics were reformulated to approximate the GRF and lower limb joint dynamics as the weighted sum of the CoM, the basic input-output relationship of the ANN. Based on this approximation, the ANN was designed to have input nodes receiving the kinematics of the CoM, with the layers being fully connected. Using this ANN, eleven lower limb joint dynamics and the GRF data at various walking speeds were predicted from a single and only IMU measurement near the sacrum. The predicted segment angles, joint torques, and GRF data showed a reasonably good match with the experimental data despite the very simple network structure ( Figure 1) with a relatively small number of training data sets (Figures 5-8). Although trained with the CoM data measured at a specific gait speed, the ANN could predict multiple joint kinetics at various gait speeds that have not been used for network training. These results imply that the CoM is important as a dynamic determinant of multi-segment kinetics during walking, and that the tradeoff between data quantity and wearable convenience can be solved by utilizing a machine learning algorithm based on the dynamic characteristics of human walking.
Based on a biomechanical analysis using 2D spring mechanics, we showed that the 2D acceleration data of the CoM could predict the GRFs and joint dynamics in 2D fairly well. Due to the point mass assumption of the CoM, no rotational motion of the CoM was assumed, and no gyroscope data were used for prediction. We may improve the quality of the prediction heuristically with the full usage of IMU sensors, such as 3D accelerations and 3D angular rates from the gyroscope. We do not, however, have sound explanations as to why the addition of 3D data would help in obtaining better prediction. We are investigating whether the 3D GRF could be well formulated using 3D CoM kinematics as a form of simple mechanics, as was the case for 2D. Despite the lack of quantitative formulation, there is considerable experimental evidence to show that full usage of the IMU data could predict the unmeasured GRFs, metabolic cost, and joint angles [9,33,34]. Three dimensional acceleration and angular rate measured from the 17 IMUs attached to the whole body have been used to successfully predict 3D GRF data [9]. Similar to our study, the 2D joint angles of lower limbs have been estimated from IMU measurement at the foot, shank, and thigh [34]. The 3D acceleration and gyroscope data were processed using general regression neural networks and the authors found that foot acceleration was sufficient to predict the other joint angles. Despite the larger number of sensors, the prediction errors were slightly larger (4~9 degrees) than ours (2~3 degrees), implying that the CoM kinematics are closely related to lower limb kinematics. The predicted quantities were not limited to the GRF or joint kinematics. Metabolic cost is one of the key determinants for human motion [35][36][37][38]. Measurements of energy expenditure, however, are limited to indoor settings, and the measurement procedures are often tedious. Researchers attempted to estimate the metabolic cost during walking from IMU data measured at the hip joint using Bayesian linear regression [33]. With the support of quantitative biomechanical analysis of the target motion to be predicted, the IMU data could be used more effectively and efficiently. There are studies that use machine learning technology to predict human motion, independent of gait biomechanics, but often with the cost of increased measurement complexity of the prediction errors. Oh et al. successfully predicted the joint torques and GRFs of human walking using a general regression neural network (GRNN), with the cost of increased measurement of the whole body joint kinematics using an optical marker system [10], which is therefore limited for on-field applications. The errors of the vertical and A-P GRF were comparable to ours (Table 8). Recently, without a great loss of prediction accuracy, the number of IMUs were reduced to two, attached at each shank of the subject [11]. To complement the reduced sensor information, a feed-forward neural network was used to produce the best estimate of the GRFs. In this recent study, however, the training datasets included part of the data from the test subject [11], which may have decreased the prediction accuracy when LOO validation was performed. Because the GRF pattern could differ between subjects [39,40], we used a conservative validation method, by excluding the test subjects' data from the training dataset, and following the LOO validation method to guarantee generalization.  [41]. The predicted joint angles of lower limbs showed relatively small errors, in the range of 3~5 degrees, compared to others, but the test data were included in the training data set, which would produce less-accurate predictions when LOO validation is performed. Joint torque predictions have been performed using joint angle data [42] as well as surface electromyography (sEMG) [43]. The errors of the joint torque predictions were larger than those of GRFs and/or joint angles, as was observed in our study (Table 10). Ardestani et al. used 14 sEMG data and showed relatively high accuracy in joint torque prediction [43], but, again, the inclusion of test data in the training dataset of the four subjects may contribute to a misleadingly high prediction accuracy. Compared to those of related works, we predicted more data at various walking speeds with fewer measurements from the complementary combination of the biomechanics and the ANN. The proposed prediction method using a combination of biomechanics and machine learning could be used to resolve the tradeoff between data richness and measurement convenience inherent in wearable technologies. The rehabilitation success of stroke patients can be assessed using the A-P Sensors 2020, 20, 130 20 of 29 direction of the GRF of the nonparametric leg [6][7][8], which would be challenging to measure if patients have to regularly visit a motion analysis laboratory and undergo almost over an hour of data collection procedures. Motion monitoring of those patients using simple wearables during their clinical visits or even at home in daily life would greatly increase patients' convenience, if the data from wearables could produce data of quality and quantity comparable to those obtained in the lab. A recent study demonstrated successful long-term monitoring of the range of motion (ROM) of the shoulders of elderly people using the IMU, which showed a fairly good match with data from optical motion capture systems [44]. Similarly, onsite monitoring of lower limb joint loading could be estimated from a single IMU sensor for patients with rheumatoid arthritis or anterior cruciate ligament recovery. As a quantitative comparison, the asymmetry of the GRFs of mild hemiparesis patients seemed to be approximated as 50 N [7]. Kesar et al. proposed minimal detectable changes of the kinematic and kinetic indicators of rehabilitation assessment of post-stroke patients to be about 3-10 degrees and 20-33 N, respectively [45]. In our study, the average prediction errors of segment angles are about (3.1 degree) thigh, (2.2 degree) shank and (3.4 degree) foot, respectively. The approximate errors of vertical and AP GRF are 58 N and 23 N, and those of hip, knee, ankle joint torques are 16.7 Nm, 11.4 Nm, and 15.3 Nm, respectively. Considering the observed changes of GRFs and joint kinematics of stroke patients, the proposed method could serve well for kinematics monitoring, whereas the GRF monitoring seemed to be feasible for A-P GRF but need improvement for vertical GRF. In addition to the rehabilitation research, the number of studies monitoring sports activity is rapidly increasing. Wearable monitoring systems are widely used for monitoring the effect of posture or loading conditions on kinetics [46]. Data from sacrum-mounted IMUs were used to analyze the effect of load carriage and fatigue on jumping motions. The change of kinematics due to fatigue was similar to that under load conditions. IMUs on the lower body were used to estimate joint kinematics and GRF during running [47,48]. Combining single IMU data from the sacrum with a convolutional neural network, Johnson et al. estimated 3D GRF with the error of 4-9% during running and side-stepping [14]. A leg-mounted IMU combined with a feed forward neural network predicted knee joint forces with an error rate of 17% during running, changing direction, and jumping [49]. Issues to be addressed for sports wearables include reducing the number of sensors and enhancing prediction accuracy. The difference in kinetics with different foot strike patterns in running is about 10%-20% of joint torque and 30% of GRF. [50]. A recent study showed that a magnitude of about 30% of GRF in response to the breaking force could serve as an indicator for injury risk [51]. By assuming reference force and length to be the average weight and height of our test subjects, the observed kinetic differences range approximately 200 N and 18 Nm, respectively. Considering the error range of our study, the proposed prediction method can be applicable for monitoring running. The limitations of our study, and technical issues to be resolved, lie in the main components of the prediction method, such as the biomechanical model, the processing of the IMU data, and the design of the ANN. First, the proposed prediction method utilized a spring-based biomechanical model of CoM dynamics during steady gait trials in 2D. Pathological gaits, however, often show increased asymmetry in the medial-lateral (M-L) as well as the A-P direction [52]. Thus, further research is required into applying the proposed method to the monitoring of unsteady or abnormal walking in 3D. Biomechanical modeling of human walking with a passive spring and asymmetric foot could well emulate the human-like ankle joint torque [29], implying that the CoM and ankle joint torque are mechanically coupled. The biomechanical model of the point mass CoM, however, could not represent the activation of hip joint torques, so there is no mechanical and kinematical correlation between the CoM motion and the hip joint torques. Since the ANN used for joint torque prediction was designed based on the biomechanics of the CoM, the larger prediction errors in hip joint torque presented in Figure 9A could be attributed to the limited representability of the hip joint torque by the CoM. Secondly, one of technical challenges of handling IMU data is the determination of appropriate and accurate segmentation. For long-term IMU data from motion monitoring, the very first step to perform in signal processing is to extract meaningful data from a whole data collection. A recent study used GPS information to specify the data segment of the action of interest out of approximately two weeks of IMU data in daily life [24]. To expand the proposed prediction method to gait monitoring during daily life, an appropriate activity recognition and extraction process must be undertaken. Once the data were roughly segmented for a specific action, such as walking in the hallway, further segmentation was carried out based on specific gait events, such as HS or TO. IMU data gathered near a ground contact such as the shank or foot had more accurate event detection results than when the IMU was mounted near the sacrum or lower back [53,54] possibly due to the distal measurement of the indicating signal around the foot [54][55][56]. Similar to previous studies, the identification of gait events and segmentation of the data were performed in a heuristic manner, but in a simpler way than those presented in other studies, such as those using wavelet transforms [57]. For example, the filtering frequency of the IMU data was selected heuristically by choosing the one that did not smooth the indicating signal of the gait event, while not including unnecessary peaks. The proposed event detection method produced errors comparable to those of previous studies, ranging from approximately 0.01-0.03 seconds [53,54]. The detection errors were reduced with gait speed due to the highlighted singular motion of the sacrum at faster speeds [28]. The gait event detection mostly follows the real gait event due to the delayed sacrum acceleration profile compared to that of the true CoM acceleration, which is obtained from the GRF, and resulted in shifted estimation of sacrum trajectory. Likewise, although small in magnitude, event detection errors could greatly affect the prediction results when erroneously segmented IMU data are used for training an ANN. An ANN is highly dependent on its input information, so poorly-estimated sacrum kinematics cause a deterioration in prediction accuracy. Thus, the prediction performance is sensitively dependent on the signal quality of the IMU and the following signal processing. Considering that supervised learning is performed between the segmented IMU data and the segmented GRF, erroneously segmented IMU data distort the coupling relationship between the CoM and the GRF, and the network will be trained inaccurately and inefficiently. Further, there was a relatively large error (~25%) in the velocity estimation ( Figure 4). The integration of the acceleration had drift removal and offset compensation issues, which were resolved in a heuristic manner. Thorough examination of the best estimation of the gait speed was not performed in this study, so the error rate could potentially be reduced by employing other estimation algorithms or machine learning techniques. The effects of erroneous velocity estimation on the prediction results, however, are reduced when the velocity is fed into the ANN after normalization. In addition, the validation conducted depends on the location of IMU. Since the trunk has rotational motion with respect to the vertical axis, the location of IMU in the horizontal axis can affect the prediction performance. Thus, the sensitivity analysis about the location of IMU should be performed for validation of the reliable wearable motion monitoring system. Lastly, the proposed ANN has much room for improvement of the prediction accuracy. The proposed ANN is a very simple network, designed using a simplified mathematical approximation between the GRFs and the CoM. A more realistic approximation could be achieved by adding complexity to the network, such as additional input nodes or hidden layers, or a change of activation function. A relatively small number of training datasets, and reasonably good test results compared to previous studies with more sensor data could support an efficiently trained ANN. However, a small training data set could result in the limited performance accuracy of the network, which was observed by an increase error in response to variations of the input. Since the GRF data are known to have high intra-subject variability, training with a larger dataset should produce improved prediction results for diverse test subjects.
Biomechanical knowledge about walking dynamics was used to design the ANN that predicts the unmeasured motion data from only a single IMU measurement. From the sacrum-mounted IMU data, the motion and force information of multiple joints were well predicted at various gait speeds by a simple ANN with a relatively small number of training datasets. This is attributed to the role of the CoM as a key dynamic descriptor of multi-segment lower limb coordination during the human walking behavior. The reliable prediction performance presented in this study implies the complementary roles of biomechanics domain knowledge and machine learning technology in predicting motion kinetics.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
To examine whether the joint kinetics could be predicted through the ANN with CoM kinematics as the input, we tried to approximate the joint kinetics as a form of weighted sum of the CoM kinematics. Noting that the CoM dynamics are demonstrated by the spring mechanics, we started with the equations of motion of the recently published spring loaded inverted pendulum (SLIP) model [29]. The position of the CoM, x = (x, y), is expressed as model parameters of the SLIP, such as the leg angle θ, the leg length l, the radius of the curvy foot R, the offset of the ankle joint with respect to the foot d, and the ground contact position c hs , with the subscript hs denoting values at the onset of the stance phase ( Figure A1), as follows: By applying the small angle approximation and binomial approximation and neglecting higher-order terms in Equation (A1), we can obtain the linearized mapping equation from the CoM position x to the augmented state z = [l, θ, lθ] T , which includes the state parameters, l and θ of the SLIP model, and its nonlinear term lθ as follows: where α = Rθ hs + d − c hs . kinematics as the input, we tried to approximate the joint kinetics as a form of weighted sum of the CoM kinematics. Noting that the CoM dynamics are demonstrated by the spring mechanics, we started with the equations of motion of the recently published spring loaded inverted pendulum (SLIP) model [29]. The position of the CoM, x = (x, y), is expressed as model parameters of the SLIP, such as the leg angle θ, the leg length l, the radius of the curvy foot R, the offset of the ankle joint with respect to the foot d, and the ground contact position chs, with the subscript hs denoting values at the onset of the stance phase ( Figure A1), as follows: Figure A1. Compliant walking model with an off-centered, curvy foot combined with a springy footankle segment [29]. (A) The model parameters and state variables of the model, and (B) the force and torque components. The position of the CoM in the sagittal plane; x and y are the positions of the horizontal and vertical axes, respectively, and the model state variables l and θ are the total spring length and angle, respectively, with a positive sign in the clockwise direction with respect to the vertical. The length and stiffness of the CoM-ankle are represented by la and ka, respectively. The offcentered curvy foot parameters d and R are the offset of the ankle from the center of the foot and the radius of the curvy foot, respectively. The forces at the ground contact point are presented by fx and fy in the horizontal and vertical directions, respectively. The constraint force, fa, acts on the ankle in the vertical direction of the spring force, and the constraint ankle torque, Ta, is generated to prevent foot rotation relative to the leg rotation, with a positive sign in the extension direction. With the kinematic constraint and the positions of the CoM and the ankle, (C) multibody segment angles are determined by inverse kinematics. Segment lengths and angles are represented by l and θ, respectively, with subscripts 1~3, and angles are positive in the clockwise direction.
By applying the small angle approximation and binomial approximation and neglecting higherorder terms in equation A1, we can obtain the linearized mapping equation from the CoM position x Figure A1. Compliant walking model with an off-centered, curvy foot combined with a springy foot-ankle segment [29]. (A) The model parameters and state variables of the model, and (B) the force and torque components. The position of the CoM in the sagittal plane; x and y are the positions of the horizontal and vertical axes, respectively, and the model state variables l and θ are the total spring length and angle, respectively, with a positive sign in the clockwise direction with respect to the vertical. The length and stiffness of the CoM-ankle are represented by l a and k a , respectively. The off-centered curvy foot parameters d and R are the offset of the ankle from the center of the foot and the radius of the curvy foot, respectively. The forces at the ground contact point are presented by f x and f y in the horizontal and vertical directions, respectively. The constraint force, f a , acts on the ankle in the vertical direction of the spring force, and the constraint ankle torque, T a , is generated to prevent foot rotation relative to the leg rotation, with a positive sign in the extension direction. With the kinematic constraint and the positions of the CoM and the ankle, (C) multibody segment angles are determined by inverse kinematics. Segment lengths and angles are represented by l and θ, respectively, with subscripts 1~3, and angles are positive in the clockwise direction.
The relationship between the augmented state z and the GRF and ankle torque are obtained from the CoM of the SLIP model as follows [29]: where the spring parameters k and l 0 are the stiffness and the rest length of the spring, respectively, and the geometric parameters of R and d are the radius of the curvy foot and the offset of the ankle joint from the center of the curved foot, respectively. l a and k a are the length and stiffness of the CoM-ankle, as shown in Figure A1. The constants c l and c k are the ratio of the rest length of l a to the rest length of l and the ratio of the total leg stiffness k to the CoM-ankle stiffness k a , respectively ( Figure A1). The ankle torque T a is expressed as f a and l a , which are the constraint force applied to the ankle joint and the distance from the ankle to the CoM, respectively. By applying a small angle approximation, a small offset d approximation compared to the foot length R, such as d/R << 1, and the small length deviation approximation of the leg length l, such as l = l 0 (1 + δ), with δ being assumed to be small, then f a is approximated as follows: and the above Equation (A3) can be rewritten as follows: If we ignore higher-order terms of small value, such as θ 2 ≈ 0, then the augmented kinetics f = [f x , f y , T a ] T can finally be obtained as follows: With the specified positions of the CoM and the ankle joint, the segment angles of the thigh and shank can be kinematically constrained as a function of the augmented spring states θ as follows [29]: sin θ 2 = 0.5l a sin θ · l 2 −1 ± cos θ 1 − (0.5l a l 2 −1 ) 2 sin θ 3 = (l a sin θ − l 2 sin θ 2 )l 2 −1 (A8) where subscripts 1, 2, and 3 denote the foot, shank, and thigh, respectively, and l 1 is the same as l 2 . By applying small angle approximation and a small ankle length compared to the leg lengths, such as l a < 2l 1 , the segment angles θ = [θ 1 , θ 2 , θ 3 ] T can be rewritten as follows: (1 − c k )l 1 In summary, the kinematics and kinetics of the lower limb can be approximated as a weighted sum of augmented states z of the spring mechanics and the bias as follows: From the approximated relationship between the CoM position x and the joint kinematics and kinetics [f T , θ T ], we tried to predict the 11 unmeasured joint kinetics during walking from a single IMU measurement near the CoM.