Evaluating Electromyography and Sonomyography Sensor Fusion to Estimate Lower-Limb Kinematics Using Gaussian Process Regression

Research on robotic lower-limb assistive devices over the past decade has generated autonomous, multiple degree-of-freedom devices to augment human performance during a variety of scenarios. However, the increase in capabilities of these devices is met with an increase in the complexity of the overall control problem and requirement for an accurate and robust sensing modality for intent recognition. Due to its ability to precede changes in motion, surface electromyography (EMG) is widely studied as a peripheral sensing modality for capturing features of muscle activity as an input for control of powered assistive devices. In order to capture features that contribute to muscle contraction and joint motion beyond muscle activity of superficial muscles, researchers have introduced sonomyography, or real-time dynamic ultrasound imaging of skeletal muscle. However, the ability of these sonomyography features to continuously predict multiple lower-limb joint kinematics during widely varying ambulation tasks, and their potential as an input for powered multiple degree-of-freedom lower-limb assistive devices is unknown. The objective of this research is to evaluate surface EMG and sonomyography, as well as the fusion of features from both sensing modalities, as inputs to Gaussian process regression models for the continuous estimation of hip, knee and ankle angle and velocity during level walking, stair ascent/descent and ramp ascent/descent ambulation. Gaussian process regression is a Bayesian nonlinear regression model that has been introduced as an alternative to musculoskeletal model-based techniques. In this study, time-intensity features of sonomyography on both the anterior and posterior thigh along with time-domain features of surface EMG from eight muscles on the lower-limb were used to train and test subject-dependent and task-invariant Gaussian process regression models for the continuous estimation of hip, knee and ankle motion. Overall, anterior sonomyography sensor fusion with surface EMG significantly improved estimation of hip, knee and ankle motion for all ambulation tasks (level ground, stair and ramp ambulation) in comparison to surface EMG alone. Additionally, anterior sonomyography alone significantly improved errors at the hip and knee for most tasks compared to surface EMG. These findings help inform the implementation and integration of volitional control strategies for robotic assistive technologies.


