A method to estimate ground reaction forces (GRFs) in a robot/prosthesis system is presented. The system includes a robot that emulates human hip and thigh motion, along with a powered (active) transfemoral prosthetic leg. We design a continuous-time extended Kalman filter (EKF) and a continuous-time unscented Kalman filter (UKF) to estimate not only the states of the robot/prosthesis system but also the GRFs that act on the foot. It is proven using stochastic Lyapunov functions that the estimation error of the EKF is exponentially bounded if the initial estimation errors and the disturbances are sufficiently small. The performance of the estimators in normal walk, fast walk, and slow walk is studied, when we use four sensors (hip displacement, thigh, knee, and ankle angles), three sensors (thigh, knee, and ankle angles), and two sensors (knee and ankle angles). Simulation results show that when using four sensors, the average root-mean-square (RMS) estimation error of the EKF is 0.0020 rad for the joint angles and 11.85 N for the GRFs. The respective numbers for the UKF are 0.0016 rad and 7.98 N, which are 20% and 33% lower than those of the EKF.

## Introduction

Recent advances in microelectronics and robotic technologies have enabled the development of powered prosthetic legs [1] that can help amputees walk up stairs and slopes because of the legs' net power contribution to gait; they are also able to adapt their behavior to various environmental conditions. The company Össur has a lower limb prosthesis called the Power Knee, which is the first commercially available knee to generate power during the gait cycle [2]. Recent developments in powered knee and ankle prostheses are reported in Refs. [3–8].

Bulky load cells are often employed in robots and prosthetic legs to capture external forces (GRFs) and moments during walking [9]. These data are used as feedback measurements to control the robot or prosthesis with force or impedance controls. Even when not used for feedback, force sensing is important for monitoring and evaluation of prosthesis performance and safety. However, there are several drawbacks to the use of load cells: (i) load cells are expensive; (ii) a 250-lbf load cell weighs about 1 lb, with a length of about 3 in, and thus does not easily fit in a prosthetic leg; (iii) load cell measurements tend to drift, need to be frequently offset, and are noisy and need significant signal conditioning; (iv) load cells can get damaged easily from overloading or off-axis loading; and (v) load cells consume electrical power, which is an important consideration in prosthetics.

The above-mentioned problems do not arise with angle sensors because high-resolution encoders are accurate, reliable, and inexpensive. Nevertheless, they do not measure velocity, and velocity calculated by numerical differentiation is challenging because of the difficult compromise between noise rejection and bandwidth.

There have been several methods to reduce the number of force sensors in robotics. In Ref. [10], a robot compliance controller based on a disturbance observer is presented, where the disturbance observer is used to estimate the external reaction forces. In Ref. [11], a method for force estimation of the end-effector of a Selective Compliance Assembly Robot Arm (SCARA) robot is presented, where servomotor currents and position information are used to estimate forces. In Ref. [12], the external force on the end-effector of a four degrees-of-freedom (DOF) robot manipulator is estimated with a combination of time delay estimation and input estimation. In Ref. [13], a model-based observer is used to estimate the external forces acting on a rigid body. All the aforementioned techniques depend on either the accuracy of the robot model or servomotor current data. Thus, if the accuracy of the robot model or motor current measurement degrades, force estimation can deviate from true force. In Ref. [12], the force estimation does not require a precise robot dynamic model, but the accuracy of estimated force may degrade in the presence of noise, since the external force is considered as an unknown input.

In this paper, we treat GRF as an unknown input with known dynamic properties and known bounds, and we use an EKF to estimate GRFs along with the states of the robot/prosthesis system. Although EKF is efficient in many applications, it has two important potential drawbacks. First, the derivation of the Jacobian matrix for linearization can be complex and can cause numerical implementation difficulties. Second, linearization can lead to cumulative errors which may affect the accuracy of the state estimation. To overcome these limitations, other nonlinear estimators could be used, such as the UKF or particle filter, where the state estimations are obtained without the need for derivatives and Jacobian calculations [14–16].

In this research, a powered (active) prosthetic leg is considered for transfemoral amputees. The prosthetic leg attached to a robotic hip/thigh emulator. The combined system includes four degrees-of-freedom: vertical hip displacement, thigh angle, knee angle, and ankle angle. This paper is an extended version of Ref. [17], where we designed an EKF for the online estimation of joint coordinates, velocities, and GRFs.

In this paper, we compare EKF and UKF performance in our 4DOF robot/prosthesis system. We study the EKF because it is the most commonly used nonlinear state estimator due to its versatility and simplicity of implementation. We study the UKF because it generally provides better performance than the EKF. One of the questions we study in this paper is whether the improved performance of the UKF relative to the EKF is worth the increased computational effort.

