# Exploration of the data set with human balance and posture evaluations¶

Marcos Duarte

Here is an exploration of the data set with human balance and posture evaluations of 49 subjects (27 young subjects and 22 elderly subjects, older than 60 years old).
The data is available at Figshare, DOI: 10.6084/m9.figshare.4525082.
All the subjects stood still in four different conditions where vision and the standing surface were manipulated: on a rigid surface with eyes open; on a rigid surface with eyes closed; on an unstable surface with eyes open; on an unstable surface with eyes closed.
The subjects had their balance and posture evaluated by a dual force platform setup (one force plate under each foot) and a motion capture system to quantify the whole body 3D kinematics.

See the companion website of the data set (http://demotu.org/datasets/) for more information.

## Python libraries¶

First, let's import the necessary Python libraries and configure the environment:

In [1]:
import numpy as np
from scipy.signal import butter, filtfilt
import pandas as pd
import matplotlib
import matplotlib.pyplot as plt
import matplotlib.gridspec as gridspec
%matplotlib inline
import seaborn as sns
sns.set_context("notebook", font_scale=1.3,
rc={'lines.linewidth': 1.5, 'lines.markersize': 14, 'axes.titlesize': 'x-large'})
matplotlib.rc('legend', numpoints=1, fontsize=14)
import glob
import sys, os
sys.path.insert(1, r'./../functions')
# IPython widgets:
from IPython.display import display
import ipywidgets
from ipywidgets import FloatProgress, interactive


Versions of the Python libraries used:

In [2]:
#!pip install version_information
%version_information numpy, scipy, pandas, matplotlib, seaborn

Out[2]:
SoftwareVersion
Python3.6.1 64bit [MSC v.1900 64 bit (AMD64)]
IPython6.1.0
OSWindows 10 10.0.15063 SP0
numpy1.13.0
scipy0.19.0
pandas0.20.2
matplotlib2.0.2
seaborn0.7.1
Sat Jun 24 01:57:35 2017 E. South America Standard Time

## Meta data¶

The file PDSinfo.txt contains meta data about the subjects and the experimental trials. The file has a header plus 588 rows per 29 columns (there are 12 rows for each of the 49 subjects).
Let's use the power of the pandas library to load and explore the meta data:

In [3]:
# GitHub URL:
path2 = 'https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/'
# local directory:
#path2 = r'./data'
fname = os.path.join(path2, 'PDSinfo.txt')
print(fname)
print("Information of %s subjects loaded (%s rows, %s columns)."
%(len(pd.unique(PDSinfo.Subject)), PDSinfo.shape[0], PDSinfo.shape[1]))

https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDSinfo.txt
Information of 49 subjects loaded (588 rows, 29 columns).


Here are the first 12 rows and first 10 columns of meta data:

In [4]:
PDSinfo.iloc[:12, :10]

Out[4]:
Trial Subject Vision Surface Rep Age AgeGroup Gender Height Mass
0 PDS01OR1 1 Open Rigid 1 25.67 Young M 1.72 74.3
1 PDS01OR2 1 Open Rigid 2 25.67 Young M 1.72 74.3
2 PDS01OR3 1 Open Rigid 3 25.67 Young M 1.72 74.3
3 PDS01OF1 1 Open Foam 1 25.67 Young M 1.72 74.3
4 PDS01OF2 1 Open Foam 2 25.67 Young M 1.72 74.3
5 PDS01OF3 1 Open Foam 3 25.67 Young M 1.72 74.3
6 PDS01CR1 1 Closed Rigid 1 25.67 Young M 1.72 74.3
7 PDS01CR2 1 Closed Rigid 2 25.67 Young M 1.72 74.3
8 PDS01CR3 1 Closed Rigid 3 25.67 Young M 1.72 74.3
9 PDS01CF1 1 Closed Foam 1 25.67 Young M 1.72 74.3
10 PDS01CF2 1 Closed Foam 2 25.67 Young M 1.72 74.3
11 PDS01CF3 1 Closed Foam 3 25.67 Young M 1.72 74.3

To analyze the subjects' characteristics we can drop the additional rows for each subject:

In [5]:
info = PDSinfo.drop_duplicates(subset='Subject', inplace=False)
info

Out[5]:
Trial Subject Vision Surface Rep Age AgeGroup Gender Height Mass ... Nmedication Medication Ortho-Prosthesis Ortho-Prosthesis2 Disability Disability2 Falls12m PhysicalActivity Sequence Date
0 PDS01OR1 1 Open Rigid 1 25.67 Young M 1.72 74.30 ... 0 No Yes Corrective lens No No 0 1 OR, OF, CF, CR 2016-08-01 11:00:17.753
12 PDS02OR1 2 Open Rigid 1 26.50 Young M 1.92 126.30 ... 0 No Yes 16 screws and 4 plates in the left knee No No 0 0 OR, OF, CR, CF 2016-08-01 15:52:23.189
24 PDS03OR1 3 Open Rigid 1 22.33 Young F 1.67 52.90 ... 0 No Yes Dental retention No No 0 0 OF, CF, OR, CR 2016-08-02 19:21:10.859
36 PDS04OR1 4 Open Rigid 1 33.33 Young M 1.79 75.85 ... 0 No Yes Dental braces No No 0 1 OR, CF, OF, CR 2016-08-03 16:25:52.869
48 PDS05OR1 5 Open Rigid 1 24.50 Young M 1.84 61.05 ... 0 No Yes Corrective lens No No 0 0 OR, CF, CR, OF 2016-08-04 11:51:59.927
60 PDS06OR1 6 Open Rigid 1 25.00 Young M 1.74 83.15 ... 0 No Yes Corrective lens, Clips in the vesicle No No 0 0 CF, OR, OF, CR 2016-08-05 11:25:51.116
72 PDS07OR1 7 Open Rigid 1 24.17 Young F 1.58 71.75 ... 2 Selective serotonin reuptake inhibitor, Antico... No No No No 0 4 CR, OR, CF, OF 2016-08-09 12:22:12.814
84 PDS08OR1 8 Open Rigid 1 36.33 Young M 1.82 64.00 ... 0 No No No No No 0 7 OF, OR, CF, CR 2016-08-10 11:49:27.437
96 PDS09OR1 9 Open Rigid 1 25.75 Young F 1.69 61.25 ... 1 Synthetic thyroid hormone No No No No 0 0 OR, CR, CF, OF 2016-08-11 12:08:24.469
108 PDS10OR1 10 Open Rigid 1 31.25 Young F 1.62 61.70 ... 0 No No No No No 0 2 CF, OF, OR, CR 2016-08-17 10:34:51.218
120 PDS11OR1 11 Open Rigid 1 32.08 Young M 1.92 77.55 ... 0 No Yes Corrective lens No No 0 6 OF, CR, CF, OR 2016-08-17 14:42:12.108
132 PDS12OR1 12 Open Rigid 1 24.75 Young F 1.58 48.80 ... 1 Oral contraceptive Yes Corrective lens No No 0 2 OR, CF, OF, CR 2016-08-22 11:12:36.410
144 PDS13OR1 13 Open Rigid 1 70.83 Old M 1.75 65.35 ... 0 No Yes Corrective lens No No 0 7 OR, CR, OF, CF 2016-08-24 10:50:33.248
156 PDS14OR1 14 Open Rigid 1 30.67 Young M 1.71 95.40 ... 0 No No No No No 2 6 OR, OF, CF, CR 2016-08-25 19:05:56.330
168 PDS15OR1 15 Open Rigid 1 31.58 Young F 1.53 64.95 ... 1 Angiotensin II receptor antagonist No No No No 0 5 CF, OR, OF, CR 2016-08-26 11:01:19.284
180 PDS16OR1 16 Open Rigid 1 23.75 Young M 1.81 89.30 ... 0 No Yes Corrective lens No No 0 0 CR, OR, CF, OF 2016-08-26 19:52:15.697
192 PDS17OR1 17 Open Rigid 1 31.25 Young M 1.73 77.90 ... 0 No Yes Corrective lens No No 0 0 OR, OF, CR, CF 2016-09-05 10:33:04.611
204 PDS18OR1 18 Open Rigid 1 61.75 Old F 1.50 70.55 ... 3 Angiotensin II receptor antagonist, Synthetic ... Yes Corrective lens, Dental implant e Dental prost... No No 0 3 OF, CR, OR, CF 2016-09-05 15:52:11.767
216 PDS19OR1 19 Open Rigid 1 63.25 Old M 1.72 72.25 ... 1 5-alpha reductase inhibitor No No No No 0 6 CF, OF, CR, OR 2016-09-05 17:27:08.478
228 PDS20OR1 20 Open Rigid 1 28.58 Young M 1.86 79.05 ... 0 No No No No No 0 3 CR, OR, OF, CF 2016-09-09 11:44:18.748
240 PDS21OR1 21 Open Rigid 1 28.33 Young F 1.68 59.30 ... 0 No No No No No 0 0 CF, OR, OF, CR 2016-09-09 15:16:23.127
252 PDS22OR1 22 Open Rigid 1 62.17 Old M 1.64 70.50 ... 0 No Yes Corrective lens No No 0 4 OF, CR, CF, OR 2016-09-09 17:35:14.397
264 PDS23OR1 23 Open Rigid 1 29.92 Young M 1.76 68.25 ... 0 No Yes Corrective lens No No 0 1 CF, CR, OR, OF 2016-09-12 15:15:35.184
276 PDS24OR1 24 Open Rigid 1 68.67 Old M 1.67 70.25 ... 2 Beta-blocker, HMG coenzyme A reductase inhibitor Yes Corrective lens, Denture No No 0 6 OR, CF, CR, OF 2016-09-13 17:23:37.870
288 PDS25OR1 25 Open Rigid 1 21.83 Young F 1.70 62.45 ... 3 Synthetic thyroid hormone, Oral contraceptive,... Yes Corrective lens No No 1 3 CR, OR, CF, OF 2016-09-19 15:44:05.889
300 PDS26OR1 26 Open Rigid 1 22.75 Young M 1.81 61.50 ... 0 No No No No No 2 0 OR, OF, CR, CF 2016-09-20 15:32:53.602
312 PDS27OR1 27 Open Rigid 1 25.42 Young F 1.70 61.15 ... 1 Oral contraceptive No No No No 0 0 CR, CF, OF, OR 2016-09-20 20:18:48.047
324 PDS28OR1 28 Open Rigid 1 62.83 Old M 1.74 73.10 ... 2 HMG-CoA reductase inhibitor, Synthetic thyroid... Yes Corrective lens, Denture No No 0 4 CR, OF, OR, CF 2016-09-29 11:01:47.040
336 PDS29OR1 29 Open Rigid 1 84.75 Old M 1.56 66.35 ... 1 Synthetic thyroid hormone No No No No 2 2 CR, OR, CF, OF 2016-10-14 17:19:31.212
348 PDS30OR1 30 Open Rigid 1 68.42 Old F 1.47 49.20 ... 1 Angiotensin II receptor antagonist Yes Corrective lens, Breast implant No No 1 5 CF, OF, CR, OR 2016-10-15 10:38:16.189
360 PDS31OR1 31 Open Rigid 1 61.08 Old M 1.62 74.15 ... 0 No Yes Corrective lens, Denture, 2 screws in the righ... No No 4 0 OR, CR, CF, OF 2016-11-08 10:03:16.972
372 PDS32OR1 32 Open Rigid 1 72.67 Old F 1.46 65.70 ... 1 Beta blocker+Diuretic Yes Denture No No 1 3 CR, OF, OR, CF 2016-11-10 15:33:10.661
384 PDS33OR1 33 Open Rigid 1 72.25 Old F 1.54 62.00 ... 6 Angiotensin II receptor antagonist, Biguanide,... Yes Corrective lens, Bridge teeh, PaceMaker Yes Hearing (Left ear) 0 2 CF, OR, CR, OF 2016-11-10 15:08:17.682
396 PDS34OR1 34 Open Rigid 1 63.42 Old F 1.61 59.85 ... 3 Synthetic thyroid hormone, Bisphosphonate, Lon... Yes Corrective lens, Denture, Anchor in the shoulders No No 1 3 OR, CR, OF, CF 2016-11-10 18:13:48.634
408 PDS35OR1 35 Open Rigid 1 67.67 Old F 1.60 66.30 ... 3 Angiotensin II receptor antagonist, Diuretic, ... Yes Corrective lens, Screws in the left foot No No 0 3 CR, OR, CF, OF 2016-11-14 10:48:02.103
420 PDS36OR1 36 Open Rigid 1 69.17 Old M 1.60 73.20 ... 4 Angiotensin II receptor antagonist, Proton-pum... Yes Corrective lens, Dental implant No No 3 0 CF, OF, OR, CR 2016-11-14 11:52:12.017
432 PDS37OR1 37 Open Rigid 1 80.50 Old F 1.51 52.70 ... 10 Insulin, Bisphosphonate, Biguanide, Anticonvul... Yes Corrective lens, Denture, 2 screws in the righ... No No 0 3 OR, CR, OF, CF 2016-11-22 11:13:01.410
444 PDS38OR1 38 Open Rigid 1 72.00 Old M 1.72 92.95 ... 5 Biguanide, Sulfonylurea, Angiotensin II recept... Yes Corrective lens, Denture, Right knee prosthesis Yes Visual (Right eye) 4 0 CF, CR, OR, OF 2016-11-22 14:58:04.771
456 PDS39OR1 39 Open Rigid 1 69.58 Old F 1.60 66.10 ... 2 Angiotensin II receptor antagonist, Diabetes m... Yes Corrective lens, Dental implant No No 0 3 CR, OF, OR, CF 2016-11-24 15:24:06.582
468 PDS40OR1 40 Open Rigid 1 28.42 Young F 1.54 44.30 ... 0 No No No No No 0 2 CF, CR, OR, OF 2016-11-24 16:03:28.700
480 PDS41OR1 41 Open Rigid 1 68.17 Old F 1.51 65.35 ... 2 Angiotensin II receptor antagonist, Antiinflam... Yes Corrective lens, Bridge teeh Yes Visual 0 0 OF, CF, OR, CR 2016-11-25 15:42:21.303
492 PDS42OR1 42 Open Rigid 1 63.75 Old F 1.61 55.60 ... 4 HMG-CoA reductase inhibitor, Antiinflammatory,... Yes Dental implant Yes Visual 0 2 OF, CF, CR, OR 2016-11-25 16:30:15.784
504 PDS43OR1 43 Open Rigid 1 63.08 Old M 1.78 96.00 ... 3 HMG-CoA reductase inhibitor, Synthetic thyroid... No No No No 0 2 OR, CR, OF, CF 2016-11-25 17:22:43.148
516 PDS44OR1 44 Open Rigid 1 61.33 Old M 1.70 82.25 ... 4 Angiotensin II receptor antagonist, Diuretic, ... Yes Corrective lens, Hearing aid Yes Hearing (Right ear and Left ear) 0 0 OR, OF, CF, CR 2016-11-29 17:09:00.174
528 PDS45OR1 45 Open Rigid 1 24.17 Young F 1.61 51.70 ... 0 No Yes Corrective lens, Extensor left knee orthosis Yes Physical (Cerebral Palsy) 5 2 CR, OR, CF, OF 2016-12-01 14:10:38.794
540 PDS46OR1 46 Open Rigid 1 27.75 Young M 1.72 101.80 ... 0 No Yes Corrective lens No No 0 4 CR, OF, CF, OR 2016-12-01 16:07:11.341
552 PDS47OR1 47 Open Rigid 1 34.17 Young F 1.56 53.45 ... 0 No Yes Corrective lens No No 0 2 CF, OR, OF, CR 2016-12-01 16:36:58.335
564 PDS48OR1 48 Open Rigid 1 37.92 Young M 1.55 69.55 ... 0 No Yes Corrective lens No No 0 2 OR, CR, CF, OF 2016-12-03 12:40:49.071
576 PDS49OR1 49 Open Rigid 1 64.92 Old F 1.58 60.75 ... 3 HMG-CoA reductase inhibitor, Synthetic thyroid... Yes Corrective lens, Dental implant Yes Hearing (Left ear) 0 0 CF, OF, CR, OR 2016-12-06 09:33:45.819

49 rows × 29 columns

## Age group and gender¶

Here are the number of subjetcs in the data set by age group and gender:

In [6]:
display(info[['Subject', 'AgeGroup', 'Gender']].groupby(['AgeGroup', 'Gender']).count())

Subject
AgeGroup Gender
Old F 11
M 11
Young F 12
M 15

## Disability¶

Of the 49 subjects (27 young adults and 22 elderly adults), 7 of them were classified as a person with disability:

In [7]:
display(info[['Subject', 'AgeGroup', 'Disability', 'Illness']].groupby(['AgeGroup', 'Disability', 'Illness']).count())

Subject
AgeGroup Disability Illness
Old No No 2
Yes 14
Yes Yes 6
Young No No 16
Yes 10
Yes No 1

Because there were more disabilities in the older group, 6 vs. 1, or you are looking for reference data on balance by healthy people, the data of these seven subjects can be excluded with the folloing command (uncomment the code below):

In [8]:
#info = info.ix[info.Disability=='No']
#print('Number of subjects without disability: %s' %len(info.Subject))


## Age, height, mass, and BMI¶

The corresponding mean and range (minimum and maximum) values for age, height, mass, and BMI are:

In [9]:
pd.set_option('precision', 2)
print('Age, height, mass, and BMI values across subjects per age group (mean and range):')
info.groupby(['AgeGroup'])['Age', 'Height', 'Mass', 'BMI'].agg([np.mean, np.min, np.max])

Age, height, mass, and BMI values across subjects per age group (mean and range):

Out[9]:
Age Height Mass BMI
mean amin amax mean amin amax mean amin amax mean amin amax
AgeGroup
Old 67.83 61.08 84.75 1.61 1.46 1.78 68.66 49.2 96.0 26.31 21.34 31.60
Young 28.08 21.83 37.92 1.71 1.53 1.92 70.32 44.3 126.3 23.94 18.03 34.37

Let's visualize the age, weight, height, and BMI values by gender of the subjects:

In [10]:
sns.pairplot(info[['Age', 'Gender', 'Height', 'Mass', 'BMI']], hue='Gender',
size=2.5, aspect=1.2, plot_kws={'s':60})
plt.show()


## COP and COG displacements¶

Let's load the force plate and kinematics data of one subject:

In [11]:
fname_grf = os.path.join(path2, PDSinfo.Trial[150] + 'grf' + '.txt')
fname_mkr = os.path.join(path2, PDSinfo.Trial[150] + 'mkr' + '.txt')
print(fname_grf, grf.shape)
print(fname_mkr, mkr.shape)

https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDS13CR1grf.txt (6000, 21)
https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDS13CR1mkr.txt (6000, 130)


Here are plots of the center of pressure (COP) data (from each force plate and the resultant COP):

In [12]:
def cop_plot(grf):
plt.figure(figsize=(12, 8))

gs0 = gridspec.GridSpec(1, 1)
gs0.update(bottom=0.77, top=0.96)
ax0 = plt.subplot(gs0[0])
ax0.plot(grf['LCOP_Z']*100, grf['LCOP_X']*100, 'r', label='Left')
ax0.plot(grf['RCOP_Z']*100, grf['RCOP_X']*100, 'b', label='Right')
ax0.plot(grf['COPNET_Z']*100, grf['COPNET_X']*100, 'k', label='NET')
#ax0.set_ylim([np.min(grf['COPNET_X'])*100-1, np.max(grf['COPNET_X'])*100+1])
ax0.margins(0.05)
ax0.set_xlabel('COP ml [cm]')
ax0.set_ylabel('COP ap [cm]')
ax0.locator_params(axis='y', nbins=4)
ax0.yaxis.set_label_coords(-.05, 0.5)
ax0.legend(fontsize=14, bbox_to_anchor=(0.65, 1), loc=2, borderaxespad=0., framealpha=.5)
#ax0.legend(fontsize=14, loc='best', framealpha=.5)
ax0.text(0.005, 0.8, 'A', transform=ax0.transAxes, fontsize=18)

gs1 = gridspec.GridSpec(3, 2)
gs1.update(bottom=0.04, top=0.68, hspace=0.08, wspace=0.15)
ax1 = plt.subplot(gs1[0,0])
ax1.plot(grf['Time'], grf['LCOP_X']*100, 'r', linewidth=1.0, label='ap')
ax1.margins(y=0.1)
ax1.set_xlim([0, 60])
ax1.set_xlabel('')
ax1.set_xticklabels('')
ax1.set_ylabel('LCOP [cm]')
ax1.locator_params(axis='y', nbins=3)
ax1.yaxis.set_label_coords(-.1, 0.5)
ax1.legend(fontsize=14, loc='lower left', framealpha=.5)
ax1.text(0.01, 0.8, 'B', transform=ax1.transAxes, fontsize=18)

ax2 = plt.subplot(gs1[0,1])
ax2.plot(grf['Time'], grf['LCOP_Z']*100, 'r', linewidth=1.0, label='ml')
ax2.margins(y=0.1)
ax2.set_xlim([0, 60])
ax2.set_xlabel('')
ax2.set_xticklabels('')
ax2.locator_params(axis='y', nbins=3)
ax2.legend(fontsize=14, loc='lower left', framealpha=.5)
ax2.text(0.01, 0.8, 'C', transform=ax2.transAxes, fontsize=18)

ax3 = plt.subplot(gs1[1,0])
ax3.plot(grf['Time'], grf['RCOP_X']*100, 'b', linewidth=1.0, label='ap')
ax3.margins(y=0.1)
ax3.set_xlim([0, 60])
ax3.set_xlabel('')
ax3.set_xticklabels('')
ax3.set_ylabel('RCOP [cm]')
ax3.locator_params(axis='y', nbins=3)
ax3.yaxis.set_label_coords(-.1, 0.5)
ax3.legend(fontsize=14, loc='lower left', framealpha=.5)
ax3.text(0.01, 0.8, 'D', transform=ax3.transAxes, fontsize=16)

ax4 = plt.subplot(gs1[1,1])
ax4.plot(grf['Time'], grf['RCOP_Z']*100, 'b', linewidth=1.0, label='ml')
ax4.margins(y=0.1)
ax4.set_xlim([0, 60])
ax4.set_xlabel('')
ax4.set_xticklabels('')
ax4.locator_params(axis='y', nbins=3)
ax4.legend(fontsize=14, loc='lower left', framealpha=.5)
ax4.text(0.01, 0.8, 'E', transform=ax4.transAxes, fontsize=18)

ax5 = plt.subplot(gs1[2,0])
ax5.plot(grf['Time'], grf['COPNET_X']*100, 'k', linewidth=1.0, label='ap')
ax5.margins(y=0.1)
ax5.set_xlim([0, 60])
ax5.set_xlabel('Time [s]')
ax5.set_ylabel('COPNET [cm]')
ax5.locator_params(axis='y', nbins=3)
ax5.yaxis.set_label_coords(-.1, 0.5)
ax5.legend(fontsize=14, loc='lower left', framealpha=.5)
ax5.text(0.01, 0.8, 'F', transform=ax5.transAxes, fontsize=18)

ax6 = plt.subplot(gs1[2,1])
ax6.plot(grf['Time'], grf['COPNET_Z']*100, 'k', linewidth=1.0, label='ml')
ax6.margins(y=0.1)
ax6.set_xlim([0, 60])
ax6.set_xlabel('Time [s]')
ax6.locator_params(axis='y', nbins=4)
ax6.legend(fontsize=14, loc='lower left', framealpha=.5)
ax6.text(0.01, 0.8, 'G', transform=ax6.transAxes, fontsize=18)
#plt.suptitle('Center of pressure (COP) displacement during quiet standing', fontsize=14, y=1)
plt.show()

In [13]:
cop_plot(grf)


Let's plot the horizontal displacement of the center of gravity (COG) measured by the motion capture system for the COP data shown above:

In [14]:
def copcog_plot(grf, mkr):
fig, ax = plt.subplots(2, 1, figsize=(12, 5), sharex=True)
ax[0].plot(grf['Time'], (grf['RCOP_X']-grf['RCOP_X'].mean())*100, 'b', label='COP')
ax[0].plot(mkr['Time'], (mkr['COG_X']-mkr['COG_X'].mean())*100, 'r', label='COG')
ax[0].set_xlim([0, 60])
ax[0].yaxis.set_label_coords(-.05, 0.5)
ax[0].locator_params(axis='y', nbins=6)
ax[0].legend(fontsize=12, loc='upper center', framealpha=.5)
ax[0].set_ylabel('ap direction [cm]')
ax[0].text(0.005, 0.85, 'A', transform=ax[0].transAxes, fontsize=16)

ax[1].plot(grf['Time'], (grf['RCOP_Z']-grf['RCOP_Z'].mean())*100, 'b', label='COP')
ax[1].plot(mkr['Time'], (mkr['COG_Z']-mkr['COG_Z'].mean())*100, 'r', label='COG')
ax[1].set_xlim([0, 60])
ax[1].yaxis.set_label_coords(-.05, 0.5)
ax[1].locator_params(axis='y', nbins=6)
ax[1].legend(fontsize=12, loc='upper center', framealpha=.5)
ax[1].set_xlabel('Time [s]')
ax[1].set_ylabel('ml direction [cm]')
ax[1].text(0.005, 0.85, 'B', transform=ax[1].transAxes, fontsize=16)
#plt.suptitle('COP and COG displacements during quiet standing', fontsize=14, y=1)

plt.show()

In [15]:
copcog_plot(grf, mkr)


## COP measurements¶

Let's load all files from the stabilography evaluation and calculate some variables to quantify the stabilography data. The notexbook http://nbviewer.jupyter.org/github/demotu/BMC/blob/master/notebooks/Stabilography.ipynb describes some of the most typical measurements employed to quantify postural sway using the COP data.
Let's calculate here the variables COP area COP mean velocity (resultant), and COP frequency (resultant).

These values are already available at the PDSinfoCOP.txt file in the data set.

In [ ]:
# import the python functions to calculate the COP variables
from psd import psd
from hyperellipsoid import hyperellipsoid

def COPmfreq(COP, freq):
N = COP.shape[0]
fp_ap, mf_ap, fmax_ap, Ptot_ap, F, P_ap = psd(COP[:, 0], fs=freq, window='hanning', nperseg=N/2,
noverlap=N/4, nfft=N/2, detrend='constant', show=False)
fp_ml, mf_ml, fmax_ml, Ptot_ml, F, P_ml = psd(COP[:, 1], fs=freq, window='hanning', nperseg=N/2,
noverlap=N/4, nfft=N/2, detrend='constant', show=False)
mf_res = (mf_ap*np.sum(P_ap) + mf_ml*np.sum(P_ml))/(np.sum(P_ap)+np.sum(P_ml))
#print(mf_ap, mf_ml, mf_res)
return mf_res

from scipy.signal import detrend

fp = FloatProgress(min=0, max=len(PDSinfo.Trial)-1)
display(fp)
freq = 100
PDSinfoCOP = PDSinfo.copy(deep=True)
for i, fname in enumerate(PDSinfoCOP.Trial):
filename = os.path.join(path2, fname + 'grf' + '.txt')
fp.description = '(%s/%s)' %(i+1, len(PDSinfoCOP.Trial))
fp.value = i
COP = grf[['COPNET_X', 'COPNET_Z']].values*100
area, axes, angles, center, R = hyperellipsoid(COP, show=False)
velo = np.sum(np.abs(np.sqrt(np.sum(np.diff(COP, axis=0)**2, axis=1))), axis=0)/(COP.shape[0]/freq)
mfreq = COPmfreq(COP, freq)
PDSinfoCOP.loc[i, 'COParea'] = area
PDSinfoCOP.loc[i, 'COPvelo'] = velo
PDSinfoCOP.loc[i, 'COPmfreq'] = mfreq

PDSinfoCOP.to_csv(os.path.join(path2, 'PDSinfoCOP.txt'), sep='\t', encoding='utf-8', index=False)
#PDSinfoCOP.to_excel(os.path.join(path2, 'PDSinfoCOP.xlsx'), index=False)
print('Data from %d files were processed.' %len(PDSinfoCOP.Trial))


Let's load PDSinfoCOP.txt file from the data set with the values of the variables COP area COP mean velocity (resultant), and COP frequency (resultant):

In [16]:
# GitHub URL:
path2 = 'https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/'
# local directory:
#path2 = r'./data'
fname = os.path.join(path2, 'PDSinfoCOP.txt')
encoding='utf-8')
print(fname, PDSinfoCOP.shape)

