Mathematical Model That Predicts Lower Leg Motion In Response To Electrical Stimulation
1 Department of Mechanical
Engineering,
2 Mechanical and Aeronautical
Engineering,
3 301 McKinly
Lab, Department of Physical Therapy,
Email: perumal@me.udel.edu
Abstract
Direct electrical activation of skeletal muscles of patients with upper motor neuron lesions can restore functional movements, such as standing or walking. Because responses to electrical stimulation are highly nonlinear and time varying, accurate control of muscles to produce functional movements is very difficult. Accurate and predictive mathematical models can facilitate the design of stimulation patterns and control strategies that will produce the desired force and motion. The purpose of this study was to validate our previously developed model during nonisometric contractions when the leg was allowed to move freely in response to electrical stimulation. Our results showed that the model could accurately predict the angular position and velocity of the lower leg when the muscle was stimulated with wide range of clinically relevant stimulation patterns and with different loads placed around the ankle joint.
1 Introduction
Functional Electrical Stimulation (FES) is
the coordinated electrical excitation of paralysed or weak muscles in patients
with upper motor neuron lesions to produce functional movements such as
sit-to-stand or walking. Despite many technical advances,
Although a number of mathematical models
capture the behaviour of muscles [1][2], to date none
of the models have succeeded in
2
Methods
The model has three differential equations:
, (1)
where
for i
= 1, and
for i >1;
(2)
where
, and
;
and for general nonisometric movements, hence
, (3a)
and for isovelocity
movements
, (3b)
. (3b)
Equation (1) represents the dynamics of the
rate-limiting step leading to the formation of the Ca2+-troponin
complex, CN. t (ms) is the
time since the beginning of stimulation, ti
(ms) is the time when the ith pulse is delivered, τc
is the time constant controlling the rise and decay of CN,
and R0 characterizes the magnitude of enhancement in CN
when there are two closely spaced pulses.
Equation (2) represents the development of the force due to stimulation,
F (N) . A40 (N/ms) is the scaling factor for
force at 40° of knee flexion, a (N/ms-deg2)
and b (N/ms-deg) are scaling factors to account for force at each knee
flexion angle, q (deg) is the knee flexion angle,
(deg/s) is the angular
velocity of the shank, V1
(N/deg2) is a constant for each subject to be determined at -200°/s, Km is the
sensitivity of strongly bound cross-bridges to CN, t1
(ms) is the time constant of force decline in the absence of strongly bound
cross-bridges, and t2 (ms) is the time constant of
force decline when there are strongly bound cross-bridges. Equation (3) is the equation of motion for the shank in response
to electrical stimulation [4]. Two separate equations (3a and 3b) were used to
represent the motion of the shank during general nonisometric and isovelocity movements.
(deg/s2) is the acceleration of the shank, L
is the distance from the knee joint centre of rotation to the location of F
above the ankle joint, I (kg-m2) is the mass moment of
inertia of the lower leg and the external load wrapped above the ankle joint, M
(N) is the resistance to knee extension that includes the visco-elastic
resistance of the knee joint and the weight of the lower leg, FKC (N) is the experimental
force measured by the KinCom dynamometer during isometric and isovelocity
experiments, FLOAD (N) is the external load wrapped above the
ankle joint, and
(deg) is the resting angle.
Six healthy subjects were recruited for this study. Isometric, constant shortening velocity (isovelocity) force, and varying velocity (nonisometric) angle data were collected of the human quadriceps femoris muscle in response to electrical stimulation of the human quadriceps femoris muscle. For the isometric and isovelocity contractions, force and angle data were collected using the KinCom. In contrast, for the nonisometric contractions the dynamometer arm was replaced by a low frictional resistance custom built arm that measured knee joint angle.

