Added some cmap functions.

Selected species-specific  colors.
Quite some progress on fig_invariance_thresh_lp_species.pdf.

python/color_functions.py

@@ -1,5 +1,6 @@
 import numpy as np
 import matplotlib.pyplot as plt
+from matplotlib.colors import ListedColormap, LinearSegmentedColormap
 from tkinter.colorchooser import askcolor
 from IPython import embed
 
@@ -81,6 +82,57 @@ def load_colors(path):
         return {k: (c.item() if c.size == 1 else c) for k, c in colors.items()}
     raise ValueError(f'Expected .npy or .npz file extension: {path}')
 
+# COLORMAPS:
+
+def create_listed_cmap(colors, name=None, n=None):
+    cmap = ListedColormap(colors)
+    if n is not None:
+        cmap.resampled(n)
+    if name is not None:
+        cmap.name = name
+        plt.colormaps.register(cmap)
+    return cmap
+
+def create_linear_cmap(colors, name=None, n=None):
+    cmap = LinearSegmentedColormap.from_list(colors)
+    if n is not None:
+        cmap.resampled(n)
+    if name is not None:
+        cmap.name = name
+        plt.colormaps.register(cmap)
+    return cmap
+
+def sample_cmap(cmap, n, low=None, high=None, segments=None, alpha=None):
+    if isinstance(cmap, str):
+        cmap = plt.get_cmap(cmap)
+    colors = cmap(np.linspace(0, 1, n))
+
+    if alpha is None:
+        colors = colors[:, :3]
+    elif 0.0 <= alpha <= 1.0:
+        colors[:, 3] = alpha
+
+    if segments is None and (low is not None or high is not None):
+        segments = [(0 if low is None else low, 1 if high is None else high)]
+    if segments is not None:
+        segment_colors = []
+        for start, end in segments:
+            start, end = int(start * n), int(end * n)
+            step = 1 if start <= end else -1
+            segment_colors.append(colors[start:end:step, :])
+        colors = np.vstack(segment_colors)
+    return colors
+
+def remake_cmap(cmap, n_in, n_out=None, name=None, low=None, high=None, segments=None,
+                alpha=None):
+    colors = sample_cmap(cmap, n_in, low, high, segments, alpha)
+    cmap_type = type(cmap).__name__
+    if cmap_type == 'ListedColormap':
+        return create_listed_cmap(colors, name, n_out)
+    elif cmap_type == 'LinearSegmentedColormap':
+        return create_linear_cmap(colors, name, n_out)
+    return None
+
 # ADVANCED FUNCTIONALITY:
 
 def shade_colors(color, factors, norm=True):

python/fig_invariance_log-hp.py

@@ -30,6 +30,7 @@ def plot_snippets(axes, time, snippets, ymin=None, ymax=None, **kwargs):
 # GENERAL SETTINGS:
 target = 'Omocestus_rufipes'
 data_paths = search_files(target, excl='noise', dir='../data/inv/log_hp/')