INTRODUCTION
The continued progression of robotic (i.e., powered) lower-limb assistive devices with enhanced sensing, control and actuation has the ability to improve personalized rehabilitation strategies, restore independence and improve quality of life for individuals with mobility disorders and lower-limb loss. Significant improvements of these devices include reduction in their overall weight, improved energy-transfer capabilities via injection and recycling, and the increase in overall number of powered joints (Price et al., 2019). However, when considering their physical human-robot interaction, the increased complexity of these devices, including the number of powered joints and their capabilities, is met with an increase in the degrees-of-freedom of the overall control problem. Therefore, sensing modalities that can accurately detect motion of multiple degrees-of-freedom for precise mapping to multiple joints' control parameters are critical for device translation and clinical adoption.
A breadth of wearable sensing modalities and control schemes have been explored for the accurate recognition of user movement intent. The primary control schemes include classifier-based control and volitional control, as well as some combination of the two, often referred to as shared or indirect volitional control. Classifier-based control typically identifies the mode of ambulation or specific task Huang et al., 2011;Young et al., 2014). Conversely, volitional control allows the user to execute voluntary movements by more directly mapping user intent to a device output, such as angle, velocity or torque (Ha et al., 2011;Kannape and Herr, 2014). Recent studies that evaluated online performance of shared volitional (i.e., semivolitional) control of powered lower-limb prostheses have proven successful for providing additional capabilities to the user, such as crossing over obstacles during level walking (Mendez et al., 2020). However, to enable users the access to the full range of capabilities of powered assistive devices, such as use during weight-bearing and non-weight-bearing, it may be best to enable a more fullyvolitional control scheme, and one that incorporates intent recognition of multiple degrees-of-freedom. When there are multiple powered joints, such as knees and ankles in a transfemoral prostheses or both hips and knees in an exoskeleton, volitional control over multiple degrees-offreedom simultaneously would be ideal.
Traditionally, researchers have relied on mechanical sensors as well as neural sensors for classifier-based control and volitional control techniques. Mechanical sensors are used for extracting joint kinematics and kinetics with great success, and are advantageous for aiding in shared volitional control (Mendez et al., 2020). However, these sensors record kinematic or kinetic signals as the movement occurs and thus do not always precede impending motion. Non-invasive neural sensors, e.g. surface electromyography (EMG), are commonly used to extract muscle activation information that precedes joint movement (Farina et al., 2004). However, surface EMG is susceptible to muscle crosstalk and accesses muscle activity information from superficial muscles (Farina et al., 2002). Furthermore, muscle activity is only one contributing factor to muscle force production and resulting joint motion (Disselhorst-Klug et al., 2009). Still, EMG remains the primary peripheral sensing modality to probe information about muscle contraction that precedes limb motion. Recent advances in high-density recording and motor unit decomposition have enabled an increased signal dimensionality of EMG and understanding of the motor unit activity within superficial muscles (Merletti et al., 2008;Stango et al., 2015). However, alternative sensing modalities could provide additional information about muscle contraction, especially the contraction of multiple muscles that are within close proximity of one another and those that are not superficial. Advances are needed to enable accurate prediction of joint kinematics for multiple degrees-of-freedom during varying forms of ambulation (e.g., level ground, stair and ramp ambulation).
In comparison to EMG, which provides information about muscle activity, dynamic ultrasound imaging of skeletal muscle tissue, i.e. sonomyography, has the ability to capture real-time information about muscle deformation that is related to both fiber length and contractile velocity, two additional contributing factors to muscle force and resulting joint motion (Hallock et al., 2018). Additionally, sonomyography can access information from both superficial and deep muscles, resulting in increased spatial resolution and overall dimensionality of the signal (Rabe et al., 2021a). Multiple researchers have achieved success using sonomyography for control of multiple degrees-of-freedom in the upper limb (Shi et al., 2010;Hettiarachchi et al., 2015;Dhawan et al., 2019). In addition, previous research demonstrated improvement to classification of discrete ambulation tasks by sonomyography compared to surface electromyography (Rabe et al., 2021a). Further, sonomyography of the lower limb has been used to continuously estimate knee kinematics (Jahanandish et al., 2019b;Rabe et al., 2020a), as well as hip, knee and ankle moments during basic walking tasks and isometric contraction (Rabe et al., 2020b;Zhang et al., 2020;Zhang et al., 2021). However, the ability of sonomyography to continuously predict multiple joint kinematics during widely varying ambulation tasks encompassing level ground, ramp incline, ramp decline, stair ascent and stair descent and how the performance compares to surface EMG as well as the fusion of surface EMG with sonomyography is unknown. Identifying these important relationships would help inform the implementation and capabilities of volitional controllers for advanced lower-limb assistive technologies.
The objective of this research is to use both surface EMG and sonomyography, as well as the fusion of features from both sensing modalities, as inputs to a Gaussian process regression model for the continuous estimation of healthy subjects' hip, knee and ankle angles and angular velocities during five ambulation tasks. Gaussian process regression is a kernel-based Bayesian technique for nonlinear regression that has been introduced as an alternative to musculoskeletal model-based solutions to inverse dynamics. In this study, regression models were trained and tested in task-invariant frameworks such that strides from all ambulation tasks were included in both the training dataset and testing dataset. We hypothesized that sonomyography would improve regression prediction of hip and knee kinematics in comparison to surface EMG, but that sensor fusion of sonomyography with surface EMG would be required for accurate performance of ankle kinematic estimation. We further hypothesized that prediction of joint angle would outperform prediction of joint angular velocity due to the increased noise and variability of the angular velocity signal-noise that is associated with the use of numerical differentiation to calculate this signal from motion capture data. Results from this research can inform strategies for implementing multiple degree-of-freedom control over powered assistive devices using features extracted from novel sensing modalities such as sonomyography and traditional sensing modalities like EMG.

MATERIALS AND METHODS
A schematic overview of the methods beginning with data collection through model prediction of joint kinematics is given in Figure 1. This study was reviewed and approved by The University of Texas at Dallas institutional review board. The participants provided their written informed consent prior to participating in this study.

