{'fraction_non_alphanumeric': 0.025693236089000777, 'fraction_numerical': 0.006515863907297027, 'mean_word_length': 5.036201117318436, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Predictive Modeling of Periodic Behavior for Human-Robot Symbiotic Walking

We propose in this paper Periodic Interaction Primitives - a probabilistic framework that can be used to learn compact models of periodic behavior. Our approach extends existing formulations of Interaction Primitives to periodic movement regimes, i.e., walking. We show that this model is particularly well-suited for learning data-driven, customized models of human walking, which can then be used for generating predictions over future states or for inferring latent, biomechanical variables. We also demonstrate how the same framework can be used to learn controllers for a robotic prosthesis using an imitation learning approach. Results in experiments with human participants indicate that Periodic Interaction Primitives efficiently generate predictions and ankle angle control signals for a robotic prosthetic ankle, with MAE of 2.21 degrees in 0.0008s per inference. Performance degrades gracefully in the presence of noise or sensor fall outs. Compared to alternatives, this algorithm functions 20 times faster and performed 4.5 times more accurately on test subjects.


I. INTRODUCTION
Walking is a critical motor skill which is at the center of human mobility and independence. Healthy human adults on average walk several thousand steps per day seemingly without any effort and with substantial grace and fluency to their movements. However, for many millions of people [1] affected by musculoskeletal disorders, amputations, neurologic pathologies, or other health-related issues walking is a daily struggle or even completely out of reach [2], [3]. Modern assistive robotics technology (e.g. an exoskeleton, orthotic or prosthetic) has the potential to change the lives of people affected by such conditions for the better, by replacing missing, or augmenting existing capabilities. However, methodologies are needed that allow robots to generate periodic actions that seamlessly blend with those of the human user [4].
In particular, such assistive robots need to be able to anticipate future kinematic states of the human partner given current sensor readings of their walking gait [5], thereby providing a window of opportunity for decisionmaking and control. In addition to predicting kinematics, it is suggested that assistive devices also take biomechanical and ergonomic ramifications on the human body into account [6]. Considering that different walking gaits result in varying magnitudes in internal stresses being applied to the human Fig. 1. PIP is used to learn a realtime, closed-loop controller for symbiotic walking by analyzing and predicting sensor values, as well as control signals. musculoskeletal system; it is important that assistive robots can generate estimates of biomechanical variables in a rapid and predictive fashion to avoid overexertion, injury, or even serious chronic diseases such as osteoarthritis (OA) [7]. Inverse dynamics [8] techniques can generate accurate estimates of internal biomechanical variables, but are generally slow and non-predictive in nature. In addition, they often require information from high-fidelity sensors, usually motion capture and an instrumented treadmill with force plates, which is typically not available outside of controlled laboratory environments. Hence, other methods are needed to be used in non-clinical settings with low-fidelity sensors and can run on devices with limited computational power.
Therefore we propose Periodic Interaction Primitives (PIP) -an extension to interaction primitives (IP) [9]. PIP is a data-driven approach to the predictive modeling of periodic behavior for human-robot symbiotic walking, i.e., effortless walking with an intelligent prosthesis, exoskeleton, or other assistive device. Once trained, a PIP can be used for anticipating the future states of the human user, such as kinematics from observed sensors, and inferring latent, unobserved variables, for instance internal biomechanical variables or external prosthesis or exoskeleton control without the further necessity of a motion capture system or additional biomechanical model. These features are achieved, by correlating at runtime, the information captured from multiple sensors with corresponding biomechanical or control parameters which are incorporated during the training process. PIP is differentiated from previous interaction primitives in three ways. (1) We tie the statistical processes of IP to a periodic domain through the use of von Mises basis functions and provide analytical formulations of partial derivatives with respect to phase. (2) In concert with periodicity modifications, PIP incorporates a novel phase estimation method which is effective in tightly-coupled or symbiotic systems which exhibit periodic phase behavior. (3) We describe a data augmentation scheme that allows PIPs to infer unobservable, biomechanical information through latent variable estimation.

II. RELATED WORK
Early work on control involving prosthetics and exoskeletons has primarily focused on the control aspects of the human robot system [10]. Controllers have been successfully developed for both transtibial [11] and transfemoral [12] applications. Extensions have also been made, which provide control parameters for a variety of common terrains such as navigating slopes [13] or stairs [14]. While these methods have been shown to produce robust and stable controllers, there are a number of downsides. First, these control systems are entirely reactive to the human subject. Second, they must be closely tuned to match each individuals walking gait.
In contrast, data-driven approaches to the probabilistic modeling of HRI often leverage hidden Markov models [15][16][17], which are well-suited for the joint inference problems commonly encountered in interaction. However, they often require a discretization of the state space and suffer from computational overhead too large to be feasibly placed onto the low-power microprocessors that prosthetic devices require.
One method that overcomes some of these limitations is called Probabilistic Movement Primitives (ProMP) [18], [19]. As a probabilistic formulation of the Movement Primitivies concept, this method uses a stochastic feedback controller to reproduce learned trajectory distributions, both stroke based and rhythmic, with predefined phase progression. By incorporating a probabilistic framework (utilizing the demonstration variance) ProMPs allow for greater flexibility in trajectory reproduction than deterministic methods.
In order to increase flexibility and accuracy even further Interaction Primitives were proposed, as a learning from demonstration method in which a joint probability distribution is created over the parameterized models of two interacting agents. By incorporating the full joint trajectory distribution IPs enable the inference of all modeled parameters based on a sample set of observed parameters. This approach has proven capable in several different scenarios [20][21][22][23][24]. However while potential for rhythmic behaviors has been proposed before, careful management of basis functions as well as modifications to phase estimation must be made to to take periodicity into account. Additionally, these methods do not yet address the problem of inferring temporal position in highly dynamic tightly coupled periodic systems; where the behavior of the human influences the behavior of the robot and vice versa. Instead they utilize methods such as dynamic time warping (DTW) to help determine the correct phase of the system for interaction, which are too computationally demanding to be used in run in real-time on assistive devices. This type of symbiosis between the human and robot, where there is constant reaction and interaction, necessitates new methods of quantifying and modeling interactions.

III. PREDICTIVE MODELING WITH PERIODIC INTERACTION PRIMITIVES
The approach proposed in this paper focuses on probabilistic modeling and inference of periodic behavior during walking. The main goal is to learn a predictive model that relates current multimodal sensor readings to future states of the human during the walking gait. In contrast to previous modeling approaches, we simultaneously predict: (a) observable sensor values, (b) unobservable biomechanical data, as well as potentially (c) control values for a prosthesis, exoskeleton or other assistive device.
As an overview of the algorithm presented in this paper; first a data set is collected which represents typical variations in sensor readings during walking. The data set is then augmented with biomechanical information from the subject via tools such as OpenSim [25] or Vicon Nexus [26], as well as control signals for which action an assistive device or prosthesis should take at each timestep. The final step of the learning process is to train the predictive model, which extracts and efficiently represents the correlations and dynamics inherent in the data. The goal of PIP is to represent the reciprocal coupling between observed data (sensor values), unobservable data (extracted biomechanical parameters), and control variables (angular control value for a prosthetic ankle), instead of an analytical formulation of the dynamics.

A. Training: Data Collection and Data Augmentation
One benefit of PIP is that variables which are not observed though the sensors such as: joint torques, muscle lengths, metabolic power consumption, or mechanical work can be inferred. We accomplish this by incorporating the internal biomechanical variables and their distributions w.r.t. time directly into our predictive models; in the following manner.
Each user performs a walking action spanning multiple gait cycles which produces a time series of sensor measurements. This time series is subdivided into a set of N individual demonstrations, each representing a single gait cycle, [y 1 , . . . , y Tn ] ∈ R Ds×Tn , such that D s is the total number of degrees of freedom from all sensors and T n is the number of observations for the n-th demonstration.
Given the set of individual demonstrations, the trajectories of unobservable variables [m 1 , . . . , m Tn ] ∈ R Dm×Tn , such that D m is total number of degrees of freedom for all unobservable variables, is calculated using a desired analysis tool. The original sensor observations are then augmented such that Y 1:t = [y t , m t ] ∈ R D forms the new state representation of the human at every time step, where D = D s + D m .
The objective of this approach is to turn the estimation of internal variables into a function approximation problem which will run significantly faster and will work even when variables are missing. By incorporating unobservable variables into the training process, we can later leverage the PIP framework to infer internal loads using observable variables alone.

B. Training: Learning a Periodic Interaction Primitive
Probabilistically, the objective of PIP is to infer estimates of the future states of both the observed and non-observed variablesŶ t+1:T , given observations of the human Y 1:t : This requires us to define a state transition model for each degree of freedom along with an observation model relating the unobservable variables to the observed variables based on our set of training demonstrations. As previously established [21], we can simplify this process significantly by first transforming our state representation into a latent space via basis function decomposition. That is, we want to find a vector of coefficients w d for each degree of freedom This decomposition accomplishes three goals: a) it produces a time-invariant state w d for each DOF, allowing inference of all past and future values from a single state value, b) it captures the dynamics of each DOF over time, eliminating the need for an explicit transition function, and c) it allows for correlation between different DOFs. A necessary component of the above decomposition is the introduction of a temporal dependence for the basis functions on a relative phase function φ(t) ∈ [0, 100], or for notational simplicity φ. Utilizing phase decouples the relative progress of an interaction from its absolute length, thus preserving the shape of a trajectory across temporal speeds and allows for the efficient calculation of the values for a DOF for an entire interaction: . P here represents the number of sampled points for the trajectory and has the additional benefit of acting as a smoothing function. We succinctly represent all DOFs in our state representation by concatenating them into a single weight vector: By its very nature walking is a periodic action and biomechanical analyses have therefore shown that gait cycles can be quantitatively described by phase plane evaluations [27]. So as to truly capture the periodic nature of our state variables, we must choose an appropriate basis function which operates over a cyclical domain. This is in contrast to prior work in IPs which employ Gaussian basis functions for trajectories with distinct start and end points. In this work, we utilize the von Mises function, an approximation of the wrapped normal distribution over the domain of the unit circle: where µ lies in the interval [0, 2π] and specifies the center of the basis function, α is 2π/100 (required to transform the phase into [0, 2π]), κ is a measure of inverse variance, and I 0 is the modified Bessel function of order zero.