+species_paths = search_files('*', incl='noise', dir='../data/inv/log_hp/')
 stages = ['env', 'log', 'inv']
 load_kwargs = dict(
     files=stages,
@@ -39,10 +40,6 @@ save_path = '../figures/fig_invariance_log_hp.pdf'
 compute_ratios = True
 show_diag = True
 show_noise = True
-if compute_ratios:
-    ref_data = load_data('../data/processed/white_noise_sd-1.npz', files=stages)[0]
-    ref_measures = {k: v.std() for k, v in ref_data.items() if not k.endswith('rate')}
-
 
 # GRAPH SETTINGS:
 fig_kwargs = dict(
@@ -221,6 +218,20 @@ noise_kwargs = dict(
     zorder=1.5,
 )
 
+# PREPARATION:
+if compute_ratios:
+    ref_data = load_data('../data/processed/white_noise_sd-1.npz', files=stages)[0]
+    ref_measures = {k: v.std() for k, v in ref_data.items() if not k.endswith('rate')}
+
+species_measures = []
+for species_path in species_paths:
+    species_data, _ = load_data(species_path, **load_kwargs)
+    species_measure = species_data['measure_inv']
+    if compute_ratios:
+        species_measure /= ref_measures['inv']
+    species_measures.append(species_measure)
+species_measures = np.array(species_measures).T
+
 # EXECUTION:
 for data_path in data_paths:
     print(f'Processing {data_path}')
@@ -340,6 +351,9 @@ for data_path in data_paths:
     big_axes[1].plot(noise_scales, noise_data['measure_log'], c=colors['log'], lw=lw_big)
     big_axes[1].plot(noise_scales, noise_data['measure_inv'], c=colors['inv'], lw=lw_big)
 
+    # Plot species measures:
+    big_axes[1].plot(noise_scales, species_measures, 'k', lw=lw_big, zorder=2.1)
+
     if show_diag:
         # Indicate diagonal:
         big_axes[0].plot(pure_scales, pure_scales, **diag_kwargs)

python/fig_invariance_log-hp_backup.py

@@ -0,0 +1,380 @@
+import plotstyle_plt
+import numpy as np
+import matplotlib.pyplot as plt
+from itertools import product
+from thunderhopper.filetools import search_files
+from thunderhopper.modeltools import load_data
+from color_functions import load_colors
+from plot_functions import hide_axis, ylimits, xlabel, ylabel, hide_ticks,\
+                           plot_line, strip_zeros, time_bar, zoom_inset,\
+                           letter_subplot, title_subplot
+from IPython import embed
+
+def add_snip_axes(fig, grid_kwargs):
+    grid = fig.add_gridspec(**grid_kwargs)
+    axes = np.zeros((grid.nrows, grid.ncols), dtype=object)
+    for i, j in product(range(grid.nrows), range(grid.ncols)):
+        axes[i, j] = fig.add_subplot(grid[i, j])
+    [hide_axis(ax, 'left') for ax in axes[:, 1:].flatten()]
+    [hide_axis(ax, 'bottom') for ax in axes.flatten()]
+    return axes
+
+def plot_snippets(axes, time, snippets, ymin=None, ymax=None, **kwargs):
+    ymin, ymax = ylimits(snippets, minval=ymin, maxval=ymax, pad=0.05)
+    handles = []
+    for ax, snippet in zip(axes, snippets.T):
+        handles.extend(plot_line(ax, time, snippet, ymin=ymin, ymax=ymax, **kwargs))
+    return handles
+
+
+# GENERAL SETTINGS:
+target = 'Omocestus_rufipes'
+data_paths = search_files(target, excl='noise', dir='../data/inv/log_hp/')
+species_paths = search_files('*', incl='noise', dir='../data/inv/log_hp/')
+stages = ['env', 'log', 'inv']
+load_kwargs = dict(
+    files=stages,
+    keywords=['scales', 'snip', 'measure']
+)
+save_path = '../figures/fig_invariance_log_hp.pdf'
+compute_ratios = True
+show_diag = True
+show_noise = True
+
+# GRAPH SETTINGS:
+fig_kwargs = dict(
+    figsize=(32/2.54, 32/2.54),
+)
+snip_rows = 1
+big_rows = 1
+super_grid_kwargs = dict(
+    nrows=2 * snip_rows + big_rows,
+    ncols=1,
+    wspace=0,
+    hspace=0,
+    left=0,
+    right=1,
+    bottom=0,
+    top=1
+)
+subfig_specs = dict(
+    pure=(slice(0, snip_rows), slice(None)),
+    noise=(slice(snip_rows, 2 * snip_rows), slice(None)),
+    big=(slice(-big_rows, None), slice(None)),
+)
+block_height = 0.8
+edge_padding = 0.08
+pure_grid_kwargs = dict(
+    nrows=len(stages),
+    ncols=None,
+    wspace=0.1,
+    hspace=0.15,
+    left=0.11,
+    right=0.95,
+    bottom=1 - block_height - edge_padding,
+    top=1 - edge_padding,
+    height_ratios=[1, 2, 1]
+)
+noise_grid_kwargs = dict(
+    nrows=len(stages),
+    ncols=None,
+    wspace=pure_grid_kwargs['wspace'],
+    hspace=pure_grid_kwargs['hspace'],
+    left=pure_grid_kwargs['left'],
+    right=pure_grid_kwargs['right'],
+    bottom=edge_padding,
+    top=edge_padding + block_height,
+    height_ratios=[1, 2, 1]
+)
+big_grid_kwargs = dict(
+    nrows=1,
+    ncols=3,
+    wspace=0.3,
+    hspace=0,
+    left=pure_grid_kwargs['left'],
+    right=pure_grid_kwargs['right'],
+    bottom=0.05,
+    top=1
+)
+anchor_kwargs = dict(
+    aspect='equal',
+    adjustable='box',
+    anchor=(0.5, 0.5)
+)
+
+# PLOT SETTINGS:
+fs = dict(
+    lab_norm=16,
+    lab_tex=20,
+    letter=22,
+    tit_norm=16,
+    tit_tex=20,
+    bar=16,
+)
+colors = load_colors('../data/stage_colors.npz')
+lw_snippets = 1
+lw_big = 3
+xlabels = dict(
+    big='scale $\\alpha$',
+)
+ylabels = dict(
+    env='$x_{\\text{env}}$',
+    log='$x_{\\text{dB}}$',
+    inv='$x_{\\text{adapt}}$',
+    big='$\\sigma_{\\alpha}\\,/\\,\\sigma_{\\eta}$',
+)
+xlab_big_kwargs = dict(
+    y=0,
+    fontsize=fs['lab_norm'],
+    ha='center',
+    va='bottom',
+)
+ylab_snip_kwargs = dict(
+    x=0,
+    fontsize=fs['lab_tex'],
+    rotation=0,
+    ha='left',
+    va='center',
+)
+ylab_big_kwargs = dict(
+    x=0,
+    fontsize=fs['lab_tex'],
+    ha='center',
+    va='top',
+)
+yloc = dict(
+    env=1000,
+    log=40,
+    inv=20
+)
+title_kwargs = dict(
+    x=0.5,
+    y=1,
+    ha='center',
+    va='bottom',
+    fontsize=fs['tit_norm'],
+)
+letter_snip_kwargs = dict(
+    x=0,
+    yref=0.5,
+    ha='left',
+    va='center',
+    fontsize=fs['letter'],
+)
+letter_big_kwargs = dict(
+    x=0,
+    y=1,
+    ha='left',
+    va='bottom',
+    fontsize=fs['letter'],
+)
+zoom_inset_bounds = [0.1, 0.2, 0.8, 0.6]
+zoom_kwargs = dict(
+    x0=0.45,
+    x1=0.55,
+    y0=0,
+    y1=0.0006,
+    low_left=True,
+    low_right=True,
+    ec='k',
+    lw=1,
+    alpha=1,
+)
+inset_tick_kwargs = dict(
+    axis='y',
+    length=3,
+    pad=1,
+    left=False,
+    labelleft=False,
+    right=True,
+    labelright=True,
+)
+bar_time = 5
+bar_kwargs = dict(
+    dur=bar_time,
+    y0=-0.25,
+    y1=-0.1,
+    xshift=1,
+    color='k',
+    lw=0,
+    clip_on=False,
+    text_pos=(-0.1, 0.5),
+    text_str=f'${bar_time}\\,\\text{{s}}$',
+    text_kwargs=dict(
+        fontsize=fs['bar'],
+        ha='right',
+        va='center',
+    )
+)
+diag_kwargs = dict(
+    c=(0.75, 0.75, 0.75),
+    lw=2,
+    ls='--',
+    zorder=1.9,
+)
+noise_rel_thresh = 0.95
+noise_kwargs = dict(
+    fc=(0.9, 0.9, 0.9),
+    ec='none',
+    lw=0,
+    zorder=1.5,
+)
+
+# PREPARATION:
+if compute_ratios:
+    ref_data = load_data('../data/processed/white_noise_sd-1.npz', files=stages)[0]
+    ref_measures = {k: v.std() for k, v in ref_data.items() if not k.endswith('rate')}
+
+species_measures = []
+for species_path in species_paths:
+    species_measure = load_data(species_path, **load_kwargs)[0]['measure_inv']
+    if compute_ratios:
+        species_measure /= ref_measures['inv']
+    species_measures.append(species_measure)
+species_measures = np.array(species_measures).T
+
+# EXECUTION:
+for data_path in data_paths:
+    print(f'Processing {data_path}')
+
+    # Load invariance data:
+    pure_data, config = load_data(data_path, **load_kwargs)
+    noise_data, _ = load_data(data_path.replace('.npz', '_noise.npz'), **load_kwargs)
+    pure_scales, noise_scales = pure_data['scales'], noise_data['scales']
+    t_full = np.arange(pure_data['snip_env'].shape[0]) / config['env_rate']
+
+    # Prepare overall graph:
+    fig = plt.figure(**fig_kwargs)
+    super_grid = fig.add_gridspec(**super_grid_kwargs)
+    fig.canvas.draw()
+
+    # Prepare pure-song snippet axes:
+    pure_grid_kwargs['ncols'] = pure_data['example_scales'].size
+    pure_subfig = fig.add_subfigure(super_grid[subfig_specs['pure']])
+    pure_axes = add_snip_axes(pure_subfig, pure_grid_kwargs)
+    for ax, stage in zip(pure_axes[:, 0], stages):
+        ax.yaxis.set_major_locator(plt.MultipleLocator(yloc[stage]))
+        ylabel(ax, ylabels[stage], **ylab_snip_kwargs,
+               transform=pure_subfig.transSubfigure)
+    for ax, scale in zip(pure_axes[0, :], pure_data['example_scales']):
+        pure_title = title_subplot(ax, f'$\\alpha={strip_zeros(scale)}$', **title_kwargs)
+    letter_subplot(pure_subfig, 'a', ref=pure_title, **letter_snip_kwargs)
+    pure_inset = pure_axes[0, 0].inset_axes(zoom_inset_bounds)
+    pure_inset.spines[:].set(visible=True, lw=zoom_kwargs['lw'])
+    pure_inset.tick_params(**inset_tick_kwargs)
+    hide_ticks(pure_inset, 'bottom', ticks=False)
+
+    # Prepare noise-song snippet axes:
+    noise_grid_kwargs['ncols'] = noise_data['example_scales'].size
+    noise_subfig = fig.add_subfigure(super_grid[subfig_specs['noise']])
+    noise_axes = add_snip_axes(noise_subfig, noise_grid_kwargs)
+    for ax, stage in zip(noise_axes[:, 0], stages):
+        ax.yaxis.set_major_locator(plt.MultipleLocator(yloc[stage]))
+        ylabel(ax, ylabels[stage], **ylab_snip_kwargs,
+               transform=noise_subfig.transSubfigure)
+    for ax, scale in zip(noise_axes[0, :], noise_data['example_scales']):
+        noise_title = title_subplot(ax, f'$\\alpha={strip_zeros(scale)}$', **title_kwargs)
+    letter_subplot(noise_subfig, 'b', ref=noise_title, **letter_snip_kwargs)
+    noise_inset = noise_axes[0, 0].inset_axes(zoom_inset_bounds)
+    noise_inset.spines[:].set(visible=True, lw=zoom_kwargs['lw'])
+    noise_inset.tick_params(**inset_tick_kwargs)
+    hide_ticks(noise_inset, 'bottom', ticks=False)
+
+    # Prepare analysis axes:
+    big_subfig = fig.add_subfigure(super_grid[subfig_specs['big']])
+    big_grid = big_subfig.add_gridspec(**big_grid_kwargs)
+    big_axes = np.zeros((big_grid.ncols,), dtype=object)
+    for i, scales in enumerate([pure_scales, noise_scales, noise_scales]):
+        ax = big_subfig.add_subplot(big_grid[0, i])
+        ax.set_xlim(scales[0], scales[-1])
+        ax.set_ylim(scales[0], scales[-1])
+        ax.set_xscale('symlog', linthresh=scales[1], linscale=0.5)
+        ax.set_yscale('symlog', linthresh=scales[1], linscale=0.5)
+        ax.set_aspect(**anchor_kwargs)
+        ylabel(ax, ylabels['big'], transform=big_subfig.transSubfigure, **ylab_big_kwargs)
+        if i == 0:
+            hide_ticks(ax, 'bottom')
+            letter_subplot(ax, 'c', **letter_big_kwargs)
+        else:
+            xlabel(ax, xlabels['big'], transform=big_subfig.transSubfigure, **xlab_big_kwargs)
+            letter_subplot(ax, 'd', **letter_big_kwargs)
+        big_axes[i] = ax
+
+    # Plot pure-song envelope snippets:
+    handle = plot_snippets(pure_axes[0, :], t_full, pure_data['snip_env'],
+                           ymin=0, c=colors['env'], lw=lw_snippets)[0]
+    zoom_inset(pure_axes[0, 0], pure_inset, handle, transform=pure_axes[0, 0].transAxes, **zoom_kwargs)
+
+    # Plot pure-song logarithmic snippets:
+    plot_snippets(pure_axes[1, :], t_full, pure_data['snip_log'],
+                  c=colors['log'], lw=lw_snippets)
+
+    # Plot pure-song invariant snippets:
+    plot_snippets(pure_axes[2, :], t_full, pure_data['snip_inv'],
+                  c=colors['inv'], lw=lw_snippets)
+
+    # Plot noise-song envelope snippets:
+    ymin, ymax = pure_axes[0, 0].get_ylim()
+    handle = plot_snippets(noise_axes[0, :], t_full, noise_data['snip_env'],
+                           ymin, ymax, c=colors['env'], lw=lw_snippets)[0]
+    zoom_inset(noise_axes[0, 0], noise_inset, handle, transform=noise_axes[0, 0].transAxes, **zoom_kwargs)
+
+    # Plot noise-song logarithmic snippets:
+    ymin, ymax = pure_axes[1, 0].get_ylim()
+    plot_snippets(noise_axes[1, :], t_full, noise_data['snip_log'],
+                  ymin, ymax, c=colors['log'], lw=lw_snippets)
+
+    # Plot noise-song invariant snippets:
+    ymin, ymax = pure_axes[2, 0].get_ylim()
+    plot_snippets(noise_axes[2, :], t_full, noise_data['snip_inv'],
+                  ymin, ymax, c=colors['inv'], lw=lw_snippets)
+
+    # Indicate time scale:
+    time_bar(noise_axes[-1, -1], **bar_kwargs)
+
+    if compute_ratios:
+        # Relate pure-song measures to zero scale:
+        pure_data['measure_env'] /= ref_measures['env']
+        pure_data['measure_log'] /= ref_measures['log']
+        pure_data['measure_inv'] /= ref_measures['inv']
+        # Relate noise-song measures to zero scale:
+        noise_data['measure_env'] /= ref_measures['env']
+        noise_data['measure_log'] /= ref_measures['log']
+        noise_data['measure_inv'] /= ref_measures['inv']
+
+    # Plot pure-song measures (ideal):
+    big_axes[0].plot(pure_scales, pure_data['measure_env'], c=colors['env'], lw=lw_big)
+    big_axes[0].plot(pure_scales, pure_data['measure_log'], c=colors['log'], lw=lw_big)
+    big_axes[0].plot(pure_scales, pure_data['measure_inv'], c=colors['inv'], lw=lw_big)
+
+    # Plot noise-song measures (limited):
+    big_axes[1].plot(noise_scales, noise_data['measure_env'], c=colors['env'], lw=lw_big)
+    big_axes[1].plot(noise_scales, noise_data['measure_log'], c=colors['log'], lw=lw_big)
+    big_axes[1].plot(noise_scales, noise_data['measure_inv'], c=colors['inv'], lw=lw_big)
+
+    # Plot species measures:
+    big_axes[2].plot(noise_scales, species_measures, 'k', lw=lw_big)
+
+    if show_diag:
+        # Indicate diagonal:
+        big_axes[0].plot(pure_scales, pure_scales, **diag_kwargs)
+        big_axes[1].plot(noise_scales, noise_scales, **diag_kwargs)
+
+    if show_noise:
+        # Indicate noise floor:
+        if compute_ratios:
+            span_measure = noise_data['measure_inv'][-1] - ref_measures['inv']
+            thresh_measure = ref_measures['inv'] + noise_rel_thresh * span_measure
+        else:     
+            span_measure = noise_data['measure_inv'][-1] - noise_data['measure_inv'][0]
+            thresh_measure = noise_data['measure_inv'][0] + noise_rel_thresh * span_measure
+        thresh_ind = np.nonzero(noise_data['measure_inv'] < thresh_measure)[0][-1]
+        thresh_scale = noise_scales[thresh_ind]
+        big_axes[1].axvspan(noise_scales[0], thresh_scale, **noise_kwargs)
+
+    if save_path is not None:
+        fig.savefig(save_path, bbox_inches='tight')
+    plt.show()
+
+print('Done.')
+embed()

python/fig_invariance_thresh-lp_species.py

@@ -6,8 +6,8 @@ from itertools import product
 from thunderhopper.filetools import search_files
 from thunderhopper.modeltools import load_data
 from thunderhopper.filtertools import find_kern_specs
-from color_functions import load_colors, shade_colors
+from color_functions import load_colors, shade_colors, create_listed_cmap
-from plot_functions import hide_axis, ylimits, xlabel, ylabel, super_ylabel,\
+from plot_functions import hide_axis, title_subplot, ylimits, xlabel, ylabel, super_ylabel,\
                            plot_line, plot_barcode, strip_zeros, time_bar,\
                            letter_subplot, letter_subplots, hide_ticks,\
                            super_xlabel, super_ylabel, assign_colors
@@ -125,73 +125,90 @@ def split_subplot(ax, side='right', size=10, pad=10):
     inputs = zip(*force_sequence(side, size, pad, equal_size=True))
     return [div.append_axes(s, f'{n}%', f'{p}%') for s, n, p in inputs]
 
+def shorten_species(name):
+    genus, species = name.split('_')
+    return genus[0] + '. ' + species
 
 # GENERAL SETTINGS:
-targets = [
+target_species = [
     'Omocestus_rufipes',
     'Chorthippus_biguttulus',
-    # 'Chorthippus_mollis',
+    'Chorthippus_mollis',
-    # 'Chrysochraon_dispar',
+    'Chrysochraon_dispar',
     'Gomphocerippus_rufus',
-    # 'Pseudochorthippus_parallelus',
+    'Pseudochorthippus_parallelus',
     ]
-pure_paths = search_files(targets, incl='subset', excl='noise', dir='../data/inv/thresh_lp/')
+n_species = len(target_species)
 load_kwargs = dict(
     keywords=['scales', 'measure', 'thresh']
 )
 save_path = '../figures/fig_invariance_thresh_lp_species.pdf'
+exclude_zero = True
+show_noise = True
 
 # SUBSET SETTINGS:
-thresh_percent = np.array([0.6, 0.75, 0.999])[0]
+thresh_rel = np.array([0.5, 1, 3])[0]
-kernels = np.array([
+kern_specs = np.array([
     [1, 0.008],
     [2, 0.004],
     [3, 0.002],
 ])[np.array([0, 1])]
+n_kernels = kern_specs.shape[0]
 
 # GRAPH SETTINGS:
 fig_kwargs = dict(
-    figsize=(32/2.54, 16/2.54),
+    figsize=(32/2.54, 20/2.54),
 )
-n_species = len(targets)
 super_grid_kwargs = dict(
-    nrows=2,
+    nrows=3,
-    ncols=n_species + 2,
+    ncols=1,
     wspace=0,
     hspace=0,
     left=0,
     right=1,
     bottom=0,
-    top=1
+    top=1,
+    height_ratios=[1, 4, 3]
 )
 subfig_specs = dict(
-    spec=(slice(None), slice(0, n_species)),
+    song=(0, 0),
-    big=(slice(None), slice(n_species, None))
+    feat=(1, 0),
+    space=(2, 0)
 )
-spec_grid_kwargs = dict(
+feat_grid_kwargs = dict(
     nrows=2,
     ncols=n_species,
     wspace=0.25,
-    hspace=0.1,
+    hspace=0.15,
-    left=0.1,
+    left=0.06,
-    right=0.97,
+    right=0.985,
     bottom=0.1,
     top=0.94
 )
-big_grid_kwargs = dict(
+song_grid_kwargs = dict(
-    nrows=2,
+    nrows=1,
-    ncols=1,
+    ncols=n_species,
-    wspace=0,
+    wspace=feat_grid_kwargs['wspace'],
-    hspace=0.2,
+    hspace=0,
-    left=0,
+    left=feat_grid_kwargs['left'],
-    right=1,
+    right=feat_grid_kwargs['right'],
-    bottom=spec_grid_kwargs['bottom'],
+    bottom=0.1,
-    top=spec_grid_kwargs['top']
+    top=0.8
+)
+space_grid_kwargs = dict(
+    nrows=1,
+    ncols=2,
+    wspace=0.2,
+    hspace=0,
+    left=feat_grid_kwargs['left'],
+    right=feat_grid_kwargs['right'],
+    bottom=0.05,
+    top=0.95
 )
 anchor_kwargs = dict(
     aspect='equal',
     adjustable='box',
-    anchor=(0.3, 0.5)
+    anchor=(0, 0.5)
 )
 inset_kwargs = dict(
     y0=0.7,
@@ -208,50 +225,56 @@ fs = dict(
     tit_tex=20,
     bar=16,
 )
-base_color = load_colors('../data/stage_colors.npz')['feat']
+species_colors = load_colors('../data/species_colors.npz')
-spec_cmaps = [
+kernel_shades = [0, 0.5]
-    'Reds',
+# scale_shades = [1, 0]
-    'Greens',
-    'Blues',
-]
 lw = dict(
-    spec=2,
+    song=0.5,
+    feat=3,
     kern=3
 )
+zorder = dict(
+    Omocestus_rufipes=2,
+    Chorthippus_biguttulus=2.5,
+    Chorthippus_mollis=2.4,
+    Chrysochraon_dispar=2,
+    Gomphocerippus_rufus=2,
+    Pseudochorthippus_parallelus=2,
+)
 space_kwargs = dict(
     s=30,
 )
 xlabels = dict(
-    spec='scale $\\alpha$',
+    feat='scale $\\alpha$',
-    big='$\\mu_{f_1}$'
+    space='$\\mu_{f_1}$'
 )
 ylabels = dict(
-    spec='$\\mu_f$',
+    feat='$\\mu_f$',
-    big='$\\mu_{f_2}$',
+    space='$\\mu_{f_2}$',
     bar='scale $\\alpha$',
 )
-xlab_spec_kwargs = dict(
+xlab_feat_kwargs = dict(
     y=0,
     fontsize=fs['lab_norm'],
     ha='center',
     va='bottom',
 )
-xlab_big_kwargs = dict(
+xlab_space_kwargs = dict(
     y=0,
     fontsize=fs['lab_tex'],
     ha='center',
     va='bottom',
 )
-ylab_spec_kwargs = dict(
+ylab_feat_kwargs = dict(
     x=0,
     fontsize=fs['lab_tex'],
     ha='left',
     va='center',
 )
-ylab_big_kwargs = dict(
+ylab_space_kwargs = dict(
-    x=0.03,
+    x=0,
     fontsize=fs['lab_tex'],
-    ha='center',
+    ha='left',
     va='center',
 )
 ylab_cbar_kwargs = dict(
@@ -261,28 +284,57 @@ ylab_cbar_kwargs = dict(
     va='bottom',
 )
 xloc = dict(
-    big=0.5,
+    space=0.5,
 )
 yloc = dict(
-    spec=0.5,
+    feat=0.5,
-    big=0.5
+    space=0.5
 )
-letter_spec_kwargs = dict(
+symlog_kwargs = dict(
+    linscale=0.5,
+)
+title_kwargs = dict(
+    x=0.5,
+    yref=1,
+    ha='center',
+    va='top',
+    fontsize=fs['tit_norm'],
+    fontstyle='italic'
+)
+letter_feat_kwargs = dict(
     x=0,
     yref=1,
     ha='center',
     va='top',
     fontsize=fs['letter'],  
 )
-letter_big_kwargs = dict(
+letter_space_kwargs = dict(
     x=0,
     yref=1,
     ha='center',
     va='top',
     fontsize=fs['letter'],
 )
-time_bar_kwargs = dict(
+song_bar_time = 1.0
-    dur=0.05,
+song_bar_kwargs = dict(
+    dur=song_bar_time,
+    y0=-0.1,
+    y1=0,
+    xshift=0,
+    color='k',
+    lw=0,
+    clip_on=False,
+    # text_pos=(-0.1, 0.5),
+    text_str=f'${int(1000 * song_bar_time)}\\,\\text{{ms}}$',
+    text_kwargs=dict(
+        fontsize=fs['bar'],
+        ha='right',
+        va='center',
+    )
+)
+kern_bar_time = 0.05
+kern_bar_kwargs = dict(
+    dur=kern_bar_time,
     y0=inset_kwargs['y0'],
     y1=inset_kwargs['y0'] + 0.03,
     color='k',
@@ -290,11 +342,16 @@ time_bar_kwargs = dict(
 )
 cbar_bounds = [
     0.05,
-    big_grid_kwargs['bottom'],
+    space_grid_kwargs['bottom'],
     0.15,
-    big_grid_kwargs['top'] - big_grid_kwargs['bottom']
+    space_grid_kwargs['top'] - space_grid_kwargs['bottom']
 ]
-shade_factors = [0.9, -0.9]
+noise_kwargs = dict(
+    fc=(0.9, 0.9, 0.9),
+    ec='none',
+    lw=0,
+    zorder=0.5,
+)
 
 # EXECUTION:
 
@@ -302,105 +359,165 @@ shade_factors = [0.9, -0.9]
 fig = plt.figure(**fig_kwargs)
 super_grid = fig.add_gridspec(**super_grid_kwargs)
 
-# Prepare species-specific axes:
+# Prepare song axes:
-spec_subfig = fig.add_subfigure(super_grid[subfig_specs['spec']])
+song_subfig = fig.add_subfigure(super_grid[subfig_specs['song']])
-spec_grid = spec_subfig.add_gridspec(**spec_grid_kwargs)
+song_grid = song_subfig.add_gridspec(**song_grid_kwargs)
-spec_axes = np.zeros((spec_grid_kwargs['nrows'], n_species), dtype=object)
+song_axes = np.zeros((n_species,), dtype=object)
-for i, j in product(range(spec_grid_kwargs['nrows']), range(n_species)):
+for i in range(n_species):
-    ax = spec_subfig.add_subplot(spec_grid[i, j])
+    ax = song_subfig.add_subplot(song_grid[i])
-    ax.set_xscale('symlog', linthresh=0.1, linscale=0.5)
+    hide_axis(ax, 'bottom')
-    ax.yaxis.set_major_locator(plt.MultipleLocator(yloc['spec']))
+    hide_axis(ax, 'left')
+    song_axes[i] = ax
+
+# Prepare feature invariance axes:
+feat_subfig = fig.add_subfigure(super_grid[subfig_specs['feat']])
+feat_grid = feat_subfig.add_gridspec(**feat_grid_kwargs)
+feat_axes = np.zeros((feat_grid_kwargs['nrows'], n_species), dtype=object)
+for i, j in product(range(feat_grid_kwargs['nrows']), range(n_species)):
+    ax = feat_subfig.add_subplot(feat_grid[i, j])
+    ax.yaxis.set_major_locator(plt.MultipleLocator(yloc['feat']))
     ax.set_ylim(0, 1)
-    spec_axes[i, j] = ax
+    feat_axes[i, j] = ax
-super_xlabel(xlabels['spec'], spec_subfig, spec_axes[-1, 0], spec_axes[-1, -1], **xlab_spec_kwargs)
+super_xlabel(xlabels['feat'], feat_subfig, feat_axes[-1, 0], feat_axes[-1, -1], **xlab_feat_kwargs)
-super_ylabel(ylabels['spec'], spec_subfig, spec_axes[-1, 0], spec_axes[0, 0], **ylab_spec_kwargs)
+super_ylabel(ylabels['feat'], feat_subfig, feat_axes[-1, 0], feat_axes[0, 0], **ylab_feat_kwargs)
-[hide_ticks(ax, side='bottom') for ax in spec_axes[0, :]]
+[hide_ticks(ax, side='bottom') for ax in feat_axes[0, :]]
-[hide_ticks(ax, side='left') for ax in spec_axes[:, 1:].ravel()]
+[hide_ticks(ax, side='left') for ax in feat_axes[:, 1:].ravel()]
-letter_subplots(spec_axes[0, :], labels='abc', ref=spec_subfig, **letter_spec_kwargs)
+letter_subplots(feat_axes[0, :], labels='abc', ref=feat_subfig, **letter_feat_kwargs)
 
 # Prepare kernel insets:
-x0 = np.linspace(0, 1, kernels.shape[0] + 1)[:-1] + 1 / kernels.shape[0] / 2
+x0 = np.linspace(0, 1, n_kernels + 1)[:-1] + 1 / n_kernels / 2
 x0 -= inset_kwargs['w'] / 2
 insets = []
-for i in range(kernels.shape[0]):
+for i in range(n_kernels):
     bounds = [x0[i], inset_kwargs['y0'], inset_kwargs['w'], inset_kwargs['h']]
-    inset = spec_axes[0, 0].inset_axes(bounds)
+    inset = feat_axes[0, 0].inset_axes(bounds)
     inset.set_title(rf'$k_{{{i+1}}}$', fontsize=20)
     inset.axis('off')
     insets.append(inset)
 
 # Prepare feature space axes:
-big_subfig = fig.add_subfigure(super_grid[subfig_specs['big']])
+space_subfig = fig.add_subfigure(super_grid[subfig_specs['space']])
-big_grid = big_subfig.add_gridspec(**big_grid_kwargs)
+space_grid = space_subfig.add_gridspec(**space_grid_kwargs)
-big_axes = np.zeros(super_grid_kwargs['nrows'], dtype=object)
+space_axes = np.zeros(space_grid_kwargs['ncols'], dtype=object)
-for i in range(big_axes.size):
+for i in range(space_axes.size):
-    ax = big_subfig.add_subplot(big_grid[i, 0])
+    ax = space_subfig.add_subplot(space_grid[i])
     ax.set_xlim(0, 1)
     ax.set_ylim(0, 1)
-    ax.xaxis.set_major_locator(plt.MultipleLocator(xloc['big']))
+    ax.xaxis.set_major_locator(plt.MultipleLocator(xloc['space']))
-    ax.yaxis.set_major_locator(plt.MultipleLocator(yloc['big']))
+    ax.yaxis.set_major_locator(plt.MultipleLocator(yloc['space']))
     ax.set_aspect(**anchor_kwargs)
-    # ax.set_ylabel(ylabels['big'], **ylab_big_kwargs)
+    # ax.set_ylabel(ylabels['space'], **ylab_space_kwargs)
-    ylabel(ax, ylabels['big'], transform=big_subfig.transSubfigure, **ylab_big_kwargs)
+    ylabel(ax, ylabels['space'], transform=space_subfig.transSubfigure, **ylab_space_kwargs)
-    big_axes[i] = ax
+    space_axes[i] = ax
-super_xlabel(xlabels['big'], big_subfig, big_axes[1], big_axes[1], **xlab_big_kwargs)
+super_xlabel(xlabels['space'], space_subfig, space_axes[1], space_axes[1], **xlab_space_kwargs)
-hide_ticks(big_axes[0], side='bottom')
+hide_ticks(space_axes[0], side='bottom')
-letter_subplot(big_axes[0], 'd', ref=big_subfig, **letter_big_kwargs)
+letter_subplot(space_axes[0], 'd', ref=space_subfig, **letter_space_kwargs)
 
 # Prepare colorbars:
-cbar_bounds[0] += big_axes[-1].get_position().x1
+cbar_bounds[0] += space_axes[-1].get_position().x1
-bar_axes = [big_subfig.add_axes(cbar_bounds)]
+bar_axes = [space_subfig.add_axes(cbar_bounds)]
-bar_axes.extend(split_subplot(bar_axes[0], side=['right', 'right'], size=100, pad=0))
+bar_axes.extend(split_subplot(bar_axes[0], side=['right'] * (n_species - 1),
+                              size=100, pad=0))
+
+# Prepare kernel-specific color shading:
+kern_factors = np.linspace(*kernel_shades, n_kernels)
+kern_colors_bw = shade_colors((0., 0., 0.), kern_factors)
 
 # Plot results per species:
-for i, pure_path in enumerate(pure_paths):
+min_feat = np.zeros((n_species, n_kernels), dtype=float)
-    print(f'Processing {pure_path}')
+for i, species in enumerate(target_species):
-    noise_path = pure_path.replace('.npz', '_noise.npz')
+    print(f'Processing {species}')
+
+    # Fetch species-specific recording file:
+    song_path = search_files(species, dir='../data/processed/')[0]
+
+    # Load song data:
+    song_data, _ = load_data(song_path, files='filt')
+    song, rate = song_data['filt'], song_data['filt_rate']
+
+    # Plot species snippet:
+    song_ax = song_axes[i]
+    time = np.arange(song.shape[0]) / rate
+    plot_line(song_ax, time, song, ypad=0.05, c='k', lw=lw['song'])
+    title_subplot(song_ax, shorten_species(species), ref=song_subfig, **title_kwargs)
+    time_bar(song_ax, **song_bar_kwargs)
+
+    # Fetch species-specific invariance files:
+    pure_path = search_files(species, incl='pure', dir='../data/inv/thresh_lp/')[0]
+    noise_path = search_files(species, incl='noise', dir='../data/inv/thresh_lp/')[0]
 
     # Load invariance data:
     pure_data, config = load_data(pure_path, **load_kwargs)
     noise_data, _ = load_data(noise_path, **load_kwargs)
     scales = pure_data['scales']
 
-    # Reduce to kernel subset and single threshold:
+    # Reduce to kernel subset and a single threshold:
-    thresh_ind = np.nonzero(pure_data['thresh_perc'] == thresh_percent)[0][0]
+    thresh_ind = np.nonzero(pure_data['thresh_rel'] == thresh_rel)[0][0]
-    kern_inds = find_kern_specs(config['k_specs'], kerns=kernels)
+    kern_inds = find_kern_specs(config['k_specs'], kerns=kern_specs)
     config['k_specs'] = config['k_specs'][kern_inds]
     config['kernels'] = config['kernels'][:, kern_inds]
     pure_measure = pure_data['measure_feat'][:, kern_inds, thresh_ind]
     noise_measure = noise_data['measure_feat'][:, kern_inds, thresh_ind]
+    if exclude_zero:
+        # Reduce to nonzero scales:
+        nonzero_inds = scales > 0
+        scales = scales[nonzero_inds]
+        pure_measure = pure_measure[nonzero_inds, :]
+        noise_measure = noise_measure[nonzero_inds, :]
+    min_feat[i, :] = noise_measure.min(axis=0)
 
-    # Plot invariance curves:
+    # Prepare species-specific colors:
-    pure_ax, noise_ax = spec_axes[:, i]
+    base_color = species_colors[species]
-    pure_ax.plot(scales, pure_measure, c=base_color, lw=lw['spec'])
+    kern_colors = shade_colors(base_color, kern_factors)
-    noise_ax.plot(scales, noise_measure, c=base_color, lw=lw['spec'])
+    scale_factors = np.linspace(1, 0, scales.size)
+    scale_cmap = create_listed_cmap(shade_colors(base_color, scale_factors))
+    scale_cmap_bw = create_listed_cmap(shade_colors((0., 0., 0.), scale_factors))
+
+    # Plot feature invariance curves:
+    pure_ax, noise_ax = feat_axes[:, i]
+    symlog_kwargs['linthresh'] = scales[scales > 0][0]
+    [ax.set_xscale('symlog', **symlog_kwargs) for ax in feat_axes[:, i]]
+    pure_ax.set_xscale('symlog', **symlog_kwargs)
+    noise_ax.set_xscale('symlog', **symlog_kwargs)
+    handles = pure_ax.plot(scales, pure_measure, lw=lw['feat'])
+    [h.set_color(c) for h, c in zip(handles, kern_colors)]
+    handles = noise_ax.plot(scales, noise_measure, lw=lw['feat'])
+    [h.set_color(c) for h, c in zip(handles, kern_colors)]
 
     if i == 0:
         # Indicate kernel waveforms:
         ylims = ylimits(config['kernels'], pad=0.05)
         xlims = (config['k_times'][0], config['k_times'][-1])
-        for j, inset in enumerate(insets):
+        for kern, inset, c in zip(config['kernels'].T, insets, kern_colors_bw):
-            inset.plot(config['k_times'], config['kernels'][:, j],
+            inset.plot(config['k_times'], kern, c=c, lw=lw['kern'])
-            c='k', lw=lw['kern'])
             inset.set_xlim(xlims)
             inset.set_ylim(ylims)
-        time_bar(insets[0], parent=spec_axes[0, 0], **time_bar_kwargs)
+        time_bar(insets[0], parent=feat_axes[0, 0], **kern_bar_kwargs)
 
     # Plot pure feature space:
-    handle = big_axes[0].scatter(pure_measure[:, 0], pure_measure[:, 1], 
+    from matplotlib.colors import LogNorm
-                                 c=scales, cmap=spec_cmaps[i], **space_kwargs)
+    norm = LogNorm(vmin=scales[scales > 0][0], vmax=scales[-1])
+    handle = space_axes[0].scatter(pure_measure[:, 0], pure_measure[:, 1], 
+                                   c=scales, cmap=scale_cmap, norm=norm,
+                                   zorder=zorder[species], **space_kwargs)
 
     # Plot noise feature space:
-    big_axes[1].scatter(noise_measure[:, 0], noise_measure[:, 1],
+    space_axes[1].scatter(noise_measure[:, 0], noise_measure[:, 1],
-                        c=scales, cmap=spec_cmaps[i], **space_kwargs)
+                          c=scales, cmap=scale_cmap, norm=norm,
+                          zorder=zorder[species], **space_kwargs)
 
     # Indicate scale color code:
-    big_subfig.colorbar(handle, cax=bar_axes[i])
+    space_subfig.colorbar(handle, cax=bar_axes[i])
-    bar_axes[i].set_yscale('symlog', linthresh=scales[1], linscale=0.2)
+    bar_axes[i].set_yscale('symlog', **symlog_kwargs)
-    if i < len(pure_paths) - 1:
+    if i < n_species - 1:
         hide_ticks(bar_axes[i], 'right', ticks=False)
     else:
-        ylabel(bar_axes[i], ylabels['bar'], transform=big_subfig.transSubfigure, **ylab_cbar_kwargs)
+        ylabel(bar_axes[i], ylabels['bar'], transform=space_subfig.transSubfigure, **ylab_cbar_kwargs)
+
+if show_noise:
+    # Indicate feature noise floor:
+    min_feat = min_feat.mean(axis=0)
+    space_axes[-1].add_patch(plt.Rectangle((0, 0), min_feat[0], min_feat[1], **noise_kwargs))
 
 if save_path is not None:
     fig.savefig(save_path)

python/save_inv_data_log-hp.py

@@ -1,19 +1,19 @@
-import glob
 import numpy as np
 from thunderhopper.modeltools import load_data, save_data
-from thunderhopper.filetools import crop_paths
+from thunderhopper.filetools import search_files, crop_paths
 from thunderhopper.filters import decibel, sosfilter
 from IPython import embed
 
 # GENERAL SETTINGS:
-target = 'Omocestus_rufipes'
+target = ['Omocestus_rufipes', '*'][0]
-data_paths = glob.glob(f'../data/processed/{target}*.npz')
+data_paths = search_files(target, excl='noise', dir='../data/processed/')
 save_path = '../data/inv/log_hp/'
 
 # ANALYSIS SETTINGS:
 add_noise = False
+save_snippets = target == 'Omocestus_rufipes'
 example_scales = np.array([0.1, 1, 10, 30, 100, 300])
-scales = np.geomspace(0.1, 10000, 1000)
+scales = np.geomspace(0.1, 10000, 500)
 scales = np.unique(np.concatenate((scales, example_scales)))
 
 # EXECUTION:
@@ -60,13 +60,16 @@ for data_path, name in zip(data_paths, crop_paths(data_paths)):
         data = dict(
             scales=scales,
             example_scales=example_scales,
-            snip_env=mix[:, save_inds],
-            snip_log=mix_log[:, save_inds],
-            snip_inv=mix_inv[:, save_inds],
             measure_env=measure_env,
             measure_log=measure_log,
             measure_inv=measure_inv,
             )
+        if save_snippets:
+            data.update(
+                snip_env=mix[:, save_inds],
+                snip_log=mix_log[:, save_inds],
+                snip_inv=mix_inv[:, save_inds],
+            )
         file_name = save_path + name
         if add_noise:
             file_name += '_noise'

python/save_inv_data_thresh-lp.py

@@ -8,14 +8,14 @@ from thunderhopper.model import convolve_kernels
 from IPython import embed
 
 # GENERAL SETTINGS:
-target = ['Omocestus_rufipes', '*'][0]
+target = ['Omocestus_rufipes', '*'][1]
-data_paths = search_files(target, dir='../data/processed/')
+data_paths = search_files(target, excl='noise', dir='../data/processed/')
 noise_path = '../data/processed/white_noise_sd-1.npz'
 save_path = '../data/inv/thresh_lp/'
 
 # ANALYSIS SETTINGS:
 add_noise = True
-save_snippets = add_noise and True
+save_snippets = add_noise and (target == 'Omocestus_rufipes')
 plot_results = False
 example_scales = np.array([0, 1, 10, 30, 100])
 scales = np.geomspace(0.01, 10000, 100)
@@ -50,11 +50,11 @@ for data_path, name in zip(data_paths, crop_paths(data_paths)):
     config['k_specs'] = config['k_specs'][kern_inds, :]
     config['k_props'] = [config['k_props'][i] for i in kern_inds]
 
-    # Normalize song component:
-    song /= song[segment].std()
-
     # Get normalized noise component:
     noise = pure_noise[:song.shape[0]]
+
+    # Normalize both components:
+    song /= song[segment].std()
     noise /= noise[segment].std()
 
     # Define kernel-specific threshold values based on pure-noise response SD:

python/save_noise_data.py

@@ -9,7 +9,7 @@ from IPython import embed
 save_path = '../data/processed/white_noise'
 stages = ['filt', 'env', 'log', 'inv', 'conv', 'bi', 'feat']
 sds = [1]
-dur = 10
+dur = 60
 
 # Interactivity:
 reload_saved = False

python/save_species_colors.py

@@ -0,0 +1,39 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from color_functions import load_colors, sample_cmap, color_selector
+from IPython import embed
+
+
+# Settings:
+species = [
+    'Omocestus_rufipes',
+    'Chorthippus_biguttulus',
+    'Chorthippus_mollis',
+    'Chrysochraon_dispar',
+    'Gomphocerippus_rufus',
+    'Pseudochorthippus_parallelus',
+]
+file_name = '../data/species_colors.npz'
+sample_kwargs = dict(
+    cmap='turbo',
+    n=len(species),
+    low=None,
+    high=None,
+    segments=None,
+)
+select_kwargs = dict(
+    n=len(species),
+    save=file_name,
+    labels=species,
+)
+new_start = True
+
+
+# Execution:
+if new_start:
+    colors = sample_cmap(**sample_kwargs)
+else:
+    colors = load_colors('../data/stage_colors.npz')
+colors = color_selector(colors=colors, **select_kwargs)
+plt.show()
+embed()

python/save_stage_colors.py

@@ -1,5 +1,5 @@
 import matplotlib.pyplot as plt
-from color_functions import load_colors, color_selector, hex_to_rgb
+from color_functions import load_colors, color_selector
 from IPython import embed
 
 # Settings: