tokens-analysis-pipeline

Imported from jamespherman/tokens-analysis-pipeline (AGENTS.md). Install upstream with npx skills add jamespherman/tokens-analysis-pipeline. Copyright stays with the author.

Agent Instructions

Execution Environment

MATLAB code can be executed using the -batch option:

matlab -batch "run('script_name.m')"

Or for single expressions:

matlab -batch "disp(pwd); disp(version)"

Notes:

• -batch runs without the GUI and exits on completion
• Errors will return a non-zero exit code
• Working directory and path setup may need explicit handling

Refactoring and Impact Analysis Protocol

When you are asked to refactor or modify code, your implementation plan must be comprehensive and account for all dependencies across the entire pipeline. Before writing any code, you must consider the impact of the change on all related files. Use the following checklist to ensure your plan is complete:

1. Definition: Which files define the relevant data structures or configurations (e.g., define_task_conditions.m, .md files)?
2. Creation: Which scripts create the data on a per-session basis (e.g., analyze_*.m functions)?
3. Orchestration: Which scripts read the configuration to execute the creation process (e.g., run_4factors_analysis.m)?
4. Aggregation: Which scripts pool the data from multiple sessions (e.g., aggregate_analysis_results.m)?
5. Consumption: Which scripts consume the final aggregated data (e.g., plot_aggregated_*.m functions)?

IMPORTANT: Data Structures

Before writing or modifying any code that interacts with session_data.mat files, you MUST consult the data dictionary.

The canonical definition for the main data structure, session_data, is located in: docs/preprocessing_docs/session_data_dictionary.md

Do not make assumptions about the fields or layout of this structure. The documentation is the single source of truth. For example:

• The number of clusters/neurons should be derived from session_data.spikes.cluster_info
• Trial-related event times are in the session_data.eventTimes struct
• Spike times are stored as a single vector (spikes.times) and mapped to clusters via spikes.clusters

Consulting the documentation first will prevent bugs and ensure your code is compatible with the project's data standards.

Before writing or modifying any code that uses the aggregated_analysis_data.mat file or the condition_defs analysis plan, you MUST consult the data dictionary.

The canonical definitions for the aggregated data structure and the overarching analysis plan are located in:

• docs/aggregated_data_dictionary.md
• docs/analysis_data_dictionary.md

When creating new MATLAB scripts in this directory, please ensure you add the utils directory to the MATLAB path. This is necessary for helper functions to be found.

You can use the following code snippet at the beginning of your script to achieve this:

%% Setup Paths
% Add the 'utils' directory to the path so that helper functions can be
% found.
[script_dir, ~, ~] = fileparts(mfilename('fullpath'));
addpath(fullfile(script_dir, 'utils'));

This will ensure that the path is added correctly, regardless of the current working directory.

Coding Style

Line Length

Each line of code should be limited to a maximum of 75 characters. For longer lines, please wrap them to a new line using ....

Commenting

It is strongly preferred that comments are added on the line above the code to which they refer, rather than to the right of it.

File Structure

All new MATLAB scripts should begin with a standard file header. This header should include a brief description of the script's purpose, the author's name, and the date of creation or last modification.

Example:

%% my_script.m
%
% A brief one-line description of what the script does.
%
% A more detailed description can follow, explaining the inputs, outputs,
% and any key assumptions or dependencies.
%
% Author: Your Name
% Date: YYYY-MM-DD
%

Code should be organized into logical sections using %% headers. For example:

%% Setup
% ... setup code here ...

%% Main Analysis
% ... main analysis code here ...

%% Plotting
% ... plotting code here ...

Conventions

Analysis Pipeline

1. Selecting Trials by Task

It is crucial to filter trials based on the specific task being analyzed. Neural responses are highly task-dependent, and calculations must be performed only on the relevant behavioral data. Use the session_data.trialInfo.taskCode field to select trials corresponding to a given task.

First, load the task code definitions, then create a logical mask to select the desired trial indices.

Example:

% Load the structure containing all task codes
codes = initCodes();