Subjects and Data Collection
Nine able-bodied subjects (five male and four female) were asked to complete two sets of ambulation tasks including: level walk, 10°i ncline walk, 10°decline walk, stair ascent and stair descent. Additional subject data can be found in Table 1. Prior to two sets of experiments, subjects were equipped with reflective markers for tracking kinematics, surface EMG for recording muscle activity and ultrasound sensors for collection of sonomyography data.
First, pre-gelled, self-adhesive electrodes (H124SG Covidien, Medtronic Inc., Dublin, Ireland) were placed over the belly of eight muscles with an inter-electrode distance of 2 cm. Electrodes were connected to a Bluetooth EMG unit (Shimmer3 EMG Unit, Shimmer, Dublin, Ireland) and EMG signals were recorded from the adductor magnus (AM), biceps femoris (BF), rectus femoris (RF), vastus lateralis (VL), tensor fascia latae (TFL), medial gastrocnemius (MG), tibialis anterior (TA) and soleus (SOL) muscles at a frequency of 1,200 Hz and streamed in real-time to the same lab computer as ultrasound.
In both sets of experiments, the 128-element linear array transducer of a portable, handheld ultrasound scanner (mSonics, Lonshine Technologies Inc., Beijing, China) was affixed to the thigh of each subject (on the same limb as surface EMG) via a custom-designed probe holder. The ultrasound transducer was placed on the anterior thigh during the first set of ambulation tasks, and on the posterior thigh during the second set of ambulation tasks. Priority was given to surface EMG to ensure optimal recording from the desired muscles when placing the ultrasound transducer on the anterior and posterior thigh. For anterior thigh ultrasound trials, the transducer was placed transversely at approximately half the distance between the anterior superior iliac spine and the proximal base of the patella (above the surface EMG electrodes) to collect grayscale images of the rectus femoris (RF), vastus medialis (VM) and vastus intermedius (VI) muscles. During the posterior thigh ultrasound trials, the transducer was placed transversely over the belly of the biceps femoris (BF) and semitendinosus muscles (below the surface EMG electrodes). Correct ultrasound placement, transversely over the belly of the desired muscles (Fukunaga et al., 1997), was confirmed visually via grayscale images prior to securing the FIGURE 1 | Schematic representation of methods. Beginning with data collection during five ambulation modes, sensing modality [sonomyography, surface electromyography (EMG) and sensor fusion] feature generation, regression model implementation, and ultimate hip, knee and ankle joint kinematic prediction. transducer in both anterior thigh and posterior thigh positions. Ultrasound transmit frequency was set to 7.5 MHz for all trials and the overall image gain was adjusted to optimize the brightness of the image and preserve image resolution at the deep muscle boundary. Forty-two reflective markers were placed over anatomical landmarks following guidelines by Delagi et al. on the bilateral feet, shanks, and thighs, as well as trunk and pelvis for collection of kinematic data via a ten-camera Vicon system (Vicon Motion Systems, Oxford, UK) (Delagi et al., 2011). Kinematic data was recorded to the same lab computer as ultrasound and surface EMG at a sampling rate of 100 Hz (Peng et al., 2016). Following sensor placement, all subjects were asked to complete two sets of five ambulation tasks: level walk, 10°incline walk, 10°decline walk, stair ascent and stair descent. All walking tasks were completed for 1 minute at a self-selected pace on a split-belt treadmill (Bertec, Columbus, OH, United States). Stair trials were completed on a four-step stair case, beginning with stair ascent followed by stair descent, with subjects walking in a reciprocal gait pattern at a self-selected pace. Stair trials were repeated five times and walk-to-stair and stair-to-walk transition strides were included in the respective stair analyses. Hip, knee and ankle kinematic data were calculated in Visual 3D via inverse kinematics, and a custom MATLAB program was created to enable time stamping of data for synchronization of sonomyography and surface EMG with kinematic data.