We also study the effect of various sensor sets on estimation error at different walking speeds. We start with the four main states of the robot/prosthesis system, which are hip displacement, thigh, knee, and ankle angles, as system observations to get the baseline estimation accuracy. Then, we investigate how much estimation accuracy degrades if fewer measurements are used, and we find a tradeoff between the estimation error accuracy and the number of measurements. In general, we prefer to use fewer measurements because fewer sensors mean a simpler system, less complexity, and lower cost.

Finally, we consider the stability of the EKF and boundedness of the estimation error, which are very important for prosthesis control. The stability of the EKF under various values of initial estimation error and noise covariances is explored. We analytically show that the estimation error remains bounded under certain conditions in our 4DOF robot/prosthesis system and confirm the analysis with simulation results.

This research involves the mechanical integration of a prosthesis and a test robot so that we can test the prosthesis without human trials. We will eventually need to use the Kalman filter on a prosthesis that is attached to an amputee. We develop the Kalman filter for the robot/prosthesis system [18] in this paper for later implementation on an amputee/prosthesis system.

The paper is organized as follows: In Sec. 2, the model of the robotic system and prosthetic leg is presented. In Sec. 3, the EKF for state estimation and GRF estimation is discussed. In Sec. 4, the EKF is analyzed by mathematically deriving its stability conditions. Section 5 compares performance between the EKF and UKF when different measurement sensors are used; also the convergence of the EKF is tested with different initial estimation errors and noise magnitudes. Section 6 concludes the paper and suggests future research.

## System Model

Robotic testing of transfemoral prostheses is presented here, where motion is limited to the sagittal plane and transverse motion is not considered. Typically, only the sagittal plane is considered in transfemoral prosthesis research. Although the transverse plane is ignored here, sagittal motion captures the essential dynamics of human walking [19]. The model of the test robot/prosthesis is based on the standard robotic framework. Figure 1 shows a diagram of the hip robot and prosthesis combination [20,21]. A general dynamic model for the system is given as follows:

where $q=[q1,q2,q3,q4]T$ is the vector of joint displacements (*q*_{1} is vertical hip displacement, *q*_{2} is thigh angle, *q*_{3} is knee angle, and *q*_{4} is ankle angle), $D(q)$ is the inertia matrix, $C(q,q\u0307)$ is a matrix accounting for centripetal and Coriolis effects, $B(q,q\u0307)$ is a nonlinear damping vector, *J _{e}* is the kinematic Jacobian relative to the point of application of the external forces

*F*, $g(q)$ is the gravity vector, and

_{e}*u*is the four-element vector of control signals [22]. The kinematic and dynamic models of the robot/prosthesis combination are given in Refs. [23–25], where a mixed tracking/impedance controller based on passivity methods is designed [22]. The control signals consist of hip force, and thigh, knee, and ankle torques. As Fig. 1 shows, a triangular foot with two points of ground contact is assumed. Horizontal and vertical GRFs are applied to contact points at the toe and heel. The GRFs are denoted as

*F*,

_{xh}*F*,

_{zh}*F*,

_{xt}*F*, which represent the horizontal and vertical GRFs at the heel and toe. Thus, the external force vector

_{zt}*F*in Eq. (1) comprises these four GRFs.

_{e}*x*-axis. A treadmill provides the walking surface of the robot/prosthesis system. The belt stiffness is modeled to calculate reaction forces during contact between the heel and toe with the treadmill [22,23]. The GRFs are entirely determined by kinematics and are given as follows [23]:

where *k _{b}* is the belt stiffness,

*s*is the treadmill standoff (the distance between the origin of the world coordinate system and the belt when the leg is fully extended),

_{z}*l*

_{2}and

*l*

_{3}are the lengths of link 2 (thigh) and link 3 (shank), respectively,

*a*and

_{H}*a*are the distance from the ankle joint to the heel and the toe, respectively,

_{T}*a*is the height of the ankle joint above the sole of the foot, and

_{h}*β*is the friction coefficient between the treadmill belt and the foot. The vertical positions of the toe and heel in the world coordinate system are shown in Fig. 1 as

*z*and

_{t}*z*, respectively. We thus have four states for the positions and four states for the velocities of the joint displacements.

_{h}*F*in Eq. (1). We, therefore, augment the external forces to the state vector. The augmented state vector includes the eight original states of the robot and the four GRFs and is given as follows:

_{e}It should be noted that the forces are not states of the original system model, but are augmented here to the state vector, since we need to estimate them with a state estimator. The 12-element vector of Eq. (8) will be estimated by the state estimator.