Figure
1: Schematic representation of the three
stimulation patterns tested. Bottom train (CFT50) is a constant-frequency train
with all interpulse intervals equal to 50 ms; middle train (VFT50) is a
variable-frequency train with an initial doublet of 5 ms and remaining pulses
equally spaced by 50 ms; and top train (DFT50) is a doublet-frequency train
with 5-ms doublets separated by interdoublet interval
of 50 ms. The train’s name is based on the duration of the longest interpulse
interval within that train. Each train has a maximum duration of 1 sec.
For each subject all data were collected in a single testing session. Each testing session consisted of three parts. First, subjects were tested isometrically at knee flexion angles of 15°, 40°, 65°, and 90°. Then, subjects were tested at an isovelocity speed of -200°/s (all shortening velocities are assigned negative values in this study). Finally, nonisometric tests were performed with weights of 0, 4.54 and 9.08 kg wrapped just above the ankle joint. A rest period of 5 min was provided between testing each angle, velocity, and load to avoid fatigue. The stimulation intensity was set to activate ~40% of the muscle. For isometric and isovelocity testing the muscle was stimulated with one-second long VFT20 and VFT80 (Figure. 1), while under nonisometric conditions the muscle was stimulated with one-second long CFTs, VFTs, and DFTs of interpulse intervals (IPIs) ranging from 10ms to 100ms. The experimental procedures for isometric and isovelocity studies are detailed in [4]. For nonisometric testing, stimulations were initiated when the subject’s leg was in a resting position and were truncated when either the leg had gone through a 45° excursion or all pulses within the one-second long train were delivered. The above trains were delivered every 10 sec to avoid fatigue.
First, parameters a, b, A40,
τ1, τ2, and KM
were identified under isometric conditions, then parameters V1
and M were identified at a constant shortening velocity of -200°/s, and finally parameters
, M, and
were identified under
nonisometric conditions. Under isometric and isovelocity conditions parameters
were identified by fitting the force response to the VFT20-VFT80 train
combination. For nonisometric contractions parameters were identified by
fitting the velocity response to VFT20 train at each load condition.

The accuracy of the predictions of the response to the
nonisometric condition were tested by comparing measured and modeled knee flexion angles and velocities in response to
the CFTs, VFTs, and DFTs for all the six subjects tested. We used the RMS angle
error measure to see how well the model predicted leg position in response the
stimulation at different loads. The RMS
angle error was defined as

(4)
where qmod,i is the modeled knee flexion angle, qmeas,i is the measured knee flexion, i is the increment time interval of 5 ms, and N is the total number of data points from the beginning of the stimulation to the measured knee flexion angle at the end of extension.
3 Results
Subjective evaluations showed that the model was able to predict the knee flexion angle and the angular velocity in response to different stimulation IPIs and patterns (Figure 2). The averaged angle errors were generally lower at the two higher loads and were less than 8° for most of the stimulation IPIs tested (Figure 3). In addition, the averaged RMS angle errors were generally higher at longer IPIs for each of the two loads tested.
4 Discussion and Conclusions
The model developed and tested here
accurately predicted the angular position and velocity of the shank when the
muscle was stimulated with wide range of stimulation frequencies and patterns,
and loads. Compared to other models [1][2], the
current model requires very few parameters to describe the nonisometric motion
of the limb in response to electrical stimulation. This makes the current model
a suitable candidate for adaptive control algorithms that can track parameter
variations online and provide a better control of motion during
Figure 2: Predicted
and measured knee flexion angles and shortening angular velocities during
extension from a typical subject at load of 9.08 kg for four of the 16
stimulation trains.
Figure 3: Bar graphs showing the averaged angle errors (+
standard error) for six subjects in response to the 16 stimulation trains at 0
and 9.08 kg different loads.
[1] Riener R, Quintern J, Schmidt G. Biomechanical model of the human knee evaluated by neuromuscular stimulation. J. Biomech. 29, 1157-1167, 1996.
[2]
Dorgan SJ, O’Malley MJ. A
mathematical model for skeletal muscle activated by N-let pulse trains. IEEE Trans. Rehab.
[3]
Perumal R,
[4]
Perumal R, Wexler AS, Ding J, et al. Mathematical Model For Electrically Elicited Isovelocity
Human Muscle Contractions —I: Development Of The Model. ASME J.Biomech.
Acknowledgements
This study was supported by the National Institutes Health Grant HD 36797.