#!/usr/bin/env python
# coding: utf-8

# # Introduction

# This notebook shows the difference between identifying $\mathbf{m}^*(\varphi)$ and $\mathbf{K}(\varphi)$ as opposed to identifying $\mathbf{m}^*(\varphi)$ alone from this control structure:
# 
# $$\mathbf{m}(t) = \mathbf{m}^*(\varphi) - \mathbf{K}(\varphi)\mathbf{s}(t)$$

# # Imports

# In[1]:


import sys
sys.path.append('../src')


# In[2]:


import numpy as np
import pandas
from gaitanalysis.controlid import SimpleControlSolver


# In[3]:


import utils
from gait_landmark_settings import settings


# In[4]:


get_ipython().run_line_magic('matplotlib', 'inline')


# In[5]:


from IPython.core.pylabtools import figsize
figsize(10, 8)


# # Load and Process Data for 1 Trial

# In[6]:


trials_dir = utils.config_paths()['raw_data_dir']


# In[7]:


trial_number = '068'


# In[8]:


trial = utils.Trial('068')


# In[9]:


trial.prep_data('Longitudinal Perturbation')


# In[10]:


gait_data = trial.gait_data_objs['Longitudinal Perturbation']
gait_cycles = gait_data.gait_cycles


# In[11]:


num_cycles = gait_cycles.shape[0]
solver = SimpleControlSolver(gait_cycles.iloc[:num_cycles * 3 / 4],
                             trial.sensors,
                             trial.controls,
                             validation_data=gait_cycles.iloc[num_cycles * 3 / 4:])


# # Identification only $\mathbf{m}^*(t)$

# In[12]:


gain_inclusion_matrix = np.zeros((len(trial.controls), len(trial.sensors))).astype(bool)


# In[13]:


result = solver.solve(gain_inclusion_matrix=gain_inclusion_matrix)


# In[14]:


no_control_vafs = utils.variance_accounted_for(result[-1], solver.validation_data, trial.controls)


# In[15]:


fig, axes = utils.plot_validation(result[-1], gait_data.data, no_control_vafs)


# # Identify $\mathbf{m}^*(t)$ and $\mathbf{K}(\varphi)$ (joint isolated structure)

# In[16]:


for i, row in enumerate(gain_inclusion_matrix):
    row[2 * i:2 * i + 2] = True


# In[17]:


result = solver.solve(gain_inclusion_matrix=gain_inclusion_matrix)


# In[18]:


joint_isolated_control_vafs = utils.variance_accounted_for(result[-1],
                                                           solver.validation_data,
                                                           trial.controls)


# In[19]:


fig, axes = utils.plot_validation(result[-1], gait_data.data, joint_isolated_control_vafs)


# # Identify $\mathbf{m}^*(t)$ and $\mathbf{K}(\varphi)$ (full gain matrix)

# Note that this solution will take 30 minutes to and hour if the `ignore_cov` flag is False. This is due to a speed bottleneck in `dtk.process.least_squares_variance`.

# In[20]:


result = solver.solve(ignore_cov=True)


# In[21]:


full_control_vafs = utils.variance_accounted_for(result[-1],
                                                 solver.validation_data,
                                                 trial.controls)


# In[22]:


fig, axes = utils.plot_validation(result[-1], gait_data.data, full_control_vafs)


# # Compare VAF for each identification

# In[23]:


vafs = no_control_vafs.copy()
for k, v in vafs.items():
    vafs[k] = [v, joint_isolated_control_vafs[k], full_control_vafs[k]]


# In[24]:


vaf_df = pandas.DataFrame(vafs, index=['No Control', 'Joint Isolated Control', 'Full Control'])
vaf_df.T


# In[25]:


mean_vaf = vaf_df.mean(axis=1)
mean_vaf


# The joint isolated gain matrix structure helps the base model account for ~3% more of the data. 

# In[26]:


100.0 * (mean_vaf.iloc[1] - mean_vaf.iloc[0])


# So the full gain matrix helps the base model account for ~6% more of the data.

# In[27]:


100.0 * (mean_vaf.iloc[2] - mean_vaf.iloc[0])


# # Footer

# In[28]:


get_ipython().system('git rev-parse HEAD')


# In[29]:


get_ipython().system('git --git-dir=/home/moorepants/src/GaitAnalysisToolKit/.git --work-tree=/home/moorepants/src/GaitAnalysisToolKit rev-parse HEAD')


# In[30]:


get_ipython().run_line_magic('install_ext', 'http://raw.github.com/jrjohansson/version_information/master/version_information.py')


# In[31]:


get_ipython().run_line_magic('load_ext', 'version_information')


# In[32]:


get_ipython().run_line_magic('version_information', 'numpy, scipy, pandas, matplotlib, tables, oct2py, dtk, gaitanalysis')


# In[33]:


get_ipython().system('conda list')


# In[34]:


get_ipython().system('pip freeze')