-
Notifications
You must be signed in to change notification settings - Fork 1
/
compute_gradients_Linear.py
108 lines (71 loc) · 2.91 KB
/
compute_gradients_Linear.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
import numpy as np
import sys
#from numpy.lib.function_base import gradient
sys.path.append('code_compare')
import pdb
def compute_gradients(z_data, A, t):
# COMPUTE_GRADIENTS gets gradients for cost at time t
# w.r.t. parameters at subnetwork i
# INPUT VARIABLES
# t: integer, is the time instant
# i: integer, is the subnetwork index
# z_data N x T array, containing the given z data
# N: number of sensors
# T: number of time instants
N, T = z_data.shape
# A: N x N x P array, containing VAR parameter matrices
N,N,P = A.shape
# alpha: N x M array, alpha parameters
# M: number of units in subnetwork
# w: N x M array, w parameters
# k: N x M array, k parameters
# b: N array, b parameters
# gradient of g w.r.t theta #!LEq(17)
# Function definitions
# Forward pass (function evaluation)
# (You can use your functions f, g, etc)
#!LEq(7) from the paper
tilde_y_tm = np.zeros((N, P+1)) #!LEq(7a) from the paper # here i_prime is used instead of i
for i_prime in range(N):
for p in range(1,P+1):
assert t-p >= 0
z_i_tmp = z_data[i_prime, t-p]
tilde_y_tm[i_prime, p] = z_i_tmp
hat_y_t = np.zeros((N)) #!LEq(7b)
for i_prime in range(N):
for p in range(1,P+1):
for j in range(N):
hat_y_t[i_prime] = hat_y_t[i_prime] + A[i_prime,j,p-1]*tilde_y_tm[j,p]
hat_z_t = np.zeros((N,T)) #!LEq(7c)
for i_prime in range(N):
hat_z_t[i_prime,t] = hat_y_t[i_prime]
#pdb.set_trace()
# computing cost #!LEq(7d) and vector S #!LEq(8)
cost_i = [0]*T
cost_i_test = [0]*T
if (t < int(T*0.8)):
cost_u = 0
cost_l = 0
S = 2*( hat_z_t[:,t] - z_data[:,t])
for i_prime in range(N):
cost_u = cost_u + np.square(S[i_prime]/2)
for i_prime in range(N):
cost_l = cost_l + np.square(z_data[i_prime,t])
cost_i[t] = cost_u/cost_l
else:
cost_u = 0
cost_l = 0
S = 2*( hat_z_t[:,t] - z_data[:,t]) #make the modification here for the train and test error
for i_prime in range(N): #do signal reconstruciton if possible for sensor one
cost_u = cost_u + np.square(S[i_prime]/2)
for i_prime in range(N):
cost_l = cost_l + np.square(z_data[i_prime,t])
cost_i_test[t] = cost_u/cost_l
dC_dA = np.zeros((N,N,P))
if (t < int(T*0.8)):
for i in range(N):
#dC_dA = np.zeros((N,N,P)) #this initialization was there previously before declaring globally
for j in range(N):
for p in range(1,P+1):
dC_dA[i,j,p-1] = S[i]*tilde_y_tm[j, p]
return dC_dA,cost_i[t],cost_i_test[t],hat_z_t[:,t]