Biped Gait Analysis based on Forward Kinematics Modeling using Quaternions Algebra

Gait is the main locomotion way for human beings as an autonomous decision. Due to the increase in people with walking disabilities, the precision in gait analysis for purposes in clinical diagnosis, sports medicine or biomechanical research for the design of assistive technologies is of special relevance. The literature reports notable contributions in technological developments with diverse applications; and in some cases, algorithms for characterization and gait analysis; however, more studies related to gait kinematics are necessary, such as the solution proposed in this work. In this paper, we focus on studying the forward kinematics of the lower limbs in human gait, using in a novel way quaternions algebra as mathematical tool and comparative analysis with classical methods is established. Gait analysis unlike other works is carried out by evaluating the rotational and tilting movements of the pelvis, flexion-extension of the hip and knee; as well as dorsiflexion and plantarflexion of the ankle. Finally, an assessment of normal, mild crouch and severe crouch gaits in the three anatomical planes is performed; and a metric based on the Euclidean norm in the cartesian space is used to evaluate these gaits.


INTRODUCTION
Gait analysis has been used to evaluate different conditions in sport sciences, biomechanics and clinical diagnosis [1] . In clinical environments, this assessment tool has been applied to diagnose: hemiplegia [2] , Achilles tendinopathy [3] , inversion sprains [4] , Parkinson's disease [5] , hip arthroplasty [6] , knee osteoarthritis [7] , idiopathic scoliosis [8] , cervical myelopathy [9] , among others. In the same way, gait analysis has allowed to estimate the progress of patients rehabilitated after a stroke [10] ; or simply to determine the joint displacement in the hip, knee and ankle [11] [12] . In this paper, we focus on studying the kinematics of the lower limbs in human gait based on quaternions algebra.
On the other hand, quaternions are useful to perform a rotation of vectors in a 3D space [13] . Today, they are widely used in computer graphics, multirotor tracking and control approaches, and kinematics and dynamics of rigid bodies [14] . Therefore, the free representation of the quaternions in the Euclidian space has also been used in navigation, computer-aided design and computer vision [15] .
Recently, interest in robotic developments has been increased. The most well-known methods for robot kinematics are the Denavit-Hartenberg convention [16] and geometric methods [17] . Therefore, most of the work done to model robot kinematics using quaternions continues to follow the D-H approach, wasting quaternions capacity [15] . Some robotic applications have been focused on gait approaches to establish the imitation of human gait. For example, a robotic platform for the kinetics and kinematics characterization during gait has been used [18] . Therefore, mathematical models that associate forward kinematics of position are required [19] . For this reason, in this work we present a theoretical approach for modeling forward kinematics of position of the lower limbs for human gait using quaternions algebra not based on the D-H convention.

Related work
Some methods have been proposed in the literature for the three-dimensional analysis of gait kinematics using wearable sensors and quaternions algebra [20] [21] .
In the first work, the initial orientations are computed by quaternions of the inertial sensors placed in pelvis, femur, tibia and foot, which are acquired from the acceleration data, while the angular displacement is defined from the angular velocity. Subsequently, the orientations of the independent body segments are obtained from the sensor's orientation and the calibration rotation matrix, to then synthesize a 3D model of the whole body concerning the global reference frame.
Similarly, the algebra-based position of quaternions of the frame of each joint is obtained from the gyroscope and the accelerometer signals [21] . Therefore, in addition to the existence of position estimation error derived from the use of inertial sensors, in both works, there is a low performance in the gait analysis, since anthropometry is not directly considered and, the joint system is handled by an independent structure and not as a serial chain.
On the other hand, in the work [22] , a local analysis of the stability of the joints during walking is carried out using a marker-based optical system. This implies that forward kinematics of position is also obtained from the position of the markers and the calibration test [5] [6] [7] . In addition to what has been mentioned in this section in the aforementioned works, only the hip, knee and ankle are analyzed, that is why, to have greater precision in the kinematic calculation of gait, it is necessary to include the rotational and tilting movements of the pelvis.

