—Objective: This study describes the development and ofﬂine validation of a heuristic algorithm for accurate prediction of ground terrain in a lower limb prosthesis. This method is based on inference of the ground terrain geometry using estimation of prosthetic limb kinematics during gait with a single integrated inertial measurement unit. Methods: We asked ﬁve subjects with below-knee amputations to traverse level ground, stairs, and ramps using a high-range-of-motion powered prosthesis while internal sensor data were remotely logged. We used these data to develop two terrain prediction algorithms. The ﬁrst employed a state-of-the-art machine learning approach, while the second was a directly tuned heuristic using thresholds on estimated prosthetic ankle joint translations and ground slope. We compared the performance of these algorithms using resubstitution error for the machine learning algorithm and overall error for the heuristic algorithm. Results: Our optimal machine learning algorithm attained a resubstitution error of 3:4%using 45 features, while our heuristic method attained an overall prediction error of 2:8% using only 5 features derived from estimation of ground slope and horizontal and vertical ankle joint displacement. Compared with pattern recognition, the heuristic performed better on each individual subject, and across both level and non-level strides. Conclusion and signiﬁcance: These results demonstrate a method for heuristic prediction of ground terrain in a powered prosthesis. The method is more accurate, more interpretable, and less computationally expensive than state-of-the-art machine learning methods, and relies only on integrated prosthesis sensors. Finally, the method provides intuitively tunable thresholds to improve performance for speciﬁc walking conditions. Index Terms —intent recognition, inertial measurement, ma- chine learning, wearable robotics. I. I NTRODUCTION POwered lower limb prostheses and exoskeletons have demonstrated the ability to recapitulate certain steady- state walking tasks for their users. The control methods used in these studies most often employ a state machine, which encodes the desired walking task into a series of consistent, predictable states, each actuating a corresponding control law and set of transition conditions to other states [1–4]. However, while dynamics within any given activity are relatively cyclic and predictable, transitions between activities such as sitting, standing, level-ground walking, or stair traversal are more ran- dom. When such transitions call for biomechanically different behaviors, it becomes necessary to develop a prosthetic control paradigm that can automatically detect them, thereby allowing All authors are in the Center for Extreme Bionics at the Massachusetts Institute of Technology. This work was supported by the United States Army Medical Research Acquisition Activity under contract number W81XWH-14-C-0111. R. Stol- yarov was supported by a National Defense Science and Engineering Graduate Fellowship and the Harvard-MIT Health Sciences and Technology program.the user to seamlessly progress form one desired walking task to another. A signiﬁcant determinant of gait dynamics is the ground terrain. Numerous studies have shown that walking on level ground, inclines, or stairs can have signiﬁcant effects on the kinetics and kinematics of lower limb joints [5–8]. To address the issue of these terrain-dependent lower limb biomechanics, the aim of this work is to develop a method for accurate, real- time prediction of terrain using only the sensors on board a powered prosthesis. Many approaches to solve this problem have relied on some form of remote terrain sensing including millimeter-wave radar [9] or optical distance sensing [10]. While these methods allow direct detection of certain terrains, they present inconvenience for users because they rely on a line of sight. Alternatively, most approaches that have not employed remote sensing rely on some form of pattern recognition to infer the upcoming terrain [10–20]. These methods typically involve some sort of supervised training protocol, whereby a classiﬁer is taught to correctly select from a set of possible walking tasks using signal features extracted from windows deﬁned by detected gait events. However, while these studies have demonstrated high-accuracy terrain prediction, they often produce classiﬁers that are both complex and difﬁcult to interpret. As a result, misclassiﬁcations are difﬁcult to account for and classiﬁer performance cannot be directly improved. Furthermore, while the hope is that the resulting classiﬁers will perform well in diverse conditions, it is not possible to directly tune their performance to novel conditions. In this work, we address the issues inherent to machine- learning based terrain prediction by designing a fully heuristic algorithm to accurately predict walking tasks. In doing so, we leverage our understanding of the gait dynamics necessarily imposed on the lower limbs by the underlying terrain geometry and make use of methods developed in our previous work [21] to measure these dynamics. In particular, our approach relies on the idea that terrain geometry enforces both lower limb joint trajectories and foot inclination, and that accurate estimation of these signals can consequently enable accurate and heuristic terrain prediction. Our resulting walking task predictor relies only on a set of intuitive if-then heuristics and associated thresholds on these kinematic signals. We demonstrate that the heuristic is not only as accurate as state-of-the-art machine learning methods, but that it employs only 5 features as compared to the 40-50 characteristic of the machine learning approaches. These characteristics minimize the risk of over- ﬁtting, enhance classiﬁer interpretability, and introduce the ability to tune classiﬁer performance by adjusting one of a 2 TABLE I SUBJECT ATTRIBUTES ID Sex Age Ht(m) Wt(kg) Affected side 1 F 35 1.67 90.6 L 2 M 59 1.72 89.3 L 3 M 36 1.83 80.7 L 4 M 39 1.90 80.5 R 5 M 35 1.75 65.2 L small number of intuitive and physically relevant thresholds. II. M ETHODS A. Overview We asked ﬁve subjects with unilateral transtibial amputa- tions to don a powered lower limb prosthesis and traverse various terrains including level ground, stairs, and ramps. While each subject was walking, we remotely logged data from internal prosthesis sensors, including ankle angle, ankle torque, and raw inertial measurements from a six degree- of-freedom (three orthogonal accelerometers and gyroscopes) IMU. Additionally, we ﬁlmed the subjects to facilitate manual labeling of terrains to provide a ground-truth terrain identity. Next, we used data from all subjects to develop a heuristic based on an direct estimation of terrain geometry during gait to accurately predict the terrain of every stride. This heuristic relied on estimation of ankle horizontal and vertical position using a method similar to that employed in our previous work [21], as well as a novel ground slope estimation method. In order to validate our heuristic, we also trained a linear- discriminant-analysis based terrain classiﬁer on the collected data set using a state-of-the-art method also described in our previous work [21]. Finally, we compared the overall accu- racy of the heuristic algorithm to the mean cross-validation accuracy of the machine learning algorithm. B. Data collection We asked a total of ﬁve subjects (see Table I) to complete a series of walking trials using the TF8 ankle-foot prosthesis, ﬁrst described in [22] and outlined in the following sections. Walking trials included at least 15 circuits up and down a 9- degree ramp (length 8.5 feet), at least 4 circuits up and down a 12 step staircase (0.171m rise, 0.279m run), and at least one minute of level-ground walking at self-selected, slow, and fast speeds. The number of circuits for each trial was chosen to ensure at least 20 samples were attained for each terrain. 1) Prosthesis: a) Mechanical subsystem: This study employed a novel ankle-foot prosthesis comprised of a reaction force series elastic actuator. This device was built to achieve biological kinetics and kinematics that enable operation over a range of terrain conditions. The device is a torque controlled powered prosthesis designed around a series elastic actuator that can provide peak torques up to 180 Nm across a 115 degree total operational range of motion. The system, shown in Figure 1, consists of a large gap radius motor (manufactured by T- Motor) modiﬁed to integrate a ball screw into the rotor. The Joint EncoderControl Unit Force SensorIMU Fig. 1. The TF8 mechatronic system architecture is a reaction force series elastic actuator with an on-board embedded control system. Fig. 2. Method used to estimate foot ﬂat by thresholding the difference between shank segment velocity and ankle angle velocity. ball screw transmission applies a linear force to an output moment arm that generates a torque about the ankle joint. b) Electrical subsystem: An axial load cell directly mea- sures the force in the ball screw. This force signal is evaluated along with the joint encoder measurements to determine the effective joint torque with an accuracy of 0:5Nm. The joint encoder is a 14-bit absolute encoder, AS5048 (manu- factured by Austria Microsystems). Inertial measurements are performed by the motion tracking MPU-9250 (InvenSense) included in the control unit printed circuit board assembly. The control unit consists of a customized embedded system platform based on the FlexSEA system designed by Dephy, Inc and ﬁrst described in [23]. c) Control subsystem: The system includes a motor driver and a separate mid-level control system that is based on the STM32F427 32-bit Cortex M-4 microcontroller operating at 180 Mhz. A closed-loop torque controller asserts a torque on the joint by converting the torque command to a desired motor current that the FlexSEA motor driver’s internal current con- troller commands at the motor. To enable terrain agnostic data collection we implemented a single state control scheme across 3 all terrains. Desired torque was set as d=k(0) +b_, wherek,0, andbwere respectively manually tuned spring stiffness, spring position set point, and damping coefﬁcient. This system was tuned for each subject to be within the range of subject reported comfort C. Ofﬂine processing 1) Initial processing: We divided data from all subjects into individual strides using thresholds on ankle torque and torque derivative and timers for swing and stance phases. Individual strides were deﬁned by the period from one foot strike to the next foot strike, so as to include enough information for prediction and to ensure predictions are made early enough to enable timely actuation of the prosthesis before stance. Within this window, all strides initially consisted of eight signals including ankle angle , ankle torque , three accelerations S~ a, and three rotational velocitiesS~ !. Finally, we manually labeled all strides using trial videos as ground truth evaluation of prediction accuracy. 2) Heuristic prediction algorithm: We were interested in developing a simple heuristic to accurately predict the terrain geometry on a stride-by-stride basis. In pursuing this aim, we reasoned that all terrain archetypes can be detected using thresholds on ankle vertical position, ankle horizontal position, and ground slope. For example, strides beginning with a large vertical ankle displacement and small horizontal displacement are likely stair ascent strides, while those with large initial hor- izontal displacements are likely level ground or ramp strides. Similarly, large ground slope magnitudes indicate strides are made on inclines. Thus, for a predictor whose goal it was to distinguish between strides made on level ground, inclines, and stairs, we reasoned it was necessary to accurately estimate ankle translations and ground slope. a) Estimation of ankle translations: We used a motion integration algorithm onS~ aandS~ !similar to that described in our previous work [21]. For reference we outline the steps of the algorithm below: 1) Map sensor-frame accelerations at the IMU to the ankle joint by assuming that the residual limb and prosthetic shank comprise a rigid body. 2) Integrate sagittal rotational velocity to update estimated shank sagittal orientation. 3) Use estimated shank orientation to project sensor-frame accelerations at the ankle joint to the global frame. 4) Offset global-frame vertical acceleration by gravity. 5) Integrate global-frame accelerations to obtain global- frame ankle velocities and positions. Due to bias and noise inS~ aandS~ !, integrated signals were subject to accumulating error. To bound this error, we reset these signals at least once during every stance period. While this update was also performed in the algorithm described in [21], in this work we used a different set of conditions for triggering the update, described below. For orientation resets we employed a threshold on the acceleration norm jS~ aj, expressed as: jS~ ajg ^a (1)whereg= 9:8m=s2and= 0:1m=s2. At each detected orientation update time, sagittal shank rotation angle R was adjusted by a manually tuned weighting factor cas: R2;1:=sin( ) =c(SaY=g) + (1c)R2;1 (2) R1;1:=cos( ) =q 1R2 2;1 (3) R2;2:=sin( ) =R2;1 (4) R1;2:=cos( ) =R1;1 (5) (6) where is the shank pitch angle and c= 0:02. This value was tuned to minimize the mean absolute orientation error over the course of one trial as compared to data extracted from a 12- camera Vicon motion capture system. For position resets we employed a threshold on the norm difference between sagittal shank rotational velocityS!Xand ankle rotational velocity _. As illustrated in Figure 2, the assumption underlying this calculation was that these values would be closest when the foot is ﬂat on the ground. The calculation and threshold we employed can be expressed as: =jS!X_j^! (7) where _was ﬁltered using a 50Hz, 2nd order low pass Butterworth ﬁlter and ^!= 1:2rad=s . This threshold was chosen empirically to allow for at least one velocity/position reset per stance period, and was likely high due to deﬂection of the prosthetic foot during roll-over. At the velocity reset timet(i) R, integrated signals were reset by modeling the shank as a vertical lever rotating in the sagittal plane about a ﬁxed hinge at the ankle joint: ~ pA(t(i) R) := [0; 0; 0] (8) ~ vA(t(i) R) := [0; 0; 0] (9) ~ pK(t(i) R) :=L[0;R2;1(t(i) R);R1;1(t(i) R)] (10) ~ vK(t(i) R) :=L!X[0;R1;1(t(i) R);R2;1(t(i) R)] (11) (12) where the components represent, in order, the frontal, anterior- posterior, and longitudinal axes. Ankle translational estimates were validated by character- izing their performance in measuring the geometry of various terrains. In particular, we investigated how closely end-of- stride position matched expected values for constrained terrain geometries (for example, two stair heights altitude change for stair ascent or descent and no net change in altitude for level- ground walking). b) Estimation of ground slope: We employed the foot- ﬂat detection algorithm described above for estimating the ground slope by calculating a running average of the difference between ankle angle and estimated shank pitch while the foot was detected to be ﬂat on the ground. This can be shown through simple geometry to yield the ground slope, a relationship illustrated in Figure 3. We validated our foot ﬂat detection algorithm by quantify- ing its accuracy in measuring three different ground slopes, including level ground and an ascending and descending 9- degree ramp. 4 Fig. 3. Method used to estimate ground slope, combining knowledge of shank pitch and ankle angle to achieve the most accurate possible estimate. Estimation involved averaging the value  during a time when the foot was detected to be ﬂat on the terrain. Fig. 4. Heuristic back-estimation of stride terrain based on horizontal and vertical ankle joint positions and estimated ground slope. c) Development of the heuristic: Once we were con- vinced that we could estimate these heuristic signals accurately enough, we developed a decision-tree heuristic with manually tuned coefﬁcients to distinguish between terrains based on a direct, dynamic measurement of terrain geometry across all subjects. The heuristic architecture is shown in Figure 4. We report the accuracy of this heuristic as it pertains to classifying strides made by all subjects. 3) Machine learning algorithm: We were interested in comparing the performance of our heuristic to a state-of-the- art approach based on machine learning. To this end, we employed data from all subjects to train a terrain classiﬁer using a method identical to that we employed our previous work [21]. In particular, we used 10-fold cross-validation of a linear discriminant analysis (LDA) classiﬁer to perform forward feature selection until an optimal mean fold accuracy was reached. Included in feature selection were six candidate features (including mean, range, maximum, minimum, ﬁrst and ﬁnal values) extracted from all IMU-derived signals in a window beginning at foot-off and ending at various times post foot-off. Candidate signals included raw sagittal-plane accelerations and rotational velocities, estimated shank pitch, and vertical and horizontal ankle accelerations, velocities, and positions. We also included signals pertaining to ankle joint dynamics, including torque, ankle angle, and the derivatives of each, extracted over the entirety of each stance period. In all cases we assumed a uniform prior probability vector.III. R ESULTS A. Initial processing The ﬁve subjects whose data was used for algorithm devel- opment took a total of 1594 strides, including 916 on level ground, 169 on a 9 degree ascending ramp, 230 on a 9 degree descending ramp, 140 stair ascent strides, and 139 stair descent strides. Walking speeds during the slow, medium-speed, and fast level-ground walking trials ranged from 0.5m/s to 1.2m/s as measured using a timer and tape markers on the ground. B. Translational motion tracking Our intended heuristic approach relied on accurate trans- lational motion tracking. The performance of the translation estimation algorithm is characterized by Figure 5, which shows data on all strides. Figure 5a shows sagittal plane ankle transla- tions of all forward strides taken by the ﬁve subjects along with all underlying terrain geometries of interest. Figures 5b and 5c show the full vertical and horizontal measurement distributions for each terrain, along with expected measurements for cases where the terrain geometry constrained the ankle position. These consisted of a vertical position constraint during level- ground walking and vertical and horizontal position constraints during stair ascent and descent. In particular, we assumed the following relationships between terrain geometry and expected ankle displacement: Level-ground walking should cause a net vertical ankle displacement of 0.0m Double stair ascent and double stair descent should cause a horizontal ankle displacement of two stair runs (2 x 0.279m = 0.560m) Double stair ascent and double stair descent should cause a net vertical ankle displacement of two stair rises (2 x 0.171m = 0.342m for stair ascent and -0.342m for stair descent) C. Slope estimation Figure 6 shows the performance of ground slope estimation on level ground and ramp strides, including mean errors and standard deviations separately for each terrain. In all cases, estimation errors are below 1.7 degrees. D. Prediction accuracy Our heuristic attained a prediction error of 2.8% on all strides, while the best pattern recognition algorithm attained a mean cross validation error of 3:52%and resubstitution error of 3.4%. For both heuristic and pattern-recognition-based prediction, strides in which predictions were made less than 100ms before the next stance period were predicted as level ground by default, since in real-time conditions this would not provide enough time to servo the prosthetic ankle to a desired position. These results are shown in more detail in Figure 7a, which shows the mean and standard deviation of cross validation errors for classiﬁers trained on features extracted from signals cut off at varying times post foot-off. Figure 7b shows the forward feature selection history involved 5 Fig. 5. (a) Estimated sagittal plane ankle translations for strides taken by all subjects, plotted along with underlying terrain geometries. Scattered points indicate end-of-stride estimated ankle positions. (b) Measured vertical displacement distributions for each terrain. For terrains which vertically constrain foot position we also include dashed black lines indicating the expected ankle displacement. In particular, we assumed that ankle displacement is vertically constrained during stair ascent, level-ground walking, and stair descent. Stride proportions are relative to the stride population for the corresponding terrain. (c) Measured horizontal displacement distributions for each terrain. For terrains which horizontally constrain foot position we also include dashed black lines indicating the expected ankle displacement. In particular, we assumed that ankle displacement is horizontally constrained during stair ascent and stair descent. Stride proportions are relative to the stride population for the corresponding terrain. Fig. 6. Normalized histograms of ground slope estimates for all ramp descent, forward level ground, and ramp ascent. Distributions for each terrain are normalized relative to the number of strides with the corresponding terrain label. 6 (a) (b) Fig. 7. (a) Mean fold prediction errors from performing forward feature selection using mean 10-fold cross validation accuracy of an LDA classiﬁer. Each bar represents the optimal error attained when classifying on features extracted from windows beginning at a foot-off event and ending at different cutoff times post foot-off. Values at the top of the error bars indicate the number of features for the best model. The dashed line indicates the composite classiﬁcation accuracy of our heuristic. (b) Feature selection history from using 10-fold cross validation of an LDA classiﬁer on features extracted between foot-off and 250ms post foot-off. This time point achieved the lowest optimal error of all tested cutoff times. . in training the best-performing feature set among these, and Table II shows the confusion matrices for both this classiﬁer and the heuristic algorithm. The confusion matrix for the machine learning algorithm was formed by consolidating the predictions from all folds into one matrix. IV. D ISCUSSION We have demonstrated a method for accurate and heuristic prediction of ground terrain in a lower limb prosthesis using only integrated sensors. As compared to a machine learning algorithm, this heuristic is not only as accurate but is also both tunable and interpretable, allowing direct adjustment of values depending on the speciﬁc terrain geometry or user require-ments. Finally, the method is computationally inexpensive and thus all calculations can be performed completely on-device. A. Discriminating terrains by ankle translations The results in Figure 5 show great promise for enabling heuristic discrimination between terrains using ankle trans- lational motion. Attributes supporting high discriminability between terrains include disjoint error bars and signiﬁcantly different distributions of end-of-stride vertical displacement for each terrain of interest, disjoint end-of-stride horizontal error bars between stair and non-stair terrains, and correspon- dence of both horizontal and vertical means with the expected displacements for their corresponding terrain. However, de- 7 TABLE II CONFUSION MATRICES FROM PATTERN -RECOGNITION -BASED AND HEURISTIC -BASED TERRAIN PREDICTION . VALUES ARE ROUNDED AND THUS EACH ROW MAY NOT TOTAL TO EXACTLY 100%. A BBREVIATIONS REPRESENT ,IN ORDER ,LEVEL -GROUND WALKING ,RAMP ASCENT ,RAMP DESCENT ,STAIR ASCENT ,AND STAIR DESCENT . Pattern recognition F RA RD SA SD F 95.1 3.5 0.85 0.0 0.56 RA 7.1 92.9 0.0 0.0 0.0 RD 4.7 0.43 94.8 0.0 0.0 SA 0.0 0.0 0.0 100.0 0.0 SD 0.0072 0.0 0.0072 0.0 98.6 Heuristic F RA RD SA SD F 98.4 0.66 0.77 0.22 0.0 RA 5.9 93.5 0.0 0.59 0.0 RD 7.4 0.43 92.2 0.0 0.0 SA 0.0 0.0 0.0 100.0 0.0 SD 0.70 0.0 1.3 0.0 99.3 spite these attributes, we found that we could not achieve perfect discrimination for any given terrain using a horizontal or vertical displacement alone. However, we reasoned that we could leverage the combination of high vertical separation between stair ascent and descent and high horizontal and vertical separation between stair and non-stair terrains to design an accurate, large-margin ankle position-based heuristic for identifying stair ascent and stair descent. B. Discriminating terrains by slope estimation Given a large overlap in both horizontal and vertical dis- tributions for non-stair terrains, we necessitated a separate parameter for this case, namely the estimated ground slope. Based on the standard deviations from Figures 6, our es- timator’s resolution was 1:7degrees regardless of ground slope or walking speed. This value is likely sufﬁcient for most control applications: for example, [6] found very small differences in both kinetics and kinematics between ramps whose inclinations differed by 3 degrees. Additionally, as shown by the histograms in Figure 6, the slope estimator attained perfect discrimination between level ground, a 9- degree ascending ramp, and a 9-degree descending ramp. C. Translational motion tracking error By analyzing strides made on terrains that constrain ankle trajectory, we were able to characterize the error of our translational motion tracking approach, which determines its ability to differentiate between stairs of different inclinations. These include stair terrains for horizontal domain and both stair terrains and level ground for the vertical domain. Ac- cording to 5b and 5c, both vertical and horizontal position error was smaller for stair ascent than for stair descent. For level-ground walking, vertical position error was larger than that for stair ascent but smaller than that for stair descent. We hypothesize that the magnitude of this error is related to both the magnitude and timing of shank rotation undergone when traversing each of these terrains. Greater and earlier shank rotation leads to greater drift in the shank orientation estimate,which in turn compounds error in both vertical and horizontal position tracking. Relevantly, in an analysis of stair and level ground strides, we found that stair ascent shank rotation is minimal, while level ground and stair descent shank rotation is comparable but tends to happen earlier in stair descent. These factors are important for considering alternative terrain geometries or use of this algorithm for modes requiring a large amount of shank rotation. D. Prediction accuracy We demonstrated that a simple heuristic predictor could be used to attain a lower prediction error than the best machine learning approach. We recognized that while the machine learning algorithm’s cross validation error represents the error on unseen data, the heuristic was tested on the same data used to train it. Therefore, in order to perform a fairer comparison we also calculated the resubstitution error of the optimal LDA model, thus also testing it on the same data used for training. While resubstitution error was lower than cross validation error, it was still higher than the heuristic error. E. Advantages of heuristic terrain prediction Besides its higher accuracy, our heuristic algorithm also employs only ﬁve features, indicating its high likelihood of generalizing to alternative walking conditions found outside of the laboratory environment. Three of these are estimated ground slope and vertical and horizontal ankle positions, while the other two are max ankle angle and ankle velocity at foot off. These were used to take care of certain special cases such as shufﬂing the feet or turning around, which tended to skew ground slope and ankle position estimates. Another advantage of our approach is that it is intuitive, based directly on the idea that certain aspects of lower limb kinematics are physically enforced by ground terrain geom- etry. This intuitiveness makes our model directly tunable. In comparison, typically machine learning based models either do not have tunable parameters or have tunable parameters that unpredictably or abstractly affect model behavior, without a clear link to performance. Tunability is a highly useful attribute in the ﬁeld because it opens the possibility for performance optimization given speciﬁc walking conditions. This can be done ofﬂine, such as tuning the device to better predict the speciﬁc stair geometries within the end user’s home, or in real-time, such as adapting prosthesis control to estimated load carriage or walking speed. F . Future work Avenues for future work include the real-time application and evaluation of this algorithm under conditions of terrain- speciﬁc prosthesis control. Additionally, it would be inter- esting to investigate the performance of this algorithm in transfemoral prostheses as well, in which terrain enforces the same requirements on lower limb kinematics, even though the user is less able to fulﬁll those requirements. Finally, it would be interesting to incorporate novel walking modes within this heuristic framework, without affecting the accuracy of 8 detecting the current modes. We anticipate that our algorithm can be extended to detecting a wide variety of alternative modes such as stepping over obstacles, walking on highly uneven or soft terrain, or even running and jogging. G. Acknowledgements We would like to acknowledge Tony Shu, Jean-Francois Duval, and Matthew Handford for experimental and technical support. This work was supported by the United States Army Medical Research Acquisition Activity under grant number W81XWH-14-C-0111, the National Defense Science and En- gineering Graduate Fellowship, and the Harvard-MIT Health Sciences and Technology program.