C. Inference: Phase Detection
With the basis space fully defined, we re-formulate our probabilistic objective from (1) as Intuitively, the goal is to infer the state of the interaction w t given t observations of the human and a prior estimate of the state, w 0 . However, without knowing the phase values associated with the observations Y 1:t we cannot apply a correction to our estimate, and so we must first estimate the phase φ for each observation y t . As part of the training process we create a low dimensional manifold, or subspace, encompassing the phase and observed variables. The created manifold allows for the efficient projection of observations onto the manifold so as to determine the corresponding phase.
Before we can construct this manifold model, we must first temporally align the demonstration gait cycles. This is required as gaits are nonlinearly warped in time through variations in walking speed, incline, and other environmental and biomechanical factors. According to existing analyses on complex joint movements, there are quantitative relationships between the positions of the knee, tibia, and ankle with respect to the phase plane. In this work, we exploit these relationships to temporally align our demonstrations through measurements of the angular position, angular velocity, and angular acceleration of the tibia. While the position and velocity are directly observed, the acceleration is not. Instead it is computed through the partial derivative of the linear combination of basis functions with respect to the phase: where d corresponds to the tibia angular velocity DOF. Temporal alignment is performed with dynamic time warping (DTW) [28], in which the distance between two demonstrations u and v is computed as a cost matrix C. The optimal alignment can be found as the minimum path through the cost matrix, where the cost function c(·) is the sum of the Euclidean distances between angular positions, velocities, and accelerations. Since we wish to align all demonstrated trajectories, the full set of trajectories are temporally aligned creating an approximately optimal time alignment between all demonstrations. We assume that all aligned demonstrations lie over the phase interval [0, 100], such that the first time step is 0 and the last is 100. Finally, we discretize the data over two dimensions -E angular position states and F angular velocity states -and assign each discrete state a phase value determined via the nearest neighbor in the set of aligned trajectories. This creates an E × F lookup table, in the subspace or manifold, that can be queried quickly to approximate phase from angular position and velocity. We denote the lookup function L(y t ) which accepts an observation from which the angular position and velocity are projected onto the table and the resulting phase φ is returned.

D. Inference: Generating Predictions and Controls
In order to evaluate the posterior probability of (3), we utilize a full state linear estimator, i.e. Kalman filter, to recursively apply Gaussian conditioning to a matrix of basis weights. Which evaluates the state prediction density as, For tractability, we make the assumption that the posterior state densities are normally distributed, that is, p(w t |Y 1:t , w 0 ) = N (µ t , Σ t ). Additionally, because our state has been transformed into a time-invariant representation, we no longer have a meaningful state transition function in time. That is, p(w t |Y 1:t−1 , w 0 ) = p(w t−1 |Y 1:t−1 , w 0 ). The observation function, h(·), is simply the linear combination of weighted basis functions as described in Sec. III-B and the observation matrix H is defined as: The standard update equations for the calculation of the posterior density follow: where K is the Kalman gain matrix which controls how heavily we update our state estimate based on the observed measurement while taking into account the measurement noise matrix R t . While the observation matrix given here includes both the observed and controlled DOFs, in practice only the observed DOFs are considered. This can be orchestrated by either setting the entries corresponding to the controlled DOFs to zero in H t , or regularizing them with artificial noise in R t . Future trajectories for any degree of freedom d, including the controlled ones, can be trivially calculated with Y d 1: The initial estimate, p(w 0 ) = N (µ 0 , Σ 0 ), is determined from the set of basis weights generated from the demonstrations, where µ 0 is the sample mean and Σ 0 is the sample covariance.
The PIP algorithm is shown in its entirety in Fig. 2. At the first time step, the previous state estimate is given as µ 0
Output:Ŷ Ds 1:P ∈ R P ×Ds : inferred states of D s observed variables,Ŷ Dm 1:P ∈ R P ×Dm : inferred states of D m unobservable variables, µ t : the updated state mean, Σ t : the updated state uncertainty. 1) Estimate phase φ using trained lookup function: 2) Compute the Kalman gain matrix: 3) Incorporate observations into the state estimate:

A. Data Collection
Data collection was performed in a human subject study approved by the Institutional Review Board at Arizona State University. Five participants were fitted with smart shoes, IMUs and retroreflective markers for a VICON motion capture system with 10 cameras. Four IMU devices, one on each shank and one on each femur output 9DOF inertial data as well as angular position from a proprietary embedded sensor fusion algorithm. The smart shoes measure ground reaction forces at four points: heel, first metatarsal joint, fourth metatarsal joint and toe, through silicone tubes which are wound into air bladders, placed in the soles of the shoes, and connected to barometric pressure sensors. To collect biomechanical data, each participant was asked to walk on an instrumented treadmill at 5 different speeds from 0.5m/s to 1.3m/s for 2 min at each speed in a single trial. Motion Capture marker positions and ground reaction force data was collected on the subjects and processed with the common commercially available biomechanics model: Vicon PlugInGait model (version 2) in Vicon Nexus software [26] to calculate joint angles (angles between skeletal links), forces (reflected force from one skeletal link to another), and moments (torque between skeletal links). The entire learning process from raw data to testable model output takes approximately 3 minutes. Subjects were divided into two groups three subjects for training and two for testing. Ten consecutive strides from each subject were removed from the full data set as holdout data for testing and model evaluation.

B. Experiment 1: Predictive Capabilities
The first experiment evaluates the inference and prediction capabilities of the introduced PIP method by assessing: (a) the accuracy of the estimated phase variable when the subjects are walking at different speeds, and (b) the prediction quality of the learned models on observed variables.  Fig. 3 shows an example of the PIP phase detection method as compared to the common method of dynamic time warping, with the phase manifold on the left plot and the detected phase on the right plot. The phase manifold illustrates how each point in the space is directly associated with a specific phase value, when individual observations occur, phase is directly estimated through the use of a 2D lookup function. This leads to a very accurate phase detection which works over a large range of walking speeds. The common phase detection method of Dynamic Time Warping is further compared against the PIP phase detection method to show that the PIP phase detection method is more accurate as it lacks the temporal disturbances commonly present in DTW. During this experiment the inference time of each phase detection method was collected and it was found  Table I shows the mean absolute error (MAE) for prediction of sensor values via a PIP in the top section, using the joint model from all training subjects data, tested on each subject individually. We can see for example that the error on predicting the shank is in the range 1 • to 3 • , hence yielding accurate predictions for all participants. The shank velocity is at approximately 15 • . This may appear substantial, however, the prediction error amounts to a deviation of approximately 1.5%. In addition, all predictions of the foot pressure are lower than 30 mBar and are therefore also of high accuracy. Fig. 4 shows the predictions for the shank angle and the smart shoe heel pressure. The red shaded area highlights the sensor readings, which have been used for conditioning, which generates the predictions starting at 0.25 seconds and ending at 1.4 seconds. For both variables, the PIP predictions (blue) accurately match the ground truth values (black). However, higher uncertainty in heel pressure persists during the swing phase, likely due to the variability in shoe pressure, resulting from air decompressing when lifting the heel. Table II shows how general the model is as well as a comparison to ProMP via. the errors present when tested with: PIP on the three training subjects, PIP two test subjects, and a ProMP model on two test subjects.