https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDSinfoCOP.txt (588, 32)

In [17]:
display(PDSinfoCOP[['Trial', 'Subject', 'Vision', 'Surface', 'COParea', 'COPvelo', 'COPmfreq']].head(15))

Trial Subject Vision Surface COParea COPvelo COPmfreq
0 PDS01OR1 1 Open Rigid 3.09 0.86 0.14
1 PDS01OR2 1 Open Rigid 2.08 0.81 0.17
2 PDS01OR3 1 Open Rigid 5.08 0.90 0.12
3 PDS01OF1 1 Open Foam 9.83 1.92 0.26
4 PDS01OF2 1 Open Foam 8.41 1.62 0.23
5 PDS01OF3 1 Open Foam 8.63 1.70 0.22
6 PDS01CR1 1 Closed Rigid 2.84 1.14 0.26
7 PDS01CR2 1 Closed Rigid 1.71 0.92 0.29
8 PDS01CR3 1 Closed Rigid 2.49 0.93 0.18
9 PDS01CF1 1 Closed Foam 9.17 2.60 0.33
10 PDS01CF2 1 Closed Foam 13.25 2.96 0.34
11 PDS01CF3 1 Closed Foam 11.54 2.42 0.25
12 PDS02OR1 2 Open Rigid 2.46 1.08 0.22
13 PDS02OR2 2 Open Rigid 7.54 1.21 0.12
14 PDS02OR3 2 Open Rigid 3.56 0.97 0.17
In [18]:
PDSinfoCOP = PDSinfoCOP.groupby(['Subject', 'Vision', 'Surface', 'AgeGroup'], as_index=False).median()
display(PDSinfoCOP[['Subject', 'Vision', 'Surface', 'COParea', 'COPvelo', 'COPmfreq']].head(12))
print('%s subjects.' %len(pd.unique(PDSinfo.Subject)))

