State Space Methods and Examples for Computational Models of Human Movement

,

Based on the state space formulations of NE, we introduce a six-sub module, 18-link three-dimensional skeletal model.Each sub module is made of a three-link system.The six sub modules represent the two lower extremities (Pope, Crowninshield, Miller, & Johnson, 1976;Ren, Jones, & Howard, 2006), the torso (Minetti, Cisotti, & Mian, 2011;Winter & Stark, 1988), the two arm-lower arm-hand limbs (Charles & Hogan, 2010;Rankin & Neptune, 2012;Reich & Daunicht, 2000;Williams, Schmidt, Disselhorse-Klug, 2005).and the head and neck three segment.This 18-link module can also be used to model each hand and each foot.The hand module has a three-link palm to which five three-link fingers are attached.Therefore, a 90-link human skeletal model becomes available.
For ease of presentation, we consider point connections and holonomic constraints at all joints.All three-link sub modules have three translation and nine rotation degrees of freedom.Each sub module has its own internal musclelike actuators that span joints and may span a segment or a sub module (Ajay & Pandy, 2007;Buchner, Hines, & Hemami, 1988;Winter & Stark, 2012).
The different sub modules are coupled by cross-module forces of connection and muscle-like actuators (Rankin & Neptune, 2012).The system has 72 degrees of freedom and five internal cross-joint forces.Two ground contact forces support the system.The degrees of freedom can be reduced, by projection, to 54 degrees of rotational freedom plus three degrees of translational freedom.
In Section Two, a single body is formulated in four different state spaces.Constraint forces are addressed in section Three.Multi rigid body systems and three-body sub modules are presented in Section Four.The six-sub module system and a particular actuation system are presented in Section Five.Computational experiments are briefly presented in Section Six.Discussions and conclusions follow in Section Seven.

Single Rigid Body
The state space formulation of the single rigid body is an important step in developing multi-rigid-body dynamics that do not use tensors (Pellionisz & Llinas, 1979, 1980) and are amenable to dimensional augmentation and dimensionality reduction by projection.All required matrices are synthesized from 3 × 3 matrices.The states are the three Euler or Bryant angles Θ and the angular velocity Ω (as expressed in BCS) for the NE method, They are Θ and Θ in K and LA formulations, Θ and the angular momenta P in H method.
For systems with a natural skeleton, the NE formulation of the dynamics is easier to implement.For combinations of natural and man-made systems, such as an athlete with skates or bicycles in the air, Kane's method seems more appropriate.For deform-able and flexible bodies, the L and H methods may be more suitable.For understanding the energy conversion machinery of a single muscle, all four methodologies may be necessary.
Transformations between these state space formulations are desirable for the development of computational software development and use.All are compact and computer-adapted formulations that facilitate human interaction with computers (Marshall, Jensen, & Wood, 1985;Mansouri & Reinbolt, 2012).The LA and H methods provide energy functions that come handy in defining Lyapunov functions (Antsaklis & Michel, 1997;Winter, Prince, Powell, & Zabjik, 1996) and stability (Honarvar & Nakashima, 2012).

Newton-Euler Formulation
The NE method (Hemami, 2002) makes use of three matrices: A, B and C. Matrix A describes the transformation of vectors from the body coordinate system (BCS) to the inertial coordinate system (ICS).Matrix B relates Θ and Ω. Matrix C describes the constraints of interaction with, connection to and contact with other objects and with the environment.For examples of NE formulations, see (Amirouche, 1992;Hemami & Wyman, 2012) and Appendix A.

Kane's Formulation
The objective is to derive Kane's equations (Kane & Levinson, 1985;O Reilly, 2008) for a rigid body, connected to the origin of the inertial coordinate system (ICS) at its center of gravity.Let the three-vector T be the generalized torques acting on the system as defined by Kane.The key to the transformation is to relate the generalized holonomic torque vector T appearing in Kane's equations to the torque vector N applied to the BCS axes of the rigid body.The power p delivered by T is the inner product The same power in the NE method is given by p From the above two equations and the equality: Let us define the transpose of B as B t , the inverse B as B i and the inverse of the transpose of B as B it .Therefore, Thus, by multiplying the NE equations by B it , one obtains Kane's equations.Conversely, by multiplying Kane's equations by B t , one obtains NE equations.
There is one computational difficulty with matrix B. It is singular at this problem can be computationally addressed by using a finite large number for the inverse of the cosine function at the critical value or change the coordinate system in that vicinity.