C. Experiment 2: Latent Variable Estimation
The second experiment aims at evaluating how accurately the method predicts latent unobservable variables, e.g., biomechanical information that was added to the collected data during the data augmentation step (see Sec. III-A). Fig. 5 depicts the ground truth knee force along with the predicted values of the knee force as inferred via the PIP. The bottom section of Table (I) shows the kinematic and biomechanical variables that are inferred using this technique, namely ankle and knee: moments, forces and angles; along with the MAE for each subject. It is critical to note here, that the ground truth was obtained from the data augmentation process and therefore not composed of direct sensor measurements. Still,  Besides prediction accuracy, a critical aspect when using machine learning components on an assistive device is safety. It is important that predictions degrade gracefully (meaning minimal deterioration of the prediction accuracy) with noisy or missing sensor data. Fig. 5 shows an example for the graceful degradation as exhibited during inference in PIP, where sensors start failing one after the other. In the first section of the figure, all sensors are observed when making predictions. As sensors fail the uncertainty of the prediction grows and that the prediction accuracy slightly decreases. However, in general, the plot of the prediction trajectory still largely follows the ground truth.

D. Experiment 3: Learning Control with PIPs
The last experiment evaluates the ability of PIPs to generate control signals for a powered lower limb prosthesis along with predictions. More specifically, we investigated if ankle angles from a healthy walking gait can be imitated and reproduced when wearing the active prosthesis using PIP. Each participant was asked to wear the SpringActive Odyssey prosthesis [29], a battery powered lower limb prosthesis 1. Imitation of the subjects walking gait is achieved by inferring in real-time the (unobservable) ankle angle and sending it to the prosthesis for execution. Fig. 6 depicts the ankle angle that resulted from this experiment. The black line is an example of ankle angle when no prosthesis is worn. Ankle angles recorded when the prosthesis was used in conjunction with both PIP in blue and the commercial control software originally shipped with the SpringActive Odyssey in yellow, with three example strides of each. It can be seen that while PIP closely tracks the expected ankle angle, the commercial control software substantially deviates from how the participant would naturally move. A ProMP controller was not compared against in this experiment, since the high ankle angle error seen in Table II would possibly make this controller unsafe.

V. CONCLUSION
In this paper, we proposed Periodic Interaction Primitives, a probabilistic framework that can be used to learn compact models of periodic behavior in symbiotic systems. PIP facilitates prediction of future states of walking, as well as generation of control signals for assistive devices through its integrated phase estimation process. Most notably this approach extends Interaction Primitives to latent variables which are either unobservable or hard to measure in practice. For example, when collecting training data, we recorded motion capture readings to calculate offline biomechanical data such as ankle angles and internal forces of the subjects; but did not require additional motion capture readings at runtime. Instead, we use other low-cost sensors to estimate the biomechanical variables at each time-step.
Results in experiments with human participants indicate that Periodic Interaction Primitives efficiently generates predictions and ankle angle control signals for a robotic prosthetic ankle, even of subjects that the model was not trained on, with MAE of 2.21 • in 0.0008s per inference. Even on devices with limited compute power it is expected that inference rates of 1000Hz can easily be achieved. Compared to ProMP, with MAE of 10.10 • in 0.0169s per inference which is too slow and inaccurate for use with a robotic prosthetic. Furthermore it was shown that performance degrades gracefully in the presence of noise or sensor fall outs. This work shows that the approach can quickly generate predictive models and controls using a purely data-driven methodology.