This notebook is a simple demonstration of how you may be able to predict the phase of the gait cycle from the states from a simple linear regression.
warning This notebook doesn't seem to run with Pandas < 0.14.0. I get an error calling results.rsquared
.
import sys
sys.path.append('../src')
import matplotlib.pyplot as plt
from pandas import concat
import statsmodels.formula.api as smf
from gaitanalysis.gait import plot_gait_cycles
import utils
from gait_landmark_settings import settings
%matplotlib inline
from IPython.core.pylabtools import figsize
figsize(14, 10)
Load the path to the directory with the experimental data.
trials_dir = utils.config_paths()['raw_data_dir']
This is the "training" data set. We simple collect all of the gait cycles from the perturbed portion of one trial from subject "A".
trial_number = '068'
trial = utils.Trial('068')
trial.prep_data('Longitudinal Perturbation')
gait_cycles = trial.gait_data_objs['Longitudinal Perturbation'].gait_cycles
Loading pre-cleaned data: /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/cleaned-data/068-longitudinal-perturbation.h5 _write_event_data_frame_to_disk took 0.04 s to execute. Loading pre-computed inverse dynamics from /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/gait-data/068-longitudinal-perturbation.h5. _write_inverse_dynamics_to_disk took 0.08 s to execute. Loading pre-computed gait cycles from /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/gait-data/068-longitudinal-perturbation.h5. _section_into_gait_cycles took 0.45 s to execute.
Now create a big data frame which has a column with the percent gait cycle phi
. It also contains columns for all of the state values.
dfs = []
for i, df in gait_cycles.iteritems():
dfs.append(df)
bigdf = concat(dfs)
bigdf['phi'] = bigdf.index.values
Now specify a patsy formula for the module. We simply want to predict phi
using all of the states as the multivariate linear regressors.
sensors, controls = utils.load_sensors_and_controls()
model = 'phi ~ Q("' + '") + Q("'.join(sensors) + '")' #+ ' + Q("FP1.ForY") + Q("FP2.ForY")'
model
'phi ~ Q("Right.Ankle.PlantarFlexion.Angle") + Q("Right.Ankle.PlantarFlexion.Rate") + Q("Right.Knee.Flexion.Angle") + Q("Right.Knee.Flexion.Rate") + Q("Right.Hip.Flexion.Angle") + Q("Right.Hip.Flexion.Rate") + Q("Left.Ankle.PlantarFlexion.Angle") + Q("Left.Ankle.PlantarFlexion.Rate") + Q("Left.Knee.Flexion.Angle") + Q("Left.Knee.Flexion.Rate") + Q("Left.Hip.Flexion.Angle") + Q("Left.Hip.Flexion.Rate")'
results = smf.ols(model, data=bigdf).fit()
The results show that we have a decent model with a high $R^2$ value.
print(results.summary())
OLS Regression Results ============================================================================== Dep. Variable: phi R-squared: 0.904 Model: OLS Adj. R-squared: 0.904 Method: Least Squares F-statistic: 6753. Date: Fri, 24 Apr 2015 Prob (F-statistic): 0.00 Time: 12:12:24 Log-Likelihood: 8567.3 No. Observations: 8580 AIC: -1.711e+04 Df Residuals: 8567 BIC: -1.702e+04 Df Model: 12 Covariance Type: nonrobust ========================================================================================================= coef std err t P>|t| [95.0% Conf. Int.] --------------------------------------------------------------------------------------------------------- Intercept 1.1771 0.036 32.944 0.000 1.107 1.247 Q("Right.Ankle.PlantarFlexion.Angle") 0.4949 0.014 34.506 0.000 0.467 0.523 Q("Right.Ankle.PlantarFlexion.Rate") 0.0037 0.001 4.858 0.000 0.002 0.005 Q("Right.Knee.Flexion.Angle") -0.0861 0.011 -7.510 0.000 -0.109 -0.064 Q("Right.Knee.Flexion.Rate") -0.1039 0.001 -110.601 0.000 -0.106 -0.102 Q("Right.Hip.Flexion.Angle") -0.6223 0.027 -22.968 0.000 -0.675 -0.569 Q("Right.Hip.Flexion.Rate") 0.0712 0.004 16.062 0.000 0.063 0.080 Q("Left.Ankle.PlantarFlexion.Angle") 0.0477 0.014 3.423 0.001 0.020 0.075 Q("Left.Ankle.PlantarFlexion.Rate") 0.0123 0.001 16.125 0.000 0.011 0.014 Q("Left.Knee.Flexion.Angle") 0.0050 0.011 0.438 0.662 -0.017 0.027 Q("Left.Knee.Flexion.Rate") 0.0500 0.001 56.917 0.000 0.048 0.052 Q("Left.Hip.Flexion.Angle") 0.8544 0.027 31.859 0.000 0.802 0.907 Q("Left.Hip.Flexion.Rate") -0.1215 0.004 -29.435 0.000 -0.130 -0.113 ============================================================================== Omnibus: 1487.314 Durbin-Watson: 2.342 Prob(Omnibus): 0.000 Jarque-Bera (JB): 5339.237 Skew: 0.849 Prob(JB): 0.00 Kurtosis: 6.472 Cond. No. 150. ============================================================================== Warnings: [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.
Now we load in data from the same subject but at a different speed to test the model.
trial_number = '069'
trial = utils.Trial('068')
trial.prep_data('Longitudinal Perturbation')
gait_cycles = trial.gait_data_objs['Longitudinal Perturbation'].gait_cycles
Loading pre-cleaned data: /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/cleaned-data/068-longitudinal-perturbation.h5 _write_event_data_frame_to_disk took 0.06 s to execute. Loading pre-computed inverse dynamics from /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/gait-data/068-longitudinal-perturbation.h5. _write_inverse_dynamics_to_disk took 0.38 s to execute. Loading pre-computed gait cycles from /home/moorepants/Manuscripts/gait-control-direct-id-paper/processed-data/gait-data/068-longitudinal-perturbation.h5. _section_into_gait_cycles took 0.15 s to execute.
dfs = []
for i, df in gait_cycles.iteritems():
dfs.append(df)
bigdf = concat(dfs)
bigdf['phi'] = bigdf.index.values
The following plot shows the actual percent gait cycle and that predicted by the linear model given the states from the new trial.
plt.plot(results.predict(bigdf)[:400])
plt.plot(bigdf['phi'].values[:400])
plt.legend(['Prediction', 'Actual'])
plt.ylabel('Percent Gait Cycle')
<matplotlib.text.Text at 0x7f7828374a50>
!git rev-parse HEAD
6af0c229bb0a9edf78f8489d5f0efe02d5685a57
!git --git-dir=/home/moorepants/src/GaitAnalysisToolKit/.git --work-tree=/home/moorepants/src/GaitAnalysisToolKit rev-parse HEAD
9e80dfdcfe0a14b44e0ebcbadb6e9e827d215c3c
%install_ext http://raw.github.com/jrjohansson/version_information/master/version_information.py
Installed version_information.py. To use it, type: %load_ext version_information
%load_ext version_information
%version_information gaitanalysis, dtk, numpy, scipy, pandas, matplotlib, tables, statsmodels, oct2py
Software | Version |
---|---|
Python | 2.7.9 64bit [GCC 4.4.7 20120313 (Red Hat 4.4.7-1)] |
IPython | 3.0.0 |
OS | Linux 3.13.0 49 generic x86_64 with debian jessie sid |
gaitanalysis | 0.2.0dev |
dtk | 0.4.0 |
numpy | 1.9.2 |
scipy | 0.15.1 |
pandas | 0.16.0 |
matplotlib | 1.4.2 |
tables | 3.1.1 |
statsmodels | 0.6.1 |
oct2py | 3.1.0 |
Fri Apr 24 12:12:26 2015 PDT |
!conda list
# packages in environment at /home/moorepants/anaconda/envs/gait-direct-id-paper: # backports.ssl-match-hostname 3.4.0.2 <pip> cairo 1.12.18 0 dateutil 2.4.1 py27_0 dynamicisttoolkit 0.4.0 py27_0 fastcache 1.0.2 py27_0 fontconfig 2.11.1 2 freetype 2.4.10 0 gaitanalysistoolkit (/home/moorepants/src/GaitAnalysisToolKit) 0.2.0.dev0 <pip> hdf5 1.8.14 0 ipython 3.0.0 py27_0 ipython-notebook 3.0.0 py27_1 jinja2 2.7.3 py27_1 jsonschema 2.4.0 py27_0 libpng 1.5.13 1 libsodium 0.4.5 0 libxml2 2.9.0 0 markupsafe 0.23 py27_0 matplotlib 1.4.2 np19py27_0 mistune 0.5.1 py27_0 nose 1.3.4 py27_1 numexpr 2.3.1 np19py27_0 numpy 1.9.2 py27_0 oct2py 3.1.0 py27_0 openssl 1.0.1k 1 pandas 0.16.0 np19py27_1 patsy 0.3.0 np19py27_0 pip 6.1.1 py27_0 pixman 0.26.2 0 ptyprocess 0.4 py27_0 py2cairo 1.10.0 py27_2 pygments 2.0.2 py27_0 pyparsing 2.0.1 py27_0 pyqt 4.10.4 py27_0 pyside 1.2.1 py27_1 pytables 3.1.1 np19py27_2 python 2.7.9 2 python-dateutil 2.4.2 py27_0 pytz 2015.2 py27_0 pyyaml 3.11 py27_0 pyzmq 14.5.0 py27_0 qt 4.8.5 0 readline 6.2 2 scipy 0.15.1 np19py27_0 seaborn 0.5.1 np19py27_0 setuptools 15.1 py27_1 shiboken 1.2.1 py27_0 sip 4.15.5 py27_0 six 1.9.0 py27_0 snakeviz 0.3.1 <pip> sqlite 3.8.4.1 1 ssl_match_hostname 3.4.0.2 py27_0 statsmodels 0.6.1 np19py27_0 sympy 0.7.6 py27_0 system 5.8 2 tables 3.1.1 <pip> terminado 0.5 py27_0 tk 8.5.18 0 tornado 4.1 py27_0 yaml 0.1.4 0 zeromq 4.0.4 0 zlib 1.2.8 0
!pip freeze
backports.ssl-match-hostname==3.4.0.2 DynamicistToolKit==0.4.0 fastcache==1.0.2 -e git+git@github.com:csu-hmc/GaitAnalysisToolKit.git@9e80dfdcfe0a14b44e0ebcbadb6e9e827d215c3c#egg=GaitAnalysisToolKit-origin_speedup-inverse-dynamics ipython==3.0.0 Jinja2==2.7.3 jsonschema==2.4.0 MarkupSafe==0.23 matplotlib==1.4.2 mistune==0.5.1 nose==1.3.4 numexpr==2.3.1 numpy==1.9.2 oct2py==3.1.0 pandas==0.16.0 patsy==0.3.0 ptyprocess==0.4 Pygments==2.0.2 pyparsing==2.0.1 PySide==1.2.1 python-dateutil==2.4.2 pytz==2015.2 PyYAML==3.11 pyzmq==14.5.0 scipy==0.15.1 seaborn==0.5.1 six==1.9.0 snakeviz==0.3.1 statsmodels==0.6.1 sympy==0.7.6 tables==3.1.1 terminado==0.5 tornado==4.1