Subject Vision Surface COParea COPvelo COPmfreq
0 1 Closed Foam 11.54 2.60 0.33
1 1 Closed Rigid 2.49 0.93 0.26
2 1 Open Foam 8.63 1.70 0.23
3 1 Open Rigid 3.09 0.86 0.14
4 2 Closed Foam 21.35 3.96 0.33
5 2 Closed Rigid 5.07 1.65 0.22
6 2 Open Foam 12.75 2.75 0.33
7 2 Open Rigid 3.56 1.08 0.17
8 3 Closed Foam 6.94 2.25 0.35
9 3 Closed Rigid 2.19 1.03 0.25
10 3 Open Foam 10.33 1.80 0.20
11 3 Open Rigid 3.14 0.95 0.16
49 subjects.


Let's plot the variables COP area, COP velocity and COP mean frequency for each subject by age group at the different conditions:

In [19]:
fig, ax = plt.subplots(3, 2, figsize=(12, 8))

sns.stripplot(x='Vision', y='COParea', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[0,0], s=6)
ax[0,0].set_xlabel('')
ax[0,0].set_xticklabels('')
ax[0,0].set_ylabel('COP area (cm$^2$)')
ax[0,0].yaxis.set_label_coords(-.1, .5)
ax[0,0].legend('')
ax[0,0].set_title('Surface = Rigid', fontsize=16)
sns.stripplot(x='Vision', y='COParea', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[0,1], s=6)
ax[0,1].set_xlabel('')
ax[0,1].set_xticklabels('')
ax[0,1].set_ylabel('')
ax[0,1].legend(title='AgeGroup', loc='upper left')
ax[0,1].set_title('Surface = Foam', fontsize=16)

