Gait analysis method based on combination of smart phone and inertial sensor

Abstract

A gait analysis method based on combination of a smart phone and an inertial sensor belongs to the technical field of gait analysis and comprises the following steps: S1, acquiring walking data of the inertial sensor, transmitting the data to the smart phone through Bluetooth and uploading the data to a server for processing; s2, calculating the attitude by using an improved Madgwick filtering algorithm; s3, displacement is calculated by converting the coordinate system; s4, detecting a gait event: performing filtering processing on sagittal plane angular velocity data to detect a heel landing event and a toe off-ground event; and S5, calculating space-time gait parameters. According to the method, the defects of a single smart phone and an inertial sensor in the aspects of sensor number, measurement data precision, wearing comfort and gait parameter calculation diversity are overcome. According to the invention, the improved Madgwick filtering algorithm is used to realize the attitude resolving algorithm, and the space-time gait parameters are accurately calculated by combining the gait event detection algorithm. The method does not need specific calibration motion, does not need magnetometer data, and has better adaptability.

Description

Technical Field

[0001] The present invention belongs to the technical field of gait analysis, and particularly relates to a gait analysis method based on the combination of a smart phone and an inertial sensor. Background Art

[0002] Gait analysis, as a non-invasive assessment tool, plays an important role in clinical diagnosis and health assessment. So far, the gold standard gait analysis techniques include optical motion capture systems, force plates, or expensive instrumented walkways. They usually require a high degree of proficiency, as well as complex system setups and specific usage environments, which limit their accessibility and feasibility in non-laboratory environments.

[0003] In recent years, most studies have discussed low-cost wearable systems to provide the feasibility of instrumented gait analysis for a wider community and non-professional users, and to promote its use independent of gait laboratories and in daily life environments. The gait analysis method based on inertial measurement units (IMUs) and smart phones has the potential to be used in this way. However, the current research has the following problems:

[0004] (1) Although smart phones are portable, they are usually limited to a single device (such as the prior art application numbers: CN202110012482.2, CN202410360943.9). Gait analysis requires the simultaneous analysis of the left and right feet. The types of gait parameters obtained only by relying on the built-in sensors of smart phones are limited, and these parameters are not comprehensive enough for clinical evaluation or judgment.

[0005] (2) Smart phones often use a unique fixing strap (such as the prior art application number: CN202411175256.6). When the smart phone is fixed to the target object through the fixing strap or tape, good measurement accuracy can usually be obtained. However, considering the volume and weight of the smart phone, the comfort and biocompatibility of wearing still need to be further discussed.

[0006] (3) The performance (accuracy and sensitivity) of the built-in sensors of smart phones cannot fully reach the level of commercial inertial measurement units. High-performance built-in sensors are only available in more advanced models (such as the prior art literature Yan W, Bastos L, Performance assessment of the android smartphone’s IMU in a GNSS / INS coupled navigation model[J]. IEEE Access, 2019, 7: 171073-171083.), the high drift and error sensitivity of low-cost sensors lead to a high degree of uncertainty in the measurement results (such as the existing technical literature Real Ehrlich C, Blankenbach J. Indoor localization for pedestrians with real-time capability using multi-sensor smartphones[J]. Geo-spatial Information Science, 2019, 22(2): 73-88.), there are performance parameter differences between the built-in sensors of different models of smartphones (such as the existing technical literature Ahmed HU, Hu L, Yang X, et al. Effects of smartphone sensor variability in road roughness evaluation[J]. International Journal of Pavement Engineering, 2022, 23(12): 4404-4409.), such as sensor measurement range, sampling frequency, sensitivity error, and these parameters vary due to the hardware and operating systems of different manufacturers (such as the existing technical literature Chen KH, Hsu YW, Yang JJ, et al. Evaluating the specifications of built-in accelerometers in smartphones on fall detection performance[J]. Instrumentation Science & Technology, 2018, 46(2): 194-206.), which may lead to inconsistent measurement data and make it difficult to adapt the algorithm.

[0007] (4) Inertial sensors generally require precise initial attitude and position calibration (such as the existing patent application numbers: CN201711398588.0, CN202210661225.6), otherwise it will lead to deviations in subsequent measurement results. In practical applications, it may be difficult for users to perform accurate initial calibration, affecting the accuracy of gait analysis. At the same time, inertial sensors are easily interfered by the external environment, such as magnetic field interference, vibration, etc. These interferences will affect the output data of the sensors, and thus affect the reliability of gait analysis.