% Find the indices of trials that belong to the 'tokens' task
tokens_trial_indices = session_data.trialInfo.taskCode == codes.uniqueTaskCode_tokens;

% Subsequent analyses should use these indices to filter the data
valid_cue_on_times = session_data.eventTimes.CUE_ON(tokens_trial_indices);

2. Condition Mask Compatibility

When creating logical masks for different experimental conditions, ensure they are compatible with the data they are intended to select. This means the masks should have the same number of elements as the trials for the specific task being analyzed.

For example, the define_task_conditions.m function is designed to generate condition masks specifically for 'tokens' trials. It first filters the session data to include only rewarded tokens trials and then creates masks that are the same length as this filtered set of trials. This ensures that the masks can be directly applied to core_data arrays, which are also filtered for tokens trials.

Interpreting Specific Task Data

The gSac_jph task has special properties that can be leveraged during analysis. Memory-guided saccade trials within this task are intentionally placed at the neuron's estimated receptive/movement field center. This experimental design allows for two major simplifications:

1. The recorded side of the SC (scSide) can be inferred directly from the target location (e.g., targets in the left visual field imply a right SC recording).
2. All of these trials can be considered "contralateral" for the purposes of analysis, as they are designed to elicit a strong response from the recorded neurons.

Project Root Directory

When writing code that needs to access files or directories relative to the project root, it is important to construct paths correctly. The project root for tokens-analysis-pipeline can be found by using the following MATLAB code:

project_root = fullfile(findOneDrive, 'Code', 'tokens-analysis-pipeline');

This will ensure that paths are resolved correctly, regardless of the directory from which the code is executed.

For example, to access the figures directory from within a script, you should use:

figures_dir = fullfile(project_root, 'figures');

This is especially important for scripts located in the code/ directory, which are executed from /code in the runtime environment.

Screening SC Neurons

When screening neurons from Superior Colliculus (SC) recordings, use the screen_sc_neurons.m function. This function has been designed to automate several key steps of the analysis for 'SC type' sessions.

Function Usage

The function is called as follows:

[selected_neurons, sig_epoch_comp, scSide] = screen_sc_neurons(session_data);

Automatic scSide Determination

The function automatically determines the recorded side of the SC ('left' or 'right'). It does this by comparing the average visual-evoked firing rates for targets presented in the left versus the right visual field. The side with the stronger contralateral response is identified as the recorded hemisphere. The determined side is returned as the scSide output variable.

Trial Selection

The function specifically uses memory-guided saccade trials from the gSac_jph task to perform neuron selection. It identifies these trials using the taskCode and the isVisSac field from trialInfo. If no such trials are found, it has a fallback mechanism to use any rewarded memory-guided saccade trials present in the session.

sig_epoch_comparison Output

A key output of this function is sig_epoch_comp. This is a boolean matrix of size [nClusters x 3] that indicates whether a neuron showed a statistically significant change in firing rate for three critical comparisons:

1. Visual epoch vs. Baseline
2. Delay epoch vs. Baseline
3. Saccade epoch vs. Baseline

This variable is crucial for the functional classification of SC neurons and should be stored for later analysis. The recommended practice is to store it in the session_data structure, for example:

session_data.metadata.scSide = scSide;
session_data.metadata.sig_epoch_comparison = sig_epoch_comp;

This ensures that the results of the screening are saved with the session's data.

Pupil Analysis Conventions

• Pupil preprocessing functions should output to session_data.analysis.pupil
• Follow the same event alignment pattern as neuronal data (use alignAndBinSpikes pattern)
• New analysis plans should be added to define_task_conditions.m

Plotting Conventions

When creating complex, multi-panel figures for data verification or analysis, please adhere to the following conventions to ensure clarity, consistency, and scientific accuracy.

1. Use barStairsFill for PSTHs

For plotting Peristimulus Time Histograms (PSTHs), use the barStairsFill.m utility function instead of a standard plot. This function creates a filled, staircase-style plot that accurately represents the binned nature of the data without implying a false sense of continuity between bins.