sns.stripplot(x='Vision', y='COPvelo', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[1,0], s=6)
ax[1,0].set_xlabel('')
ax[1,0].set_xticklabels('')
ax[1,0].set_ylabel('COP velocity (cm/s)')
ax[1,0].yaxis.set_label_coords(-.1, .5)
ax[1,0].legend('')
sns.stripplot(x='Vision', y='COPvelo', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[1,1], s=6)
ax[1,1].set_xlabel('')
ax[1,1].set_xticklabels('')
ax[1,1].set_ylabel('')
ax[1,1].legend('')

sns.stripplot(x='Vision', y='COPmfreq', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[2,0], s=6)
ax[2,0].set_ylabel('COP mean freq. (Hz)')
ax[2,0].yaxis.set_label_coords(-.1, .5)
ax[2,0].legend('')

sns.stripplot(x='Vision', y='COPmfreq', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[2,1], s=6)
ax[2,1].set_ylabel('')
ax[2,1].legend('')

plt.tight_layout()
plt.show()


The mean and standard deviation values for the variables COP area, COP velocity and COP mean frequency are:

In [20]:
pd.set_option('precision', 2)
PDSinfoCOP.groupby(['AgeGroup', 'Vision', 'Surface'])['COParea', 'COPvelo', 'COPmfreq'].agg([np.mean, np.std])