Lagrange's Formulation
The derivation of the Lagrange'e equations of motion (McCuskey, 1959;O Reilly, 2008) is carried out in (Hemami, 2002) starting with the Lagrangian (Lam, 1998;Wittenburg, 1989).We provide an alternative derivation here which is more general and may also be easily applicable to multi body systems.The equation for Θ is differentiated with respect to time.
By multiplying both side of Eq. 5 by B it J(B i ), and recalling one obtains the Lagrangian equations: The above development points to an indirect method of deriving the Lagrangian angular acceleration variables: (( Θ) from available NE (or K) state space equations.One adds the angular accelerations as outputs to the state space formulation.These outputs are functions of the states and inputs, i.e., the structure is that of a Mealy machine as opposed to Moore machines where the output is a function of the state alone.The system has a feed through component (Bay, 1999).

Hamilton Formulation
In the Hamiltonian formulation, the state variables are defined as Θ and P, where It follows that H = 0.5P t J L P.
The latter equation appears to be compact and simple, but it involves derivation of the tensor ∂B/∂Θ or the simpler form ∂P T B/∂Θ.
The relation between P and Ω is

Forces of Constraint
Briefly, forces of constraint either appear in the equations of motion (Pope, Crowninshield, Miller, & Johnson, 1976;Hemami & Wyman, 1979;Popovich, Goswami, & Herr, 2005) or are not part of the equations of motion (Kane & Levinson, 1985).According to Kane's method, additional degrees of freedom are added such that the forces are computable.
When the forces appear in the dynamics, two methods exist to derive the forces.The physical constraints may be differentiated and combined with the system equations.The latter system can be simultaneously solved for the accelerations and the forces (Smith, 1973) Alternatively, the accelerations can be eliminated and the forces derived as functions of the state of the system and the inputs (Hemami & Wyman, 1979).
For the single rigid body system above, the procedure is as follows: According to Kane's method, one adds three translation degrees of freedom X to the hinge.Let X and Θ be the six generalized coordinates.Let Γ be the vector force of constraint at the center of gravity (COG).Let T be the sum of all torques acting on the body -some being independent couples and some, due to external forces, represented here by their resultant force acting at COG: Similarly let Kane's equations of motion, for rotation, are F * r + F r = 0.For the translation, we have Now, the constraints are imposed, i.e., Ẍ = 0.
From which it follows that Γ is equal to the negative of all applied forces.

Multi-segment Rigid Body Systems
With the above development, the Lagrangian equations of motion for a multi rigid body system can be derived.However, use of tensors for derivation and partial derivations may be necessary as the dimensions increase.
An alternative is to use K or NE dynamics in state space form, and add Θ as Mealy outputs of the latter as stated earlier.
Co-activation amounts to position and velocity feedback (Hemami, Barin, & Pai, 2006).Lyapunov stability can be used to show that the system is stable in the regions of interest.We further assume the CNS produces the trajectories.

The Sub modules Dynamics
A skeletal model of a human or animal could be constructed from connecting a large number of rigid bodies (Amirouche, 1992;Wittenburg, 1989) While the large system could be constructed with individual rigid bodies, a better compromise would be to construct three rigid body sub modules first.The larger system would then be constructed from a number of such sub modules.The movements of the sub modules may be conceived, programmed and/or implemented by parts of the CNS.Therefore, the physical model of the skeleton may be in certain ways closer to the natural structures of the CNS.
There are three different sub module types in this system.The simplest one is that of the upper extremities and the neck and head.These three modules each have two internal forces of connection (constraint) and one external force.The two lower extremities have two internal and two external forces as shown in Figure 1.

The Torso Sub Module
The torso sub module, as shown in Figure 2, has five external forces and two internal forces.Initially, the BCS of the three links coincide with the ICS.The sequence of angles is roll, pitch, and yaw.The rotational space is considered first, and it is assumed that the torque vector N, expressed in the BCS operates on the system.The other two simpler sub modules (Hemami & Wyman, 2012;Hemami & Zheng, 2012) can be inferred from the torso sub module.
The 36-dimensional states of positions and angular and translational velocities are The velocity vector V is A projection operator P eliminates the internal forces and maintains the five external forces.

Inter Sub Module Dynamics
The six sub module system is illustrated in Figure 3: • Module 1 is the torso, composed of three segments.
• Module 2 is the right leg.
• Module 3 is the left leg • Module 4 is the right hand.
• Module 5 is the left hand • Module 6 is the three-segment neck and head.
The external forces are Γ 02 and Γ 0 3 at the support points with the ground.The external forces on the torso are Γ 12 and Γ 13 for the left foot, Γ 14 and Γ 15 between the torso and the hands, and finally Γ 16 between the torso and the neck Two alternative models of the dynamics are considered.In the first model, the system remains of dimension 72, and the formulation allows the joint forces and the forces of support to be available in the dynamics and continuously observable (Figure 4).The alternative is to eliminate the five joint forces and reduce the dimension to 57.Only the two support forces remain in the equations.For the 57-dimensional case, the torso sub module is of dimension 12 and the remaining five sub modules each are of dimension nine.Both support forces go to zero when the system is in the air.The system can also stand on one leg.For simplicity, we only consider the 57-dimensional case in this paper.Because the neck, the hip and the shoulder joints all have three degrees of freedom, 15 holonomic constraints allow a projection, similar to the one defined in the Appendix A, to be developed as sketched below.
The equations for the six sub modules are symbolically written: Initially all segments are vertical and all are aligned with the inertial coordinate system (ICS).The three segments in every module are numbered, respectively, 1, 2 and 3, starting from the bottom.Every module has its own system of inputs: N i where i = 1....6.
There are inter module inputs N 12 , N 13 , N 24 ,N 25 and N 26 The translational motion of the torso module will be kept in the system equations and represents the translation component of the system: The translation of the system could be represented by the translation of one point of the six-sub module system.We consider the movement of the center of gravity of the mid-torso segment for this purpose.
Suppose the dynamics of the six-sub module system are programmed by the same subroutine that provides the 12 equations of motion from the 36-dimensional equations of motion of the three-link.The 24 position and velocity variables in the equations of motion must be augmented to 36 equations as input to the subroutine.Therefore, a more accurate block diagram is the one presented in Figure 5 where the augmentation and reduction of dimension are marked and allow for inclusion of passive structures (Hemami & Hemami, 2014).
. The six-sub module system where the expansion and reduction of the dimensinality are also marked to accommodate passive structures.
Let the individual segments of each sub module be numbered by 1, 2 and 3 upwards.The variables in the system have three subscripts: sub module number, segment number and attitude angle (1 for roll, 2 for pitch and 3 for yaw.)With the subscripts ij, i refers to the sub module and j refers to the segment.The torso segments are labeled 11,12 and 13.The right foot, leg and thigh are labeled 21, 22 and 23.The left foot, leg and thigh are labeled 31, 32 and 33.The right upper arm, lower arm and hand are labeled 41, 42 and 43 -the arm is assumed to be stretched upward vertically.The left upper arm, lower arm and hand are labeled 51, 52 and 53.The neck and head segments are labeled 61, 62 and 63.For every one of the eighteen rigid bodies above, a 3 × 3 matrix A i j defines the transformation from of the rigid body's BCS to ICS.We briefly introduce the contact constraints that relate the individual sub module's translation to that of the torso.The holonomic connection constraints that define the projection operator, as is done in appendix A for the torso module, are described here for the whole module.Let the translation of the torso module be X m and let the remaining five translational vectors for the other five modules be X22, X32, X42 , X52 and X62.It follows that the translations can be eliminated by the following projection: In this projection, the translational degrees of freedom of the torso module are kept and the other 15 translation degrees of freedom as well as the forces between the segments are eliminated.Only the contact forces of the feet remain in the 57-dimensional dynamics.The first derivative of Eq.13 defines the projection operator for reducing the dimension of the system from 72 to 57.The second derivative helps carry out the reduction and derivation of the reduced order system equations.Let Ha be the latter projection operator.This operator is parsed into six horizontal strips: the first 12 rows pertain to module 1, and the next successive nine rows sequentially are related to modules 2, 3, 4, 5 and 6.These parsed strips are named, respectively, Ha1, Ha2, Ha3, Ha4, Ha5, Ha6.
The inner product of each of these strips with H, as defined in Eq. 12, defines the dynamics of the individual sub modules and their couplings.
To complete the model, an additional module is necessary to compute the forces of support while the contact with the ground is in effect for either individual foot or for both feet.This step is well understood (Hemami & Wyman, 1979, 1980).It is also to be noted that, up to this point ( namely, the development of the dynamic sub-modules) all operations are of the inner product type and summation.Only integration of the dynamics and computation of the support forces requires inversion of appropriate matrices.

Control and Support
A set of N moments of force act on the module.These inputs will be specified later.Each moment of force is brought about by a pair of muscle-like actuators (Ajay & Pandy, 2007;Brooks, 1986;Kandel, Schwartz, & Jessell, 1991).Three pairs of muscles (Ajay & Pandy, 2007) control each of the three degrees of freedom, respectively, at every joint.We assume the system is on tiptoe and the toe extensors and flexors that control the sagittal motion of the toes can be replicated in the frontal and medial planes such that there is one pair of actuators for every degree of freedom.Similarly we assume three pairs at the ankle and three pairs at the knee joint.The structures of the nine pairs of muscles may not be anatomically correct for humans, but they simplify the modeling and simulations here.
The model here can be made more anatomically correct by adding a three-segment foot and more representative musculature (Humphrey, 2009).For bi-articular (such as the gastrocnemius) and multi-articular (such as the rectus femoris) muscles, modeling difficulties arise partly because the function of many such intertwined structures are not well understood or articulated.Slings and caps can be added here.
The precise understanding of the function of the cap is necessary in order to formulate the action of arbitrary muscles in the dynamics of any musculo-skeletal system.Briefly, if one uses Lagrangian dynamics, the incremental work dW performed by an actuator on the system is given by F.dl where F and dl are, respectively, the muscle force and the incremental change in the length of the muscle.Let Θ be the vector of the totality of the angles of the skeletal elements.It follows that the effect of the force F can be entered in the dynamics as: (∂l/∂Θ)F.
Therefore, it is important to be able to specify the dependence of the length of the muscle on θ.To clarify this point, we consider a sagittal plane approximation to a sling where an actuator spans the toes, the bottom of the foot and around the heel.It is connected to the tendon of another actuator as shown in Figure 6.The incremental change of the length of the muscle is where ∆l 3 could potentially depend on several more angles of the system.
The first term in the above expression can summarize one function of the five flexors: flexor digitorum brevis, quadratus plantae, flexor hallucis brevis, flexor digiti minimi brevis and lumbricales (Warfel, 1978).When the toes are firmly in touch with the ground and the limb is unloaded, the action of lifting the heel is facilitated and can be fast due to the collective action of these five muscles.This means there will be a moment of force acting on the foot to plantar-flex it, and the heel lifts easily due to the simultaneous effort of these five flexors.
This action amounts to having available feedback and/or feed-forward by these five muscles for the purpose of lifting the heel when the lower limb is not loaded.These feed-backs amount to state feedback in the system theoretic sense: it is shown below that linear state feedback stabilizes the attitude of the rigid body system in physical space.

Computational Experiment
The objective is to carry the arm through a point to point movement, such as pointing (Soechting & Flanders, 1989) and directing music, and then carry the reverse movement to get the arm back to its original position.
The arm has nine degrees of freedom and the motion is observed imagined or intended in the ICS system.A neural network (Arbib, Erdi, & Szentagothai, 1998;Bullock, Grossberg, & Gyenther, 1996;Kandel, Schwartz, & Jessell, 1991) with two populations of neurons is designed (Houk, 1995;Ajemian & Bullock, 2000;Hemami & Dariush, 2012) to transform the ICS trajectory to the trajectories of four of the nine Bryant angles.All other Bryant angles and velocities are zero.The fingers position in ICS is described by Suppose this motion is, arbitrarily, carried out by pitch and yaw of the arm, yaw of the lower arm and yaw of the hand which is assumed to be the same as the yaw of the arm.This reduces the computational problem to three unknown angles that are the solution of the following three equations: Given a desired X p (t), expressed in ICS, which is received for example from visual observation, the task of the central controller is to solve for the three Bryant angles in Eq. 15 above as functions of time.We make use of the three neural population circuits described in Hemami and Dariush (2012) to construct the Bryant angles from X p (t). Suppose the Bryant angles, in a simpler form and for two seconds, are parameterized as: The arm motion was simulated to demonstrate the feasibility of the three-link sub modules and the control strategy of co activation (Hogan, 1984).All the physical parameters are listed in Appendix B. The four Bryant angle inputs and the actual four Bryant outputs of the system are shown in Figure 7.The desired trajectory of the tip of the fingers, and the actual trajectory of the fingers in the ICS system are shown in Figure 8.

Discussion and Conclusions
We have contrasted different state space formulations for a rigid body.These state space formulations can be extended to multiple rigid body systems where dimensionality reduction procedures are easy to implement by projection operators.With the addition of Mealy outputs, the K or NE simulation programs can be adapted to compute Lagrangian outputs.
We have developed a six sub module computational tool for studying the machinery, control and coordination of human movement.For simplicity, all joints in the modules and between modules were assumed to have three degrees of freedom.Anatomically accurate systems of muscles can be added to this model.One may vary the number of muscles, disable them, and try different synergies.One can compute internal and support forces to compare with available experimental measurements.One may add more modules to have a more accurate spinal column, add modular tail structures for modeling animals as well as hand, foot, digit and toe modules.For testing the six-sub modules system, a pair of agonist-antagonist actuators were added for each of the 54 rotational degrees of freedom.The system can be supplied with a variety of sensory modalities: vision, vestibular, and tactile for more comprehensive and instructive models.Parallel computations (Hebb, 1949;Adeli, 1992;Singer, 2000) can be implemented to make the system more physiologically realistic.
The current model has two points of contact, with the ground providing support.This means all the movements are on tiptoe and mass-less toes provide the support.The human foot is very complex.One of the many challenges is to identify, define and precisely describe a three-link three-dimensional module for each foot (Humphrey & Hemami, 2010).Specifically, one of the many problems in formulating the foot is to articulate the dynamics of the joint between the toes and the metatarsals, sometimes refereed to as a sling.Briefly stated, one of the simpler functions of the sling is to change the direction of the force of the muscle or several muscles.This is similar to the function of the kneecap for changing the direction of action of the rectus femoris' heads at the knee (Rasch & Burke, 1978;Warfel, 1978).The mechanism is needed for almost all bi-articular muscles that span a joint or segment.
Philosophically speaking, one may cite some advantages for this computational method.The ideas presented here may be, pedagogically, superior in teaching rigid-body dynamics.They yield to more efficient and labor-saving methods of simulation and system integration.Finally, many internal states can be added and are candidates to be traced in the CNS of natural systems, and, therefore, help to unravel some of the magnificence of the brain mysteries in signal processing.
segments.Let G j be the vector of gravity for body j.Let N j be the moment of input forces acting on body j.Let the internal constraint forces between the mid-torso and the lower and upper torso be, respectively, Γ 2 and Γ 3 .The Newton-Euler equations of motion for the three segment torso, as free bodies, are: The second body represents middle torso The equations of the upper torso with the three connections to the arms and the neck are With these definitions, the 18 Newton-Euler equations can be summarily written as where V is defined in Eq. 11, and where the coefficient of Γ is an 18 × 21 matrix H1 c and where N 1 , f 1 and G 1 are 18-vectors.
Let Z be the 18-dimensional vector of translational and rotational positions states, defined by Eq. 10: Let the center of gravity of the mid-torso segment be X m = X 2 and define vector W: The projection from Z to W Z = P2(W).
allows one to derive the equation of motion for the torso module.By this projection, six degrees of translational freedom as well as Γ 2 and Γ 3 are eliminated.The projection P is derived from the six holonomic connection constraint equations We differentiate Eq. 22 with respect to time to derive the projection operator.The result of the differentiation is as follows: For convenience, let The 18 × 12 matrix H1 a is Differentiating Eq. 22 once more, with respect to time, renders The 18-vector H1 b is

Figure 1 .
Figure1.The three-link sub-module with two internal and two external forces.The lengths of the three segments are, respectively, l 1 , l 2 and l 3 .

Figure 2 .
Figure 2. The three-link torso sub module with two internal and five external forces

Figure 3 .
Figure 3.The six-sub module biped with five internal and two external support forces.

Figure 4 .
Figure 4.The block diagram of the six-sub module system

Figure 6 .
Figure 6.An actuator spanning two joints and, at one end, connected to a tendon

Figure 7 .
Figure 7.The four Bryant input angles (arm pitch and yaw and lower arm and hand yaw angles) are shown on the left side.The same angles as outputs of the three-link system are shown on the right side.

Figure 8 .
Figure 8.The desired trajectory of the hand and the actual trajectory in the ICS are shown.