# code for simulating size tuning curves using the edog model # import import matplotlib.pyplot as plt from operator import itemgetter from edog.tools import* import spatint_utils from importlib import reload import os # reload modules reload(spatint_utils) # load params spatint_utils.plot_params() _, _, optocolor = spatint_utils.get_colors() class SimEdog: def __init__(self, start=1, stop=41): """Init class Parameters ----- start: int starting width for inhibitory feedback kernel stop: int end width for inhibitory feedback kernel width + 1 """ self.start = start self.stop = stop self.step = 1 self.lowbounds = None self.highbounds = None self.sids = None self.rfcsds = None self.smallrateds = None self.largerateds = None self.siss = None self.rfcsss = None self.wavenumber = None self.patch_diameter = None self.w_rc_mix = None def describe(self): """Get parameter boundaries in which experimental results are qualitatively replicated. Creates two dictionaries, one for the low boundaries, one for the high boundaries. Each dict contains one key for each parameter (smallrateds: response to preferred size; largerateds: response to largest stimulus; sids: surround suppression; rfcsds: preferred size) plus one key containing the intersection of all parameter boundaries (all).""" if self.sids is None: # run simulations if this has not been done self.iterate() # init boundary dicts lowbounds = {} highbounds = {} # check for which kernel size CT feedback increases responses to optimal stimulus size smallrateds_increaseis = (np.where(np.array(self.smallrateds) > 0))[0] + self.start lowbounds['smallrateds'] = np.min(smallrateds_increaseis) highbounds['smallrateds'] = np.max(smallrateds_increaseis) print('Responses to optimal stimulus diameter inreases with feedback for kernel ' 'widths between %d and %d degree.' % (lowbounds['smallrateds'], highbounds['smallrateds'])) # check for which kernel size CT feedback decreases responses to largest stimulus size largerateds_decreaseis = (np.where(np.array(self.largerateds) < 0))[0] + self.start lowbounds['largerateds'] = np.min(largerateds_decreaseis) highbounds['largerateds'] = np.max(largerateds_decreaseis) print('Responses to largest stimulus diameter decreases with feedback for kernel ' 'widths between %d and %d degree.' % (lowbounds['largerateds'], highbounds['largerateds'])) # check for which kernel size CT feedback increases surround suppression sids_increaseis = (np.where(np.array(self.sids) > 0))[0] + self.start lowbounds['sids'] = np.min(sids_increaseis) highbounds['sids'] = np.max(sids_increaseis) print('SI increases with feedback for kernel widths between %d and %d degree.' % (lowbounds['sids'], highbounds['sids'])) # check for which kernel size CT feedback decreases receptive field center size rfcsds_decreaseis = (np.where(np.array(self.rfcsds) < 0))[0] + self.start lowbounds['rfcsds'] = np.min(rfcsds_decreaseis) highbounds['rfcsds'] = np.max(rfcsds_decreaseis) print('RFCS decreases with feedback for kernel widths between %d and %d degree.' % (lowbounds['rfcsds'], highbounds['rfcsds'])) # check range for all paramters all_lowi = max(lowbounds.items(), key=itemgetter(1))[0] all_low = lowbounds[all_lowi] all_highi = min(highbounds.items(), key=itemgetter(1))[0] all_high = highbounds[all_highi] lowbounds['all'] = all_low highbounds['all'] = all_high print('Qualitatively matching range between %d and %d degree.' % (lowbounds['all'], highbounds['all'])) # store dicts as attributes self.lowbounds = lowbounds self.highbounds = highbounds def iterate(self): """Run size tuning simulation for range of inhibitory FB kernel widths specified for conditions in which feedback is intact (weight: 1) and feedback is abolished ( weight: 0). For each stimulation store percent change in surround suppression ( sids); preferred size (rfcsds), response to preferred size (smallrateds), and response to the largest stimulus size (largerateds), as well as absolute values for suppression index (siss) and preferred size (rfcsss) in both conditions""" # init lists to store values sids = [] rfcsds = [] smallrateds = [] largerateds = [] siss = [] rfcsss = [] # create vector of inhibitory feedback widths inhib_fbs = np.arange(self.start, self.stop, self.step) # iterate over widths for inhib_fb in inhib_fbs: # run simulation curves = self.simulate(inhib_fb=inhib_fb) # compute target values sid, rfcsd, smallrated, largerated, sis, rfcss = self.characterize(curves=curves) # store values in list sids.append(sid) rfcsds.append(rfcsd) smallrateds.append(smallrated) largerateds.append(largerated) siss.append(sis) rfcsss.append(rfcss) # store as attributes self.sids = sids self.rfcsds = rfcsds self.smallrateds = smallrateds self.largerateds = largerateds self.siss = siss self.rfcsss = rfcsss def simulate(self, inhib_fb=None): """Run the edog model and return the two tuning curves Parameters ---- inhib_fb: int width of inhibitory feedback kernel """ # get parameters parentdir = os. path. dirname(os. getcwd()) filename = parentdir + '/data/params_mouse.yaml' params = parse_parameters(filename) nt, nr, dt, dr = itemgetter("nt", "nr", "dt", "dr")(params["grid"]) k_id, w_id, patch_diameter = itemgetter("k_id", "w_id", "patch_diameter")( params["stimulus"]) A_g, a_g, B_g, b_g = itemgetter("A", "a", "B", "b")(params["ganglion"]) w_rg, A_rg, a_rg = itemgetter("w", "A", "a")(params["relay"]["Krg"]) w_rig, A_rig, a_rig = itemgetter("w", "A", "a")(params["relay"]["Krig"]) w_rc_mix = itemgetter("w")(params["relay"]["Krc_mix"]) A_rc_mix_in, a_rc_mix_in = itemgetter("A", "a")(params["relay"]["Krc_mix"]["Krc_in"]) A_rc_mix_ex, a_rc_mix_ex = itemgetter("A", "a")(params["relay"]["Krc_mix"]["Krc_ex"]) if inhib_fb is not None: # if input to inhib_fb is not none replace default value of inhibitory FB kernel # width a_rc_mix_in = inhib_fb * pq.deg # initiate variables tuning_curve = np.zeros([len(w_rc_mix), len(patch_diameter)]) cen_size = np.zeros(len(w_rc_mix)) supp_index = np.zeros(len(w_rc_mix)) tuning_curve[:] = np.nan cen_size[:] = np.nan supp_index[:] = np.nan # run model for i, w in enumerate(w_rc_mix): network = create_spatial_network(nt=nt, nr=nr, dt=dt, dr=dr, A_g=A_g, a_g=a_g, B_g=B_g, b_g=b_g, w_rg=w_rg, A_rg=A_rg, a_rg=a_rg, w_rig=w_rig, A_rig=A_rig, a_rig=a_rig, w_rc_in=w, A_rc_in=A_rc_mix_in, a_rc_in=a_rc_mix_in, w_rc_ex=w, A_rc_ex=A_rc_mix_ex, a_rc_ex=a_rc_mix_ex) angular_freq = network.integrator.temporal_angular_freqs[int(w_id)] wavenumber = network.integrator.spatial_angular_freqs[int(k_id)] spatiotemporal_tuning = spatiotemporal_size_tuning(network=network, angular_freq=angular_freq, wavenumber=wavenumber, patch_diameter=patch_diameter) # store tuning curve, preferred size and suppression index tuning_curve[i, :] = spatiotemporal_tuning[0, :] # normalize by no feedback condition curves = tuning_curve / tuning_curve[0].max() # store anglular frequencies, FB weighs and patch diameters self.wavenumber = network.integrator.spatial_angular_freqs[int(k_id)] self.patch_diameter = patch_diameter self.w_rc_mix = w_rc_mix return curves def characterize(self, curves): """Compute differences between size tuning curves with and without feedback. Return percent change in surround suppression (sid); preferred size (rfcsd), response to preferred size (smallrated), and response to the largest stimulus size (largerated), as well as absolute values for suppression index (sis) and preferred size (rfcss) Parameters ---- curves: ndarray 2D array containing size tuning curves in both conditions Returns ---- sid: float percent change in suppression index rfcsd: float percent change in preferred size smallrated: float percent change in response to small stimulus largerated: float percent change in response to large stimulus sis: list suppression indices under both condtitions rfcss: list preferred sizes under both conditions """ # get patch_diameters patch_diameter = self.patch_diameter # init lists to store target values sis = [] rfcss = [] smallrates = [] largerates = [] # preferred size is defined by stimulus where an increase in one degree does not # increase response by 0.05%. This is analogous to how visTRN size tuning is analyzed perc_thres = 0.05 # iterate over the two conditions for icurve in range(curves.shape[0]): # get x and y y = curves[1-icurve, :] # start with control condition x = patch_diameter # compute RF center size: check where 1deg in increase of visual stim fails to # increase resp by perc_thres (output 2 to deal with divide by 0 -> will set # first increase to 100%) rel_change = np.divide(y[1:], y[:-1], out=np.ones_like(y[1:])*2, where=y[:-1] != 0) # compute percent change perc_change = (rel_change - 1) * 100 # check where percent change is lower than perc_thres max_resp = np.where(perc_change < perc_thres)[0] # receptive field center size is either first max_resp or last stimulus size if any(max_resp): rfcs = max_resp[0] else: rfcs = x[-1] if icurve == 0: # store rfcs of control curve cont_rfcs = rfcs # compute surround size rfss = int(x[-1].item()) # determine modeled firing rates smallrate_cont = y[cont_rfcs] smallrate = y[rfcs] largerate = y[rfss] # compute suppression index, if si is smaller 0 because of small increase for large # stimuli, set it to 0 si = max(((smallrate - largerate) / smallrate), 0) # store values sis.append(si) rfcss.append(rfcs) smallrates.append(smallrate_cont) largerates.append(largerate) # compute percent change sid = ((sis[0] / sis[1]) - 1) * 100 rfcsd = ((rfcss[0] / rfcss[1]) - 1) * 100 smallrated = ((smallrates[0] / smallrates[1]) - 1) * 100 largerated = ((largerates[0] / largerates[1]) - 1) * 100 return sid, rfcsd, smallrated, largerated, sis, rfcss def plot(self, figsize=(6, 6), inhib_fb=None, quant=False, ax=None): """Plot the two tuning curves for inhibitory feedback kernel width Parameters ---- figsize: tuple figure size (width, height) inhib_fb: int inhibitory feedback kernel width quant: bool if true print info regarding ax: mpl axis axis for plotting """ # simulate data curves = self.simulate(inhib_fb=inhib_fb) # get xvalues patch_diameter = self.patch_diameter # set colors colors = [optocolor, 'k'] if ax is None: # init fig if doesnt exist fig, ax = plt.subplots(figsize=figsize) # plot curves for curve, color in zip(curves, colors): ax.plot(patch_diameter, curve, '-', color=color) # layout ax.set_ylabel("Normalized response") ax.set_xlabel("Diameter ($\degree$)") ax.set_xticks((0, 25, 50, 75)) ax.set_xticklabels((0, 25, 50, 75)) ax.set_xlim([0, 75]) ax.set_ylim([0, 1.5]) ax.spines['bottom'].set_bounds(0, 75) if quant: # plot information about differences # init variables title = "Mixed feedback" wavenumber = self.wavenumber # get quant measures sid, rfcsd, smallrated, largerated, sis, rfcss = self.characterize(curves=curves) # insert model results ax.text(50, 0.1, 'opto si = %0.2f' % sis[1]) ax.text(50, 0.14, 'control si = %0.2f' % sis[0]) ax.text(50, 0.18, 'perc change si = %0.2f' % sid) ax.text(50, 0.22, 'opto rfcs = %0.2f' % rfcss[1]) ax.text(50, 0.26, 'control rfcs = %0.2f' % rfcss[0]) ax.text(50, 0.3, 'perc change rfcs = %0.2f' % rfcsd) ax.text(50, 0.34, 'perc change smallrate = %0.2f' % smallrated) ax.text(50, 0.38, 'perc change largerate = %0.2f' % largerated) ax.text(50, 0.42, 'model') # add legends and titles ax.set_title(title) fig.text(0.39, 0.97, "Patch grating ({})".format(round(wavenumber.item(), 2)), fontsize=12) # plot preferred size for rfcs, color in zip(rfcss[::-1], colors): ax.axvline(rfcs, color=color)