Problem statement
Generally, the normal gait pattern is established in the joint space [23] and independently from each anatomical movement of each joint structure. Which limits the global perception of the performance of the joint system and the use of well-known metrics such as speed, cadence and stride length, which are established in the operational space (gait space) and are related to anthropometry of the individual. Therefore, the calculation of direct kinematics is necessary to determine the gait pattern in cartesian space as a function of joint variables.

Proposed solution
In this work, a method to calculate the forward kinematics of position of the articular system of the lower limbs during gait, using quaternion algebra as a mathematical tool is proposed. To this end, the modeling of the articular system of the lower extremities is established as an open serial chain, which allows obtaining a global performance of the system and a reference pattern of normal gait in the operational space. Also, the analysis in the gait space allows the evaluation of the metrics of that space and some not common metrics for this purpose such as: the Euclidean distance, the areas and the centroid between each reference.
Finally, to evaluate abnormal performances in the cartesian space for reference pattern, a comparative analysis is performed in the 3 anatomical planes for a normal gait and 2 types of crouched gaits, of which their forward position kinematics is calculated from the articular anatomical positions related to the rotation and inclination of the pelvis, dorsiflexion of the hips and knees, as well as dorsiflexion and plantarflexion of the ankles.

Paper organization
In this work, the anthropometry and the joint parameters considered as the starting point are adapted to [20] and [24] , respectively. Subsequently, mathematical modeling is proposed to calculate the forward kine- To describe the work done, this paper is organized as follows: Related work, problem statement and proposed solution are described in the introduction section. In the materials and methods section, the mathematical modeling for the calculation of the forward kinematics of position of joint systems and lower limbs using the methods: i) geometric, ii) Denavit-Hartenberg and iii) quaternions algebra, is presented.
The analysis and discussion of results is shown in the corresponding section and finally, the last section focuses on the conclusions.

MATERIALS AND METHODS
Kinematics is the study of the motion of mechanical systems without regard to the cause of the motion.
The most well-known methods for forward kinematics of position of robots are the Denavit-Hartenberg convention [16] and geometric [17] methods. In this work, the lower limb joint system is considered as a set of rigid links connected together at various joints during the swing phase of the gait cycle. Therefore, the methodologies used in robotics such as: i) geometric method ii) D-H) and iii) quaternions [25] , can be used to calculate the forward kinematics of position of the lower limbs during gait, which is described as follows.

Hip -knee -ankle system (HAK)
To develop the kinematic analysis in this work, the anatomical terms describing the relationships of the different parts of the body are based on the anatomical positions of sagittal, frontal and transverse planes; and their main directions [23] . Firstly, we develop a synthesis of forward kinematics of position of 3 DoF system as a planar robot, whose home position in ( Figure   2) of the first link is on the y 0R (-) axis of the origin of the base frame OΣ 0R , and the origins of the orthonormal frames of the robot correspond with the reference points of the joints of the human's lower right limb model (Table 1). In this case, the anatomical sagittal plane corresponds to the x-y plane of the robot, while that z 0R (+), x 0R (+) and y 0R (+) axes correspond to the right, front and top directions, respectively [23] . The model features lower limb right joint as 3 rigid-body segment 1) femur, 2) tibia and 3) foot and the anthropometry is adopted from [20] . The relative motion of these segments is defined successively by quaternions algebra.
In a 3 DoF planar robot, the width of each servomotor and the thickness of the joint bar are determined by β 1R , β 2R and β 3R , for simplicity in this case, these values are equal to 0. The placement of the z iR axes coincide with the axes of rotation of the joints, while the x iR axis is assigned in the direction of the link and the y iR axis according to the right-hand rule. Meanwhile, q 1R , q 2R and q 3R are the hip and knee flexo-extension, and ankle dorsiflexion and plantarflexion, with respect to the z 0R , z 1R and z 2R axis, respectively. (1)

Geometric method
The position of the right big toe [x 3R , y 3R , z 3R ] T ϵ R 3x1 without orientation, could be obtained geometrically adding a link in the x 0R axis direction, whose joint is arranged in the Σ 2R frame origin ( Figure 2) [17] , such that, the coordinates of such position are determined as (2) . a) Open kinematic chain that represents the frames that define joint movements in the sagittal plane. b) Correspondence of the frames to the musculoskeletal reference points in the female model [27] .
where c qiR and s qiR are the cosine and sine transcendental functions that depend on q iR , respectively.
Similarly, the position of the left big toe could be obtained by the D-H convention.

Denavit -Hartenberg convention
If we assume the same configuration ( Figure 2), by means of the D-H convention [16] the parameters of the joints and the links of the right lower limb are shown in (Table 2).

Gait type Movement model
Severe crouch (SC) crouch4 Thus, the homogeneous transformation matrices for each link can be written as Then, applying trigonometric identities D-H procedure [16] , the generalized homogeneous transformation matrix results as where x 3R , y 3R and z 3R are equivalents to the Equation (1), Equation (2) and Equation (3) of the geometric method. Therefore, the end position of the point OΣ 3R (x 3R , y 3R , z 3R ) is the same using the geometric method and the D-H convention.

Quaternions-based forward kinematics
The four-dimensional space H is formed by the real axis and three orthogonal axis, spanned by the principal imaginaries vectors i= (1,0,0), j=(0,1,0) and k= (0,0,1), which obey Hamilton rules [28] : i 2 = j 2 = k 2 = ijk= -1. Where multiplication of these imaginaries resembles cross product, such that ij= k, jk= i, ki= j, ji= -k, kj= A quaternion Q= r + xi + y j + z k consists of a real part r and a pure part v= xi + yj + zk [13] . Let Q 1 = a 1 + v 1 and Q 2 = a 2 + v 2 two quaternions, then, their product is calculated using the dot product and the cross product as: The quaternion Q= a+v also decomposes into a+bu, which resembles a complex number, where the imaginary u is a unit-three vector Thus, |(|u|)|= 1 and x, y, z are the cartesian coordinates of v. Let Q= a + bu, so its conjugate is = a -bu [29] .
On the other hand, a rotation of θ around axis u is represented as the unit quaternion.
Given a unit quaternion Q that represents a rotation, then, the rotation around an arbitrary pure vector v ϵ R 3 is where = cos ( ) -u sin ( ) is the conjugate.
Then, if we represent the HAK system (Figure 2), by quaternions algebra, from Equation (10) the quaternion representation of the rotation around the z 0R axis is Where u 1R = [0, 0, 0, k] is a unit vector, q 1R is the rotation angle around the z 0R axis and OΣ is the rotation angle around the z 0R axis and OΣ 1R = [0, 0, -l 1R j, 0] is the representation of the home position of the frame Σ 1R origin. Then, the rotated final from OΣ 1R using Equation Let us       coordinates of M. Let * = 0 + NO, so its conjugate is 290 * U = 0 − NO [29]. 291 On the other hand, a rotation of 9 around axis O is 292 represented as the unit quaternion: 293 Given a unit quaternion * that represents a rotation, 294 then, the rotation around an arbitrary pure vector M 5 + % 295 is 296 Then, if we represent the HAK system (Figure 2), by 298 quaternions algebra, from Equation (10) the quaternion 299 representation of the rotation around the & !" axis is 300 freedom increase [16] as in this model. While the quaternion-based method presented as following, represents a flexible and precise tool for this approach.
In system PHAK (Figure 3), the model features the lower limbs as 7 rigid-body segments 1) pelvis, 2) right femur, 3) left femur, 4) right tibia, 5) left tibia, 6) right foot and 7) left foot. The relative motion of these segments is defined successively by quaternions algebra.
The anthropometry is adopted from [20] . The cartesian frames Σ 0 , Σ 1 , Σ 2 , Σ 3 and Σ 4 correspond to the references  toe is calculated from Equation (11)   gait cycle. These values were adopted from de 382 gait2392_simbody model [24] for each gait (Table 4). 383 This model is a three-dimensional, 23-degree-of-384 freedom computer model of the human musculoskeletal 385 system. While, ( 1-, ( 1> , ( 2-, ( 2> , ( 3-, ( 3> , ( 4-, ( 4> , pelvis, 386 femur, tibia and foot anthropometry was adopted from 387 [20]. in that corresponding sample can be made. The data of 400 the cartesian positions of each joint are stored to be 401 compared later. 402 Table 3. Movement models used from 403 gait2392_simbody.osim [24] . toe is calculated from Equation (11)   as, dorsiflexion and plantarflexion of the ankle, for both 380 extremities, respectively. ] is the number of samples per 381 gait cycle. These values were adopted from de 382 gait2392_simbody model [24] for each gait (Table 4). 383 This model is a three-dimensional, 23-degree-of-384 freedom computer model of the human musculoskeletal 385 femur, tibia and foot anthropometry was adopted from 387 [20]. in that corresponding sample can be made. The data of 400 the cartesian positions of each joint are stored to be 401 compared later. 402 Table 3. Movement models used from 403 gait2392_simbody.osim [24] . 404  (Figure 4). In the same way.
the algorithm is applicable for the HAK system. Where, the input data q 1R , q 1l , q 2R , q 2L , q 3R , q 3L , q 4R , q 4l , q 5R , q 5l are the values of the angles for rotation and tilt of the pelvis, flexo-extension of hip and knee, as well as, dorsiflexion and plantarflexion of the ankle, for both extremities, respectively. N is the number of samples per gait cycle. These values were adopted from de gait2392_simbody model [24] for each gait (Table 4).
This model is a three-dimensional, 23-degree-of-freedom computer model of the human musculoskeletal system. While, l 1R , l 1L , l 2R , l 2L , l 3R , l 3L , l 4R , l 4L , pelvis, femur, tibia and foot anthropometry was adopted from [20] . The proposed method in this work allows a flexible and transparent coupling of the joint angles

RESULTS AND DISCUSSION
In this section, the most relevant results of this work are presented and described. Then, in order to present the movement pattern of the gait model [24] TABLE 3. Movement models used from gait2392_simbody.osim [24] .

Tabla 1
Orthonormal cartesian frames Joint references of The right lower limb Σ !" ( !" , !" , !" ) Hip Σ #" ( #" , #" , #" ) Knee  In the 3 figures, the gait pattern generated by the mathematical quaternion-based model is similar to the OpenSim® biomechanical model [24] , which is very important, since a more practical and less complex implementation can be done with the method based on quaternions than with the D-H convention.
Furthermore, as it can be seen, it is possible to approximate the gait pattern of a 23 DoF model [24] , with an 8 pattern is generally a problem presented by people with cerebral palsy [30] [31] , as well as spastic diplegia and quadriplegia [32] [33] . Disability as a result of a cerebrovascular accident or other types of accident has an evolution that dominates from the stroke phase, immediately after the accident, to a phase of low to severe spasticity, due to the inattention of the patient The position of the joints is useful to improve the anatomical models or attachments that help to walk.
Likewise, obtaining the position by forward kinematics, the cadence can be obtained from the calculation of the angular velocity, which cannot be obtained from the joint space. Also, through forward kinematics it is possible to obtain the distance between the joints during the gait cycle and from there determine the dependency that exists between them. Therefore, using metrics such as these, through these comparative metrics, it is possible to contribute to diagnosis and decision-making. Generally, these metrics are used to determine age, gender, pathologies [34] or physiotherapeutic progress [35] .
The most common metrics used in gait analysis are speed, cadence, stride length, toe angle, number of daily steps [23] . There are some less used ones based on the calculation of the area of the silhouette to determine gender [36] . This represents an area of opportunity to study and propose new metrics evaluated in the operational space. While the Euclidean distance is the

ETHICAL STATEMENT
The authors declare that they have no potential conflicts of interest regarding the research, authorship and / or publication of the article.

SUPPLEMENTARY MATERIAL
Visualizations of the positions of the lower limbs during the gait cycle in the cartesian space using the proposed method in this work and OpenSim®, in the three anatomical planes of normal, mild crouch and, severe crouch gaits are shown at: https://youtu.be/frCdRDTfLz4