Equations (2)–(7) are a preliminary approximation to GRF and are accurate only if the amputee (or robot in our case) walks in a highly controlled environment with a known walking surface stiffness. In future work, it will be important to study the robustness of the Kalman filter to modeling errors and particularly to errors in these GRF equations.

We should mention that there are at least two other possible approaches to obtain GRF estimates. In the first alternative approach, GRF could be considered as an unknown input with no prior assumptions about its model. In that case, the EKF would have difficulty estimating GRF, since its effect on the system would be indistinguishable from process noise. In the second alternative approach, state estimation of the robot/prosthesis system could be performed by a nonlinear observer, such as the EKF, and then the estimated states of the robot could be substituted into the GRF Eqs. (2)–(7) to obtain the GRF estimates. However, this approach would be less flexible than the approach that we are using. Augmenting the GRF states onto the original state as in Eq. (8) allows us additional tuning flexibility in the artificial process noise that is included in the augmented states, as we will see in Sec. 5.

## Extended Kalman Filtering for Robot/Prosthesis State Estimation

*x*(

*t*), the input

*u*(

*t*), and the output

*y*(

*t*) are in $\mathbb{R}q,\u2009\mathbb{R}p$ and $\mathbb{R}m$, respectively. Moreover, the noise terms

*w*(

*t*) and

*v*(

*t*) are in $\mathbb{R}l$ and $\mathbb{R}k$ and are uncorrelated, zero-mean white noise processes with identity covariances.

*G*(

*t*) and

*D*(

*t*) are time-varying matrices of size

*q*×

*l*and

*m*×

*k*. If the nonlinear functions

*f*and

*h*are sufficiently smooth in

*x*, Taylor series expansions can be performed. The following Jacobian matrices are used to linearize the system:

*x*(0) is a random vector with covariance

*P*(0)

*x*(0) is assumed to be uncorrelated with

*w*(

*t*) and

*v*(

*t*). The EKF equations are given as

*K*(

*t*) and

*P*(

*t*) are the Kalman gain and estimation error covariance matrix. Covariance matrices

*Q*(

*t*) and

*R*(

*t*) are

*q*×

*q*and

*m*×

*m*, respectively, and are defined as follows:

## The Stability of the Extended Kalman Filter

The analysis of the stability properties of continuous-time extended Kalman filters is complex and has been treated only for especial cases. Safonov and Athans [28] considered the stability of the constant and modified gain EKF. The boundedness of the EKF estimates were investigated when used as a parameter estimator for linear systems [29,30]. Reif et al. [31,32] determined that the estimation error in both the discrete and continuous-time EKF is exponentially bounded under certain conditions. They showed that the estimation error remains bounded if the initial estimation error and disturbances are sufficiently small and the nonlinear system satisfies a detectability rank condition [33,34]. Here, we extend Reif's work [31] to obtain the stability conditions of the state estimator of the robot/prosthesis system. The main purpose of this section is to derive the stability of the estimator in terms of initial condition errors and disturbances.

*X*are the remaining terms of the Taylor series. The estimation error $\xi (t)=x(t)\u2212x\u0302(t)$ can be written recursively as follows:

In stability theory, there are two commonly used methods to prove boundedness: supermartingales for stochastic differential equations, and Lyapunov functions for deterministic differential equations. Supermartingales can be regarded as the stochastic equivalent of Lyapunov functions in some cases [35]. In order to analyze the error dynamics of the EKF, we first present some preliminary results [31,35]. Note that $\Vert .\Vert $ represents the spectral norm of a matrix or the Euclidian norm of a vector.

*Suppose there is a stochastic process*$V[\xi (t),t]$

*with*$\xi (t)\u2208\mathbb{R}q$

*governed by Eq.*

*(24)*

*and there exist real numbers*$v1,v2,\mu ,\u03d1\u22650$

*such that*

*Proof of Lemma 4.1*. See Appendix A.

Theorem 4.1. *Consider the continuous-time nonlinear dynamic system Eqs.**(9)**and**(10)**and the extended Kalman filter Eqs.**(13)**–**(15)**. Suppose that the following conditions hold:*

