A Method to Accurately Estimate the Muscular Torques of Human Wearing Exoskeletons by Torque Sensors

In exoskeletal robots, the quantification of the user’s muscular effort is important to recognize the user’s motion intentions and evaluate motor abilities. In this paper, we attempt to estimate users’ muscular efforts accurately using joint torque sensor which contains the measurements of dynamic effect of human body such as the inertial, Coriolis, and gravitational torques as well as torque by active muscular effort. It is important to extract the dynamic effects of the user’s limb accurately from the measured torque. The user’s limb dynamics are formulated and a convenient method of identifying user-specific parameters is suggested for estimating the user’s muscular torque in robotic exoskeletons. Experiments were carried out on a wheelchair-integrated lower limb exoskeleton, EXOwheel, which was equipped with torque sensors in the hip and knee joints. The proposed methods were evaluated by 10 healthy participants during body weight-supported gait training. The experimental results show that the torque sensors are to estimate the muscular torque accurately in cases of relaxed and activated muscle conditions.


Introduction
In recent years, there has been increasing interest in using robotic devices to assist in the rehabilitative training of people with motion impairments. Most of the initially developed OPEN ACCESS rehabilitation robots only provide passive-mode training, which moves the user's limbs along a predefined fixed trajectory. In recent years, many researchers have insisted that robotic assistance should be adaptive according to the user's contribution for more effective and optimal training [1,2]. In this robotic training paradigm, the quantification of the user's muscular effort is important to make the robot's behavior adaptive and to inform the user of their contribution to the training [1]. For example, in rehabilitative training for neurological disorders (e.g., after stroke or spinal cord injury), the patient's motor performance can be measured and evaluated by the muscular effort estimation.
There are two widely used methods for quantifying the user's muscular effort in a rehabilitation robot: by measuring electromyography (EMG) using surface electrodes attached to the user's skin; and by estimating muscular torque based on inverse-dynamics analysis. Although measuring EMG signals has advantages in terms of detecting user intentions with accurate timing, it has practical limitations. For example, the attachment of the electrodes is time-consuming, and complex signal processing is required [3,4]. This paper considers inverse-dynamics-based muscular torque estimation for practical use.
The muscular torque of the human user can be estimated by measuring the applied external torque at each joint of the exoskeletal robot, and by removing the inertial, Coriolis, and gravitational torques of the user's limb (referred to as "passive torque" throughout the paper to distinguish it from the torque generated by muscle). Computation of the passive torque requires accurate estimates of anthropometric and inertial characteristics of the limb segment, such as mass, center of mass location, and moment of inertia (often referred to as body segment inertial parameters; BSIPs). In inverse-dynamics analyses of human movement, BSIPs are typically estimated from anthropometric models [1,[5][6][7]. Although such anthropometric data provide simple solutions for the researcher, they may not match that of the actual user because they cannot provide comprehensive solutions for variations in the gender, race, age, and body type of users [8,9].
In this paper, we present a method for estimating the user's muscular torque using joint torque sensors and its implementation in an actual exoskeletal robot. In particular, we focus on the identification of user-specific inertial parameters rather than using typical anthropometric models. This approach is important because the isolation of active muscular effort from joint torque measurements critically relies on the accuracy of the dynamic model of the user's limb. The wheelchair-integrated lower limb exoskeleton robot EXOwheel was used as a test bed, and 10 subjects participated in the experiments. This paper provides a mathematical formulation of the joint torque resulting from the user in the exoskeleton and experimental procedures for identifying user-specific parameters. The performance of the proposed method was verified by experiments on body-weight-supported gait training.