[0008] (5) Related research on smartphones and inertial sensors in gait analysis is usually limited to one plane or direction and only focuses on certain events in the gait cycle (heel strike and toe off) (such as the prior art application numbers: CN201811122096.3, CN202011367253.4). Summary of the Invention

[0009] To solve the above problems, the present invention proposes a gait analysis method based on the combination of a smartphone and an inertial sensor, including the following steps:

[0010] S1. Collect the walking data of the inertial sensor, and the data is transmitted to the smartphone via Bluetooth and uploaded to the server for processing;

[0011] S2. Use the improved Madgwick filtering algorithm to solve the attitude: The improved Madgwick filtering algorithm fuses the accelerometer and gyroscope data to calculate the attitude quaternion; the algorithm removes the magnetometer term and introduces an adaptive step size to improve the accuracy;

[0012] S3. Calculate the displacement by converting the coordinate system: The acceleration data is rotated from the device coordinate system to the earth coordinate system through quaternion rotation and the gravity component is removed; combined with the zero velocity update method, the three-dimensional displacement is calculated by integration;

[0013] S4. Detect gait events: Detect the heel strike and toe off events through the filtered sagittal plane angular velocity data;

[0014] S5. Calculate the spatio-temporal gait parameters: Calculate the step frequency, step length and angle parameters based on the gait events and displacement data; the results are displayed in real time through the smartphone and can be exported as a structured file.

[0015] Further, in the step S1, the inertial sensor is fixed at the back of the shoe to collect the three-axis acceleration and angular velocity data.

[0016] Further, in the step S2, it includes the following steps:

[0017] S21. Under the original quaternion update framework based on gradient descent, remove the terms related to the magnetometer and only retain the accelerometer and gyroscope data;

[0018] S22. Optimize the construction of the objective function and the Jacobian matrix to make it only related to the acceleration measurement;

[0019] S23. Introduce an adaptive step size adjustment strategy to dynamically adjust the gradient descent step size according to the measurement noise characteristics of the accelerometer, and improve the convergence speed and robustness.

[0020] Further, in the step S21, remove the terms related to the magnetometer:

[0021] Original objective function:

[0022]