- (i)There exist
*p*_{1},*p*_{2},*c*_{1},*c*_{2},*q*_{1},*q*_{2},*r*_{1},*r*_{2}such that the following bounds hold for every $t\u22650$:$p1I\u2264P(t)\u2264p2I$(31)$c1\u2264\Vert C(t)\Vert \u2264c2$(32)$q1I\u2264Q(t)\u2264q2I$(33)$r1I\u2264R(t)\u2264r2I$(34) - (ii)There exist $\epsilon \phi $,
*ε*, $K\phi ,\u2009Kx\u22650$ such that the nonlinear functions in Eq. (26) are bounded by_{x}$\Vert \phi (x,x\u0302,u)\Vert \u2264K\phi \Vert x\u2212x\u0302\Vert 2$(35)for$\Vert X(x,x\u0302)\Vert \u2264Kx\Vert x\u2212x\u0302\Vert 2$(36)*x*, $x\u0302\u2208\mathbb{R}q$ and $u\u2208\mathbb{R}p$ with $\Vert x\u2212x\u0302\Vert \u2264\epsilon \phi $ and $\Vert x\u2212x\u0302\Vert \u2264\epsilon x$, respectively.

*Proof of Theorem* 4.1. See Appendix B.

*K*, we can consider a compact subset $N\u2208\mathbb{R}q$; the constants $K\phi $,

_{x}*K*can be obtained by

_{x}If *f* and *h* are twice differentiable with respect to *x* for every $x\u2208N$, the spectral norm of the Hessian matrices of *f _{i}* and

*h*are bounded, where

_{i}*f*and

_{i}*h*are the components of

_{i}*f*and

*h*, respectively (see Refs. [31] and [38, Chap. 7]). In the given robot/prosthesis system, the Hessian matrix of

*f*is not a function of control input

*u*.

Note that the nonlinear system has to be uniformly detectable in order to satisfy the boundedness condition of Eq. (31) [39,40]. The detectability of the nonlinear system results in the boundedness of the estimation error bounds for the solution *P*(*t*) of the Riccati differential equation (15) based on Eq. (31). We present the following definition for the detectability of nonlinear stochastic systems [39].

*The pair*$[(\u2202f/\u2202x)(x,u),(\u2202h/\u2202x)(x)],\u2009x\u2208\mathbb{R}q,\u2009u\u2208\mathbb{R}p$

*is uniformly detectable if there exists a bounded matrix valued function*$\Lambda (x)$

*such that*

holds for all *x*, $\eta \u2208\mathbb{R}q$, where *λ _{m}* is the maximum eigenvalue of $[(\u2202f/\u2202x)(x,u)+\Lambda (x)(\u2202h/\u2202x)(x)]$.

## Simulation Results

*A*matrix in Eq. (11) is 12 × 12. For our first simulation, we use four states as the measurements: vertical hip position, thigh angle, knee angle, and ankle angle. Therefore, the dimension of the measurement matrix

*C*in Eq. (12) is 4 × 12. The initial value of the state vector and estimate, and the diagonal covariance matrices

*Q*and

*R*for the process and measurement noise are given as

where $x\u03020$ differs from *x*(0) due to random initial state estimation errors, $\Delta t$ denotes the discretization step size with $\Delta t=0.5ms$, *Q* and *R* represent continuous-time noise covariances, and *Q _{d}* and

*R*represent discrete-time noise covariances [27, Chap. 8]. The process noise represents small values of unmodeled dynamics and parameter uncertainties. The

_{d}*Q*matrix is chosen based on prior experience with the accuracy of system dynamic modeling; the last four elements of the diagonal matrix

*Q*are zero since the GRF model is assumed to be perfectly known. The

*R*matrix is usually straightforward to determine on the basis of our knowledge of the accuracy of the measurement system.

### Reduced Sensor Sets and Unscented Kalman Filtering.

In a real-world scenario, it would not be practical to use a hip displacement sensor on an amputee because of its invasiveness. The other sensors—thigh angle, knee angle, and ankle angle—can easily be used in the real world because they can be mounted on the prosthesis. In this paper, we use a hip displacement sensor to provide a baseline scenario for the reduced sensor sets, which we study later in this section.

More observations result in more accurate state estimation. However, our goal is to use the fewest possible number of sensors for estimation, which result in a reduction of complexity and cost in the prosthesis hardware. Here, we remove the hip displacement sensor and use only three measurements: the thigh, knee, and ankle angles. For our final test in this section, we will use only two measurements: the knee and ankle angles. Moreover, we will also compare the performance of the EKF and UKF for the state and GRF estimation. We use the same initial values for state and estimated state as used in the EKF for normal walking. We also test the performance of the filters with different gait speeds, and we will see that no extra filter tuning is needed for fast and slow walking.

The original discrete-time form of the UKF has been widely used for state estimation of discrete-time systems. However, the discrete-time UKF cannot be directly applied to continuous-time filtering problems, in which the process and measurement equations are modeled as continuous-time stochastic systems. In Ref. [15], the continuous-time UKF was derived from the discrete-time UKF in matrix form, and this is the algorithm that we implement here. The tuning parameters in the unscented transforms were chosen to be $\alpha =0.9$, *β* = 2, and $\kappa =\u22123$.

Table 1 compares the accuracy of EKF state estimation in terms of the RMSE in three different gait modes (fast walk, normal walk, slow walk). The EKF works well for different walking speeds with a reduced measurement set. It can be seen, as expected, that for all three walking speeds, the estimation errors with four measurements are smaller than with three or two measurements.

Gait modes | $No.$ measurements | $x1$ (m) | $x2$ (rad) | $x3$ (rad) | $x4$ (rad) | $x5$ (m/s) | $x6$ (rad/s) | $x7$ (rad/s) | $x8$ (rad/s) | $x9$ (N) | $x10$ (N) | $x11$ (N) | $x12$ (N) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Normal walking | 4 | 0.003 | 0.002 | 0.003 | 0.005 | 0.022 | 0.051 | 0.084 | 0.322 | 3.145 | 9.758 | 4.247 | 21.015 |

3 | 0.008 | 0.003 | 0.003 | 0.005 | 0.057 | 0.091 | 0.214 | 0.524 | 5.015 | 16.247 | 9.365 | 38.146 | |

2 | 0.036 | 0.009 | 0.004 | 0.004 | 0.085 | 0.352 | 0.314 | 0.558 | 9.524 | 23.287 | 15.174 | 49.561 | |

Fast walking | 4 | 0.005 | 0.001 | 0.002 | 0.002 | 0.013 | 0.055 | 0.0946 | 0.218 | 3.124 | 10.015 | 4.647 | 22.144 |

3 | 0.009 | 0.002 | 0.004 | 0.005 | 0.048 | 0.137 | 0.376 | 0.933 | 8.715 | 22.014 | 15.247 | 46.984 | |

2 | 0.042 | 0.007 | 0.004 | 0.006 | 0.092 | 0.445 | 0.765 | 1.234 | 14.845 | 28.956 | 21.568 | 53.175 | |

Slow walking | 4 | 0.007 | 0.001 | 0.001 | 0.001 | 0.015 | 0.036 | 0.065 | 0.159 | 2.399 | 8.995 | 3.512 | 18.621 |

3 | 0.008 | 0.004 | 0.003 | 0.006 | 0.058 | 0.143 | 0.275 | 0.856 | 5.067 | 19.324 | 13.189 | 42.872 | |

2 | 0.035 | 0.007 | 0.005 | 0.006 | 0.058 | 0.355 | 0.586 | 1.058 | 10.188 | 25.357 | 18.957 | 50.143 |

Gait modes | $No.$ measurements | $x1$ (m) | $x2$ (rad) | $x3$ (rad) | $x4$ (rad) | $x5$ (m/s) | $x6$ (rad/s) | $x7$ (rad/s) | $x8$ (rad/s) | $x9$ (N) | $x10$ (N) | $x11$ (N) | $x12$ (N) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Normal walking | 4 | 0.003 | 0.002 | 0.003 | 0.005 | 0.022 | 0.051 | 0.084 | 0.322 | 3.145 | 9.758 | 4.247 | 21.015 |

3 | 0.008 | 0.003 | 0.003 | 0.005 | 0.057 | 0.091 | 0.214 | 0.524 | 5.015 | 16.247 | 9.365 | 38.146 | |

2 | 0.036 | 0.009 | 0.004 | 0.004 | 0.085 | 0.352 | 0.314 | 0.558 | 9.524 | 23.287 | 15.174 | 49.561 | |

Fast walking | 4 | 0.005 | 0.001 | 0.002 | 0.002 | 0.013 | 0.055 | 0.0946 | 0.218 | 3.124 | 10.015 | 4.647 | 22.144 |

3 | 0.009 | 0.002 | 0.004 | 0.005 | 0.048 | 0.137 | 0.376 | 0.933 | 8.715 | 22.014 | 15.247 | 46.984 | |

2 | 0.042 | 0.007 | 0.004 | 0.006 | 0.092 | 0.445 | 0.765 | 1.234 | 14.845 | 28.956 | 21.568 | 53.175 | |

Slow walking | 4 | 0.007 | 0.001 | 0.001 | 0.001 | 0.015 | 0.036 | 0.065 | 0.159 | 2.399 | 8.995 | 3.512 | 18.621 |

3 | 0.008 | 0.004 | 0.003 | 0.006 | 0.058 | 0.143 | 0.275 | 0.856 | 5.067 | 19.324 | 13.189 | 42.872 | |

2 | 0.035 | 0.007 | 0.005 | 0.006 | 0.058 | 0.355 | 0.586 | 1.058 | 10.188 | 25.357 | 18.957 | 50.143 |

Table 2 shows that the UKF achieves smaller RMS estimation errors than the EKF for almost all measurement sets and walking speeds. Table 3 summarizes the average performance of the EKF and UKF. We see that the UKF achieves 30% improvement in the average RMSE of GRFs with four measurements. Although the UKF performs better than the EKF, its computational effort is between two and three times that of the EKF in terms of the number of multiplications and additions [15].

Gait modes | $No.$ measurements | $x1$ (m) | $x2$ (rad) | $x3$ (rad) | $x4$ (rad) | $x5$ (m/s) | $x6$ (rad/s) | $x7$ (rad/s) | $x8$ (rad/s) | $x9$ (N) | $x10$ (N) | $x11$ (N) | $x12$ (N) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Normal walking | 4 | 0.003 | 0.002 | 0.002 | 0.002 | 0.015 | 0.025 | 0.08 | 0.184 | 2.547 | 7.985 | 3.521 | 18.687 |

3 | 0.006 | 0.003 | 0.002 | 0.005 | 0.066 | 0.089 | 0.252 | 0.604 | 4.170 | 14.178 | 8.115 | 36.027 | |

2 | 0.035 | 0.009 | 0.003 | 0.004 | 0.079 | 0.362 | 0.264 | 0.452 | 10.015 | 23.057 | 16.185 | 47.111 | |

Fast walking | 4 | 0.004 | 0.001 | 0.002 | 0.002 | 0.011 | 0.015 | 0.0614 | 0.205 | 2.552 | 8.095 | 4.012 | 19.574 |

3 | 0.007 | 0.002 | 0.004 | 0.005 | 0.033 | 0.139 | 0.289 | 0.721 | 7.075 | 20.584 | 13.449 | 44.159 | |

2 | 0.037 | 0.007 | 0.003 | 0.006 | 0.078 | 0.415 | 0.611 | 1.002 | 11.899 | 24.152 | 20.222 | 49.788 | |

Slow walking | 4 | 0.005 | 0.001 | 0.001 | 0.001 | 0.008 | 0.022 | 0.054 | 0.160 | 2.005 | 6.778 | 3.444 | 16.550 |

3 | 0.007 | 0.004 | 0.003 | 0.005 | 0.044 | 0.186 | 0.201 | 0.699 | 4.997 | 18.875 | 12.854 | 39.921 | |

2 | 0.028 | 0.007 | 0.005 | 0.006 | 0.048 | 0.368 | 0.502 | 0.917 | 9.187 | 24.179 | 17.045 | 47.854 |

Gait modes | $No.$ measurements | $x1$ (m) | $x2$ (rad) | $x3$ (rad) | $x4$ (rad) | $x5$ (m/s) | $x6$ (rad/s) | $x7$ (rad/s) | $x8$ (rad/s) | $x9$ (N) | $x10$ (N) | $x11$ (N) | $x12$ (N) |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|

Normal walking | 4 | 0.003 | 0.002 | 0.002 | 0.002 | 0.015 | 0.025 | 0.08 | 0.184 | 2.547 | 7.985 | 3.521 | 18.687 |

3 | 0.006 | 0.003 | 0.002 | 0.005 | 0.066 | 0.089 | 0.252 | 0.604 | 4.170 | 14.178 | 8.115 | 36.027 | |

2 | 0.035 | 0.009 | 0.003 | 0.004 | 0.079 | 0.362 | 0.264 | 0.452 | 10.015 | 23.057 | 16.185 | 47.111 | |

Fast walking | 4 | 0.004 | 0.001 | 0.002 | 0.002 | 0.011 | 0.015 | 0.0614 | 0.205 | 2.552 | 8.095 | 4.012 | 19.574 |

3 | 0.007 | 0.002 | 0.004 | 0.005 | 0.033 | 0.139 | 0.289 | 0.721 | 7.075 | 20.584 | 13.449 | 44.159 | |

2 | 0.037 | 0.007 | 0.003 | 0.006 | 0.078 | 0.415 | 0.611 | 1.002 | 11.899 | 24.152 | 20.222 | 49.788 | |

Slow walking | 4 | 0.005 | 0.001 | 0.001 | 0.001 | 0.008 | 0.022 | 0.054 | 0.160 | 2.005 | 6.778 | 3.444 | 16.550 |

3 | 0.007 | 0.004 | 0.003 | 0.005 | 0.044 | 0.186 | 0.201 | 0.699 | 4.997 | 18.875 | 12.854 | 39.921 | |

