published on eLife

Author: Eloy Parra Barrero

analyze/batch_analyses.py

@@ -12,10 +12,10 @@
 
 def analyze(session, group_name):
 
-    if group_name == experimental_group_name:
-        tracking = initialize((Tracking,), session, group_name)['Tracking']
-        run_tracking(tracking)
-        del tracking
+    # if group_name == experimental_group_name:
+    #     tracking = initialize((Tracking,), session, group_name)['Tracking']
+    #     run_tracking(tracking)
+    #     del tracking
 
     # place field sizes and skews
     firing_fields = initialize((FiringFields,), session, group_name)['FiringFields']

analyze/batch_config.py

@@ -4,20 +4,20 @@
 
 # select which sessions to analyze
 sessions = [
-    "ec013.156",
-    "ec013.206",
-    "ec013.374",
-    "ec013.386",
-    "ec013.388",
-    "ec013.395",
-    "ec013.412",
-    "ec013.413",
-    "ec013.440",
-    "ec013.454",
-    "ec013.466",
-    "ec013.469",
-    "ec013.502",
-    "ec013.531",
+    # "ec013.156",
+    # "ec013.206",
+    # "ec013.374",
+    # "ec013.386",
+    # "ec013.388",
+    # "ec013.395",
+    # "ec013.412",
+    # "ec013.413",
+    # "ec013.440",
+    # "ec013.454",
+    # "ec013.466",
+    # "ec013.469",
+    # "ec013.502",
+    # "ec013.531",
     # only these three sessions were used in the models for ec013:
     "ec013.555",
     "ec013.556",
@@ -56,11 +56,15 @@
 
 # select what to analyze: experimental data and/or models
 group_names = [
-    experimental_group_name,  # experimental data
+    # experimental_group_name,  # experimental data
     # 'Time',  # temporal sweep model
     # 'Position',  # spatial sweep model
-    # 'SpeedDep57',  # behavior-dependent sweep model
-    # 'VariableNoise'  # variable noise model
+    # 'SpeedDepVanilla'  # behavior-dependent sweep model
+
+    # 'VariableNoise'  # variable noise model with extra variance
+    'VariableNoiseVanilla'  #
+
+    # 'SpeedDep57',  # behavior-dependent sweep model with extra variance
 ]
 
 figures_path = f"{figures_path}/ALL"

analyze/batch_plots.py

@@ -58,7 +58,7 @@
     font_size = 7
     line_width = 1
 
-    thin_line_width = 0.5
+    thin_line_width = 0.75  # 0.5
 
     plt.rcParams.update({'font.family': 'sans-serif', 'font.sans-serif': 'Arial',
                          'font.size': font_size, 'mathtext.default': 'regular', 'axes.linewidth': thin_line_width,

analyze/config.py

@@ -44,10 +44,12 @@
 fig, ax = plt.subplots(1, 2, sharey='row', gridspec_kw={'width_ratios': [1, 0.2]}, figsize=(6.5*cm, 6*cm))
 ax[0].set_visible(False)
 # ax[1].invert_xaxis()
-ax[1].axis('off')
-ax[1].set_xlim([-0.25, 1.05])
-rect = plt.Rectangle((-0.2, bottom), 0.2, top - bottom, facecolor='whitesmoke')
-ax[1].add_patch(rect)
+# ax[1].axis('off')
+for spine in ['top', 'bottom', 'right']:
+    ax[1].spines[spine].set_visible(False)
+ax[1].set_xlim([0, 1.2])
+# rect = plt.Rectangle((-0.2, bottom), 0.2, top - bottom, facecolor='whitesmoke')
+# ax[1].add_patch(rect)
 ax[1].set_ylim([bottom, top])
 ax[1].set_xlabel("Firing\nrate")
 
@@ -122,8 +124,8 @@ def add_inset(left, bottom, width, height, x_label=None, y_label=None):
 
 # FIRING FIELDS FROM THE MODEL
 super_group_name = "BDSweepDemo"
-group_name = "BDSweepDemo"
-t_sim = 300
+group_name = "BDSweepDemoRevision"
+t_sim = 1400
 
 # LFP
 p = general_parameters['LFP']
@@ -158,9 +160,9 @@ def add_inset(left, bottom, width, height, x_label=None, y_label=None):
 
 # Model Spikes
 p = general_parameters[f'SpeedSpikes|{group_name}']
-num_cells = 6
+num_cells = 2
 # field_centers = [c - bottom + 100 for c in [138, 150, 162, 238, 250, 262]]
-field_centers = [138, 150, 162, 238, 250, 262]
+field_centers = [150, 250]
 
 model_spikes = SpeedSpikes(super_group_name, group_name, "Spikes", lfp, tracking, num_cells, p["ds"], p["dt"],
                            p["phase_range"], p["phase_current"], p["firing_rate_0"], p["firing_rate_slope"],
@@ -173,13 +175,31 @@ def add_inset(left, bottom, width, height, x_label=None, y_label=None):
 p = general_parameters['FiringFields']
 firing_fields = FiringFields(super_group_name, group_name, "FiringFields", model_spikes.spikes, tracking,
                              p['firing_rate_sigma'], p['consecutive_nans_max'], save_figures=False)
+firing_fields.find_fields_candidates(p['min_spikes'], p['min_peak_firing_rate'], p['firing_rate_threshold'],
+                                     p['peak_prominence_threshold'])
+
+path_arrow_color = "#D65F5F"
+field_arrow_color = "#D65F5F"
 
 for field_num, firing_rates in enumerate(firing_fields.smooth_rate_maps):
     normalized_rates = firing_rates[0] / np.nanmax(firing_rates[0])
-    ax[1].plot(normalized_rates, firing_fields.positions, color=colors[field_num])
-
-ax[1].axvline(-0.2, color='gray')
-ax[1].axvline(0, color='gray')
+    ax[1].plot(normalized_rates, firing_fields.positions, color='C0')
+    field_bounds = firing_fields.cand_bounds[field_num*2]
+    ax[1].axhline(field_bounds[0], color='C7', linestyle='dotted', zorder=0)
+    ax[1].axhline(field_bounds[1], color='C7', linestyle='dotted', zorder=0)
+
+    arrow_x = 1.1
+    ax[1].annotate('', (arrow_x, field_bounds[0]), xytext=(arrow_x, field_bounds[1]),
+                   arrowprops=dict(arrowstyle='<|-|>', shrinkA=0, shrinkB=0,
+                                   facecolor=field_arrow_color, edgecolor=field_arrow_color))
+    ax[1].annotate(f"{field_bounds[1] - field_bounds[0]:.0f}",
+                   (1.2, (field_bounds[0] + field_bounds[1]) / 2),
+                   color=field_arrow_color, verticalalignment='center')
+
+# ax[1].axvline(-0.2, color='gray')
+# ax[1].axvline(0, color='gray')
+ax[1].set_xticks([])
+ax[1].set_yticks([125, 150, 175, 200, 225, 250, 275])
 
 # fig.tight_layout(pad=2, w_pad=5)
 fig.savefig(f"figures/BD_sweep {datetime.now().strftime('%Y-%m-%d %H:%M:%S')}.pdf", bbox_inches='tight')

analyze/demo_BD_sweep.py

@@ -200,6 +200,29 @@
     "exponential_factor": -0.025
   },
 
+    "VariableSpikes|VariableNoiseVanilla":
+  {
+    "num_cells": 20,
+    "ds": 1.0,
+    "dt": 0.002,
+    "phase_range": [0, 360],
+    "phase_current": 180,
+    "firing_rate_0": 15,
+    "firing_rate_slope": 0.2,
+    "theta_modulation": 0.35,
+    "fields_from": "characteristic_speed",
+    "theta_d_slope": 0,
+    "theta_d_offset": 30,
+    "theta_d_min": 0,
+    "theta_d_add_sigma": 0,
+    "theta_d_mul_sigma": 0,
+    "size_to_theta_d": 0.23,
+    "size_min": 3,
+    "size_add_sigma": 0,
+    "phase_sigma_0": 120,
+    "exponential_factor": -0.025
+  },
+
   "SpeedSpikes|SpeedDep57":
   {
     "num_cells": 20,
@@ -218,6 +241,24 @@
     "shift_sigma": 0
   },
 
+    "SpeedSpikes|SpeedDepVanilla":
+  {
+    "num_cells": 20,
+    "ds": 1.0,
+    "dt": 0.002,
+    "phase_range": [0, 360],
+    "phase_current": 180,
+    "firing_rate_0": 15,
+    "firing_rate_slope": 0.2,
+    "theta_modulation": 0.35,
+    "theta_time": 0.57,
+    "multiplicative_sigma": 0.0,
+    "additive_sigma": 0.0,
+    "size_to_theta_d": 0.3,
+    "size_min": 4,
+    "shift_sigma": 0
+  },
+
   "SpeedSpikes|BDSweepDemo":
   {
     "num_cells": 6,
@@ -234,5 +275,23 @@
     "size_to_theta_d": 0.25,
     "size_min": 3,
     "shift_sigma": 0
+  },
+
+  "SpeedSpikes|BDSweepDemoRevision":
+  {
+    "num_cells": 6,
+    "ds": 1.0,
+    "dt": 0.002,
+    "phase_range": [0, 360],
+    "phase_current": 180,
+    "firing_rate_0": 15,
+    "firing_rate_slope": 0.2,
+    "theta_modulation": 0.35,
+    "theta_time": 0.5,
+    "multiplicative_sigma": 0,
+    "additive_sigma": 0,
+    "size_to_theta_d": 0,
+    "size_min": 6,
+    "shift_sigma": 0
   }
 }

analyze/fields/Demo.Position.pkl

@@ -7,7 +7,7 @@ def run_cell_coordination(cell_coordination):
     print("Analyzing cell coordination...")
     p = general_parameters['CellCoordination']
 
-    cell_coordination.full_histograms()
+    # cell_coordination.full_histograms()
     cell_coordination.coordination_by_speed(speed_groups, num_phase_bins=p['num_phase_bins'],
                                             sample_size=p['sample_size'])
 

analyze/fields/Demo.PositionDemo.pkl

@@ -31,9 +31,8 @@ def run_path_lengths(path_lengths):
 
     path_lengths.single_cycles(min_peak_prob=p['min_peak_prob'], min_phase_coverage=p['min_phase_coverage'],
                                min_phase_extent=p['min_phase_extent'], radon_fit_params=p['radon_fit_params'],
-                               max_cycles_to_plot=128,
-                               # max_cycles_to_plot=0,
-                               from_run_types=(1, 1), cycles_per_figure=(8, 8), cycles_fig_size=(9.7*cm, 9*cm))
+                               max_cycles_to_plot=128, from_run_types=(1, 1), cycles_per_figure=(8, 8),
+                               cycles_fig_size=(9.7*cm, 9*cm))
 
     # if path_lengths.group_name == experimental_group_name:
     #     # print("NO displacement compensation:")

analyze/fields/Demo.Time.pkl

@@ -7,17 +7,15 @@ def run_firing_fields(firing_fields):
     print("Calculating place fields...")
     p = {**general_parameters['ALL'], **general_parameters['FiringFields']}
 
-    # firing_fields.plot_traces(field_nums=(0, 1, 2))
+    # # firing_fields.plot_traces(field_nums=(0, 1, 2))
     firing_fields.plot_heatmap()
 
     firing_fields.field_sizes_by_speed(speed_groups=speed_groups, min_peak_firing_rate=p['min_peak_firing_rate'],
                                        threshold=p['firing_rate_threshold'],
                                        peak_prominence_threshold=p['peak_prominence_threshold'],
-                                       min_occupancy=p['min_occupancy'],
-                                       min_spread=p['min_spread'],
-                                       plot_fields=False,
+                                       min_occupancy=p['min_occupancy'], min_spread=p['min_spread'],
                                        # plot_fields=True, field_nums=(22, 9, 10), fields_per_plot=3,
-                                       # fig_size=(6*cm, 10*cm), constrained_layout=True
+                                       fig_size=(6.5*cm, 10*cm), constrained_layout=True
                                        )
 
     firing_fields.field_sizes_vs_stuff(plot=True)

analyze/fields/Demo.TimeDemo.pkl

@@ -22,11 +22,8 @@ def run_phase_vs_position(phase_vs_position):
     print(0)
     phase_vs_position.pool_all(full_speed_groups=speed_groups, fit_type=p['fit_type'], min_spikes=p['pooled_min_spikes'],
                                pool_by_pass_speed=False, spike_speed_threshold=False, min_occupancy=p['min_occupancy'],
-                               min_spread=p['min_spread'], plot_fits=True, plot_occupancy=True,
-                               fig_size=(12, 9),
+                               min_spread=p['min_spread'], plot_fits=True, plot_occupancy=True, fig_size=(18*cm, 13*cm),
                                # fields_per_plot=3, field_nums=(8, 20, 22), fig_size=(7.8*cm, 10.1*cm), constrained_layout=0
-                               # fields_per_plot=4, field_nums=(1, 6, 11, 16), fig_size=(8.5*cm, 10.1*cm), constrained_layout=0
-                               # fields_per_plot=4, field_nums=(5, 22, 24, 26), fig_size=(8.5 * cm, 10.1 * cm), constrained_layout=0
                                )
 
     # single passes for one field
@@ -50,7 +47,7 @@ def run_phase_vs_position(phase_vs_position):
 
 if __name__ == '__main__':
     # session = "ec013.156"
-    session = "ec013.454"  # fields 4 and 16 nice for single passes
+    session = "ec013.454"  # fields 4 and 16 nice for single passes; also nice for comparing to variable noise model
     # session = "ec013.395"
     # session = "ec014.468"  # fields 8, 20, 22 nice for pooled at different speed bins
     # session = "ec014.639"
@@ -60,10 +57,11 @@ def run_phase_vs_position(phase_vs_position):
     # session = "Buddy.06272013"
     # session = "Cicero.09172014"
 
-    group_name = experimental_group_name
+    # group_name = experimental_group_name
     # group_name = "SpeedDepSharp"
     # group_name = "Time"
     # group_name = "VariableNoise3"
+    group_name = "VariableNoiseVanilla"
 
     # with open(f"fields/{session}.{group_name}.pkl", 'rb') as fields_file:
     #     fields = pickle.load(fields_file)['complete_fields']

analyze/fields/ec013.454.VariableNoise3.pkl

@@ -20,14 +20,15 @@ def run_tracking(tracking):
 
 if __name__ == '__main__':
 
+    session = "ec013.156"
     # session = "ec013.388"
     # session = "ec013.395"
     # session = "ec013.412"
     # session = "ec014.468"
     # session = "ec014.639"
     # session = "ec016.233"
 
-    session = "vvp01.1815"
+    # session = "vvp01.1815"
     # session = "gor01.1216"
 
     # session = "Achilles.10252013_alt"

analyze/fields/ec013.469.EXPERIMENTAL.pkl

@@ -9,8 +9,8 @@
 # session = "ec014.468"
 # session = "ec014.639"
 # session = "ec016.233"
-session = "Buddy.06272013"
-# session = "Achilles.10252013"
+# session = "Buddy.06272013"
+session = "Achilles.10252013"
 # session = "Gatsby.08022013"
 
 p = general_parameters['LFP']
@@ -36,8 +36,8 @@
 lfp.calculate_phase_locking(data_path, s['dataset'], s['session_set'], s['session'], s['discarded_intervals'], plot=True)
 lfp.comp_instantaneous_frequency()
 
-lfp.plot(time_interval=(0, 500), channels=[0])
+# lfp.plot(time_interval=(0, 500), channels=[0])
 # lfp.plot(time_interval=(17079.2, 17300), channels=[24])  # Gatsby
-# lfp.plot(time_interval=(18080, 182000), channels=[0])  # Achilles
+lfp.plot(time_interval=(18080, 182000), channels=[0])  # Achilles
 
 plt.show()