[0023] Among them, the first 3 lines correspond to the gravitational acceleration error term respectively represent the differences between the estimated attitude and the measured acceleration in the x, y, and z directions. The last 3 lines correspond to the geomagnetic field direction error term respectively represent the differences between the estimated attitude and the measured geomagnetic field direction in the x, y, and z directions; is the normalized accelerometer measurement value, is the normalized magnetometer measurement value, b x , b z is the component of the geomagnetic field direction in the sensor coordinate system, q = [q 0 , q 1 , q 2 , q 3 T is the currently estimated attitude quaternion;

[0024] Objective function after removing the magnetometer term:

[0025]

[0026] Original Jacobian matrix:

[0027]

[0028] Jacobian matrix after removing the magnetometer term:

[0029]

[0030] Furthermore, in the step S23, an adaptive step size adjustment is introduced:

[0031] Dynamically adjust the step size according to the characteristics of the accelerometer measurement noise. Adaptive step size adjustment formula:

[0032]

[0033] Among them, β k is the gradient descent step size at the k-th moment, β 0 is the initial step size, γ is the attenuation coefficient, is the variance of the accelerometer measurement at the k-th moment;

[0034] Calculation formula for the accelerometer measurement variance:

[0035]

[0036] ​Among them, N is the window length for variance calculation, which is the average value measured by the accelerometer within the window.

[0037] Furthermore, according to the simplified objective function and the adaptive step size, modify the quaternion update formula:

[0038] Original quaternion update formula:

[0039]

[0040] Improved quaternion update formula:

[0041]

[0042] Among them, Δq is the quaternion increment obtained by integrating the angular velocity: Δq = is the gradient of the objective function at q k-1 : ω k = [ω k,x , ω k,y , ω k,z T represents the angular velocity vector measured by the gyroscope at time k; T is the sampling time interval;

[0043] Quaternion normalization formula:

[0044] Among them, q k is the modulus of the quaternion:

[0045] Furthermore, in the step S3, the coordinate system transformation process is as follows:

[0046] According to the acceleration data A sensor = [A x , A y , A z T collected by the inertial sensor during walking, and solve the attitude data quaternion Q = [Q 0 , Q 1 , Q 2 , Q 3 T , convert the acceleration data into a pure quaternion Rotate the acceleration in the inertial sensor coordinate to the earth coordinate system; the coordinate system transformation formula is as follows:

[0047]

[0048] Among them, Q * is the conjugate of the quaternion Q; represents quaternion multiplication; ​​​

[0049] Remove the gravitational acceleration, and the formula is as follows:

[0050] A earth =[A earth x,A earth y,A earth z]-[0,0,g]#(11)

[0051] where, A earth x, A earth y, A earth z are the triaxial accelerations in the Earth coordinate system; where, x is the east direction, y is the north direction, z is the vertical direction; g is the gravitational acceleration constant;

[0052] The formula for calculating the velocity is as follows:

[0053]

[0054] where, stationary is the motion state (0: relative motion; 1: relative rest); f s is the sampling rate;

[0055] During the non-stationary period, the velocity drift V drift is calculated by linear interpolation; let the start and end of the non-stationary period be t start and t end , the drift rate is the drift amount is V drift (t), and the final velocity is obtained by subtracting the integral drift amount;

[0056]

[0057] V = V - V drift #(14)

[0058] The displacement P(t) is obtained by integrating the final velocity:

[0059]

[0060] Furthermore, in the step S4, the gait events are detected by the sagittal plane angular velocity of the inertial sensor and the vertical direction position change information of the reflective marker points of the Vicon motion capture system respectively.

[0061] Furthermore, in the step S5, the gait parameter calculation

[0062] Based on the heel strike point HS i and the toe off point TO i of the current gait cycle i, the following gait parameters are calculated:

[0063] (1) Step frequency: The number of steps taken per minute. One step is considered to be between the landing points of two consecutive heel strikes.

[0064]

[0065] (2) Stride time: The duration between the heel strike of one leg and the next heel strike of the same leg.

[0066] Stride time=HS i+1 -HS i #(17)

[0067] (3) Stance phase: This phase starts from the heel strike of one foot and ends when the toes of that foot leave the ground.

[0068]

[0069] (4) Swing phase: This phase starts from the toe-off of one foot and ends when the same foot strikes the ground again.

[0070]

[0071] (5) Stance time: The period during which the foot is in contact with the ground in a gait cycle.

[0072] Stance time=TO i -HS i #(20)

[0073] (6) Swing time: The period during which the foot is off the ground in a gait cycle.

[0074] Swing time=HS i+1 -TO i #(21)

[0075] (7) Double stance phase: Refers to the period during which both feet are in contact with the ground in a gait cycle. The calculation formula for the double stance phase is usually based on the gait cycle and the time of the single-leg stance phase, approximately 20 - 30%.

[0076] Dual stance phase=Left(Right)stance phase - Right(left)swingphase#(22)

[0077] (8) Double stance time: The total time during which both feet are in contact with the ground in a gait cycle.

[0078] Dual stance time=Stride time * Dual stance phase#(23)

[0079] (9) Step length: The distance between the times when the sole touches the ground in the horizontal plane (X and Y) displacement information of the sensor in the Earth coordinate system;

[0080]

[0081] where p x is the displacement along the x-axis of the Earth coordinate system; p y is the displacement along the y-axis of the Earth coordinate system; t is the time point when the sole touches the ground detected;

[0082] (10) Toe-off height: Defined as the maximum gap in the vertical direction between the foot and the ground during walking;

[0083]

[0084] where P max|z is the local maximum value of the displacement along the z-axis of the Earth coordinate system; P footflat|z is the displacement value corresponding to the sole touching the ground of the displacement along the z-axis of the Earth coordinate system;

[0085] (11) Heel strike angle and toe-off angle

[0086] According to the acceleration and angular velocity, combined with the Madgwick algorithm to achieve attitude solution to calculate the quaternion, convert the quaternion to Euler angles, and calculate the heel strike angle and toe-off angle according to the pitch angle in the Euler angles.

[0087] The beneficial effects of the present invention are:

[0088] In view of the above situation, the present invention aims to explore a gait analysis method using a combination of a smart phone and a wearable device, to make up for the defects of a single smart phone and an inertial sensor in terms of the number of sensors, the accuracy of measurement data, wearing comfort, and the diversity of gait parameter calculation. The advantages of the present invention are:

[0089] (1) The present invention proposes a gait analysis method combining a smart phone and an inertial sensor, using the smart phone to connect an external inertial sensor, obtaining comprehensive and accurate gait analysis results through multiple sensors, and using the smart phone to achieve real-time visualization display.

[0090] (2) The present invention proposes an attitude solution algorithm using an improved Madgwick filtering algorithm, and accurately calculates spatio-temporal gait parameters by combining a gait event detection algorithm. This method does not require a specific calibration movement, nor does it require magnetometer data, and has better adaptability.

[0091] (3) In order to understand the motion pattern more deeply, the study obtains more gait parameter information from the sagittal plane and the transverse (horizontal) plane. Description of the Drawings

[0092] Figure 1 Schematic diagram of data acquisition for the present invention;

[0093] Figure 2 Schematic diagram of the quaternion calculation result for the present invention;

[0094] Figure 3 Schematic diagram of the three-axis acceleration in the device coordinate system for the present invention;

[0095] Figure 4 Schematic diagram of the three-axis acceleration in the earth coordinate system for the present invention;

[0096] Figure 5 Flowchart of coordinate system transformation and displacement calculation for the present invention;

[0097] Figure 6 Schematic diagram of sagittal plane angular velocity and gait event detection for the present invention;

[0098] Figure 7 Schematic diagram of the vertical displacement of the reflective marker points for the present invention;

[0099] Figure 8 Schematic diagram of the gait cycle for the present invention;

[0100] Figure 9 Schematic diagram of the horizontal plane displacement for the present invention;

[0101] Figure 10 Schematic diagram of the vertical direction displacement for the present invention;

[0102] Figure 11 Schematic diagram of the pitch angle for the present invention;

[0103] Figure 12 Schematic diagram of the wearing of the inertial sensor and the reflective marker points for the present invention. Detailed implementation manners

[0104] To make the technical means and the achieved purpose of the present invention easy to understand, the present invention will be further described below in conjunction with the detailed implementation manners. A gait analysis method based on the combination of a smart phone and an inertial sensor, the main steps are as follows:

[0105] S1. Collect the walking data of the inertial sensor: The inertial sensor is fixed at the heel of the shoe to collect the three-axis acceleration and angular velocity data. The data is transmitted to the smart phone via Bluetooth and uploaded to the server for processing.

[0106] S2. Calculate the attitude using the improved Madgwick filtering algorithm: The improved Madgwick filtering algorithm fuses the accelerometer and gyroscope data to calculate the attitude quaternion. The algorithm removes the magnetometer term and introduces an adaptive step size to improve the accuracy.

[0107] S3. Transform the coordinate system and calculate displacement: The acceleration data is rotated from the device coordinate system to the earth coordinate system by quaternion and the gravity component is removed. Combining with the zero velocity update method, the three-dimensional displacement is calculated by integration.

[0108] S4. Detect gait events: Detect the heel strike and toe off events by filtering the sagittal plane angular velocity data.

[0109] S5. Calculate spatio-temporal gait parameters: Calculate the step frequency, step length and angular parameters based on the gait events and displacement data. The results are displayed in real time on the smartphone and can be exported as a structured file.

[0110] S6. Verify the reliability and effectiveness of the method: Evaluate the measurement consistency of the gait parameters measured by the IMU and the Vicon optical system.

[0111] In the step S1, data acquisition:

[0112] Current studies basically choose to collect gait data with shoes on or barefoot. Since people usually wear shoes in daily life, in order to enable the wearable device based on the inertial measurement unit to perform gait analysis in the non-laboratory environment of daily life and at the same time facilitate the fixation and alignment of the sensors, the present invention chooses to fix the sensors at the back of the shoe instead of directly on the bare foot.

[0113] As Figure 1 shown, the inertial sensors are fixed at the back of the subject's shoe by (medical) tape, the smartphone is connected to the inertial sensors via Bluetooth, the sensors collect the acceleration and angular velocity data of walking in real time and upload them to the server via HTTP requests, the server calls the gait parameter parsing algorithm to calculate the gait parameters, responds to the smartphone, displays the gait parameter results, and at the same time can export the results as an Excel file to view the detailed information (gait symmetry, gait variability, single gait cycle parameters) for convenient secondary processing and analysis.

[0114] In the step S2, attitude solution is realized based on the improved Madgwick filtering algorithm:

[0115] Although the existing Madgwick filtering algorithm can fuse the data of the accelerometer, gyroscope and magnetometer for attitude estimation, in some application scenarios, the magnetometer may be interfered by the external magnetic field, resulting in a decrease in the attitude estimation accuracy. Therefore, an improved algorithm is needed to achieve robust attitude estimation only using the accelerometer and gyroscope data without magnetometer data.

[0116] In the absence of magnetometer data, attitude estimation mainly depends on the gravity direction information provided by the accelerometer and the angular velocity integration information provided by the gyroscope. The improvement idea is:

[0117] (1) Under the original quaternion update framework based on gradient descent, remove the terms related to the magnetometer and only retain the accelerometer and gyroscope data.

[0118] (2) Optimize the construction of the objective function and the Jacobian matrix so that they are only related to the acceleration measurement.

[0119] (3) Introduce an adaptive step size adjustment strategy to dynamically adjust the gradient descent step size according to the measurement noise characteristics of the accelerometer, improving the convergence speed and robustness.

[0120] The Madgwick filter algorithm is an attitude estimation algorithm based on gradient descent optimization. By fusing the measurement data of the accelerometer, gyroscope, and magnetometer, it estimates the attitude quaternion of the sensor. Its basic principle is to construct an objective function to make the estimated attitude conform to the sensor measurement data, and then minimize the objective function through the gradient descent method to iteratively update the attitude quaternion.

[0121] (1) Remove the terms related to the magnetometer

[0122] Original objective function:

[0123]

[0124] Among them, the first 3 lines correspond to the gravity acceleration error terms respectively represent the differences between the estimated attitude and the measured acceleration in the x, y, and z directions. The last 3 lines correspond to the geomagnetic field direction error terms respectively represent the differences between the estimated attitude and the measured geomagnetic field direction in the x, y, and z directions. is the normalized accelerometer measurement value. is the normalized magnetometer measurement value. b x ,b z is the component of the geomagnetic field direction in the sensor coordinate system. q = [q 0 ,q 1 ,q 2 ,q 3 T is the currently estimated attitude quaternion.

[0125] Objective function after removing the magnetometer terms:

[0126]

[0127] Original Jacobian matrix:

[0128]

[0129] Jacobian matrix after removing the magnetometer terms:​

[0130]

[0131] (2) Introduction of adaptive step size adjustment

[0132] In the original algorithm, the gradient descent step size β is a fixed hyperparameter. To improve the adaptability and robustness of the algorithm, the step size can be dynamically adjusted according to the characteristics of the accelerometer measurement noise. Adaptive step size adjustment formula:

[0133]

[0134] where β k is the gradient descent step size at the k-th moment, β 0 is the initial step size, γ is the attenuation coefficient, is the variance of the accelerometer measurement at the k-th moment.

[0135] Calculation formula for the variance of the accelerometer measurement:

[0136]

[0137] where N is the window length for variance calculation, is the average value of the accelerometer measurements within the window.

[0138] (3) Quaternion update

[0139] Modify the quaternion update formula according to the simplified objective function and the adaptive step size.

[0140] Original quaternion update formula:

[0141]

[0142] Improved quaternion update formula:

[0143]

[0144] where Δq is the quaternion increment obtained by integrating the angular velocity: is the gradient of the objective function at q k-1 : ω k = [ω k,x , ω k,y , ω k,z T represents the angular velocity vector measured by the gyroscope at time k; T is the sampling time interval.

[0145] (4) Quaternion normalization

[0146] ​To ensure the unit property of the quaternion, it is necessary to normalize the quaternion after each update. The quaternion normalization formula:

[0147]

[0148] where q k is the modulus of the quaternion: The quaternion calculation result is as Figure 2 shown.

[0149] In the step S3, coordinate system transformation and displacement calculation

[0150] (1) Coordinate system transformation

[0151] According to the acceleration data A sensor =[A x ,A y ,A z T collected by the inertial sensor during walking, and solve the attitude data quaternion Q=[Q 0 ,Q 1 ,Q 2 ,Q 3 T . Convert the acceleration data into a pure quaternion Rotate the acceleration in the inertial sensor coordinate to the Earth coordinate system. Figure 3 is the schematic diagram of the original acceleration in the device coordinate system, Figure 4 is the schematic diagram of the acceleration in the Earth coordinate system, Figure 5 is the flow chart of coordinate system transformation and displacement calculation. The coordinate system transformation formula is as follows:

[0152]

[0153] where, Q * is the conjugate of the quaternion Q; represents quaternion multiplication.

[0154] Remove the gravitational acceleration, the formula is as follows

[0155] A earth =[A earth x,A earth y,A earth z]-[0,0,g]#(11)

[0156] where, A earth x, A earth y, A earth z are the three-axis accelerations in the Earth coordinate system respectively; where, x is the east direction, y is the north direction, z is the vertical direction; g is the gravitational acceleration constant.

[0157] ​​(2) Calculation speed

[0158] During walking, the speed of the human body is theoretically zero at the moment of each foot touching the ground (i.e., the moment of relative rest). By detecting these stationary moments, these points can be used as the reference points for speed update, thereby correcting the integration error and improving the positioning accuracy, the zero speed update method. The formula is as follows:

[0159]

[0160] Among them, stationary is the motion state (0: relative motion; 1: relative rest); f s is the sampling rate.

[0161] (3) Calculate the integration drift during non-stationary periods

[0162] During non-stationary periods, the drift of the speed V drift is calculated by linear interpolation. Let the start and end of the non-stationary period be t start and t end , the drift rate is the drift amount is V drift (t), and the final speed is obtained by subtracting the integration drift amount.

[0163]

[0164] V = V - V drift #(14)

[0165] (4) Calculate displacement

[0166] The displacement P(t) is obtained by integrating the final speed.

[0167]

[0168] In the step S4, gait event detection

[0169] Gait events (heel strike and toe off) are detected by the sagittal plane angular velocity of the inertial sensor and the vertical direction position change information of the reflective marker points of the Vicon motion capture system, as Figure 6 and Figure 7 shown.

[0170] The core of gait analysis lies in gait event detection. Key events in the gait cycle are detected by the sensitive axis data of the sensor wearing position. The sensitive axis data is the data with the most periodic change characteristics of the sensor. For the sensor worn on the lower limb position, the sagittal plane angular velocity is usually selected as the sensitive axis data to detect gait events.

[0171] Such as Figure 6As shown, gait event detection starts from the zero-crossing point of the sagittal plane angular velocity. The zero-crossing point is used as the heel strike point, and the gait cycle is divided. The second local minimum in the gait cycle is used as the toe-off point. The sagittal plane angular velocity is filtered to calculate the motion state (stationary: 1 for relatively stationary; 0 for relatively moving) to determine foot contact and heel lift-off. The pseudocode in Table 1 illustrates the specific process of the gait event detection algorithm, where the threshold is an empirical value obtained by observing experimental data.

[0172] Table 1 Gait Event Detection Pseudocode

[0173]

[0174] In the step S5 described above, gait parameter calculation

[0175] As Figure 8 shown, based on the heel strike point HS i and the toe-off point TO i of the current gait cycle i, the following gait parameters are calculated:

[0176] (1) Step frequency: The number of steps taken per minute. The distance between two heel strike points is regarded as one step.

[0177]

[0178] (2) Stride time: The duration between heel strike of the same leg and the next heel strike.

[0179] Stride time=HS i+1 -HS i #(17)

[0180] (3) Stance phase: This phase starts from the heel strike of one foot and ends until the toe-off of the same foot.

[0181]

[0182] (4) Swing phase: This phase starts from the toe-off of one foot and ends until the next heel strike of the same foot.

[0183]

[0184] (5) Stance time: The time period during which the foot is in contact with the ground in the gait cycle.

[0185] Stance time=TO i -HS i #(20)

[0186] (6) Swing time: The time period during which the foot is off the ground in the gait cycle.

[0187] Swing time=HS i+1 -TO i #(21)

[0188] (7) Double support phase: Refers to the time period during a gait cycle when both feet are in contact with the ground simultaneously. The calculation formula for the double support phase is usually based on the gait cycle and the time of the single-leg support phase, approximately 20 - 30%.

[0189] Dual stancephase=Left(Right)stancephase - Right(left)swingphase#(22)

[0190] (8) Double support time: The total time during a gait cycle when both feet are in contact with the ground simultaneously.

[0191] Dual stance time=Stride time * Dual stance phase#(23)

[0192] (9) Stride length: In the horizontal plane (X and Y) displacement information of the sensor in the Earth coordinate system, the distance between the moments when the sole touches the ground. The displacement information of the horizontal plane is as Figure 9 shown.

[0193]

[0194] Among them, p x is the displacement along the x-axis of the Earth coordinate system; p y is the displacement along the y-axis of the Earth coordinate system; t is the time point when the sole touches the ground detected.

[0195] (10) Lifting height of the foot: Defined as the maximum gap in the vertical direction between the foot and the ground during walking. The vertical displacement of the foot is as Figure 10 shown.

[0196]

[0197] Among them, P max|z is the local maximum value of the displacement along the z-axis of the Earth coordinate system; P footflat|z is the displacement value corresponding to the sole touching the ground of the displacement along the z-axis of the Earth coordinate system.

[0198] (11) Heel strike angle and toe off angle

[0199] According to the acceleration and angular velocity, combined with the Madgwick algorithm to achieve attitude solution and calculate the quaternion, and convert the quaternion into Euler angles. Calculate the heel strike angle and toe off angle according to the pitch angle in the Euler angles, as Figure 11 shown.

[0200] In step S6, reliability and validity verification

[0201] To verify the reliability and validity of the gait analysis system proposed by the present invention, 20 subjects (11 males, 9 females; age 26.40±1.55; height 172.33±8.31 cm; weight 65.37±11.36 kg) were recruited.

[0202] The IMU inertial measurement unit (IM948, Chenyi Electronics, sampling rate: 100 Hz) was fixed at the heel (Achilles tendon) position of the subject. At the same time, the Vicon optical motion capture system (sampling rate: 100 Hz) was used as the gold standard device, and reflective markers were fixed at the heel and toe positions of the subject respectively, as Figure 12 shown. Each participant was required to walk a distance of at least 7 meters at slow (Low), normal (Normal), and fast (High) speeds of their own choice.

[0203] Table 2 shows the data distribution and standard error of the spatio-temporal gait parameters measured by the IMU and Vicon for the participants. The results show that at three different walking speeds, the gait parameters calculated by the IMU all have good reliability (ICCs>0.75). At slow and normal speeds, they have excellent reliability (ICCs>0.9). There are no significant differences in the gait parameters between the IMU and Vicon at three different walking speeds (P>0.05). Compared with slow and normal speeds, the ICC value is lower and the deviation value is higher at fast speed.

[0204] Table 2 Differences in gait parameters between IMU and Vicon

[0205]

[0206]

[0207] As described above, the above are only the preferred specific embodiments of the present invention, but the protection scope of the present invention is not limited thereto. Any person skilled in the art within the technical scope disclosed by the present invention, according to the technical solution and its concept of the present invention, making equivalent substitutions or changes should be covered by the protection scope of the present invention.

Claims

1. A gait analysis method based on a combination of a smart phone and an inertial sensor, characterized in that: The steps include: S1, collect walking data from the inertial sensor, transmit the data to the smartphone via Bluetooth and upload it to the server for processing; S2. Calculate attitude using the improved Madgwick filter algorithm: The improved Madgwick filter algorithm fuses the accelerometer and gyroscope data to calculate the attitude quaternion; the algorithm removes the magnetometer term and introduces an adaptive step size to improve accuracy; S3. Calculate displacement by converting the coordinate system: The acceleration data is converted from the device coordinate system to the earth coordinate system by quaternion rotation and the gravity component is removed; Combined with the zero velocity update method, the three-dimensional displacement is calculated by integration; S4, detecting gait events: detecting heel strike and toe-off events through filtering of sagittal plane angular velocity data; S5. Calculate spatiotemporal gait parameters: Calculate the frequency, stride length and angle parameters based on gait events and displacement data; the results are displayed in real time on a smartphone and can be exported as a structured file.

2. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 1, characterized in that: In step S1, an inertial sensor is fixed at the heel of the shoe to collect three-axis acceleration and angular velocity data.

3. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 1, characterized in that: The step S2 includes the following steps: S21. In the original quaternion update framework based on gradient descent, the magnetometer related items are removed, and only the accelerometer and gyroscope data are retained; S22, optimization objective function and Jacobian matrix construction, so that it is only related to acceleration measurement; S23. Introduce an adaptive step size adjustment strategy to dynamically adjust the gradient descent step size according to the measurement noise characteristics of the accelerometer to improve the convergence speed and robustness.

4. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 3, characterized in that: In step S21, the magnetometer related items are removed: Original objective function: Among them, the first 3 rows correspond to the gravity acceleration error term They represent the difference between the estimated attitude and the measured acceleration in the x, y, and z directions respectively. The last three rows correspond to the error terms in the direction of the geomagnetic field. They represent the differences between the estimated attitude and the measured geomagnetic field direction in the x, y, and z directions respectively; is the normalized accelerometer measurement, is the normalized magnetometer measurement value, b x ,b z is the component of the geomagnetic field direction in the sensor coordinate system, q = [q0, q1, q2, q3] T is the current estimated attitude quaternion; Objective function after removing the magnetometer term: Original Jacobian matrix: Jacobian matrix after removing the magnetometer term:

5. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 4, characterized in that: In step S23, an adaptive step length adjustment is introduced: The step size is dynamically adjusted according to the characteristics of the accelerometer measurement noise. The adaptive step size adjustment formula is: Among them, β k is the gradient descent step size at the kth moment, β0 is the initial step size, γ is the attenuation coefficient, is the variance of the accelerometer measurement at the kth moment; The calculation formula of accelerometer measurement variance is: Where N is the window length for variance calculation, is the average of the accelerometer measurements within the window.

6. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 5, characterized in that: According to the simplified objective function and adaptive step size, modify the quaternion update formula: Original quaternion update formula: Improved quaternion update formula: Among them, Δq is the quaternion increment obtained by integrating the angular velocity: The objective function is k-1 The gradient at: ω k =[ω k,x ,ω k,y ,ω k,z ] T represents the angular velocity vector measured by the gyroscope at time k; T is the sampling time interval; Quaternion normalization formula: where q k is the modulus of the quaternion:

7. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 1, characterized in that: In step S3, the coordinate system transformation process is as follows: According to the inertial sensor, the walking acceleration data A is collected sensor =[A x ,A y ,A z ] T , and solve the attitude data quaternion Q = [Q0, Q1, Q2, Q3] T , convert acceleration data into pure quaternion Rotate the acceleration of the inertial sensor coordinates to the earth coordinate system; the coordinate system transformation formula is as follows: Among them, Q * is the conjugate of the quaternion Q; Represents quaternion multiplication; Remove the acceleration due to gravity, the formula is as follows: A earth =[A earth x,A earth y,A earth z]-[0,0,g]#(11) Among them, A earth x,A earth y, A earth z are the three-axis accelerations in the earth coordinate system; x is the east direction, y is the north direction, and z is the vertical direction; g is the gravitational acceleration constant; The formula for calculating speed is as follows: Among them, stationary is the motion state (0: relative motion; 1: relative stillness); f s is the sampling rate; During the non-stationary period, the velocity drift V drift Calculated by linear interpolation; assuming that the start and end of the non-stationary period are t start and t end , the drift rate is The drift is V drift (t), the final velocity is obtained by subtracting the integrated drift; V=VV drift #(14 The displacement P(t) is obtained by integrating the final velocity:

8. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 1, characterized in that: In step S4, gait events are detected respectively by using the sagittal plane angular velocity of the inertial sensor and the vertical position change information of the reflective marker point of the Vicon motion capture system.

9. The gait analysis method based on the combination of a smart phone and an inertial sensor as claimed in claim 1, characterized in that: In step S5, gait parameter calculation Take the heel position HS of the current gait cycle i i and toe-off point TO i Based on this, the following gait parameters are calculated: (1) Cadence: the number of steps taken per minute, with the distance between the two heel points being considered a step; (2) stride time: the duration from heel strike to the next heel strike of the same leg; Stride time=HS i+1 -HS i #(17) (3) Stance phase: This phase starts when the heel of one foot strikes the ground and ends when the toe of that foot leaves the ground; (4) Swing phase: This phase starts when the toe of one foot leaves the ground and ends when the heel of the foot strikes the ground again; (5) Stance phase time: the time period during the gait cycle when the foot is in contact with the ground; Stance time=TO i -HS i #(20) (6) Swing phase time: the time period during the gait cycle when the foot leaves the ground; Swing time=HS i+1 -TO i #(21) (7) Double support phase: refers to the period of time when both feet are in contact with the ground during a gait cycle. The calculation formula of the double support phase is usually based on the time of the gait cycle and the single-leg support phase, which is about 20-30%; Dualstance phase=Left(Right)stance phase-Right(left)swingphase#(22 (8) Double support time: the total time that both feet are in contact with the ground during a gait cycle; Dual stance time=Stride time*Dual stance phase#(23) (9) Stride length: the distance between the moment when the foot touches the ground in the horizontal plane (X and Y) displacement information of the sensor in the earth coordinate system; Among them, p x is the x-axis displacement of the earth coordinate system; p y is the y-axis displacement of the earth coordinate system; t is the time point when the sole of the foot touches the ground; (10) Foot lift height: defined as the maximum vertical gap between the foot and the ground during walking; Among them, P max|z is the local maximum displacement of the z-axis of the earth coordinate system; P footflat|z The displacement value corresponding to the foot landing on the ground in the z-axis displacement of the earth coordinate system; (11) Heel strike angle and toe lift angle According to the acceleration and angular velocity, the Madgwick algorithm is combined to realize the posture solution to calculate the quaternion, which is converted into Euler angles. The heel landing angle and toe-off angle are calculated according to the pitch angle in the Euler angle.

Images