Sensor Feature Generation
Raw EMG signals were first processed to remove any baseline offset, then a fourth-order Butterworth bandpass filter was applied with a low cut-off at 20 Hz and a high cut-off at 450 Hz to remove motion artifact and high frequency noise, respectively (De Luca et al., 2010;Chowdhury et al., 2013). FIGURE 2 | Hip, knee and ankle joint angle and angular velocity during ten strides of each ambulation task (level walk, incline walk, decline walk, stair ascent and stair descent) from a single representative subject.
Frontiers in Robotics and AI | www.frontiersin.org March 2022 | Volume 9 | Article 716545 Then, a sliding analysis window method was used for extracting six features from processed surface EMG data and has been described in detail previously (Huang et al., 2011;Rabe et al., 2021a). The six features included four time-domain features (mean absolute value, number of slope sign-changes, number of zero-crossings and waveform length), as well as the first and second coefficient of a fourth-order autoregressive model. Sliding windows of length 200 ms were used with a 50 ms overlap, generating a 20 Hz signal of six EMG features from each of the eight muscles, resulting in 48 features per time-point (Smith et al., 2011).
The changes in grayscale ultrasound image intensity are correlated to changes in muscle density, as well as muscle architecture features such as muscle thickness and fascicle dynamics. Many researchers have evaluated these ultrasoundbased features of muscle morphology for correlations with muscle force production, muscle contraction, and joint motion as well as overall muscle strength and muscle fatigue (Kurokawa et al., 2001;Muraoka et al., 2005;Blazevich et al., 2006;Han et al., 2013;Panizzolo et al., 2013;Li et al., 2020). Based on previous research demonstrating these features are useful for estimation of knee kinematics from sonomyography, mean intensity and temporal intensity features were extracted from each ultrasound imaging frame, as described in detail previously (Jahanandish et al., 2019a;Jahanandish et al., 2019b;Rabe et al., 2020a;Rabe et al., 2021a;Rabe et al., 2021b). The image sequence from each trial was split by heel strikes to create an ultrasound image sequence for each stride. A spatial filter with a block size of 3 × 3 mm was used to extract mean intensity of each 3 × 3 mm block. Then, this 2-dimensional array of mean image intensity features was rearranged into 1-dimension by horizontally concatenating rows of features ranging from superficial to deep image features. The temporal features were created by taking the time derivative of each feature set between consecutive frames. Finally, the mean intensity and temporal intensity features were combined to create a single sonomyography feature set consisting of 520 features per frame.

Regression Model
The Gaussian process regression model is a Bayesian nonparametric kernel-based model that has been introduced as a real-time approximation of inverse dynamics solutions as an alternative for musculoskeletal model-based calculations (Rasmussen, 2004;Nguyen-Tuong et al., 2008;Wang et al., 2008). These models offered accurate function approximation given high-dimensional inputs at a relatively low computational demand to their inverse dynamics counterpart. In this study, we make the assumption that the distribution of joint kinematics within a subject, across multiple strides of a given task (e.g., Figure 2), can be characterized by a Gaussian distribution. Based on these preliminary analyses and previous work, we chose to implement a Gaussian process regression model with a rational quadratic kernel for solving the joint kinematics given the sonomyography, surface EMG, and sensor fusion feature sets. This model has been described in detail previously (Rasmussen, 2004;Zhang et al., 2018), but we assume the general function f (x) is distributed as a Gaussian function 2 | Mean (SD) root mean square error (RMSE) and range-normalized RMSE (nRMSE) of hip angle and angular velocity during level walking, incline walking, decline walking, stair ascent and stair descent. Joint kinematics were predicted by Gaussian process regression models trained and tested on features from five sensing modalities: 1) surface electromyography (EMG), 2) anterior sonomyography (Ant. SMG), 3) posterior sonomyography (Pos. SMG), 4) sensor fusion of Ant. SMG with EMG (Ant. Fusion), and 5) sensor fusion of Pos. SMG with EMG (Pos. Fusion). Overall average values are mean across all ambulation tasks. with a mean function m(x) and covariance, or kernel, The mean function reflects the expected function value at an input, x, and the rational quadratic kernel is described as θ is the maximum a posteriori estimates, σ is the signal standard deviation or amplitude parameter, and α is the non-negative parameter of the covariance or scale mixture, and l is the lengthscale (Schulz et al., 2018). A block coordinate descent method was used to solve for model parameters during training with a block size of 1,000 observations, gradient tolerance of 0.001, step tolerance of 0.001 and a maximum of 100,000 iterations to reach convergence. Each of the three feature sets: sonomyography, surface EMG, and sensor fusion were used as inputs to separate Gaussian process regression models for subject-dependent training and testing. Leave-one-stride-out cross-validation was used to prevent over-fitting of the model. In other words, a pooled dataset containing strides from all ambulation tasks was created and one stride of each ambulation task was removed for testing, then looped through the strides such that each stride of each ambulation task was the test stride once.

Statistical Analysis
Root mean square errors (RMSE) of each ambulation task and overall root mean square errors across ambulation tasks were calculated, and descriptive statistics (mean and standard deviations (SD)) were reported. A one-way ANOVA was completed to compare RMSE across sensing modalities within ambulation tasks at the hip, knee and ankle joint (α = 0.05). When a group significant difference was found, subsequent Tukey-Kramer multiple comparisons tests were completed to determine significant differences between sensing modalities. Additionally, RMSE were normalized to the range of the measured kinematics for comparison across joints and tasks. Additionally, an adjusted coefficient of determination (R 2 ) for nonlinear regression was calculated as a goodness-of-fit measure for the predicted kinematic trajectories.
where R 2 1 − Residual Sum of Squares Total Sum of Squares , N is the size of the training set, and p is the number of variables in the regression model.

RESULTS
Fusion of features from anterior sonomyography with surface EMG resulted in the lowest overall RMSE averaged across all ambulation tasks at the hip, knee and ankle ( Table 2, Table 3, Table 4). At the hip, there was a significant difference in overall average RMSE when comparing hip angles predicted from features of EMG with anterior sonomyography, as well as both anterior and posterior sensor fusion. Additionally, there was a significant difference in overall average RMSE when comparing hip angles from features of posterior sonomyography with anterior sensor fusion. There was a significant difference in overall average RMSE when comparing the hip angular velocities from features of EMG with anterior sensor fusion. At the knee, there was a significant difference in overall knee angle RMSE from features of EMG with all other sensing modalities, in addition to features from posterior sonomyography with anterior sensor fusion. Significant differences were found between overall RMSE of knee angular velocity prediction from features of EMG with anterior sonomyography and anterior sensor fusion. There was an additional significant difference in overall RMSE of knee angular velocity prediction between posterior sonomyography and anterior sensor fusion features. Lastly, at the ankle, there were significant differences between overall RMSE of ankle angle between EMG and anterior sonomyography, anterior sensor fusion and posterior sensor fusion. There was a significant difference between RMSE of ankle angle prediction based on posterior sonomyography features and anterior sensor fusion features. There were no significant differences between RMSE of all sensing modalities prediction of ankle angular velocity. Comparing across joints, the lowest RMSEs for both angle and angular velocity are observed at the hip and ankle, where the overall magnitude of the joint angle and angular velocity are lower (Figure 3). The normalized RMSE revealed that the hip angle and angular velocity were the lowest, there was an increase in normalized RMSE of knee angle and angular velocity, and the greatest normalized RMSE was observed when predicting ankle angle and angular velocity across most sensing modalities. Within the hip, knee and ankle joints individually, the normalized RMSE reveals that there was increased difficulty predicting joint angular velocity in comparison to joint angle for all sensing modalities.

Hip Estimation Performance
Anterior sensor fusion resulted in the best predictive performance of the Gaussian process regression model for both hip angle and angular velocity during all ambulation tasks, as evidenced by RMSE, normalized RMSE and the adjusted R 2 values (Table 2; Figure 4). There were no significant differences in RMSE of hip angle and angular velocity between anterior sonomyography alone and anterior sensor fusion, although anterior sensor fusion consistently outperformed anterior sonomyography. Surface EMG resulted in the worst predictive performance for all walking and stair tasks, followed by posterior sonomyography. The regression model predicted hip angle and angular velocity best during incline walking strides, followed by level walking strides for all sensing modalities. The lowest performance 4 | Mean (SD) root mean square error (RMSE) and range-normalized RMSE of ankle angle and angular velocity during level walking, incline walking, decline walking, stair ascent and stair descent. Joint kinematics were predicted by Gaussian process regression models trained and tested on features from five sensing modalities: 1) surface electromyography (EMG), 2) anterior sonomyography (Ant. SMG), 3) posterior sonomyography (Pos. SMG), 4) sensor fusion of Ant. SMG with EMG (Ant. Fusion), and 5) sensor fusion of Pos. SMG with EMG (Pos. Fusion). Overall average values are mean across all ambulation tasks. (greatest RMSE and normalized RMSE, lowest adjusted R 2 ) of regression model prediction of hip angle and angular velocity was observed during stair strides for all sensing modalities, with stair descent resulting in the greatest error.

Knee Estimation Performance
Similar to the hip, anterior sensor fusion resulted in the best predictive performance of the Gaussian process regression model for both knee angle and angular velocity during all ambulation tasks (Table 3; Figure 5). There were no significant differences in RMSE of knee angle and angular velocity between anterior sonomyography alone and anterior sensor fusion. Surface EMG resulted in the worst predictive performance for all ambulation tasks, followed by posterior sonomyography. The regression model predicted knee angle during level and incline walking with the lowest error for all sensing modalities compared to decline walking and stair strides. The lowest error for knee angular velocity prediction was observed during incline walking for regression models trained from all sensing modalities, followed by decline walking, level walking, and stair strides, as evidenced by the normalized RMSE and adjusted R 2 values.