Mathematical Formulation
Experiments were performed with an EXOwheel robot [10], shown in Figure 1a. The EXOwheel is designed to support exercise and rehabilitative training in the daily lives of individuals with disabilities. The exoskeleton provides assistive joint torques via electric motors in the sagittal plane at the hip and knee joints. The user is connected to the exoskeleton through three attachment points: the thigh, shank, and foot. The length of the thigh and shank in the exoskeleton can be manually adjusted to fit the user's leg length. Figure 1b shows a schematic diagram of an exoskeleton joint. The exoskeleton is equipped at each joint with an encoder for the motor's position and a sensor for the joint torque which is located between the motor and link ("torque sensor" in Figure 1b). The exoskeleton's hip and knee joints have the same configuration. The technical specifications of the torque sensor are: sensing range: ±120 Nm, resolution: 0.015 Nm, non-linearity: 0.03% full scale and repeatability: 0.02% full scale. For a human lower extremity wearing an exoskeletal robot, we considered a two-segmental model in the sagittal plane as illustrated in Figure 2. Several assumptions are made to simplify the calculations: (1) The human leg consists of rigid segments, and each segment is connected with a fixed hinge joint.
(2) The human leg is rigidly linked to the exoskeletal robot, and both systems have the same kinematics.
(3) The model only considers motion in the sagittal plane. (4) The shank and foot are treated as one rigid segment (i.e., lower leg). The model consists of two rigid segments (thigh and lower leg) and two pin joints (hip and knee). Each segment of the model is defined by five parameters: length (L), mass (m), position of center of mass in the directions parallel and perpendicular to the link (a and b, respectively), and moment of inertia (Iz). θ is the joint angle, and positive values of θ represent counter-clockwise rotation. Subscript 1 refers to variables of the thigh segment and hip joint, whereas subscript 2 refers to the lower leg segment and knee joint.
The equation of motion for the human lower limb model is expressed as [11]: ) and "active" torque (muscular torque; i.e., τM). The external torque, τEXT, is the applied torque from the environment.
In the EXOwheel, the external torque can be expressed by subtracting the torque needed to move the exoskeleton from the torque generated by the robot's actuator: Indexes H, R, and HR are omitted in Equation (4) where g is gravitational acceleration, θ12 = θ1 + θ2, JHR = JH + JR, XHR = XH + XR , and YHR = YH + YR. From the location of the torque sensor (between the actuator and link; see Figure 1b), the measured torque can be described as: where τS is the torque measured by the exoskeleton's torque sensor. The active muscular torque of the human user can be estimated from the measured torque τS as follows: where ( ) and "hats" are placed on the parameters to denote the estimated values. As shown in Equations (9) and (10), it is important to extract the dynamic effects of the human body and exoskeleton accurately from the measured torque for estimating user's active muscular torque. The user's BISP values in Equation (4) (mH,i, aH,i, and IzH,i) are typically estimated from the literature. Table 1 shows three widely used BSIP estimation models: two models derived from cadaver studies [12,13] and one model derived from the gamma-ray scanning of living subjects [14]. In the table, each segment's mass (M) is described as a ratio of the total weight, whereas the center of mass location (CM) and radius of gyration (RG) are described as segment length ratios: where W is the subject's weight. The BSIPs estimated from Table 1 are not identical to those of actual users; in this study, the BISPs are not only estimated from the literature but also measured and validated with actual human subjects in the following sections. Passive elastic torque, P(θ), is a torque generated by the mechanisms of the joint surface, ligaments, and connective tissue. This torque is weak relative to the gravitational and inertial torques, but it becomes significant at the end range of motion [15]. A model for estimating passive elastic torque is based on Riener's double-exponential equations [16].
The unit of the angle is radians and the unit of torque is Nm.

Experimental Procedure
To improve the accuracy of the BSIPs, we identified certain parameters from the series of experiments. All experiments were performed on 10 healthy subjects ( Table 2). The experimental procedure was approved by the Ethics Committee of Sogang University (approval number: Sogang-IRB-2014-08), and written informed consent was obtained from all participants. Table 2. Subject characteristics.  Figure 3 shows the experimental setup for parameter identification. The subjects were placed in the EXOwheel, and the exoskeleton was connected to their legs. Two experimental configurations were used for parameter identification.