Example:

% Good: Using barStairsFill for a PSTH WITH ONE PLOT ON A SET OF AXES:
mean_psth = mean(rates, 1);
hBS = barStairsFill(time_vector, zeros(size(mean_psth)), mean_psth);
delete(hBS(2))
set(hBS(1), 'FaceColor', 'k');
set(hBS(3), 'Color', 'k')

% Good: Using barStairsFill for a PSTH WITH MORE THAN ONE PLOT ON A SET OF AXES:
mean_psth1 = mean(rates1, 1);
mean_psth2 = mean(rates2, 1);
hold on;
hBS1 = barStairsFill(time_vector, zeros(size(mean_psth1)), mean_psth1);
delete(hBS1(1:2))
set(hBS1(3), 'Color', 'b');
hBS2 = barStairsFill(time_vector, zeros(size(mean_psth2)), mean_psth2);
delete(hBS2(1:2))
set(hBS2(3), 'Color', 'r');

2. De-clutter Axes in Multi-Panel Figures

To improve readability and reduce visual clutter in figures with multiple subplots arranged in a grid, remove redundant axis labels.

• For any given column of plots, only the bottom-most plot should have X-axis tick labels.
• For any given row of plots, only the left-most plot should have Y-axis tick labels.

Example:

% h is an array of handles to the axes in a 2x3 grid
% Remove x-labels from the top row
set(h(1:3), 'XTickLabel', []);
% Remove y-labels from the middle and right columns
set(h([2,3,5,6]), 'YTickLabel', []);

3. Proportional Subplot Widths for Time-Series Data

When comparing different time epochs in adjacent subplots, the width of each subplot column should be proportional to the duration of the time window it represents. This provides an intuitive visual comparison of the temporal dynamics.

Do not use automated subplot tools like subplot. Instead, ONLY use mySubPlot. Unlike a call to 'subplot' (which shouldn't be used) which is formatted as 'subplot(a,b,c), a call to 'mySubPlot' shoud be formatted like: 'mySubPlot([a,b,c]).

4. Use a Single, Informative sgtitle

Use a single, comprehensive main title for the entire figure using sgtitle(). This title should provide the key takeaway or context for the figure as a whole.

While sgtitle is primary, individual title() calls on subplots are permitted when they are necessary to provide specific context for each panel. For example, in a figure comparing multiple different analyses, each subplot's title can clarify the specific comparison it displays. The plot_aggregated_roc_comparison.m script serves as a good example of this convention. Column and row labels should still be used for axis identification (i.e., ylabel on the leftmost plots and xlabel on the bottom plots).

5. Set Interpreter to 'none' for Titles with File Paths

When including filenames, unique IDs, or other text with special characters (like underscores _) in a title, always set the 'Interpreter' property to 'none'. This prevents MATLAB from treating the underscores as subscript indicators and ensures the text is displayed literally.

Example:

unique_id = 'My_Session_2025_09_09';
sgtitle(sprintf('Data Verification for: %s', unique_id), 'Interpreter', 'none');

6. Save figures using 'pdfSave.m'

Example: With a figure handle ('fig_handle') and a filename defined for saving the pdf (fig_filename):

pdfSave(fig_filename, fig_handle.Position(3:4)/72, fig_handle);

Coordinate System and Angle Conventions

All angular measurements adhere to the standard mathematical convention for polar coordinates.

• trialInfo.targetTheta: This specific variable stores the target angle in degrees, multiplied by 10 for strobing purposes. It must be divided by 10 before use in any geometric calculations.
• Angle Reference: 0° corresponds to the positive x-axis (the 3 o'clock position).
• Angle Direction: Angles increase in a counter-clockwise direction.
• Angle Range: After division by 10, angles are expressed in degrees, typically in the range [0, 360).
• Visual Field Definitions:
  • Right Visual Field: (theta >= 0 AND theta < 90) OR (theta > 270 AND theta < 360)
  • Left Visual Field: (theta > 90 AND theta < 270)