2 | 0.028 | 0.007 | 0.005 | 0.006 | 0.048 | 0.368 | 0.502 | 0.917 | 9.187 | 24.179 | 17.045 | 47.854 |

$No.$ measurements | Filter | Ave. RMSE $angular$ positions (rad) | Ave. RMSE$angular$ velocities (rad/s) | Ave. RMSE GRFs (N) |
---|---|---|---|---|

4 | EKF | 0.0020 | 0.1205 | 11.8525 |

UKF | 0.0016 | 0.0896 | 7.9792 | |

3 | EKF | 0.0039 | 0.3943 | 19.8650 |

UKF | 0.0037 | 0.3533 | 18.7003 | |

2 | EKF | 0.0058 | 0.6297 | 26.7279 |

UKF | 0.0055 | 0.5437 | 25.0579 |

$No.$ measurements | Filter | Ave. RMSE $angular$ positions (rad) | Ave. RMSE$angular$ velocities (rad/s) | Ave. RMSE GRFs (N) |
---|---|---|---|---|

4 | EKF | 0.0020 | 0.1205 | 11.8525 |

UKF | 0.0016 | 0.0896 | 7.9792 | |

3 | EKF | 0.0039 | 0.3943 | 19.8650 |

UKF | 0.0037 | 0.3533 | 18.7003 | |

2 | EKF | 0.0058 | 0.6297 | 26.7279 |

UKF | 0.0055 | 0.5437 | 25.0579 |

### Stability of the Extended Kalman Filter.

It has been shown that the estimation error in the continuous-time Kalman filter remains bounded under certain conditions. Small initial estimation errors and small noise magnitudes are required conditions to obtain error bounds. In this section, we study the stability of the state estimates of the robot/prosthesis system. First, we need to check that the Jacobian matrices in the robot/prosthesis system satisfy the uniform detectability condition of Definition 4.1. To find an appropriate $\Lambda (x)$, we assume that $[(\u2202f/\u2202x)(x,u)+\Lambda (x)(\u2202h/\u2202x)(x)=\u2212Iq\xd7q]$ and solve for $\Lambda (x)=(\u2212I\u2212(\u2202f/\u2202x))(\u2202h/\u2202x)T((\u2202h/\u2202x)(\u2202h/\u2202x)T)\u22121$. Note that with our measurement system, $((\u2202h/\u2202x)(\u2202h/\u2202x)T)\u22121$ is an identity matrix and is thus invertible. Therefore, the uniform detectability condition of Definition 4.1 is satisfied in the robotic/prosthesis system. To test the stability of the robotic/prosthesis EKF, the initial value $x\u0302(0)$ and the measurement and process noise covariances *R*, *Q* are chosen as shown in Table 4, where

Small initial error and small noise | Small initial error and large noise | Large initial error and small noise | |
---|---|---|---|

$x\u0302(0)$ | $x\u0302good(0)$ | $x\u0302good(0)$ | $x\u0302poor(0)$ |

Q | $10\u22125I12\xd712$ | $10\u22122I12\xd712$ | $10\u22125I12\xd712$ |

R | $10\u22124I4\xd74$ | $I4\xd74$ | $10\u22124I4\xd74$ |

Error behavior | Bounded | Divergent | Divergent |

Fig. 5 | Fig. 6 | Fig. 7 |

Small initial error and small noise | Small initial error and large noise | Large initial error and small noise | |
---|---|---|---|

$x\u0302(0)$ | $x\u0302good(0)$ | $x\u0302good(0)$ | $x\u0302poor(0)$ |

Q | $10\u22125I12\xd712$ | $10\u22122I12\xd712$ | $10\u22125I12\xd712$ |

R | $10\u22124I4\xd74$ | $I4\xd74$ | $10\u22124I4\xd74$ |

Error behavior | Bounded | Divergent | Divergent |

Fig. 5 | Fig. 6 | Fig. 7 |

$x\u0302good(0)$ is a reasonably accurate initial condition for the state estimate, and we will see that it is good enough to ensure convergence of the EKF. $x\u0302bad(0)$ is a significantly worse initial condition for the state estimate, and we will see that it is not good enough to ensure convergence of the EKF. The *R* and *Q* values in Table 4 have been chosen to demonstrate conditions that, respectively, ensure, or do not ensure, EKF convergence. These values for $x\u0302good(0),\u2009x\u0302bad(0)$, and *R* and *Q* were chosen by trial and error to demonstrate their effect on EKF convergence.

