1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798991001011021031041051061071081091101111121131141151161171181191201211221231241251261271281291301311321331341351361371381391401411421431441451461471481491501511521531541551561571581591601611621631641651661671681691701711721731741751761771781791801811821831841851861871881891901911921931941951961971981992002012022032042052062072082092102112122132142152162172182192202212222232242252262272282292302312322332342352362372382392402412422432442452462472482492502512522532542552562572582592602612622632642652662672682692702712722732742752762772782792802812822832842852862872882892902912922932942952962972982993003013023033043053063073083093103113123133143153163173183193203213223233243253263273283293303313323333343353363373383393403413423433443453463473483493503513523533543553563573583593603613623633643653663673683693703713723733743753763773783793803813823833843853863873883893903913923933943953963973983994004014024034044054064074084094104114124134144154164174184194204214224234244254264274284294304314324334344354364374384394404414424434444454464474484494504514524534544554564574584594604614624634644654664674684694704714724734744754764774784794804814824834844854864874884894904914924934944954964974984995005015025035045055065075085095105115125135145155165175185195205215225235245255265275285295305315325335345355365375385395405415425435445455465475485495505515525535545555565575585595605615625635645655665675685695705715725735745755765775785795805815825835845855865875885895905915925935945955965975985996006016026036046056066076086096106116126136146156166176186196206216226236246256266276286296306316326336346356366376386396406416426436446456466476486496506516526536546556566576586596606616626636646656666676686696706716726736746756766776786796806816826836846856866876886896906916926936946956966976986997007017027037047057067077087097107117127137147157167177187197207217227237247257267277287297307317327337347357367377387397407417427437447457467477487497507517527537547557567577587597607617627637647657667677687697707717727737747757767777787797807817827837847857867877887897907917927937947957967977987998008018028038048058068078088098108118128138148158168178188198208218228238248258268278288298308318328338348358368378388398408418428438448458468478488498508518528538548558568578588598608618628638648658668678688698708718728738748758768778788798808818828838848858868878888898908918928938948958968978988999009019029039049059069079089099109119129139149159169179189199209219229239249259269279289299309319329339349359369379389399409419429439449459469479489499509519529539549559569579589599609619629639649659669679689699709719729739749759769779789799809819829839849859869879889899909919929939949959969979989991000100110021003100410051006100710081009101010111012101310141015101610171018101910201021102210231024102510261027102810291030103110321033103410351036103710381039104010411042104310441045104610471048104910501051105210531054105510561057105810591060106110621063106410651066106710681069107010711072107310741075107610771078107910801081108210831084108510861087108810891090109110921093109410951096109710981099110011011102110311041105110611071108110911101111111211131114111511161117111811191120112111221123112411251126112711281129113011311132113311341135113611371138113911401141114211431144114511461147114811491150115111521153115411551156115711581159116011611162116311641165116611671168116911701171117211731174117511761177117811791180118111821183118411851186118711881189119011911192119311941195119611971198119912001201120212031204120512061207120812091210121112121213121412151216121712181219122012211222122312241225122612271228122912301231123212331234123512361237123812391240124112421243124412451246124712481249125012511252125312541255125612571258125912601261126212631264126512661267126812691270127112721273127412751276127712781279128012811282128312841285128612871288128912901291129212931294129512961297129812991300130113021303130413051306130713081309131013111312131313141315131613171318131913201321132213231324132513261327 |
- import itertools
- import os
- import matplotlib.legend as mlegend
- import matplotlib.pyplot as plt
- import numpy as np
- import pandas as pd
- from matplotlib.patches import Ellipse
- from matplotlib.patches import Polygon
- from matplotlib.patches import Rectangle
- from mpl_toolkits.axes_grid1 import ImageGrid
- from pypet import Trajectory
- from pypet.brian2 import Brian2MonitorResult
- from scripts.spatial_maps.spatial_network_layout import Interneuron, get_position_mesh
- from scripts.spatial_maps.correlation_length_fit.correlation_length_fit import \
- correlation_length_fit_dict
- from scripts.model_figure.figure_utils import remove_frame, remove_ticks, add_length_scale, cm_per_inch, \
- panel_size, head_direction_input_colormap
- from scripts.spatial_network.perlin_map.run_simulation_perlin_map import DATA_FOLDER, TRAJ_NAME, \
- get_input_head_directions, POLARIZED, CIRCULAR, NO_SYNAPSES
- plt.style.use('../../model_figure/figures.mplstyle')
- FIGURE_SAVE_PATH = '../../../figures/' + TRAJ_NAME + '/'
- FIGURE_SAVE_FORMAT = '.pdf'
- save_figs = False
- def tablelegend(ax, col_labels=None, row_labels=None, title_label="", *args, **kwargs):
- """
- Place a table legend on the axes.
- Creates a legend where the labels are not directly placed with the artists,
- but are used as row and column headers, looking like this:
- title_label | col_labels[1] | col_labels[2] | col_labels[3]
- -------------------------------------------------------------
- row_labels[1] |
- row_labels[2] | <artists go there>
- row_labels[3] |
- Parameters
- ----------
- ax : `matplotlib.axes.Axes`
- The artist that contains the legend table, i.e. current axes instant.
- col_labels : list of str, optional
- A list of labels to be used as column headers in the legend table.
- `len(col_labels)` needs to match `ncol`.
- row_labels : list of str, optional
- A list of labels to be used as row headers in the legend table.
- `len(row_labels)` needs to match `len(handles) // ncol`.
- title_label : str, optional
- Label for the top left corner in the legend table.
- ncol : int
- Number of columns.
- Other Parameters
- ----------------
- Refer to `matplotlib.legend.Legend` for other parameters.
- """
- #################### same as `matplotlib.axes.Axes.legend` #####################
- handles, labels, extra_args, kwargs = mlegend._parse_legend_args([ax], *args, **kwargs)
- if len(extra_args):
- raise TypeError('legend only accepts two non-keyword arguments')
- if col_labels is None and row_labels is None:
- ax.legend_ = mlegend.Legend(ax, handles, labels, **kwargs)
- ax.legend_._remove_method = ax._remove_legend
- return ax.legend_
- #################### modifications for table legend ############################
- else:
- ncol = kwargs.pop('ncol')
- handletextpad = kwargs.pop('handletextpad', 0 if col_labels is None else -2)
- title_label = [title_label]
- # blank rectangle handle
- extra = [Rectangle((0, 0), 1, 1, fc="w", fill=False, edgecolor='none', linewidth=0)]
- # empty label
- empty = [""]
- # number of rows infered from number of handles and desired number of columns
- nrow = len(handles) // ncol
- # organise the list of handles and labels for table construction
- if col_labels is None:
- assert nrow == len(
- row_labels), "nrow = len(handles) // ncol = %s, but should be equal to len(row_labels) = %s." % (
- nrow, len(row_labels))
- leg_handles = extra * nrow
- leg_labels = row_labels
- elif row_labels is None:
- assert ncol == len(col_labels), "ncol = %s, but should be equal to len(col_labels) = %s." % (
- ncol, len(col_labels))
- leg_handles = []
- leg_labels = []
- else:
- assert nrow == len(
- row_labels), "nrow = len(handles) // ncol = %s, but should be equal to len(row_labels) = %s." % (
- nrow, len(row_labels))
- assert ncol == len(col_labels), "ncol = %s, but should be equal to len(col_labels) = %s." % (
- ncol, len(col_labels))
- leg_handles = extra + extra * nrow
- leg_labels = title_label + row_labels
- for col in range(ncol):
- if col_labels is not None:
- leg_handles += extra
- leg_labels += [col_labels[col]]
- leg_handles += handles[col * nrow:(col + 1) * nrow]
- leg_labels += empty * nrow
- # Create legend
- ax.legend_ = mlegend.Legend(ax, leg_handles, leg_labels, ncol=ncol + int(row_labels is not None),
- handletextpad=handletextpad, **kwargs)
- ax.legend_._remove_method = ax._remove_legend
- return ax.legend_
- def get_closest_scale(traj, scale):
- available_lengths = sorted(list(set(traj.f_get("scale").f_get_range())))
- closest_length = available_lengths[np.argmin(np.abs(np.array(available_lengths) - scale))]
- if closest_length != scale:
- print("Warning: desired correlation length {:.1f} not available. Taking {:.1f} instead".format(
- scale, closest_length))
- corr_len = closest_length
- return corr_len
- def get_closest_fit_correlation_length(traj, fit_corr_len):
- corr_len_fit_dict = correlation_length_fit_dict(traj, map_type='perlin_map', load=True)
- available_lengths = list(corr_len_fit_dict.values())
- closest_length = available_lengths[np.argmin(np.abs(np.array(available_lengths) - fit_corr_len))]
- closest_corresponding_scale = list(corr_len_fit_dict.keys())[list(corr_len_fit_dict.values()).index(closest_length)]
- if closest_length != fit_corr_len:
- print("Warning: desired fit correlation length {:.1f} not available. Taking {:.1f} instead".format(
- fit_corr_len, closest_length))
- return closest_corresponding_scale
- def colorbar(mappable):
- from mpl_toolkits.axes_grid1 import make_axes_locatable
- import matplotlib.pyplot as plt
- last_axes = plt.gca()
- ax = mappable.axes
- fig = ax.figure
- divider = make_axes_locatable(ax)
- cax = divider.append_axes("right", size="5%", pad=0.05)
- cbar = fig.colorbar(mappable, cax=cax)
- plt.sca(last_axes)
- return cbar
- def plot_firing_rate_map_excitatory(traj, direction_idx, plot_run_names, exemplary_excitatory_neuron_id):
- max_val = 0
- for run_name in plot_run_names:
- fr_array = traj.results.runs[run_name].firing_rate_array
- f_rates = fr_array[:, direction_idx]
- run_max_val = np.max(f_rates)
- if run_max_val > max_val:
- # if traj.derived_parameters.runs[run_name].morphology.morph_label == POLARIZED:
- # n_id_max_rate = np.argmax(f_rates)
- max_val = run_max_val
- n_id_polar_plot = exemplary_excitatory_neuron_id
- # Mark the neuron that is shown in polar plot
- ex_positions = traj.results.runs[plot_run_names[0]].ex_positions
- polar_plot_x, polar_plot_y = ex_positions[n_id_polar_plot]
- # Vertices for the plotted triangle
- tr_scale = 30.
- tr_x = tr_scale * np.cos(2. * np.pi / 3. + np.pi / 2.)
- tr_y = tr_scale * np.sin(2. * np.pi / 3. + np.pi / 2.) + polar_plot_y
- tr_vertices = np.array(
- [[polar_plot_x, polar_plot_y + tr_scale], [tr_x + polar_plot_x, tr_y], [-tr_x + polar_plot_x, tr_y]])
- height = 1 * panel_size
- width = 3 * panel_size
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, (0.05, 0.15, 0.85, 0.7), axes_pad=panel_size / 6.0, cbar_location="right",
- cbar_mode="single",
- cbar_size="7%", cbar_pad=panel_size / 10.0,
- nrows_ncols=(1, 3))
- for ax, run_name in zip(axes, [plot_run_names[i] for i in [2, 0, 1]]):
- label = traj.derived_parameters.runs[run_name].morphology.morph_label
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- firing_rate_array = traj.results[run_name].firing_rate_array
- number_of_excitatory_neurons_per_row = int(np.sqrt(traj.N_E))
- c = ax.pcolor(X, Y, np.reshape(firing_rate_array[:, direction_idx], (number_of_excitatory_neurons_per_row,
- number_of_excitatory_neurons_per_row)),
- vmin=0, vmax=max_val, cmap='Reds')
- # ax.set_title(adjust_label(label))
- # ax.add_artist(Ellipse((polar_plot_x, polar_plot_y), 20., 20., color='k', fill=False, lw=2.))
- # ax.add_artist(Ellipse((polar_plot_x, polar_plot_y), 20., 20., color='w', fill=False, lw=1.))
- ax.add_artist(Polygon(tr_vertices, closed=True, fill=False, lw=2.5, color='k'))
- ax.add_artist(Polygon(tr_vertices, closed=True, fill=False, lw=1.5, color='w'))
- for spine in ax.spines.values():
- spine.set_edgecolor("grey")
- spine.set_linewidth(1)
- remove_ticks(ax)
- # fig.suptitle('spatial firing rate map', fontsize=16)
- ax.cax.colorbar(c)
- ax.cax.annotate("fr (Hz)", xy=(1, 1), xytext=(3, 3), xycoords="axes fraction", textcoords="offset points")
- # fig.tight_layout()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'C_firing_rate_map_excitatory' + FIGURE_SAVE_FORMAT, dpi=300)
- plt.close(fig)
- def plot_firing_rate_map_inhibitory(traj, direction_idx, plot_run_names, selected_inhibitory_neuron):
- max_val = 0
- for run_name in plot_run_names:
- fr_array = traj.results.runs[run_name].inh_firing_rate_array
- f_rates = fr_array[:, direction_idx]
- run_max_val = np.max(f_rates)
- if run_max_val > max_val:
- max_val = run_max_val
- n_id_polar_plot = selected_inhibitory_neuron
- # Mark the neuron that is shown in Polar plot
- inhibitory_axonal_cloud_array = traj.results.runs[plot_run_names[1]].inhibitory_axonal_cloud_array
- polar_plot_x = inhibitory_axonal_cloud_array[n_id_polar_plot, 0]
- polar_plot_y = inhibitory_axonal_cloud_array[n_id_polar_plot, 1]
- width = 3 * panel_size
- height = panel_size
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, (0.05, 0.15, 0.85, 0.7), axes_pad=panel_size / 6.0, cbar_location="right",
- cbar_mode="single",
- cbar_size="7%", cbar_pad=panel_size / 10.0,
- nrows_ncols=(1, 3))
- for ax, run_name in zip(axes, [plot_run_names[i] for i in [2, 0, 1]]):
- label = traj.derived_parameters.runs[run_name].morphology.morph_label
- inhibitory_axonal_cloud_array = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- inh_positions = [[p[0], p[1]] for p in inhibitory_axonal_cloud_array]
- X, Y = get_position_mesh(inh_positions)
- inh_firing_rate_array = traj.results[run_name].inh_firing_rate_array
- number_of_inhibitory_neurons_per_row = int(np.sqrt(traj.N_I))
- c = ax.pcolor(X, Y, np.reshape(inh_firing_rate_array[:, direction_idx], (number_of_inhibitory_neurons_per_row,
- number_of_inhibitory_neurons_per_row)),
- vmin=0, vmax=max_val, cmap='Blues')
- # ax.set_title(adjust_label(label))
- circle_r = 40.
- ax.add_artist(Ellipse((polar_plot_x, polar_plot_y), circle_r, circle_r, color='k', fill=False, lw=4.5))
- ax.add_artist(Ellipse((polar_plot_x, polar_plot_y), circle_r, circle_r, color='w', fill=False, lw=3))
- for spine in ax.spines.values():
- spine.set_edgecolor("grey")
- spine.set_linewidth(0.5)
- remove_ticks(ax)
- # fig.colorbar(c, ax=ax, label="f (Hz)")
- # fig.suptitle('spatial firing rate map', fontsize=16)
- cbar = ax.cax.colorbar(c)
- cbar_ticks = range(0,151,30)
- cbar.ax.set_yticks(cbar_ticks)
- ax.cax.annotate("fr (Hz)", xy=(1, 1), xytext=(3, 3), xycoords="axes fraction", textcoords="offset points")
- # fig.tight_layout()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'C_firing_rate_map_inhibitory' + FIGURE_SAVE_FORMAT, dpi=300)
- plt.close(fig)
- return max_val
- def normal_labels(label):
- if label == POLARIZED:
- label = 'polarized'
- elif label == NO_SYNAPSES:
- label = 'no interneurons'
- return label
- def short_labels(label):
- if label == POLARIZED:
- label = 'pol.'
- elif label == CIRCULAR:
- label = 'cir.'
- elif label == "no conn":
- label = "no inh."
- return label
- def plot_input_map(traj, run_name, figsize=(panel_size, panel_size), figname='input_map' + FIGURE_SAVE_FORMAT):
- n_ex = int(np.sqrt(traj.N_E))
- width, height = figsize
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, (0.15, 0.1, 0.7, 0.8), axes_pad=panel_size / 3.0, cbar_location="right", cbar_mode=None,
- nrows_ncols=(1, 1))
- ax = axes[0]
- traj.f_set_crun(run_name)
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- head_dir_preference = np.array(traj.results.runs[run_name].ex_tunings).reshape((n_ex, n_ex))
- scale_length = traj.morphology.long_axis * 2
- traj.f_restore_default()
- ax.pcolor(X, Y, head_dir_preference, vmin=-np.pi, vmax=np.pi, cmap=head_direction_input_colormap)
- for spine in ax.spines.values():
- spine.set_edgecolor("grey")
- spine.set_linewidth(0.5)
- remove_ticks(ax)
- start_scale_x = 100
- end_scale_x = start_scale_x + scale_length
- start_scale_y = -70
- end_scale_y = start_scale_y
- add_length_scale(ax, scale_length, start_scale_x, end_scale_x, start_scale_y, end_scale_y)
- ax.cax.set_visible(False)
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + figname)
- plt.close(fig)
- def plot_example_input_maps(traj, figsize=(panel_size, panel_size), figname='perlin_example_input_maps' + FIGURE_SAVE_FORMAT):
- n_ex = int(np.sqrt(traj.N_E))
- width, height = figsize
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, (0.1, 0.1, 0.85, 0.85), axes_pad=panel_size / 12.0, nrows_ncols=(3, 3))
- fit_corr_len_list = [20, 50, 100]
- seed_list = [1, 2, 3]
- corr_and_seed = itertools.product(fit_corr_len_list, seed_list)
- corr_len_fit_dict = correlation_length_fit_dict(traj, map_type='perlin_map', load=True)
- for idx, (ax, (fit_corr_len, seed)) in enumerate(zip(axes, corr_and_seed)):
- corresp_scale = get_closest_fit_correlation_length(traj, fit_corr_len)
- par_dict = {'seed': seed, 'correlation_length': corresp_scale}
- run_name = filter_run_names_by_par_dict(traj, par_dict)[0]
- traj.f_set_crun(run_name)
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- head_dir_preference = np.array(traj.results.runs[run_name].ex_tunings).reshape((n_ex, n_ex))
- scale_length = traj.morphology.long_axis * 2
- traj.f_restore_default()
- ax.pcolor(X, Y, head_dir_preference, vmin=-np.pi, vmax=np.pi, cmap=head_direction_input_colormap)
- for spine in ax.spines.values():
- spine.set_edgecolor("grey")
- spine.set_linewidth(0.5)
- remove_ticks(ax)
- ax.set_ylabel('{:3.0f} um'.format(corr_len_fit_dict[corresp_scale]))
- inhibitory_axonal_cloud_array = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- plt_polar_id = 255
- plt_polar_coordinates = inhibitory_axonal_cloud_array[plt_polar_id]
- plt_polar_neuron = Interneuron(plt_polar_coordinates[0],
- plt_polar_coordinates[1],
- traj.morphology.long_axis,
- traj.morphology.short_axis,
- plt_polar_coordinates[2])
- circ_radius = np.sqrt(traj.morphology.long_axis * traj.morphology.short_axis)
- plt_circ_id = 65
- plt_circ_coordinates = inhibitory_axonal_cloud_array[plt_circ_id]
- plt_circ_neuron = Interneuron(plt_circ_coordinates[0],
- plt_circ_coordinates[1],
- circ_radius,
- circ_radius,
- plt_circ_coordinates[2])
- plt_interneurons = [plt_polar_neuron, plt_circ_neuron]
- for p in plt_interneurons:
- ell = p.get_ellipse()
- edgecolor = 'black'
- alpha = 1
- zorder = 10
- linewidth = 1.
- ell.set_edgecolor(edgecolor)
- ell.set_alpha(alpha)
- ell.set_zorder(zorder)
- ell.set_linewidth(linewidth)
- ax.add_artist(ell)
- if idx == 8:
- start_scale_x = 100
- end_scale_x = start_scale_x + scale_length
- start_scale_y = -70
- end_scale_y = start_scale_y
- add_length_scale(ax, scale_length, start_scale_x, end_scale_x, start_scale_y, end_scale_y)
- ax.cax.set_visible(False)
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + figname)
- plt.close(fig)
- def plot_axonal_clouds(traj, plot_run_names):
- n_ex = int(np.sqrt(traj.N_E))
- height = 1 * panel_size
- width = 3 * panel_size
- cluster_positions = [(250, 250), (750, 200), (450, 600)]
- cluster_sizes = [4, 4, 4]
- selected_neurons = []
- inhibitory_positions = traj.results.runs[plot_run_names[0]].inhibitory_axonal_cloud_array[:, :2]
- for cluster_position, number_of_neurons_in_cluster in zip(cluster_positions, cluster_sizes):
- selection = get_neurons_close_to_given_position(cluster_position, number_of_neurons_in_cluster,
- inhibitory_positions)
- selected_neurons.extend(selection)
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, 111, axes_pad=0.15, cbar_location="right", cbar_mode="single", cbar_size="7%",
- nrows_ncols=(1, 3))
- for ax, run_name in zip(axes, [plot_run_names[i] for i in [2, 0, 1]]):
- traj.f_set_crun(run_name)
- label = traj.derived_parameters.runs[run_name].morphology.morph_label
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- inhibitory_axonal_cloud_array = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- axonal_clouds = [Interneuron(p[0], p[1], traj.morphology.long_axis, traj.morphology.short_axis, p[2]) for p in
- inhibitory_axonal_cloud_array]
- head_dir_preference = np.array(traj.results.runs[run_name].ex_tunings).reshape((n_ex, n_ex))
- c = ax.pcolor(X, Y, head_dir_preference, vmin=-np.pi, vmax=np.pi, cmap=head_direction_input_colormap)
- # ax.set_title(normal_labels(label))
- ax.set_aspect('equal')
- remove_frame(ax)
- remove_ticks(ax)
- # fig.colorbar(c, ax=ax, label="Tuning")
- if label != NO_SYNAPSES and axonal_clouds is not None:
- for i, p in enumerate(axonal_clouds):
- ell = p.get_ellipse()
- if i in selected_neurons:
- edgecolor = 'black'
- alpha = 1
- zorder = 10
- linewidth = 1
- else:
- edgecolor = 'gray'
- alpha = 0.5
- zorder = 1
- linewidth = 0.3
- ell.set_edgecolor(edgecolor)
- ell.set_alpha(alpha)
- ell.set_zorder(zorder)
- ell.set_linewidth(linewidth)
- ax.add_artist(ell)
- traj.f_set_crun(plot_run_names[0])
- scale_length = traj.morphology.long_axis * 2
- traj.f_restore_default()
- start_scale_x = 100
- end_scale_x = start_scale_x + scale_length
- start_scale_y = -70
- end_scale_y = start_scale_y
- add_length_scale(axes[0], scale_length, start_scale_x, end_scale_x, start_scale_y, end_scale_y)
- axes[1].set_yticklabels([])
- axes[2].set_yticklabels([])
- axes[0].cax.set_visible(False)
- traj.f_restore_default()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'B_i_axonal_clouds' + FIGURE_SAVE_FORMAT)
- plt.close(fig)
- def get_neurons_close_to_given_position(cluster_position, number_of_neurons_in_cluster, positions):
- position = np.array(cluster_position)
- distance_vectors = positions - np.expand_dims(position, 0).repeat(positions.shape[0],
- axis=0)
- distances = np.linalg.norm(distance_vectors, axis=1)
- selection = list(np.argpartition(distances, number_of_neurons_in_cluster)[:number_of_neurons_in_cluster])
- return selection
- def plot_orientation_maps_diff_scales(traj):
- n_ex = int(np.sqrt(traj.N_E))
- scale_run_names = []
- plot_scales = [0.0, 100.0, 200.0, 300.0]
- for scale in plot_scales:
- par_dict = {'seed': 1, 'correlation_length': get_closest_scale(traj, scale), 'long_axis': 100.}
- scale_run_names.append(*filter_run_names_by_par_dict(traj, par_dict))
- fig, axes = plt.subplots(1, 4, figsize=(18., 4.5))
- for ax, run_name, scale in zip(axes, scale_run_names, plot_scales):
- traj.f_set_crun(run_name)
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- head_dir_preference = np.array(traj.results.runs[run_name].ex_tunings).reshape((n_ex, n_ex))
- # TODO: Why was this transposed for plotting? (now changed)
- c = ax.pcolor(X, Y, head_dir_preference, vmin=-np.pi, vmax=np.pi, cmap='twilight')
- ax.set_title('Correlation length: {}'.format(scale))
- fig.colorbar(c, ax=ax, label="Tuning")
- # fig.suptitle('axonal cloud', fontsize=16)
- traj.f_restore_default()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'orientation_maps_diff_scales' + FIGURE_SAVE_FORMAT)
- plt.close(fig)
- def plot_orientation_maps_diff_scales_with_ellipse(traj):
- n_ex = int(np.sqrt(traj.N_E))
- scale_run_names = []
- plot_scales = [0.0, 400.0, 800.0]
- real_scales = [get_closest_scale(traj, scale) for scale in plot_scales]
- for scale in plot_scales:
- par_dict = {'seed': 1, 'correlation_length': get_closest_scale(traj, scale), 'long_axis': 100.}
- scale_run_names.append(*filter_run_names_by_par_dict(traj, par_dict))
- inhibitory_positions = traj.results.runs[scale_run_names[1]].inhibitory_axonal_cloud_array[:, :2]
- selected_polar_neuron = get_neurons_close_to_given_position((300, 300), 1, inhibitory_positions)[0]
- selected_circular_neuron = get_neurons_close_to_given_position((600, 600), 1, inhibitory_positions)[0]
- width = panel_size
- height = 1.5 * panel_size
- fig = plt.figure(figsize=(width, height))
- axes = ImageGrid(fig, 111, axes_pad=0.15,
- nrows_ncols=(len(plot_scales), 1))
- for ax, run_name, scale in zip(axes, scale_run_names, real_scales):
- traj.f_set_crun(run_name)
- X, Y = get_position_mesh(traj.results.runs[run_name].ex_positions)
- inhibitory_axonal_cloud_array = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- axonal_clouds = [Interneuron(p[0], p[1], traj.morphology.long_axis, traj.morphology.short_axis, p[2]) for p in
- inhibitory_axonal_cloud_array]
- head_dir_preference = np.array(traj.results.runs[run_name].ex_tunings).reshape((n_ex, n_ex))
- # TODO: Why was this transposed for plotting? (now changed)
- c = ax.pcolor(X, Y, head_dir_preference, vmin=-np.pi, vmax=np.pi, cmap=head_direction_input_colormap)
- # ax.set_title('Correlation length: {}'.format(scale))
- # fig.colorbar(c, ax=ax, label="Tuning")
- ax.set_xticks([])
- ax.set_yticks([])
- p1 = axonal_clouds[selected_polar_neuron]
- ell = p1.get_ellipse()
- ell._linewidth = 2.
- ax.add_artist(ell)
- p2 = axonal_clouds[selected_circular_neuron]
- circ_r = 2 * np.sqrt(p2.a * p2.b)
- circ = Ellipse((p2.x, p2.y), circ_r, circ_r, fill=False, zorder=2, edgecolor='k')
- circ._linewidth = 2.
- ax.add_artist(circ)
- ax.annotate("{:.0f}".format(scale), xy=(1.0, 0.5), xytext=(2, 0), xycoords="axes fraction", textcoords="offset "
- "points",
- va="center", ha="left")
- remove_frame(ax)
- # fig.suptitle('axonal cloud', fontsize=16)
- traj.f_restore_default()
- add_length_scale(axes[1], 200, -300, -100, 50, 50)
- axes[0].annotate("input maps", xy=(-0.5, 1.0), xycoords="axes fraction", rotation=90, ha="left", va="top")
- fig.tight_layout()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'F_orientation_maps_diff_scales_with_ellipse' + FIGURE_SAVE_FORMAT)
- plt.close(fig)
- def plot_polar_plot_excitatory(traj, plot_run_names, selected_neuron_idx,
- figname=FIGURE_SAVE_PATH + 'D_polar_plot_excitatory' + FIGURE_SAVE_FORMAT):
- directions = np.linspace(-np.pi, np.pi, traj.input.number_of_directions, endpoint=False)
- directions_plt = list(directions)
- directions_plt.append(directions[0])
- height = panel_size
- width = panel_size
- fig, ax = plt.subplots(1, 1, figsize=(height, width), subplot_kw=dict(projection='polar'))
- # head_direction_indices = traj.results.runs[plot_run_names[0]].head_direction_indices
- # sorted_ids = np.argsort(head_direction_indices)
- # plot_n_idx = sorted_ids[-75]
- plot_n_idx = selected_neuron_idx
- line_styles = ['dotted', 'solid', 'dashed']
- colors = ['r', 'lightsalmon', 'grey']
- labels = ['pol. ', 'cir. ', 'no inh.']
- line_widths = [1.5, 1.5, 1]
- zorders = [10, 2, 1]
- max_rate = 0.0
- ax.plot([], [], color='white',label=' ')
- hdis = []
- for run_idx, run_name in enumerate(plot_run_names):
- # label = traj.derived_parameters.runs[run_name].morphology.morph_label
- label = labels[run_idx]
- hdi = traj.results.runs[run_name].head_direction_indices[selected_neuron_idx]
- hdis.append(hdi)
- tuning_vectors = traj.results.runs[run_name].tuning_vectors
- rate_plot = [np.linalg.norm(v) for v in tuning_vectors[plot_n_idx]]
- run_max_rate = np.max(rate_plot)
- if run_max_rate > max_rate:
- max_rate = run_max_rate
- rate_plot.append(rate_plot[0])
- # plt_label = '{:s} {:.2f}'.format(short_labels(label), hdi)
- plt_label = '{:s}'.format(short_labels(label))
- ax.plot(directions_plt, rate_plot, linewidth=line_widths[run_idx],
- label=plt_label, color=colors[run_idx], linestyle=line_styles[run_idx], zorder=zorders[run_idx])
- # ax.set_title('Firing Rate')
- # ax.plot([0.0, 0.0], [0.0, 1.05 * max_rate], color='red', alpha=0.25, linewidth=4.)
- # TODO: Set ticks for polar
- ticks = [40., 80.]
- ax.set_rgrids(ticks, labels=["{:.0f} Hz".format(ticklabel) if idx == len(ticks) - 1 else "" for idx, ticklabel in
- enumerate(ticks)], angle=60)
- ax.set_thetagrids([0, 90, 180, 270], labels=[])
- ax.xaxis.grid(linewidth=0.4)
- ax.yaxis.grid(linewidth=0.4)
- leg = ax.legend(loc="lower right", bbox_to_anchor=(1.15, -0.2), handlelength=1, fontsize="medium")
- leg.get_frame().set_linewidth(0.0)
- hdi_box_x, hdi_box_y = (0.86, -0.09)
- hdi_box_dy = 0.14
- hdi_box = ax.text(hdi_box_x, hdi_box_y + 3 * hdi_box_dy, 'HDI', fontsize='medium', transform=ax.transAxes,zorder=9.)
- hdi_box = ax.text(hdi_box_x, hdi_box_y + 2 * hdi_box_dy, '{:.2f}'.format(hdis[0]), fontsize='medium',transform=ax.transAxes, zorder=9.)
- hdi_box = ax.text(hdi_box_x, hdi_box_y + hdi_box_dy, '{:.2f}'.format(hdis[1]), fontsize='medium', transform=ax.transAxes,zorder=9.)
- hdi_box = ax.text(hdi_box_x, hdi_box_y, '{:.2f}'.format(hdis[2]), fontsize='medium', transform=ax.transAxes,zorder=9.)
- ax.axes.spines["polar"].set_visible(False)
- if save_figs:
- plt.savefig(figname)
- plt.close(fig)
- def plot_polar_plot_inhibitory(traj, plot_run_names, selected_neuron_idx, figname=FIGURE_SAVE_PATH +
- 'D_polar_plot_inhibitory' + FIGURE_SAVE_FORMAT):
- directions = np.linspace(-np.pi, np.pi, traj.input.number_of_directions, endpoint=False)
- directions_plt = list(directions)
- directions_plt.append(directions[0])
- height = panel_size
- width = panel_size
- fig, ax = plt.subplots(1, 1, figsize=(height, width), subplot_kw=dict(projection='polar'))
- # head_direction_indices = traj.results.runs[plot_run_names[0]].inh_head_direction_indices
- # sorted_ids = np.argsort(head_direction_indices)
- # plot_n_idx = sorted_ids[-75]
- plot_n_idx = selected_neuron_idx
- line_styles = ['dotted', 'solid']
- colors = ['b', 'lightblue']
- labels = ['pol. ', 'cir. ']
- line_widths = [1.5, 1.5]
- zorders = [10, 2]
- ax.plot([], [], color='white',label=' ')
- hdis = []
- for run_idx, run_name in enumerate(plot_run_names[:2]):
- # ax = axes[max_hdi_idx, run_idx]
- # label = traj.derived_parameters.runs[run_name].morphology.morph_label
- label = labels[run_idx]
- hdi = traj.results.runs[run_name].inh_head_direction_indices[selected_neuron_idx]
- hdis.append(hdi)
- tuning_vectors = traj.results.runs[run_name].inh_tuning_vectors
- rate_plot = [np.linalg.norm(v) for v in tuning_vectors[plot_n_idx]]
- rate_plot.append(rate_plot[0])
- # plt_label = '{:s} {:.2f}'.format(short_labels(label), hdi)
- plt_label = '{:s}'.format(short_labels(label))
- ax.plot(directions_plt, rate_plot, linewidth=line_widths[run_idx],
- label=plt_label,
- color=colors[run_idx], linestyle=line_styles[run_idx], zorder=zorders[run_idx])
- # ax.set_title('Inh. Firing Rate')
- # TODO: Set ticks for polar
- # ticks = [np.round(max_rate / 3.), np.round(max_rate * 2. / 3.), np.round(max_rate)]
- ticks = [40., 80., 120.]
- ax.set_rgrids(ticks, labels=["{:.0f} Hz".format(ticklabel) if idx == len(ticks) - 1 else "" for idx, ticklabel in
- enumerate(ticks)], angle=60)
- ax.set_thetagrids([0, 90, 180, 270], labels=[])
- ax.xaxis.grid(linewidth=0.4)
- ax.yaxis.grid(linewidth=0.4)
- leg = ax.legend(loc="lower right", bbox_to_anchor=(1.15, -0.15), handlelength=1, fontsize="medium")
- leg.get_frame().set_linewidth(0.0)
- hdi_box_x, hdi_box_y = (0.86, -0.04)
- hdi_box_dy = 0.14
- hdi_box = ax.text(hdi_box_x, hdi_box_y+2*hdi_box_dy, 'HDI', fontsize='medium', transform=ax.transAxes, zorder=9.)
- hdi_box = ax.text(hdi_box_x, hdi_box_y+hdi_box_dy, '{:.2f}'.format(hdis[0]), fontsize='medium', transform=ax.transAxes, zorder=9.)
- hdi_box = ax.text(hdi_box_x, hdi_box_y, '{:.2f}'.format(hdis[1]), fontsize='medium', transform=ax.transAxes, zorder=9.)
- ax.axes.spines["polar"].set_visible(False)
- if save_figs:
- plt.savefig(figname)
- plt.close(fig)
- def plot_hdi_histogram_combined_and_overlayed(traj, plot_run_names, ex_polar_plot_id, in_polar_plot_id, cut_off_dist):
- labels = []
- inh_hdis = []
- exc_hdis = []
- no_conn_hdi = 0.
- for run_idx, run_name in enumerate(plot_run_names):
- label = traj.derived_parameters.runs[run_name].morphology.morph_label
- if label != NO_SYNAPSES:
- labels.append(normal_labels(label))
- inh_head_direction_indices = traj.results.runs[run_name].inh_head_direction_indices
- inh_axonal_cloud = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- inh_cut_off_ids = (inh_axonal_cloud[:, 0] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (inh_axonal_cloud[:, 1] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- inh_hdis.append(sorted(inh_head_direction_indices[inh_cut_off_ids]))
- exc_head_direction_indices = traj.results.runs[run_name].head_direction_indices
- ex_positions = traj.results.runs[run_name].ex_positions
- exc_cut_off_ids = (ex_positions[:, 0] >= cut_off_dist) & \
- (ex_positions[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (ex_positions[:, 1] >= cut_off_dist) & \
- (ex_positions[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- exc_hdis.append(sorted(exc_head_direction_indices[exc_cut_off_ids]))
- else:
- exc_head_direction_indices = traj.results.runs[run_name].head_direction_indices
- no_conn_hdi = np.mean(exc_head_direction_indices)
- # Look for a representative excitatory neuron
- hdi_mean_dict = {}
- excitatory_hdi_means = [np.mean(hdis) for hdis in exc_hdis]
- hdi_mean_dict["polar_exc"] = excitatory_hdi_means[0]
- hdi_mean_dict["circular_exc"] = excitatory_hdi_means[1]
- inhibitory_hdi_means = [np.mean(hdis) for hdis in inh_hdis]
- hdi_mean_dict["polar_inh"] = inhibitory_hdi_means[0]
- hdi_mean_dict["circular_inh"] = inhibitory_hdi_means[1]
- width = 2 * panel_size
- height = 1.2 * panel_size
- fig, axes = plt.subplots(2, 1, figsize=(width, height))
- plt.subplots_adjust(wspace=0, hspace=0.1)
- bins = np.linspace(0.0, 1.0, 21, endpoint=True)
- max_density = 0
- for i in range(2):
- # i = i + 1 # grid spec indexes from 0
- # ax = fig.add_subplot(gs[i])
- ax = axes[i]
- density_e, _, _ = ax.hist(exc_hdis[i], color='r', edgecolor='r', alpha=0.3, bins=bins, density=True)
- density_i, _, _ = ax.hist(inh_hdis[i], color='b', edgecolor='b', alpha=0.3, bins=bins, density=True)
- max_density = np.max([max_density, np.max(density_e), np.max(density_i)])
- ax.axvline(np.mean(exc_hdis[i]), color='r')
- ax.axvline(np.mean(inh_hdis[i]), color='b')
- ax.axvline(no_conn_hdi, color='dimgrey', linestyle='--', linewidth=1.5)
- ax.set_ylabel(labels[i], rotation='vertical')
- # plt.axis('on')
- if i == 0:
- ax.set_xticklabels([])
- else:
- ax.set_xlabel('HDI')
- remove_frame(ax, ["top", "right", "bottom"])
- max_density = 1.2 * max_density
- fig.subplots_adjust(left=0.2, right=0.95, bottom=0.2)
- axes[0].annotate('% cells', (0, 1.0), xycoords='axes fraction', va="bottom", ha="right")
- axes[1].annotate("no ihn.\n{:.2f}".format(no_conn_hdi), xy=(no_conn_hdi, max_density),
- xytext=(-2, 0), xycoords="data",
- textcoords="offset points",
- va="top", ha="right", color="dimgrey")
- for i, ax in enumerate(axes):
- ax.annotate("{:.2f}".format(np.mean(exc_hdis[i])), xy=(np.mean(exc_hdis[i]), max_density),
- xytext=(2, 0), xycoords="data",
- textcoords="offset points",
- va="top", ha="left", color="r")
- # i_ha = "left" if i == 1 else "right"
- # i_offset = 2 if i == 1 else -2
- i_ha = "right"
- i_offset = -1
- ax.annotate("{:.2f}".format(np.mean(inh_hdis[i])), xy=(np.mean(inh_hdis[i]), max_density),
- xytext=(i_offset, 0), xycoords="data",
- textcoords="offset points",
- va="top", ha=i_ha, color="b")
- for ax in axes:
- ax.set_ylim(0, max_density)
- # plt.annotate('probability density', (-0.2,1.5), xycoords='axes fraction', rotation=90, fontsize=18)
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'E_hdi_histogram_combined_and_overlayed' + FIGURE_SAVE_FORMAT.format(cut_off_dist))
- plt.close(fig)
- return hdi_mean_dict
- def get_neurons_with_given_hdi(polar_hdi, circular_hdi, max_number_of_suggestions, plot_run_names, traj, type):
- polar_run_name = plot_run_names[0]
- circular_run_name = plot_run_names[1]
- polar_ex_hdis = traj.results.runs[polar_run_name].head_direction_indices if type == "ex" else traj.results.runs[
- polar_run_name].inh_head_direction_indices
- circular_ex_hdis = traj.results.runs[circular_run_name].head_direction_indices if type == "ex" else \
- traj.results.runs[
- polar_run_name].inh_head_direction_indices
- neuron_indices = get_indices_of_closest_values(polar_ex_hdis, polar_hdi,
- circular_ex_hdis,
- circular_hdi, 0.1 * np.abs(
- polar_hdi - circular_hdi), max_number_of_suggestions)
- return neuron_indices
- def get_indices_of_closest_values(first_list, first_value, second_list, second_value, absolute_tolerance_list_one,
- number_of_indices):
- is_close_in_list_one = np.abs(first_list - first_value) < absolute_tolerance_list_one
- indices_close_in_list_one = np.where(is_close_in_list_one)[0]
- indices_closest_in_list_two = indices_close_in_list_one[np.argpartition(np.abs(second_list[
- indices_close_in_list_one] - second_value),
- number_of_indices)]
- return indices_closest_in_list_two[:number_of_indices]
- def filter_run_names_by_par_dict(traj, par_dict):
- run_name_list = []
- for run_idx, run_name in enumerate(traj.f_get_run_names()):
- traj.f_set_crun(run_name)
- paramters_equal = True
- for key, val in par_dict.items():
- if (traj.par[key] != val):
- paramters_equal = False
- if paramters_equal:
- run_name_list.append(run_name)
- traj.f_restore_default()
- return run_name_list
- def plot_exc_and_inh_hdi_over_simplex_grid_scale(traj, plot_run_names, cut_off_dist):
- corr_len_expl = traj.f_get('scale').f_get_range()
- seed_expl = traj.f_get('seed').f_get_range()
- label_expl = [traj.derived_parameters.runs[run_name].morphology.morph_label for run_name in traj.f_get_run_names()]
- label_range = set(label_expl)
- exc_hdi_frame = pd.Series(index=[corr_len_expl, seed_expl, label_expl])
- exc_hdi_frame.index.names = ["corr_len", "seed", "label"]
- inh_hdi_frame = pd.Series(index=[corr_len_expl, seed_expl, label_expl])
- inh_hdi_frame.index.names = ["corr_len", "seed", "label"]
- for run_name, corr_len, seed, label in zip(plot_run_names, corr_len_expl, seed_expl, label_expl):
- ex_tunings = traj.results.runs[run_name].ex_tunings
- inh_hdis = []
- exc_hdis = []
- inh_head_direction_indices = traj.results.runs[run_name].inh_head_direction_indices
- inh_axonal_cloud = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- inh_cut_off_ids = (inh_axonal_cloud[:, 0] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (inh_axonal_cloud[:, 1] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- inh_hdis.append(sorted(inh_head_direction_indices[inh_cut_off_ids]))
- exc_head_direction_indices = traj.results.runs[run_name].head_direction_indices
- ex_positions = traj.results.runs[run_name].ex_positions
- exc_cut_off_ids = (ex_positions[:, 0] >= cut_off_dist) & \
- (ex_positions[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (ex_positions[:, 1] >= cut_off_dist) & \
- (ex_positions[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- exc_hdis.append(sorted(exc_head_direction_indices[exc_cut_off_ids]))
- exc_hdi_frame[corr_len, seed, label] = np.mean(exc_hdis)
- inh_hdi_frame[corr_len, seed, label] = np.mean(inh_hdis)
- # TODO: Standard deviation also for the population
- exc_hdi_n_and_seed_mean = exc_hdi_frame.groupby(level=[0, 2]).mean()
- exc_hdi_n_and_seed_std_dev = exc_hdi_frame.groupby(level=[0, 2]).std()
- inh_hdi_n_and_seed_mean = inh_hdi_frame.groupby(level=[0, 2]).mean()
- inh_hdi_n_and_seed_std_dev = inh_hdi_frame.groupby(level=[0, 2]).std()
- markersize = 4.
- exc_style_dict = {
- NO_SYNAPSES: ['dimgrey', 'dashed', '', 0],
- POLARIZED: ['red', 'solid', '^', markersize],
- CIRCULAR: ['lightsalmon', 'solid', '^', markersize]
- }
- inh_style_dict = {
- NO_SYNAPSES: ['dimgrey', 'dashed', '', 0],
- POLARIZED: ['blue', 'solid', 'o', markersize],
- CIRCULAR: ['lightblue', 'solid', 'o', markersize]
- }
- width = 2 * panel_size
- height = 1.2 * panel_size
- fig, ax = plt.subplots(1, 1, figsize=(width, height))
- for label in sorted(label_range, reverse=True):
- if label == NO_SYNAPSES:
- no_conn_hdi = exc_hdi_n_and_seed_mean[get_closest_scale(traj, 200), label]
- ax.axhline(no_conn_hdi, color='grey', linestyle='--')
- ax.annotate(short_labels(label), xy=(1.0, no_conn_hdi), xytext=(0, -2), xycoords='axes fraction',
- textcoords="offset points",
- va="top", \
- ha="right",
- color="dimgrey")
- continue
- exc_hdi_mean = exc_hdi_n_and_seed_mean[:, label]
- exc_hdi_std = exc_hdi_n_and_seed_std_dev[:, label]
- inh_hdi_mean = inh_hdi_n_and_seed_mean[:, label]
- inh_hdi_std = inh_hdi_n_and_seed_std_dev[:, label]
- corr_len_range = exc_hdi_mean.keys().to_numpy()
- exc_col, exc_lin, exc_mar, exc_mar_size = exc_style_dict[label]
- inh_col, inh_lin, inh_mar, inh_mar_size = inh_style_dict[label]
- simplex_grid_scale = corr_len_range * np.sqrt(2)
- ax.plot(simplex_grid_scale, exc_hdi_mean, label='exc., ' + label, marker=exc_mar, color=exc_col, linestyle=exc_lin,
- markersize=exc_mar_size, alpha=0.5)
- plt.fill_between(simplex_grid_scale, exc_hdi_mean - exc_hdi_std,
- exc_hdi_mean + exc_hdi_std, alpha=0.3, color=exc_col)
- ax.plot(simplex_grid_scale, inh_hdi_mean, label='inh., ' + label, marker=inh_mar, color=inh_col, linestyle=inh_lin,
- markersize=inh_mar_size, alpha=0.5)
- plt.fill_between(simplex_grid_scale, inh_hdi_mean - inh_hdi_std,
- inh_hdi_mean + inh_hdi_std, alpha=0.3, color=inh_col)
- ax.set_xlabel('simplex grid scale')
- ax.set_ylabel('HDI')
- ax.axvline(get_closest_scale(traj, 200.0) * np.sqrt(2), color='k', linewidth=0.5, zorder=0)
- ax.set_ylim(0.0, 1.0)
- # ax.set_xlim(0.0, np.max(corr_len_range))
- remove_frame(ax, ["right", "top"])
- tablelegend(ax, ncol=2, bbox_to_anchor=(1.1, 1.1), loc="upper right",
- row_labels=None,
- col_labels=[short_labels(label) for label in sorted(label_range - {"no conn"}, reverse=True)],
- title_label='', borderaxespad=0, handlelength=2, edgecolor='white')
- fig.subplots_adjust(bottom=0.2, left=0.2)
- # plt.legend()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'F_hdi_over_grid_scale' + FIGURE_SAVE_FORMAT)
- plt.close(fig)
- def plot_exc_and_inh_hdi_over_fit_corr_len(traj, plot_run_names, cut_off_dist):
- corr_len_expl = traj.f_get('scale').f_get_range()
- seed_expl = traj.f_get('seed').f_get_range()
- label_expl = [traj.derived_parameters.runs[run_name].morphology.morph_label for run_name in traj.f_get_run_names()]
- label_range = set(label_expl)
- exc_hdi_frame = pd.Series(index=[corr_len_expl, seed_expl, label_expl])
- exc_hdi_frame.index.names = ["corr_len", "seed", "label"]
- inh_hdi_frame = pd.Series(index=[corr_len_expl, seed_expl, label_expl])
- inh_hdi_frame.index.names = ["corr_len", "seed", "label"]
- for run_name, corr_len, seed, label in zip(plot_run_names, corr_len_expl, seed_expl, label_expl):
- ex_tunings = traj.results.runs[run_name].ex_tunings
- inh_hdis = []
- exc_hdis = []
- inh_head_direction_indices = traj.results.runs[run_name].inh_head_direction_indices
- inh_axonal_cloud = traj.results.runs[run_name].inhibitory_axonal_cloud_array
- inh_cut_off_ids = (inh_axonal_cloud[:, 0] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (inh_axonal_cloud[:, 1] >= cut_off_dist) & \
- (inh_axonal_cloud[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- inh_hdis.append(sorted(inh_head_direction_indices[inh_cut_off_ids]))
- exc_head_direction_indices = traj.results.runs[run_name].head_direction_indices
- ex_positions = traj.results.runs[run_name].ex_positions
- exc_cut_off_ids = (ex_positions[:, 0] >= cut_off_dist) & \
- (ex_positions[:, 0] <= traj.parameters.input_map.sheet_size - cut_off_dist) & \
- (ex_positions[:, 1] >= cut_off_dist) & \
- (ex_positions[:, 1] <= traj.parameters.input_map.sheet_size - cut_off_dist)
- exc_hdis.append(sorted(exc_head_direction_indices[exc_cut_off_ids]))
- exc_hdi_frame[corr_len, seed, label] = np.mean(exc_hdis)
- inh_hdi_frame[corr_len, seed, label] = np.mean(inh_hdis)
- # TODO: Standard deviation also for the population
- exc_hdi_n_and_seed_mean = exc_hdi_frame.groupby(level=[0, 2]).mean()
- exc_hdi_n_and_seed_std_dev = exc_hdi_frame.groupby(level=[0, 2]).std()
- inh_hdi_n_and_seed_mean = inh_hdi_frame.groupby(level=[0, 2]).mean()
- inh_hdi_n_and_seed_std_dev = inh_hdi_frame.groupby(level=[0, 2]).std()
- markersize = 4.
- exc_style_dict = {
- NO_SYNAPSES: ['dimgrey', 'dashed', '', 0],
- POLARIZED: ['red', 'solid', '^', markersize],
- CIRCULAR: ['lightsalmon', 'solid', '^', markersize]
- }
- inh_style_dict = {
- NO_SYNAPSES: ['dimgrey', 'dashed', '', 0],
- POLARIZED: ['blue', 'solid', 'o', markersize],
- CIRCULAR: ['lightblue', 'solid', 'o', markersize]
- }
- # colors = ['blue', 'grey', 'lightblue']
- # linestyles = ['solid', 'dashed', 'solid']
- # markers = [verts, '', 'o']
- corr_len_fit_dict = correlation_length_fit_dict(traj, map_type='perlin_map', load=True)
- width = 2 * panel_size
- height = 1.2 * panel_size
- fig, ax = plt.subplots(1, 1, figsize=(width, height))
- for label in sorted(label_range, reverse=True):
- if label == NO_SYNAPSES:
- no_conn_hdi = exc_hdi_n_and_seed_mean[get_closest_scale(traj, 200), label]
- ax.axhline(no_conn_hdi, color='grey', linestyle='--')
- ax.annotate(short_labels(label), xy=(1.0, no_conn_hdi), xytext=(0, -2), xycoords='axes fraction',
- textcoords="offset points",
- va="top", \
- ha="right",
- color="dimgrey")
- continue
- exc_hdi_mean = exc_hdi_n_and_seed_mean[:, label]
- exc_hdi_std = exc_hdi_n_and_seed_std_dev[:, label]
- inh_hdi_mean = inh_hdi_n_and_seed_mean[:, label]
- inh_hdi_std = inh_hdi_n_and_seed_std_dev[:, label]
- corr_len_range = exc_hdi_mean.keys().to_numpy()
- exc_col, exc_lin, exc_mar, exc_mar_size = exc_style_dict[label]
- inh_col, inh_lin, inh_mar, inh_mar_size = inh_style_dict[label]
- fit_corr_len = [corr_len_fit_dict[corr_len] for corr_len in corr_len_range]
- ax.plot(fit_corr_len, exc_hdi_mean, label='exc., ' + label, marker=exc_mar, color=exc_col, linestyle=exc_lin,
- markersize=exc_mar_size, alpha=0.5)
- plt.fill_between(fit_corr_len, exc_hdi_mean - exc_hdi_std,
- exc_hdi_mean + exc_hdi_std, alpha=0.3, color=exc_col)
- ax.plot(fit_corr_len, inh_hdi_mean, label='inh., ' + label, marker=inh_mar, color=inh_col, linestyle=inh_lin,
- markersize=inh_mar_size, alpha=0.5)
- plt.fill_between(fit_corr_len, inh_hdi_mean - inh_hdi_std,
- inh_hdi_mean + inh_hdi_std, alpha=0.3, color=inh_col)
- ax.set_xlabel('correlation length (um)')
- ax.set_ylabel('HDI')
- ax.axvline(corr_len_fit_dict[get_closest_scale(traj, 200.0)], color='k', linewidth=0.5, zorder=0)
- ax.set_ylim(0.0, 1.0)
- # ax.set_xlim(0.0, np.max(corr_len_range))
- remove_frame(ax, ["right", "top"])
- tablelegend(ax, ncol=2, bbox_to_anchor=(1.1, 1.1), loc="upper right",
- row_labels=None,
- col_labels=[short_labels(label) for label in sorted(label_range - {"no conn"}, reverse=True)],
- title_label='', borderaxespad=0, handlelength=2, edgecolor='white')
- fig.subplots_adjust(bottom=0.2, left=0.2)
- # plt.legend()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'F_hdi_over_corr_len_scaled' + FIGURE_SAVE_FORMAT)
- plt.close(fig)
- def get_phase_difference(total_difference):
- """
- Map accumulated phase difference to shortest possible difference
- :param total_difference:
- :return: relative_difference
- """
- return (total_difference + np.pi) % (2 * np.pi) - np.pi
- def plot_firing_rate_similar_vs_diff_tuning(traj, plot_run_names, figsize=(9, 9)):
- # The plot that Imre wanted
- n_bins = traj.parameters.input.number_of_directions
- fig, ax = plt.subplots(1, 1, figsize=figsize)
- dir_bins = np.linspace(-np.pi, np.pi, n_bins + 1)
- plot_fr_array = []
- labels = []
- similarity_threshold = np.pi / 6.
- directions = np.linspace(-np.pi, np.pi, traj.input.number_of_directions, endpoint=False)
- for run_idx, run_name in enumerate(plot_run_names):
- label = traj.derived_parameters.runs[run_name].morphology.morph_label
- labels.append(short_labels(label))
- fr_similar_tunings = []
- fr_different_tunings = []
- ex_tunings = traj.results.runs[run_name].ex_tunings
- firing_rate_array = traj.results[run_name].firing_rate_array
- for tuning, firing_rates in zip(ex_tunings, firing_rate_array):
- for idx, dir in enumerate(directions):
- if np.abs(get_phase_difference(tuning - dir)) <= similarity_threshold:
- fr_similar_tunings.append(firing_rates[idx])
- elif np.abs(get_phase_difference(tuning + np.pi - dir)) <= similarity_threshold:
- fr_different_tunings.append(firing_rates[idx])
- plot_fr_array.append([np.mean(fr_similar_tunings), np.mean(fr_different_tunings)])
- x = np.arange(3) # the label locations
- width = 0.35 # the width of the bars
- plot_fr_array = np.array(plot_fr_array)
- # rects1 = ax.bar(x - width / 2, plot_fr_array[:, 0], width,
- # label=r'$\theta_{pref} \pm 30°$')
- # rects2 = ax.bar(x + width / 2, plot_fr_array[:, 1], width,
- # label=r'$\theta_{opp} \pm 30°$')
- rects1 = ax.bar(x - width / 2, plot_fr_array[:, 0], width,
- label=r'sim.')
- rects2 = ax.bar(x + width / 2, plot_fr_array[:, 1], width,
- label=r'diff.')
- ax.set_xticks(x)
- ax.set_xticklabels(labels)
- ax.spines['right'].set_visible(False)
- ax.spines['top'].set_visible(False)
- # ax.set_title('Mean firing rate for tunings similar and different to input')
- # ax.set_ylabel('Mean firing rate')
- ax.annotate(r'$\overline{\mathrm{fr}}$ (Hz)', (0.05, 1.0), xycoords='axes fraction', va="bottom", ha="right")
- leg = ax.legend(loc="upper left", bbox_to_anchor=(0.0, 1.2), handlelength=1, fontsize="medium")
- leg.get_frame().set_linewidth(0.0)
- def autolabel(rects):
- """Attach a text label above each bar in *rects*, displaying its height."""
- for rect in rects:
- height = rect.get_height()
- ax.annotate('{}'.format(np.round(height)),
- xy=(rect.get_x() + rect.get_width() / 2, height),
- xytext=(0, 3), # 3 points vertical offset
- textcoords="offset points",
- ha='center', va='bottom')
- autolabel(rects1)
- autolabel(rects2)
- fig.tight_layout()
- if save_figs:
- plt.savefig(FIGURE_SAVE_PATH + 'SUPPLEMENT_B_firing_rate_similar_vs_diff_tuning' + FIGURE_SAVE_FORMAT, dpi=200)
- plt.close(fig)
- def get_firing_rates_along_preferred_axis(traj, run_name, neuron_idx):
- firing_rates = traj.results[run_name].firing_rate_array[neuron_idx, :]
- tuning = traj.results[run_name].ex_tunings[neuron_idx]
- anti_tuning = tuning + np.pi if tuning + np.pi < np.pi else tuning - np.pi
- tuning_idx = np.argmin(np.abs(directions - tuning))
- anti_tuning_idx = np.argmin(np.abs(directions - anti_tuning))
- firing_at_the_preferred_direction = firing_rates[tuning_idx]
- firing_at_the_opposite_direction = firing_rates[anti_tuning_idx]
- return firing_at_the_preferred_direction, firing_at_the_opposite_direction
- def get_hdi(traj, run_name, neuron_idx, type):
- return traj.results.runs[run_name].head_direction_indices[neuron_idx] if type=="ex" else traj.results.runs[
- run_name].inh_head_direction_indices[neuron_idx]
- def plot_colorbar(figsize=(2, 2), figname=None):
- azimuth_no = 360
- zenith_no = 15
- azimuths = np.linspace(-180, 180, azimuth_no)
- zeniths = np.linspace(0.85, 1, zenith_no)
- values = azimuths * np.ones((zenith_no, azimuth_no))
- fig, ax = plt.subplots(subplot_kw=dict(projection='polar'), figsize=figsize)
- ax.pcolormesh(azimuths * np.pi / 180.0, zeniths, values, cmap=head_direction_input_colormap)
- # ax.set_yticks([])
- ax.set_thetagrids([0, 90, 180, 270])
- ax.tick_params(pad=-2)
- ax.set_ylim(0, 1)
- ax.grid(True)
- y_tick_labels = []
- ax.set_yticklabels(y_tick_labels)
- gridlines = ax.yaxis.get_gridlines()
- [line.set_linewidth(0.0) for line in gridlines]
- gridlines = ax.xaxis.get_gridlines()
- [line.set_linewidth(0.0) for line in gridlines]
- ax.axes.spines["polar"].set_visible(False)
- plt.subplots_adjust(left=0.25, right=0.75, bottom=0.25, top=0.75)
- # plt.show()
- if figname is not None:
- plt.savefig(figname, transparent=True)
- plt.close(fig)
- if __name__ == "__main__":
- traj = Trajectory(TRAJ_NAME, add_time=False, dynamic_imports=Brian2MonitorResult)
- NO_LOADING = 0
- FULL_LOAD = 2
- traj.f_load(filename=os.path.join(DATA_FOLDER, TRAJ_NAME + ".hdf5"), load_parameters=FULL_LOAD,
- load_results=NO_LOADING)
- traj.v_auto_load = True
- save_figs = True
- print("# Plotting script polarized interneurons")
- print()
- map_length_scale = 200.0
- map_seed = 1
- exemplary_head_direction = 0
- print("## Map specifications")
- print("\tinput map scale: {:.1f} um".format(map_length_scale))
- print("\tmap seed: {:d}".format(map_seed))
- print()
- print("## Input specification")
- print("\tselected head direction: {:.0f}°".format(exemplary_head_direction))
- print()
- print("## Selected simulations")
- plot_scale = get_closest_scale(traj, map_length_scale)
- par_dict = {'seed': map_seed, 'scale': plot_scale}
- plot_run_names = filter_run_names_by_par_dict(traj, par_dict)
- run_name_dict = {}
- for run_name in plot_run_names:
- traj.f_set_crun(run_name)
- run_name_dict[traj.derived_parameters.runs[run_name].morphology.morph_label] = run_name
- for network_type, run_name in run_name_dict.items():
- print("{:s}: {:s}".format(network_type, run_name))
- directions = get_input_head_directions(traj)
- direction_idx = np.argmin(np.abs(np.array(directions) - np.deg2rad(exemplary_head_direction)))
- selected_neuron_excitatory = 1052
- selected_inhibitory_neuron = 28
- print("## Figure specification")
- print("\tpanel size: {:.2f} cm".format(panel_size * cm_per_inch))
- print()
- plot_colorbar(figsize=(0.8 * panel_size, 0.8 * panel_size), figname=FIGURE_SAVE_PATH + "A_i_colormap.svg")
- plot_input_map(traj, run_name_dict[POLARIZED], figname="A_i_exemplary_input_map"+FIGURE_SAVE_FORMAT,
- figsize=(panel_size, panel_size))
- # plot_example_input_maps(traj, figsize=(2 * panel_size, 2 * panel_size))
- plot_axonal_clouds(traj, plot_run_names)
- #
- plot_firing_rate_map_excitatory(traj, direction_idx, plot_run_names, selected_neuron_excitatory)
- in_max_rate = plot_firing_rate_map_inhibitory(traj, direction_idx, plot_run_names, selected_inhibitory_neuron)
- #
- hdi_means = plot_hdi_histogram_combined_and_overlayed(
- traj, plot_run_names,
- selected_neuron_excitatory,
- selected_inhibitory_neuron,
- cut_off_dist=100.)
- #
- plot_polar_plot_excitatory(traj, plot_run_names, selected_neuron_excitatory)
- plot_polar_plot_inhibitory(traj, plot_run_names, selected_inhibitory_neuron)
- plot_firing_rate_similar_vs_diff_tuning(traj, plot_run_names, figsize=(1.2*panel_size, 1.2*panel_size))
- plot_exc_and_inh_hdi_over_simplex_grid_scale(traj, traj.f_get_run_names(), cut_off_dist=100.)
- plot_exc_and_inh_hdi_over_fit_corr_len(traj, traj.f_get_run_names(), cut_off_dist=100.)
- if not save_figs:
- plt.show()
- traj.f_restore_default()
|