Out[20]:
COParea COPvelo COPmfreq
mean std mean std mean std
AgeGroup Vision Surface
Old Closed Foam 22.00 11.57 4.10 1.16 0.39 0.09
Rigid 3.31 1.90 1.32 0.27 0.30 0.08
Open Foam 15.77 8.10 2.79 0.63 0.30 0.08
Rigid 2.98 1.48 1.10 0.19 0.26 0.09
Young Closed Foam 13.04 6.48 2.73 0.51 0.32 0.06
Rigid 2.73 1.27 1.03 0.21 0.22 0.06
Open Foam 10.41 4.86 1.85 0.33 0.21 0.07
Rigid 2.55 1.24 0.90 0.13 0.19 0.05

And the corresponding plots for the mean and standard deviation values for the variables COP area, COP velocity and COP mean frequency are:

In [21]:
fig, ax = plt.subplots(3, 2, figsize=(12, 8))

sns.pointplot(x='Vision', y='COParea', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'], estimator=np.mean, ci=95, kind='point',
split=True, dodge=True, size=4, ax=ax[0,0], s=5)
ax[0,0].set_xlabel('')
ax[0,0].set_xticklabels('')
ax[0,0].set_ylabel('COPNET area (cm$^2$)')
ax[0,0].yaxis.set_label_coords(-.1, .5)
ax[0,0].legend(title='AgeGroup', loc='upper center')
ax[0,0].set_title('Surface = Rigid', fontsize=16)
ax[0,0].text(0.01, 0.87, 'A', transform=ax[0,0].transAxes, fontsize=16)