The simulation results are presented in Figs. 5–7, where the unknown state $x6(t)$ in Eq. (8) (angular velocity of thigh), the estimated state $x\u03026(t)$, and its estimation error $\xi 6(t)$ are plotted. We can see in Fig. 5 that for small initial estimation error and small noise, the estimation error is bounded. However, Figs. 6 and 7 show that for large initial estimation errors or large noise magnitudes, the estimation error is no longer bounded.

Extensive simulation tests show that boundedness of the estimation error is obtained in the numerical simulation if $\Vert \xi (0)\Vert \u2264\epsilon =0.35,\u2009\u2009G(t)GT(t)\u2264\delta =2\xd710\u22122$ and $D(t)DT(t)\u2264\delta =2\xd710\u22127$. The theoretical results for *ε* and *δ* via Eqs. (B21) and (B24) from Appendix B yield smaller bounds for stability: $\epsilon \u22646.5\xd710\u22125$ and $\delta \u22648.7\xd710\u221213$ with $N={x\u2208\mathbb{R}q,\Vert x\Vert \u22641000}$. Since the estimation error does not diverge in practice even with larger initial condition errors and noise terms than those given by the theorem, we conclude that the bounds are very conservative in this system. However, we note that it is not possible to test all possible conditions that exceed the theoretical error bounds; therefore, we naturally expect the theoretical results to yield more conservative stability bounds than the simulation results.

We want the state to be estimated as accurately as possible for eventual implementation in a state-feedback controller. The controller that we used in this research is robust to estimation errors [22], although more accurate state estimates are always desirable in order to reduce controller errors. We have not explored the relationship between estimation error and controller error in this paper but leave it as an important area for future research.

## Conclusion and Future Work

We designed an EKF and an UKF to estimate not only the states of a robot/prosthesis system but also the external forces acting on the prosthetic foot. This approach removes the need for heavy and bulky load cells that are otherwise needed for GRF estimation. We achieved satisfactory estimation errors in various gait speeds for the robot/prosthesis system using four, three, and two measurements. The average RMS estimation errors of the EKF for the thigh, knee, and ankle angles in normal walking with four measurements is 0.0033 rad, in comparison with that the UKF which is 0.0020 rad. Although the UKF outperforms the EKF, it requires more computational effort than the EKF, which will be a consideration for real-time implementation.

We proved mathematically that the estimation error in the EKF is exponentially bounded if the initial estimation errors and disturbances are sufficiently small and if the nonlinear system satisfies a detectability rank condition. In simulation tests, we verified that the estimation errors remained bounded for small initial estimation errors and small disturbances. However, the filter is unstable if the initial estimation errors or disturbances are too large.

As far as we know, the present research is the first time that GRF has been estimated for prostheses using state estimation techniques. As we noted in the introduction, some other methods have been used for external force estimation in robotics. For instance, [12] achieved a 1.72% estimation error for a 4DOF robot, compared to our results, which show a 1.85% error when using the EKF with four sensors during normal walking. Our results are, therefore, quantitatively similar to Ref. [12], although the comparison may not be fair since the robotic systems in the two approaches are much different. Other publications in the area of external force estimation do not include quantitative results [10,11,13].

Future work will include increasing and verifying the robustness of the filters, and experimental implementation and verification of the EKF and UKF, first on robotic hardware [21] and then in human trials.

• NSF Grant No. 1344954.

### Appendix A

*Proof of Lemma 4.1*. Let $\xi (t)$ be the unique solution to the stochastic differential in Eq. (24) and let $V[\xi (t),t]$ be a stochastic process which is a non-negative and twice differentiable. Suppose that $V[\xi (t),t]$ satisfies the conditions in Eqs. (27) and (28). By applying Ito's formula [42] we can derive

*t*yields

### Appendix B

*Proof of Theorem 4.1*. To prove this theorem, we need to construct a stochastic Lyapunov function to satisfy the conditions of Lemma 4.1. The stochastic Lyapunov function is chosen as follows:

*P*(

*t*) is positive definite with probability one. Therefore, from Eq. (31), it follows that

*K*is defined as

_{d}*ε*as

In summary, the assumptions of Lemma 4.1 Eqs. (27) and (28) are satisfied by Eqs. (B3) and (B22), where $(1/p2)=v1,\u2009(1/p1)=v2,\u2009((q1p1r2+c12p1p22)/2p22r2)=\u03d1$ and $Kd\delta =\mu $. We conclude that the estimation error is exponentially bounded in mean square under the conditions in Eqs. (37)–(39).

*ε*from Eq. (B21) and calculate

*δ*as follows. We need to take care that $\epsilon \u0303\u2264\Vert \xi (t)\Vert \u2264\epsilon $ with some $\epsilon \u0303\u2264\epsilon $ satisfies the inequality in Eq. (B22)