Skip to content

Commit

Permalink
Merge pull request #517 from PrincetonUniversity/fix/gilzenrat/match-…
Browse files Browse the repository at this point in the history
…matlab

Fix/gilzenrat/match matlab
  • Loading branch information
KristenManning authored Nov 8, 2017
2 parents 0bfbd34 + 2211ad7 commit 22ecfbc
Show file tree
Hide file tree
Showing 3 changed files with 76 additions and 44 deletions.
56 changes: 31 additions & 25 deletions Scripts/Models/GilzenratModel.py
Original file line number Diff line number Diff line change
Expand Up @@ -19,17 +19,21 @@

# --------------------------------- Global Variables ----------------------------------------

# Mode ("coherence")
C = 0.95
high_C = True
if high_C:
C = 0.95 # Mode ("coherence")
initial_h_of_v = 0.07
initial_w=0.14
else:
C = 0.55
initial_h_of_v = 0.2
initial_w = 0.2

# Uncorrelated Activity
d = 0.5

# Initial values
initial_h_of_v = 0.07
# initial_h_of_v = 0.07
# get v from h of v
initial_v = (initial_h_of_v - (1-C)*d)/C
# initial_w = 0.14
initial_w=0.14

# g(t) = G + k*w(t)

Expand All @@ -40,11 +44,11 @@

# numerical integration
time_step_size = 0.02
# number_of_trials = int(20/time_step_size)
number_of_trials = 1
number_of_trials = int(20/time_step_size)
# number_of_trials = 1

# noise
standard_deviation = 0.22*(time_step_size**0.5)
standard_deviation = 0.22

# --------------------------------------------------------------------------------------------

Expand All @@ -58,8 +62,8 @@