sns.pointplot(x='Vision', y='COParea', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
split=True, dodge=True, size=4, ax=ax[0,1], s=5)
ax[0,1].set_xlabel('')
ax[0,1].set_xticklabels('')
ax[0,1].set_ylabel('')
ax[0,1].legend('')
ax[0,1].set_title('Surface = Foam', fontsize=16)
ax[0,1].text(0.01, 0.87, 'B', transform=ax[0,1].transAxes, fontsize=16)

sns.pointplot(x='Vision', y='COPvelo', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'],
split=True, dodge=True, size=4, ax=ax[1,0], s=5)
ax[1,0].set_xlabel('')
ax[1,0].set_xticklabels('')
ax[1,0].set_ylabel('COPNET velocity (cm/s)')
ax[1,0].yaxis.set_label_coords(-.1, .5)
ax[1,0].legend('')
ax[1,0].text(0.01, 0.87, 'C', transform=ax[1,0].transAxes, fontsize=16)

sns.pointplot(x='Vision', y='COPvelo', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
split=True, dodge=True, size=4, ax=ax[1,1], s=5)
ax[1,1].set_xlabel('')
ax[1,1].set_xticklabels('')
ax[1,1].set_ylabel('')
ax[1,1].legend('')
ax[1,1].text(0.01, 0.87, 'D', transform=ax[1,1].transAxes, fontsize=16)