Subjects Age (years) Height (m) Weight (kg) BMI a (kg/m 2 ) L 1 (m) L 2 (m)
(1) Identification for the lower leg ( Figure 3a): The subject sat on the wheelchair seat, which allowed the lower leg to swing. During the experiment, the hip was fixed at 90° flexion, and the knee was moved with a sinusoidal trajectory. (2) Identification for the thigh (Figure 3b): The subject was in the standing position, which allowed the entire leg to swing. During the experiment, the hip was moved with a sinusoidal trajectory, and the knee was fixed at 90° flexion.
In each experiment, the hip and knee angles were imposed by the exoskeleton, and the subject was asked to fully relax his leg to allow the leg to move passively against the torque imposed by the exoskeleton. Each experiment was carried out five times for 40 s. A conventional proportional derivative (PD) controller was applied at each joint of the exoskeleton to track the desired position. The exciting trajectories were selected as a sum of harmonic sine and cosine functions within the frequency range during the gait training (0.3-1.5 Hz). After the experiments on the subjects were completed, experiments on the exoskeleton only (i.e., exoskeleton not worn by a human subject) were performed in the same configuration to separate the parameters of the exoskeleton and human limb.

Data Analysis and Parameter Identification
The measured joint angles and torques were recorded during the experiment. The angular velocity and acceleration were obtained off-line from the joint angle using a zero-phase digital low-pass filter (Butterworth, second-order, cut-off frequency of 20 Hz). To detect the user's actual muscular effort during the experiment, EMG data were measured by surface electrodes attached to five major muscles of the right leg in all subjects ( Figure 3b): gluteus maximus (GMAX), vastus medialis (VM), rectus femoris (RF), semitendinosus (ST), and gastrocnemius (GAS). The raw EMG signals collected from the electrodes were amplified by a Bagnoli TM 8-channel system (Delsys Inc.) to a gain of 1000 and then band-pass filtered (Butterworth, fourth-order, cut-off frequency of 20-450 Hz). All signals were collected at a sample rate of 1 kHz.
Based on the measured EMG data, muscle activation, Ach, is classified into two states: "activated" and "inactivated". The single-threshold method is used as the classification rule: where χch(k) is the amplified, band-pass filtered, and full-wave rectified EMG signal of each channel (e.g., GMAX, VL, RF, ST, and GAS) at a discrete time instant k. Zch is the threshold value, which was set to mean plus three standard deviations (Mean + 3SD) of χch when the muscle is relaxed.
The parameters J, X, and Y of the exoskeleton and human limb were estimated using the off-line least-squares method. For the human limb parameters, data that satisfy passive conditions (i.e., Ach = 0 for all EMG channels) were used. The estimation procedures are shown schematically in Figure 4. The procedures were split into two parts: identification of the exoskeleton parameters (Step 1: R2 segment; Step 2: R1 segment) and identification of the human limb parameters with the exoskeleton's parameters assumed to be known (Step 3: H2 segment; Step 4: H1 segment).  7), the relation between the measured signals, input u and output y, are expressed as: Step Step 2 (R1 segment; Step 3 (H2 segment; Step 4 (H1 segment; Lower leg of the HM+EXO Whole leg of the HM+EXO where * j u ∈ℜ is the modified input signal, 3 j ∈ℜ w is the regression vector as a function of geometry and , , y y y  , and 3 j φ ∈ℜ is the vector of the segment parameters. The components of the vectors uj*, wj, and ϕj in Equation (18) can be found in Appendix.
When N measurements are used, Equation (18) can be represented as: A least-squares method [17] can be used to estimate the parameters in Equation (23). Figure 5 shows the input and output data recorded in the identification experiments. The exciting trajectories were selected as the sum of harmonic sine and cosine functions considering the joint angle range and frequency in typical gait training in the EXOwheel system. As shown in the second row of Figure 5, the range of the joint angle was 5 to 40 deg for the hip joint and −5 to −70 deg for the knee joint. The range of angular velocity was ±120 deg/s for the hip joint and ±165 deg/s for the knee joint. The kinematic data of the exoskeleton and the combined human-exoskeleton system were identical for each joint (i.e., y1 = y3, y2 = y4). Table 3 and 4 represent the results of the parameter identification in the exoskeletal robot and human body, respectively. In Table 3, the identified parameters of exoskeletal robot vary since the link length of the thigh and lower leg in exoskeleton are adjusted to fit the subject's leg length. The parameters of human body, however, vary significantly though the segment length and weight of individuals are similar (Table 4). For example, S3 and S5 have highly similar weights (W = 74.3 kg for S3; W = 75.1 kg for S5) and segment lengths (L1 = 0.43, L2 = 0.55 for both subjects), but the identified human limb parameters show significant discrepancy between the two subjects (e.g., XH1 = 4.40, XH2 = 1.50 for S3; XH1 = 3.56, XH2 = 1.27 for S5). In particular, parameter YH, which represents the distribution of the mass in the direction perpendicular to the link, deviated significantly among all of the subjects regardless of the weight and segment length.

Identification Results
In Table 4, the identified parameters of the human body are compared from that of the two anthropometric models: the De Leva model and the Dempster model. The observed results indicate that the anthropometric data are not sufficient to estimate user-specific parameters. There are two possible reasons for the discrepancies between the parameters identified from the experiments and the anthropometric models. First, normalizing the parameters with respect to weight and segment length alone might cause poor estimation because the identified parameters are sensitive to the body shape of the subject (i.e., spatial distribution of the segment's mass). Second, the assumptions applied in the anthropometric-model-based estimation (e.g., the radius of gyration of the each segment is constant and the anthropometric segment boundaries coincide with that of the exoskeleton) might increase the discrepancies. Step 2 (hip joint data; exoskeleton); (c) Step 3 (knee joint data; combined human-exoskeleton system); (d) Step 4 (hip joint data; combined human-exoskeleton system).

Validation of the Identification Results
The quality of the identified parameters was evaluated in body-weight supported gait training. During the experiment, the exoskeleton was controlled to move along a given trajectory while the subject was asked to fully relax his leg muscles. The trajectories of the hip and knee joint angles were determined from the physiological gait pattern [18], and the gait speed was selected to be 2 km/h (typical range of gait speed in rehabilitation training; approximate stride period of 2 s). During the training, the subject's body weight was fully supported by the EXOwheel's electrical lifter for no ground contact. The experiment was carried out for 10 stride cycles on each subject.
As the subject's leg moves passively, the measured joint torque only contains motion-dependent "passive" torque (i.e., in Equation (10)). In the relaxed muscle condition, the estimated torque should coincide with the measured torque if the identified parameters are accurate. Figure 6 shows the results of the validation experiment for subject 1, which include measured kinematics (Figure 6a), measured EMG data (Figure 6b), and a comparison between measured and estimated torque (Figure 6c). Figure 6b shows that all of the muscles were in the "inactivate" state, which indicates the subject was fully relaxed throughout the experiment. As shown in Figure 6c, for both the hip and knee joints, there was good agreement between the estimated and measured torques, validating the accuracy of the identified parameters. The performance of the identified model was quantitatively evaluated by inspecting the root-mean-square (RMS) error of the torque differences between the estimated and measured torque: In calculating the RMS error, the estimated torque was calculated using the three different parameter sets: experimentally identified parameters and the parameters estimated from the De Leva and Dempster models. According to the results shown in Table 5 for all subjects, the RMS errors range 1.5-2.2 Nm for the hip joint and 0.7-1.2 Nm for the knee joint when torques are computed using the experimentally identified parameters. In contrast, large RMS errors are observed when using the parameters estimated from the De Leva model (3.3-6.5 Nm for the hip joint; 2.0-3.9 Nm for the knee joint) and the Dempster model (2.8-7.3 Nm for the hip joint; 2.3-4.5 Nm for the knee joint). These results demonstrate that the experimentally identified parameter set is acceptable and significantly increases the accuracy of the inverse-dynamics analysis relative to the parameter set estimated from typical anthropometric models. 18 19

Experimental Validation of Active Muscular Torque Estimation
As indicated in Equation (9), the reliability of the estimated active muscular torque is guaranteed if the passive torque is accurately removed from the measured torque . Although the accuracy of the estimated passive torque was analyzed in the previous section, it remains unclear whether the amplitude of the estimation error is sufficiently small to recognize a user's muscular effort. To validate the estimation of active muscular torque in Equation (9), we investigated the correlation between the estimated muscular torque and EMG data. Prior to the experiment on gait training, the relationship between the EMG data and joint torque was established with isometric contractions. Figure 7 shows the experimental setup for the isometric calibration procedures. The subjects were placed in the EXOwheel, and the exoskeleton was connected to their legs. The subject stood upright during the isometric hip flexion and extension (Figure 7a), whereas the subject sat on the wheelchair seat during isometric knee flexion and extension (Figure 7b). The subjects were asked to perform five to eight isometric flexion and extensions for each joint. The raw rectified EMG data (i.e., χch) were low-pass filtered (Butterworth, fourth-order, 2 Hz) to obtain linear envelopes. The resulting EMG envelope of each muscle was normalized to the respective muscle's maximal voluntary contraction (MVC). The subjects were instructed to perform MVCs against manual isometric resistance with specific test positions for each muscle of interest [19]. Three trials were performed with 2 min of resting time between trials; the average of the three highest EMG peaks was used as the MVC value for normalization.
To calculate muscular torque from the EMG signals, the model reported by Olney and Winter [4] was used in this study. The muscular torques at the hip and knee joint during isometric contraction were calculated from EMG data as follows: where ̂ , and ̂ , are, respectively, the muscular torques at the hip and knee joints calculated using EMG under isometric contraction, is the normalized linear envelope of EMG, and is a constant relating the amplitude of the to the joint torque. is determined by a least-squares curve-fitting procedure as follows: where α1 = {α RF,1 , αGMAX,1, αST,1}, α2 = {αVM,2, αRF,2, αST,2, αGAS,2}, and ̂ , ( ) = , ( ). Note that the measured torque at the joint torque sensor only contains the user's muscular torque in the experimental setup shown in Figure 7, which indicates that ̂ , is reliable during the isometric calibration procedures regardless of the accuracy of the inverse dynamics model. The optimized values of α computed by the least-squares method are (Mean ± SD): α1 = {1.88 ± 0.76, 1.08 ± 0.42, 1.70 ± 0.83}, α2 = {2.53 ± 0.77, 0.89 ± 0.37, 2.91 ± 1.12, 1.51 ± 0.35}. The results of the isometric calibration for the hip and knee joint of subject 1 are presented in Figure 8a,b, respectively. In each figure, a time plots of the EMG signals recorded from GMAX, VM, RF, ST, GAS, and the corresponding muscular torques are presented. In the EMG graphs of Figure 8, thin grey lines represent the raw-rectified EMG signals, and thick blue lines represent EMG linear envelopes. As shown in the torque graphs of Figure 8, good agreement was observed between the muscular torque measured from the torque sensor (solid red line in the figure) and the torque calculated using EMG (dashed black line in the figure). The R 2 values and normalized RMS errors between and , for all 10 subjects are presented in Table 6. The average R 2 values were high (ranging between 0.951 and 0.976), and the normalized RMS errors were low (ranged between 5.03% and 8.12%) for all isometric contractions; this indicates that the EMG to torque processing model produced accurate estimates of the joint torque.  After the calibration was completed, the experiment was conducted on gait training. The experimental setup was the same as that used in the model-validation experiment, i.e., body weight supported gait training with pre-defined gait pattern (gait speed: 2 km/h). In this case, the subjects were asked to generate more force than required to follow the pre-defined gait trajectory (i.e., gait with exaggerated flexion/extension of hip and knee joint). The experiment was performed for 10 stride cycles on each subject.
To calculate muscular torque from the EMG signals during dynamic contractions, the model reported by Olney and Winter [4] was used in this study. The muscular torques at the hip and knee joint during dynamic contraction were calculated as follows: where ̂ , and ̂ , are, respectively, the muscular torques at the hip and knee joints calculated using EMG under dynamic contraction. β is a constant, in deg −1 , depending on the difference between the joint angle and the angle at which the isometric calibration trials were conducted ( = 0, = − 90 deg). γ is a constant, in (deg/s) −1 , accounting for the variations in angular velocities. The optimized values of β and γ were determined by a least-squares curve-fitting procedure. The optimization problem was formulated as follows: where ̂ , ( ) = , ( ) − ̂ , ( ).
The optimized values of β and γ computed by the least-squares method are presented in Table 7. The values of α were derived from the isometric calibration procedure. Figure 9 shows the results of the gait experiment for subject 1, which include the measured EMG, measured joint angle, measured joint torque, and estimated muscular torque. The subject was fully relaxed up to 5 s, and then generated muscular force with exaggerated flexion/extension of the hip and knee joint. During the periods when the subject was passive (i.e., 0-5 s), the measured torque exhibited a repetitive pattern (Figure 9b). This behavior can be explained by the fact that most of the measured torque was induced by the passive torque, such as the inertial, Coriolis/centrifugal, and gravitational torque of the subject's limb. It can be observed that the active muscular torque estimated using inverse dynamics (solid red lines in Figure 9d) was close to zero over the period 0-5 s, indicating the passive torque was accurately removed from the measured torque.   During the periods when the subject was active (i.e., after 5 s), the pattern observed for the measured torque was not repetitive because it featured the active muscular torque in addition to the passive torque (Figure 9b). Figure 9d shows a comparison of the muscular torque estimated using inverse dynamics, , and the muscular torque calculated using EMG, , . As expected, good agreement was observed between (solid red line in the figure) and , (dashed black line in the figure). The R 2 values and normalized RMS errors between and , for all 10 subjects are presented in Table 7. As seen in Figure 9d and Table 7, the estimated muscular torque and EMG data are closer when using the experimentally identified parameter set than using the parameter set estimated from typical anthropometric models. The average R 2 is 0.935 for the hip joint and 0.924 for the knee joint when using the experimentally identified parameter set and it shows that the estimated muscular torque during gait experiment is highly correlated with the EMG data. This indicates that most of the estimated muscular torque was induced by the subject's muscle activation (i.e., neural command) and the proposed method can be effectively used to estimate the user's muscular effort.

Conclusions
This study presents a method for estimating users' active muscular torque measured by sensor systems typically used in exoskeletal rehabilitation robots (i.e., encoder and torque sensors). The key step in this method was to accurately identify the inertial parameters of the user's limb. Parameters experimentally derived from actual users are significantly more accurate than those obtained by widely used anthropometric models. Experimental results on gait training validate the proposed muscular torque estimation method.
The method proposed in this study is valid only when the user moves in the air with no ground contact because ground reaction force (GRF) cannot be measured with the current EXOWheel system. This method can also be extended to over-ground walking with measurement of the GRF. When the foot is in contact with the ground, however, the accuracy of the GRF measurements will likely become the critical factor in estimating muscular torque rather than the accuracy of the user's limb inertial parameters. A future study will examine the performance of this method during ground contact.
The advantages of the proposed method can be summarized as follows: (1) no additional sensors for measuring bio-signals are necessary because this method only uses common rehabilitation robot sensors and (2) the method provides user's muscular effort in terms of joint torque, which is adequate as a feedback signal for the controller in joint space; this method can be used for the shared control or adaptive control of exoskeletons.