# Implement self-excitatory (auto) and mutually inhibitory (hetero) connections within the decision layer
decision_layer = GilzenratTransferMechanism(size=2,
initial_value=np.array([[1,0]]),
matrix=np.matrix([[1,0],[0,-1]]),
initial_value=np.array([[0.0, 0.0]]),
matrix=np.matrix([[1,-1],[-1,1]]),
#auto=1.0,
#hetero=-1.0,
time_step_size=time_step_size,
Expand All @@ -74,15 +78,15 @@
#To do Markus: specify recurrent self-connrection weight for response unit to 2.00
response = GilzenratTransferMechanism(size=1,
initial_value=np.array([[2.0]]),
matrix=np.matrix([[0.5]]),
function=Logistic(bias=2),
matrix=np.matrix([[2.0]]),
function=Logistic(bias=2.0),
time_step_size=time_step_size,
noise=NormalDist(mean=0.0,standard_dev=standard_deviation).function,
name='RESPONSE')

# Implement response layer with input_state for ObjectiveMechanism that has a single value
# and a MappingProjection to it that zeros the contribution of the decision unit in the decision layer
LC = LCControlMechanism(
LC = LCControlMechanism(integration_method="EULER",
time_step_size_FHN=time_step_size, # integrating step size
mode_FHN=C, # coherence: set to either .95 or .55
uncorrelated_activity_FHN=d, # Baseline level of intrinsic, uncorrelated LC activity
Expand All @@ -103,7 +107,7 @@
threshold_FHN=0.5, #Parameter describing shape of the FitzHugh–Nagumo cubic nullcline for the fast excitation variable v
objective_mechanism=ObjectiveMechanism(
function=Linear,
# monitored_output_states=[(decision_layer, None, None, np.array([[0.3],[0.0]]))],
monitored_output_states=[(decision_layer, None, None, np.array([[0.335],[0.0]]))],
# monitored_output_states=[{PROJECTION_TYPE: MappingProjection,
# SENDER: decision_layer,
# MATRIX: np.array([[0.3],[0.0]])}],
Expand All @@ -114,7 +118,6 @@

for signal in LC._control_signals:
signal._intensity = k*initial_w + G

# ELICITS WARNING:
decision_process = Process(pathway=[input_layer,
input_weights,
Expand All @@ -123,14 +126,17 @@
response],
name='DECISION PROCESS')


lc_process = Process(pathway=[decision_layer,
# CAUSES ERROR:
# np.array([[1,0],[0,0]]),
LC],
name='LC PROCESS')

task = System(processes=[decision_process, lc_process])
#
# lc_process = Process(pathway=[decision_layer,
# # CAUSES ERROR:
# # np.array([[1,0],[0,0]]),
# LC],
# name='LC PROCESS')

task = System(processes=[decision_process,
# lc_process
]
)

# stimulus
stim_list_dict = {input_layer: np.repeat(np.array([[0,0],[1,0]]),10/time_step_size,axis=0)}
Expand Down
59 changes: 41 additions & 18 deletions psyneulink/components/functions/function.py
Original file line number Diff line number Diff line change
Expand Up @@ -4201,7 +4201,7 @@ def function(self,

# Compute function based on integration_type param
# Gilzenrat: previous_value + (-previous_value + variable)*self.time_step_size + noise --> rate = -1
value = previous_value + (rate*previous_value + new_value)*self.time_step_size + noise
value = previous_value + (rate*previous_value + new_value)*self.time_step_size + noise*(self.time_step_size**0.5)

adjusted_value = value + offset
# If this NOT an initialization run, update the old value
Expand Down Expand Up @@ -5549,6 +5549,7 @@ def __init__(self,
time_constant_w=12.5,
mode=1.0,
uncorrelated_activity=0.0,
integration_method="RK4",
params: tc.optional(dict)=None,
owner=None,
prefs: is_pref_set = None,
Expand All @@ -5575,6 +5576,7 @@ def __init__(self,
threshold=threshold,
mode=mode,
uncorrelated_activity=uncorrelated_activity,
integration_method=integration_method,
time_constant_w=time_constant_w,
params=params)

Expand All @@ -5590,7 +5592,28 @@ def __init__(self,

self.auto_dependent = True

def _validate_params(self, request_set, target_set=None, context=None):
super()._validate_params(request_set=request_set,
target_set=target_set,
context=context)
if self.integration_method not in {"RK4", "EULER"}:
raise FunctionError("Invalid integration method ({}) selected for {}".
format(self.integration_method, self.name))

def _euler_FHN(self, previous_value_v, previous_value_w, previous_time, slope_v, slope_w, time_step_size):

slope_v_approx = slope_v(previous_time,
previous_value_v,
previous_value_w)

slope_w_approx = slope_w(previous_time,
previous_value_w,
previous_value_v)

new_v = previous_value_v + time_step_size*slope_v_approx
new_w = previous_value_w + time_step_size*slope_w_approx

return new_v, new_w

def _runge_kutta_4_FHN(self, previous_value_v, previous_value_w, previous_time, slope_v, slope_w, time_step_size):

Expand Down Expand Up @@ -5705,23 +5728,23 @@ def dw_dt(time, w, v):
# val = (v - 0.5*w)

return val

# new_v = self._runge_kutta_4(previous_time=self.previous_t,
# previous_value=self.previous_v,
# slope=dv_dt,
# time_step_size=self.time_step_size)*self.scale + self.offset
#
# new_w = self._runge_kutta_4(previous_time=self.previous_t,
# previous_value=self.previous_w,
# slope=dw_dt,
# time_step_size=self.time_step_size)*self.scale + self.offset

approximate_values = self._runge_kutta_4_FHN(self.previous_v,
self.previous_w,
self.previous_t,
dv_dt,
dw_dt,
self.time_step_size)
if self.integration_method == "RK4":
approximate_values = self._runge_kutta_4_FHN(self.previous_v,
self.previous_w,
self.previous_t,
dv_dt,
dw_dt,
self.time_step_size)
elif self.integration_method == "EULER":
approximate_values = self._euler_FHN(self.previous_v,
self.previous_w,
self.previous_t,
dv_dt,
dw_dt,
self.time_step_size)
else:
raise FunctionError("Invalid integration method ({}) selected for {}".
format(self.integration_method, self.name))

if not context or INITIALIZING not in context:
self.previous_v = approximate_values[0]
Expand Down
5 changes: 4 additions & 1 deletion psyneulink/library/subsystems/agt/lccontrolmechanism.py
Original file line number Diff line number Diff line change
Expand Up @@ -503,6 +503,7 @@ def __init__(self,
# modulated_mechanisms:tc.optional(tc.any(list,str)) = None,
modulated_mechanisms=None,
modulation:tc.optional(_is_modulation_param)=ModulationParam.MULTIPLICATIVE,
integration_method="RK4",
initial_w_FHN=0.0,
initial_v_FHN=0.0,
time_step_size_FHN=0.1,
Expand Down Expand Up @@ -531,6 +532,7 @@ def __init__(self,
objective_mechanism=objective_mechanism,
modulated_mechanisms=modulated_mechanisms,
modulation=modulation,
integration_method=integration_method,
initial_v_FHN=initial_v_FHN,
initial_w_FHN=initial_w_FHN,
time_step_size_FHN=time_step_size_FHN,
Expand All @@ -553,7 +555,8 @@ def __init__(self,

super().__init__(system=system,
objective_mechanism=objective_mechanism,
function=FHNIntegrator(initial_v=initial_v_FHN,
function=FHNIntegrator( integration_method=integration_method,
initial_v=initial_v_FHN,
initial_w=initial_w_FHN,
time_step_size=time_step_size_FHN,
t_0=t_0_FHN,
Expand Down

0 comments on commit 22ecfbc

Please sign in to comment.