sns.pointplot(x='Vision', y='COPmfreq', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Rigid'],
split=True, dodge=True, size=4, ax=ax[2,0], s=5)
ax[2,0].set_ylabel('COPNET mean freq. (Hz)')
ax[2,0].yaxis.set_label_coords(-.1, .5)
ax[2,0].legend('')
ax[2,0].text(0.01, 0.87, 'E', transform=ax[2,0].transAxes, fontsize=16)

sns.pointplot(x='Vision', y='COPmfreq', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoCOP[PDSinfoCOP.Surface=='Foam'],
split=True, dodge=True, size=4, ax=ax[2,1], s=5)
ax[2,1].set_ylabel('')
ax[2,1].legend('')
ax[2,1].text(0.01, 0.87, 'F', transform=ax[2,1].transAxes, fontsize=16)

plt.show()


## Kinematics¶

Let's load and visualize the kinematic data for one trial:

In [22]:
sns.set_style("whitegrid")

fname_mkr = os.path.join(path2, PDSinfo.Trial[150] + 'mkr' + '.txt')
mkrm = mkr.iloc[:, 1::].mean().values
[4,5],                                                            # shoulder
[6,7],[7,8],[8,9],[9,10],[10,11],                                 # spine
[12,13],                                                          # front trunk
[14,15],[14,16],[15,17],[16,17],[14,19],[16,19],[15,18],[17,18],  # pelvis
[20,22],[22,23],[23,20],                                          # right tigh
[26,28],[26,30],[30,31],[31,26],                                  # right leg
[21,24],[24,25],[25,21],                                          # left tigh
[27,29],[27,32],[32,33],[33,27],                                  # left leg
[34,36],[36,37],[37,38],[38,34],                                  # right foot
[35,39],[39,40],[40,41],[41,35]]                                  # left foot

from mpl_toolkits.mplot3d import Axes3D
%matplotlib inline
fig = plt.figure(figsize=(8, 10))
ax.view_init(15, 20)

for s in seg:
ax.plot(mkrm[[s[0]*3,s[1]*3]], mkrm[[s[0]*3+2,s[1]*3+2]], mkrm[[s[0]*3+1,s[1]*3+1]],
linewidth=2, c='b')

ax.scatter(mkrm[0:-3:3], mkrm[2:-1:3], mkrm[1:-2:3], c='r', s=30, depthshade=False)
ax.scatter(mkrm[-3], mkrm[-1], mkrm[-2], c='g', s=60, depthshade=False)

ax.set_xlim3d([np.nanmin(mkrm[0::3])-.4, np.nanmax(mkrm[0::3])+.4])
ax.set_ylim3d([np.nanmin(mkrm[2::3])-.4, np.nanmax(mkrm[2::3])+.4])
ax.set_zlim3d([np.nanmin(mkrm[1::3]), np.nanmax(mkrm[1::3])])
ax.set_xlabel('\n' + 'X [m]', linespacing=2)
ax.set_ylabel('\n' + 'Z [m]', linespacing=2)
ax.set_zlabel('\n' + 'Y [m]', linespacing=2)
ax.invert_yaxis()

plt.show()


## Joint angular displacement at sagital plane¶

Lert's plot the joint angles of a subject during one trial:

In [23]:
sns.set_style("darkgrid")

angz_l = ['LHe_Th_Angle_Z', 'RHe_Th_Angle_Z',
'LTh_Pel_Angle_Z', 'RTh_Pel_Angle_Z',
'LHip_Angle_Z', 'RHip_Angle_Z',
'LKnee_Angle_Z', 'RKnee_Angle_Z',
'LAnkle_Angle_Z', 'RAnkle_Angle_Z']
'L Hip', 'R Hip', 'L Knee', 'R Knee', 'L Ankle', 'R Ankle']

def angz_plot(angles, angz_l):
fig, ax = plt.subplots(5, 2, figsize=(12, 7), sharex=True)
axs = ax.flatten()
color = ['r', 'b']
letters = 'ABCDEFGHIJ'
for i, ax in enumerate(axs):
ax.plot(angles['Time'], angles[[angz_l[i]]], color[i%2], label=label[int(i)])
ax.locator_params(axis='y', nbins=4)
ax.set_xlim([0, 60])
ax.text(0.01, 0.82, letters[i], transform=ax.transAxes, fontsize=18)
ax.legend(fontsize=14, loc='best', framealpha=.5)
if i > 7:
ax.set_xlabel('Time [s]')
if not i%2:
ax.set_ylabel('Angle [$^o$]')
ax.yaxis.set_label_coords(-.1, .5)

