forked from DiegoAE/HMMLFM
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmocap_ICM.py
110 lines (90 loc) · 3.4 KB
/
mocap_ICM.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
import GPy
from hmm.continuous.LFMHMMcontinuous import LFMHMMcontinuous
from hmm.continuous.ICMHMMcontinuousMO import ICMHMMcontinuousMO
from matplotlib import pyplot as plt
import numpy as np
seed = np.random.random_integers(10000)
np.random.seed(seed)
print "Using GPy version: ", GPy.__version__
data = GPy.util.datasets.cmu_mocap('05', ['20'], sample_every=1)
print data['info']
Y = data['Y']
nsamples, nfeatures = Y.shape
print "Y's shape ", Y.shape
channel_ids = [9, 16] # ltibia and rtibia
for idx in channel_ids:
plt.plot(np.arange(nsamples), Y[:, idx])
plt.show()
number_hidden_states = 3
outputs = len(channel_ids)
start_t = 0.1
end_t = 5.1
locations_per_segment = 20
n_latent_forces = 1 # TODO: currently not passing this argument to the model.
lfm_hmm = ICMHMMcontinuousMO(outputs, number_hidden_states, locations_per_segment, start_t,
end_t, verbose=True)
number_training_sequences = 1
obs = []
for s in xrange(number_training_sequences):
number_segments = 55 # fixed for now.
c_obs = np.zeros((number_segments, locations_per_segment * outputs))
for output_id in xrange(outputs):
signal = Y[:, channel_ids[output_id]]
idx = 0
for i in xrange(number_segments):
c_obs[i, output_id::outputs] = signal[idx:idx + locations_per_segment]
idx = idx + locations_per_segment - 1
obs.append(c_obs)
lfm_hmm.set_observations(obs)
print "before training"
print lfm_hmm.pi
print lfm_hmm.A
print lfm_hmm.ICMparams
train_flag = False
file = "second-ICM-MOCAP-MO"
if train_flag:
lfm_hmm.train()
lfm_hmm.save_params("/home/diego/tmp/Parameters/MOCAP", file)
else:
lfm_hmm.read_params("/home/diego/tmp/Parameters/MOCAP", file)
print "after training"
print lfm_hmm.pi
print lfm_hmm.A
print lfm_hmm.ICMparams
obs_1 = obs[0]
hidden_states_1 = lfm_hmm._viterbi()[0]
considered_segments = len(obs_1)
number_testing_points = 100
# prediction
last_value = 0
plt.axvline(x=last_value, color='red', linestyle='--')
for i in xrange(considered_segments):
c_hidden_state = hidden_states_1[i]
c_obv = obs_1[i]
# predicting more time steps
t_test = np.linspace(start_t, end_t, number_testing_points)
mean_pred, cov_pred = lfm_hmm.predict(t_test, c_hidden_state, c_obv)
mean_pred = mean_pred.flatten()
cov_pred = np.diag(cov_pred)
current_outputs = np.zeros((number_testing_points, outputs))
current_covariances = np.zeros((number_testing_points, outputs))
# separating the outputs accordingly.
for j in xrange(outputs):
current_outputs[:, j] = mean_pred[j::outputs]
current_covariances[:, j] = cov_pred[j::outputs]
sl = lfm_hmm.sample_locations
for j in xrange(outputs):
plt.scatter(last_value + sl - sl[0], c_obv[j::outputs],
facecolors='none', label=[None, 'observations'][i == 0])
plt.plot(last_value + t_test - t_test[0], current_outputs, color='green',
label=[None, 'predicted mean'][i == 0])
diag_cov = np.diag(cov_pred)
plt.plot(last_value + t_test - t_test[0],
current_outputs - 2 * np.sqrt(current_covariances), 'k--')
plt.plot(last_value + t_test - t_test[0],
current_outputs + 2 * np.sqrt(current_covariances), 'k--')
last_value = last_value + end_t - start_t
plt.axvline(x=last_value, color='red', linestyle='--')
plt.show()
print "Inferred hidden states ", hidden_states_1
print "USED SEED", seed