Ankle Estimation Performance
In comparison to the hip and knee, the Gaussian process regression model performed worst at the ankle for all sensing modalities. However, comparing sensing modalities, the results remain consistent that anterior sensor fusion resulted in the lowest error of ankle angle and angular velocity, followed by anterior sonomyography alone (Table 4; Figure 6). Surface EMG resulted in the greatest error for ankle angle prediction during all ambulation tasks, with significant improvement between all other sensing modalities versus EMG during incline and decline walking. For both anterior sensor fusion and anterior sonomyography-based estimates of ankle angle, the highest performance (lowest normalized RMSE and greatest adjusted R 2 ) was observed during decline walking, followed by stair descent, incline walking, level walking, and lastly, stair ascent. There were no significant differences between RMSE of ankle angular velocity prediction by any sensing modality during any of the five ambulation tasks.

Computational Expense
All models were trained and tested offline on a single CPU (Intel(R) Core i7-7700 at 3.60 GHz). Average time for training and testing a single fold of the leave-one-stride-out framework for the rational quadratic Gaussian process regression model are given in Table 5.

DISCUSSION
The objective of this research was to evaluate features of surface EMG and sonomyography, as well as fusion of these features, for the  Predicted kinematics from Gaussian process regression models trained and tested on features from electromyography (EMG), sonomyography (SMG) and sensor fusion (Fusion) are displayed with respective standard deviations. Adjusted R 2 given as a goodness-of-fit metric for each sensing modality compared to the measured kinematics.
Frontiers in Robotics and AI | www.frontiersin.org March 2022 | Volume 9 | Article 716545 9 significant differences between RMSE of angle and angular velocity prediction based on the respective anterior or posterior sonomyography and sensor fusion features, and anterior sonomyography and sensor fusion resulted in relatively lower error compared to posterior sonomyography and sensor fusion. These results are surprising at the ankle, where the addition of FIGURE 5 | Knee angle and angular velocity as a function of the gait cycle. Measured kinematics are displayed in gray with standard deviations in shaded regions. Predicted kinematics from Gaussian process regression models trained and tested on features from electromyography (EMG), sonomyography (SMG) and sensor fusion (Fusion) are displayed with respective standard deviations. Adjusted R 2 given as a goodness-of-fit metric for each sensing modality compared to the measured kinematics.
Frontiers in Robotics and AI | www.frontiersin.org March 2022 | Volume 9 | Article 716545 information of muscle activity from muscles that span that joint did not significantly improve the regression performance. The increased performance, or reduced error, of anterior sonomyography in comparison to posterior sonomyography, could possibly be attributed to capturing features from three muscles on the anterior thigh (rectus femoris, vastus intermedius, and vastus Predicted kinematics from Gaussian process regression models trained and tested on features from electromyography (EMG), sonomyography (SMG) and sensor fusion (Fusion) are displayed with respective standard deviations. Adjusted R 2 given as a goodness-of-fit metric for each sensing modality compared to the measured kinematics.
Frontiers in Robotics and AI | www.frontiersin.org March 2022 | Volume 9 | Article 716545 11 medialis), as opposed to only two muscles on the posterior thigh (biceps femoris and semitendinosus). This gives further justification to the increased resolution, as well as ability to access deep muscle tissue, as probable explanations for the improved regression performance of sonomyography in comparison to surface EMG. These results are in agreement with our previous work demonstrating increased performance of sonomyography in comparison to surface EMG for ambulation mode classification, as well as high performance of sonomyography-based knee angular velocity prediction and hip, knee and ankle joint moment prediction (Rabe et al., 2020a;Rabe et al., 2021a;Rabe et al., 2021b).
For all sensing modalities, model error increased when comparing prediction of joint angle and angular velocity during walking tasks to stair ambulation tasks. This can likely be explained by the decrease in the number of stair strides in the training dataset compared to the level, incline and decline walking strides, or perhaps the increased variability of the stair ambulation kinematic trajectories. Another possible explanation is the need for additional high resolution (i.e., sonomyography) information from additional thigh muscles during the stair ambulation tasks to account for the increased variability during these tasks. Previous research showed that muscle activity of the individual quadriceps and hamstrings muscles is reduced during stair ascent and stair descent compared to level walking in healthy individuals; however, there is an opposite relationship in muscle coactivity where there is an increase in hamstrings/quadriceps coactivity during stair ambulation as opposed to level walking (Bae et al., 2009). These results point to the potential benefit of sonomyography from the anterior thigh muscles and posterior thigh muscles simultaneously to improve the prediction of hip, knee and ankle joint kinematics.
Accurate prediction of joint output kinematics from a reliable sensing modality is vital to the success of powered assistive devices. While many researchers evaluate control over multiple degrees-offreedom in the upper-limb, research evaluating simultaneous control strategies over the hip, knee and ankle degrees-of-freedom of the lower limb is less widely reported. The first powered lower-limb assistive devices focused only on powering a single joint, such as the knee or ankle alone (Au et al., 2007;Fite et al., 2007). However, the technological advancement of powered lower-limb assistive devices has enabled the plausibility for multiple powered joints within a single device. Minimizing the complexity of the control strategy for these devices will expediate their entrance into a rehabilitation clinic setting and the hands of users. To date, researchers have relied on pre-programmed trajectories for various ambulatory tasks, such as stair ascent or ramp walking (Lawson et al., 2013;Simon et al., 2014;Lenzi et al., 2018). These strategies have proven to be successful, but require users to adapt their mechanics to accommodate the device, rather than allowing the device mechanics to adapt to the user. This can result in further injury or long-term wear to the unaffected joints (Grabowski and D'Andrea, 2013). The results from the present study are promising for continuous prediction of hip, knee and ankle kinematics that may be able to adapt to the user.
The Gaussian process regression model was chosen based on its recent success for use as a real-time approximation of inverse dynamics solutions to musculoskeletal model-based control. We believe the Gaussian process regression model is an appropriate model for the non-linear relationship between the sonomyography and surface EMG features from multiple muscles with the output of the hip, knee and ankle joints. As an alternative to this approach, Zhang et al. evaluated various ultrasound features of tibialis anterior muscle contraction and surface EMG in a musculoskeletal modelbased algorithm for continuous estimation of isometric ankle dorsiflexion moment and compared this model to linear regression as well as a neural network approach (Zhang et al., 2021). While the model-based approach has many advantages, including establishing a functional, "white-box" relationship between the muscle features with the desired output, the results demonstrated there was no significant difference between the musculoskeletal model-based approach and two regression approaches for ankle moment estimation. However, these results indicated that surface EMG features, muscle fascicle length and muscle pennation angle resulted in higher accuracy of isometric ankle moment prediction in comparison to overall mean image echogenicity (of longitudinal images along the tibialis anterior muscle). In a previous study, Zhang et al. combined surface EMG muscle activation features with tibialis anterior pennation angle from ultrasound as inputs to a Hill-type neuromusculoskeletal model for prediction of isometric ankle dorsiflexion moment and found that the fusion of the two sensor data improved prediction accuracy in comparison to sole surface EMG or ultrasound alone (Zhang et al., 2020). Additionally, researchers evaluating contributions of sonomyography and electromyography in both the lower-and upper-limb, have found that both signals contribute to accurate prediction of joint and muscle mechanics (Shi et al., 2007;Ruiz-Muñoz and Cuesta-Vargas, 2014;Boyd and Liu, 2020). The present results further support fusion of sonomyography and surface EMG as an input to regression-based approaches for continuous prediction of joint-level output.
The predictive performance of the Gaussian process regression model was greatest when trained on sensor fusion features from anterior sonomyography and surface EMG. For almost all ambulation modes of all joints, anterior sensor fusion resulted in significantly improved hip, knee and ankle angel, as well as hip and knee angular velocity prediction in comparison to surface EMG. There were fewer instances of significant improvement of the prediction of hip, knee and ankle kinematics when comparing anterior sonomyography alone in comparison to surface EMG. However, there was no significant differences between the anterior sonomyography alone and the anterior sensor fusionbased predictions of hip, knee and ankle kinematics. Importantly, the lack of significant differences between sonomyography alone and sensor fusion of sonomyography with surface EMG indicates that the addition of separate features from surface EMG may not greatly improve the performance of the regression model at the hip, knee or ankle joints. This may be evident when considering the more traditional measurement of surface EMG we used in this study, and not true when considering other newer strategies for recording and decomposition of the surface EMG signal (Stango et al., 2015) Additional work is required to further evaluate the comparison of sensor fusion with sonomyography alone to determine if the burden of collecting data from a separate sensing modality is required for accurate prediction of lower-limb kinematics. All sonomyography and surface EMG data were recorded in real-time and synchronized to the joint kinematics during each ambulation task. However, feature generation and Gaussian process regression model training and testing was completed offline. The average time to train the regression model on the anterior sonomyography fusion dataset was 23.45 s, and the average time to test the models on strides from the hip, knee and ankle was 5.4 ms. The prediction time is within the optimal reported controller delay window of 100-125 ms for myoelectric prostheses (Farrell and Weir, 2007), as well as faster than the update rates of common lower-limb device controllers (Fluit et al., 2020). However, future work is required to determine if the hip, knee and ankle angle and angular velocity prediction times and accuracies are maintained during online implementation. Given previous success in the upper-limb for online sonomyographic control of dexterous hands , along with an established correlation between online and offline accuracy of a powered prostheses (Hargrove et al., 2015;Lu et al., 2017), we expect the accuracy of the offline joint kinematic prediction to be comparable, and at a minimum correlate with future online results.
Limitations of the present work exist beyond the offline implementation of the regression models. The subject population included only healthy individuals, without any mobility disorders or limb loss. Additional work is required to determine how the sonomyography features will translate to individuals with gait impairments as well as limb loss. Additionally, the ability of these features and the Gaussian process regression model to predict kinematics of unknown ambulation tasks has not been addressed. Future work evaluating the ability of these features to predict hip, knee and ankle kinematics during unknown tasks will improve understanding of their feasibility for implementation in various activities of daily living, such as stepping over an obstacle. Furthermore, alternative features, such as muscle cross-sectional area and aspect ratio (Guo et al., 2010), as well as alternative models, such as more complex neural networks or musculoskeletal models, may improve performance of the surface EMG and sonomyography features for predicting kinematics of unknown (or "untrained") tasks. To aid in the evaluation of online prediction of joint kinematics from surface EMG and sonomyographic features of muscle contraction of various subject populations including individuals with mobility disorders, steps should be taken to integrate the sensors into the assistive devices and their control schemes. Researchers have begun to develop electrodes for simultaneous surface EMG and grayscale ultrasound recording within a single system at the same muscle position. (Botter et al., 2019;Yeon et al., 2021). Furthermore, a reduction in the number of ultrasound elements, and perhaps a reduction to single amplitude ultrasound elements, could be possible. Sikdar et al. has demonstrated success in predicting dexterous finger movements with a single element ultrasound transducer in the upper limb (Sikdar et al., 2014), and sparsity analyses of grayscale ultrasound indicated a reduction to equally spaced single scanlines did not introduce significant error (Akhlaghi et al., 2020). Additionally, the transition to single element amplitude ultrasound sensors could reduce the ambiguity of the transfer function between amplitude signals and grayscale images provided by brightness modulated ultrasound. Theoretically, this will support improvement of the evaluation of task-and user-independent prediction from sonomyography. This work evaluated a task-invariant approach for continuous prediction of hip, knee and ankle angle and angular velocity during five common, yet widely-varying, ambulation tasks using Gaussian process regression models trained and tested on features from surface EMG, sonomyography, and sensor fusion. Anterior sonomyography sensor fusion with surface EMG significantly improved almost all joint kinematics in comparison to surface EMG alone, while anterior sonomyography alone significantly improved only some of the joint kinematics (mostly at the hip and knee) in comparison to surface EMG alone. However, there were no significant differences between sonomyography alone and sensor fusion. Additionally, the results revealed that anterior sonomyography and anterior sonomyography sensor fusion gives more accurate predictive performance in comparison to posterior sonomyography and posterior sonomyography sensor fusion, respectively. Sensor fusion and the Gaussian process regression models predicted level ground and ramp ambulation kinematics with greater accuracy in comparison to stair kinematics. Larger training datasets including additional stair strides, or perhaps sonomyography data from muscles on the anterior and posterior thigh simultaneously could possibly improve stair kinematic prediction. We believe this able-bodied study is a fundamental contribution to the evaluation and potential integration of sonomyography sensors into powered assistive devices. Improved control over multiple degrees-of-freedom of powered assistive devices is critical for translation of these devices into the daily lives of users and improving the quality of life for individuals with gait impairments or limb loss.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

ETHICS STATEMENT
The studies involving human participants were reviewed and approved by the Institutional Review Board at The University of