#plt.suptitle('Joint angular displacement at sagital plane', fontsize=16, y=1)
plt.tight_layout()
plt.show()

In [24]:
fname_ang = os.path.join(path2, PDSinfo.Trial[150] + 'ang' + '.txt')
print(fname_ang, angles.shape)
angz_plot(angles, angz_l)

https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDS13CR1ang.txt (6000, 74)


We can calculate the range of joint displacement for all subjects to characterize the patterns of joint motion for the different standing conditions. Let's use the standard deviation as measure of the amount of joint displacement instead of maximum minus minimum angles because the former is a more robust measurement.

These values are already available at the PDSinfoANG.txt file in the data set. We can skip the next cell and load this file instead.

In [ ]:
# range calculation for the angles
fp = FloatProgress(min=0, max=len(PDSinfo.Trial)-1)
display(fp)
PDSinfoANG = PDSinfo.copy(deep=True)
for i, fname in enumerate(PDSinfo.Trial):
filename = os.path.join(path2, fname + 'ang' + '.txt')
fp.description = '(%s/%s)' %(i+1, len(PDSinfoANG.Trial))
fp.value = i
ang2 = angs[angz_l].std().values.reshape(5, 2).mean(axis=1)
PDSinfoANG.loc[i, 'He_Th_Angle_Z'] = ang2[0]
PDSinfoANG.loc[i, 'Th_Pel_Angle_Z'] = ang2[1]
PDSinfoANG.loc[i, 'Hip_Angle_Z'] = ang2[2]
PDSinfoANG.loc[i, 'Knee_Angle_Z'] = ang2[3]
PDSinfoANG.loc[i, 'Ankle_Angle_Z'] = ang2[4]

PDSinfoANG.to_csv(os.path.join(path2, 'PDSinfoANG.txt'), sep='\t', encoding='utf-8', index=False)
#print('Data from %d files were processed.' %len(PDSinfoANG.Trial))

In [25]:
fname = os.path.join(path2, 'PDSinfoANG.txt')
PDSinfoANG = PDSinfoANG.groupby(['Subject', 'Vision', 'Surface', 'AgeGroup'], as_index=False).median()
print(fname, PDSinfoANG.shape)

https://raw.githubusercontent.com/demotu/datasets/master/PDS/data/PDSinfoANG.txt (196, 19)

In [26]:
fig, ax = plt.subplots(5, 2, figsize=(12, 8))

sns.stripplot(x='Vision', y='He_Th_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[0,0], s=5)
ax[0,0].set_xlabel('')
ax[0,0].set_xticklabels('')
ax[0,0].yaxis.set_label_coords(-.1, .5)
ax[0,0].legend('')
ax[0,0].set_title('Surface = Rigid', fontsize=16)
ax[0,0].locator_params(axis='y', nbins=4)
ax[0,0].text(0.01, 0.85, 'A', transform=ax[0,0].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='He_Th_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[0,1], s=5)
ax[0,1].set_xlabel('')
ax[0,1].set_xticklabels('')
ax[0,1].set_ylabel('')
ax[0,1].legend(title='AgeGroup', loc='best')
ax[0,1].set_title('Surface = Foam', fontsize=16)
ax[0,1].locator_params(axis='y', nbins=4)
ax[0,1].text(0.01, 0.85, 'B', transform=ax[0,1].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Th_Pel_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[1,0], s=5)
ax[1,0].set_xlabel('')
ax[1,0].set_xticklabels('')
ax[1,0].set_ylabel('Pelvis/Trunk')
ax[1,0].yaxis.set_label_coords(-.1, .5)
ax[1,0].legend('')
ax[1,0].locator_params(axis='y', nbins=4)
ax[1,0].text(0.01, 0.85, 'C', transform=ax[1,0].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Th_Pel_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[1,1], s=5)
ax[1,1].set_xlabel('')
ax[1,1].set_xticklabels('')
ax[1,1].set_ylabel('')
ax[1,1].legend('')
ax[1,1].locator_params(axis='y', nbins=4)
ax[1,1].text(0.01, 0.85, 'D', transform=ax[1,1].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Hip_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[2,0], s=5)
ax[2,0].set_xlabel('')
ax[2,0].set_xticklabels('')
ax[2,0].set_ylabel('Hip')
ax[2,0].yaxis.set_label_coords(-.1, .5)
ax[2,0].legend('')
ax[2,0].locator_params(axis='y', nbins=4)
ax[2,0].text(0.01, 0.85, 'E', transform=ax[2,0].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Hip_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=ax[2,1], s=5)
ax[2,1].set_xlabel('')
ax[2,1].set_xticklabels('')
ax[2,1].set_ylabel('')
ax[2,1].legend('')
ax[2,1].locator_params(axis='y', nbins=4)
ax[2,1].text(0.01, 0.85, 'F', transform=ax[2,1].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Knee_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Rigid'],
jitter=.35, split=True, size=4, ax=ax[3,0], s=5)
ax[3,0].set_xlabel('')
ax[3,0].set_xticklabels('')
ax[3,0].set_ylabel('Knee')
ax[3,0].yaxis.set_label_coords(-.1, .5)
ax[3,0].legend('')
ax[3,0].locator_params(axis='y', nbins=4)
ax[3,0].text(0.01, 0.85, 'G', transform=ax[3,0].transAxes, fontsize=16)

sns.stripplot(x='Vision', y='Knee_Angle_Z', hue='AgeGroup', order=['Open', 'Closed'],
data=PDSinfoANG[PDSinfoANG.Surface=='Foam'],
jitter=.35, split=True, size=4, ax=