diff --git a/CHANGELOG.md b/CHANGELOG.md index 0d73798..c5bf583 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,31 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] +### Added + +- **Lever push durations.** `visiomode-analysis session` picks up a sibling + `sub-_exp-_ses-_lever-durations.csv` (or one passed with + `--lever-durations`) and adds a `lever_duration` column (ms) to the trials CSV, + matching pushes to precued, hit and false alarm trials in order. If the file can't + be matched (empty or unparseable, missing `push_id`/`duration` columns, non-numeric + durations, `push_id` not running 0..n-1, or a push count that differs from the + number of lever push trials) a warning is printed and the column is left empty; no + durations file means no column. Session and + subject summaries report the mean, median, SD and IQR of lever durations for all + pushes, cued pushes, hits, false alarms and precued trials, each with a `_wc` + variant, and the session report gains a median lever push duration plot. +- **Summaries and reports from an existing trials dataframe.** `session.summary()` now + also accepts a trials dataframe from `session.get_trials()`, and + `session.generate_report()` takes an optional `trials` dataframe, so a session is + parsed once rather than re-read for every metric and plot. + +### Fixed + +- **Reaction time IQR error bars.** The median reaction time plots in the session + report passed the 25th/75th percentiles to Plotly as error bar lengths, which + Plotly treats as offsets from the median, so bars spanned `median - q25` to + `median + q75`. They now span the IQR (`q25` to `q75`). + ## [0.3.0] - 2026-09-23 ### Added diff --git a/README.md b/README.md index cf38e7d..cb56b15 100644 --- a/README.md +++ b/README.md @@ -52,6 +52,8 @@ visiomode-analysis session path/to/session.json -o output/ --no-report visiomode-analysis session path/to/session.json -o output/ --with-regressors --regressor-timestamps frame_times.csv ``` +**Lever push durations.** If a `sub-_exp-_ses-_lever-durations.csv` file (columns `push_id`, `duration` in ms) sits next to the session JSON, `session` matches each push to its lever push trial (precued trials, hits and false alarms, in trial order) and adds a `lever_duration` column to the trials CSV. Point it at a file elsewhere with `--lever-durations path/to/durations.csv`. Trials without a lever push have an empty `lever_duration`, and sessions without a durations file have no `lever_duration` column at all. If the durations can't be matched (the file is empty or malformed, `push_id` doesn't run 0..n-1, or the number of pushes doesn't match the number of lever push trials), a warning is printed and the column is left empty for that session. + **Generate regressors** for an already-processed session, aligned to an external timestamp series: ```bash @@ -72,6 +74,8 @@ visiomode-analysis session-start-time path/to/session.json visiomode-analysis subject path/to/subject_dir/ -o output/ ``` +When the trials CSVs have a `lever_duration` column, the summary also reports `lever_duration_[_][_wc]` (in ms) for all lever pushes and for the `cued`, `hits`, `false_alarms` and `precued` subsets, with `_wc` variants that include correction trials. These columns are omitted if no session has lever durations. + Run `visiomode-analysis --help` or `visiomode-analysis --help` for full option details. ### Python API @@ -85,7 +89,7 @@ trials = session.get_trials("path/to/session.json") metadata = session.get_metadata("path/to/session.json") summary = session.summary(trials) -session.generate_report(trials, metadata, output_dir="output/") +session.generate_report("path/to/session.json", output_dir="output/", trials=trials) ``` ### Input files and naming convention diff --git a/src/visiomode_analysis/reports/templates/session.html b/src/visiomode_analysis/reports/templates/session.html index 439a4b6..53203f4 100644 --- a/src/visiomode_analysis/reports/templates/session.html +++ b/src/visiomode_analysis/reports/templates/session.html @@ -121,6 +121,13 @@
Cued:uncued response breakdown
Median reaction time
{{ fig_rt_median | safe }} + + {% if fig_lever_duration_median %} +
+
Median lever push duration
+ {{ fig_lever_duration_median | safe }} +
+ {% endif %} {% if corrections_enabled %} @@ -143,6 +150,13 @@
Cued:uncued response breakdown (all trials)
Median reaction time (all trials)
{{ fig_rt_median_wc | safe }} + + {% if fig_lever_duration_median_wc %} +
+
Median lever push duration (all trials)
+ {{ fig_lever_duration_median_wc | safe }} +
+ {% endif %}
diff --git a/src/visiomode_analysis/session/__init__.py b/src/visiomode_analysis/session/__init__.py index b2854d3..4d56fad 100644 --- a/src/visiomode_analysis/session/__init__.py +++ b/src/visiomode_analysis/session/__init__.py @@ -22,6 +22,7 @@ import os import json +import warnings import click import datetime import h5py @@ -56,6 +57,11 @@ # outcome vocabulary so downstream logic only has to deal with one set of labels. LEGACY_OUTCOME_LABELS = {"hit": "correct", "false_alarm": "incorrect", "miss": "no_response"} +# Companion file holding the duration (ms) of every lever push in a session, one row per push in trial order. +LEVER_DURATIONS_SUFFIX = "_lever-durations.csv" +LEVER_DURATION = "lever_duration" +LEVER_DURATION_STATS = ("mean", "median", "sd", "iqr") + env = Environment(loader=PackageLoader("visiomode_analysis.reports", "templates"), autoescape=select_autoescape()) @@ -89,6 +95,12 @@ type=click.Path(exists=True, dir_okay=False), help="Path to a CSV/TXT file of ISO timestamps, or a mesoscopy H5 with /timestamps_aligned, for regressor generation.", ) +@click.option( + "--lever-durations", + type=click.Path(exists=True, dir_okay=False), + help=f"Path to a lever push durations CSV (push_id, duration in ms). Defaults to the `*{LEVER_DURATIONS_SUFFIX}` " + "file next to the session JSON, if there is one.", +) def session_cmd(**kwargs): """Generate a session report and extract trials from a Visiomode JSON file.""" out_dir = preprocess_session(**kwargs) @@ -137,6 +149,7 @@ def preprocess_session( no_report: bool = False, with_regressors: bool = False, regressor_timestamps: str | None = None, + lever_durations: str | None = None, ) -> str: """Generate a session summary report and trials file from a raw Visiomode JSON. @@ -146,14 +159,15 @@ def preprocess_session( no_report (bool, optional): Whether to generate a session report. Defaults to False. with_regressors (bool, optional): Whether to generate regressors for the session. Defaults to False. regressor_timestamps (str, optional): Path to a CSV/TXT file of ISO timestamps, or a mesoscopy H5, for regressor generation. Defaults to None. + lever_durations (str, optional): Path to a lever push durations CSV. Defaults to the sibling file found by `find_lever_durations`. Returns: str: Returns directory under which files were saved """ - trials_df = get_trials(path, to_csv=True, output_dir=output_dir) + trials_df = get_trials(path, to_csv=True, output_dir=output_dir, lever_durations=lever_durations) meta = get_metadata(path) if not no_report: - generate_report(path, output_dir=output_dir) + generate_report(path, output_dir=output_dir, trials=trials_df) if with_regressors: if not regressor_timestamps: @@ -244,13 +258,17 @@ def get_metadata(path: str) -> dict: } -def get_trials(path: str, to_csv: bool = False, output_dir: str = ".") -> pd.DataFrame: +def get_trials( + path: str, to_csv: bool = False, output_dir: str = ".", lever_durations: str | None = None +) -> pd.DataFrame: """Parse a Visiomode JSON file and return a trials dataframe. Optionally save to CSV. Args: path (str): Path to Visiomode JSON file. to_csv (bool, optional): Save trials dataframe as a CSV. If True, indicate which directory to save in via `output_dir`. Defaults to False. output_dir (str, optional): Output directory for trials CSV file. Only used if `to_csv` is True. Defaults to current directory. + lever_durations (str, optional): Path to a lever push durations CSV. Defaults to the sibling file found by + `find_lever_durations`; without one, the dataframe has no `lever_duration` column. Returns: pd.DataFrame: Session trials dataframe, where each row is a trial. @@ -280,6 +298,10 @@ def get_trials(path: str, to_csv: bool = False, output_dir: str = ".") -> pd.Dat # built any other way is treated the same. df["outcome"] = df["outcome"].replace(LEGACY_OUTCOME_LABELS) + lever_durations = lever_durations or find_lever_durations(path) + if lever_durations: + df = _add_lever_durations(df, lever_durations) + if to_csv: out_path = f"{output_dir}{os.sep}sub-{metadata.get('animal_id')}_exp-{metadata.get('experiment')}_ses-{str(metadata.get('session_date')).replace('-', '')}_behaviour-{metadata.get('protocol')}_trials.csv" df.to_csv(out_path) @@ -287,9 +309,100 @@ def get_trials(path: str, to_csv: bool = False, output_dir: str = ".") -> pd.Dat return df +def find_lever_durations(path: str) -> str | None: + """Return the lever push durations CSV next to a session JSON, if there is one. + + `sub-X_exp-Y_ses-Z_behaviour-P.json` pairs with `sub-X_exp-Y_ses-Z_lever-durations.csv`; a JSON without a + `_behaviour-` token pairs with `_lever-durations.csv`. + """ + json_path = Path(path) + prefix = json_path.stem.split("_behaviour-")[0] + candidate = json_path.with_name(f"{prefix}{LEVER_DURATIONS_SUFFIX}") + return str(candidate) if candidate.is_file() else None + + +def _is_lever_push(df: pd.DataFrame) -> pd.Series: + """Trials in which the lever was pushed: precued trials, hits and false alarms.""" + return (df["outcome"] == "precued") | df["sdt_type"].isin([HIT, FALSE_ALARM]) + + +def _add_lever_durations(df: pd.DataFrame, lever_durations: str) -> pd.DataFrame: + """Add a `lever_duration` (ms) column, matching the Nth push in `lever_durations` to the Nth lever push trial. + + Non-push trials are NaN. If the durations can't be matched to the trials (see `_read_lever_durations`), a warning + is raised and the whole column is left NaN. Warning rather than raising keeps a bad durations file from breaking + callers that don't need durations (e.g. `regressors`). + """ + is_push = _is_lever_push(df) + + df = df.copy() + df[LEVER_DURATION] = np.nan + try: + df.loc[is_push, LEVER_DURATION] = _read_lever_durations(lever_durations, num_push_trials=int(is_push.sum())) + except ValueError as e: + warnings.warn(f"Can't match {lever_durations} to trials: {e}; leaving {LEVER_DURATION} empty.", stacklevel=2) + return df + + +def _read_lever_durations(path: str, num_push_trials: int) -> npt.NDArray: + """Read lever push durations (ms) in push order, checking they can be matched to the session's lever push trials. + + Raises: + ValueError: If the file is empty or malformed, lacks a `push_id` or `duration` column, has non-numeric + durations, `push_id` doesn't run 0..n-1 (gaps or duplicates), or the number of pushes doesn't match + `num_push_trials`. + """ + pushes = pd.read_csv(path) # pandas' EmptyDataError and ParserError are ValueErrors + missing_columns = {"push_id", "duration"} - set(pushes.columns) + if missing_columns: + raise ValueError(f"missing column(s) {', '.join(sorted(missing_columns))}") + if sorted(pushes["push_id"].tolist()) != list(range(len(pushes))): + raise ValueError("push_id does not run from 0 without gaps or duplicates") + if len(pushes) != num_push_trials: + raise ValueError( + f"{len(pushes)} lever pushes but {num_push_trials} lever push trials (precued, hits and false alarms)" + ) + return pd.to_numeric(pushes.sort_values("push_id")["duration"]).to_numpy(dtype=float) # ValueError if non-numeric + + +def _lever_duration_subsets(df: pd.DataFrame, include_corrections: bool = True) -> dict[str, pd.Series]: + """Lever durations per trial subset, keyed by the suffix used in summary field names ("" for all pushes).""" + if not include_corrections: + df = df[df.correction == False] # noqa: E712 + durations = df[LEVER_DURATION] + return { + "": durations[_is_lever_push(df)].dropna(), + "cued": durations[df.sdt_type.isin([HIT, FALSE_ALARM])].dropna(), + "hits": durations[df.sdt_type == HIT].dropna(), + "false_alarms": durations[df.sdt_type == FALSE_ALARM].dropna(), + "precued": durations[df.outcome == "precued"].dropna(), + } + + +def _lever_duration_summary(df: pd.DataFrame) -> dict[str, float]: + """Mean, median, SD and IQR of lever durations for every subset in `_lever_duration_subsets`, with and without + correction trials, as `lever_duration_[_][_wc]` fields. Empty subsets yield NaN.""" + fields = {} + for wc_suffix, include_corrections in (("", False), ("_wc", True)): + for subset, durations in _lever_duration_subsets(df, include_corrections=include_corrections).items(): + subset_suffix = f"_{subset}" if subset else "" + stats = { + "mean": durations.mean(), + "median": durations.median(), + "sd": durations.std(ddof=1), + "iqr": durations.quantile(0.75) - durations.quantile(0.25), + } + for stat in LEVER_DURATION_STATS: + fields[f"{LEVER_DURATION}_{stat}{subset_suffix}{wc_suffix}"] = float(stats[stat]) + return fields + + def get_rts(path: str, sdt_type=None, include_corrections=True) -> npt.NDArray: - trials = get_trials(path=path) + return _select_rts(get_trials(path=path), sdt_type=sdt_type, include_corrections=include_corrections) + +def _select_rts(trials: pd.DataFrame, sdt_type=None, include_corrections=True) -> npt.NDArray: + """Reaction times of the cued trials with a response, optionally of one SDT type only.""" if not include_corrections: trials = trials[trials.correction == False] # noqa: E712 @@ -299,26 +412,28 @@ def get_rts(path: str, sdt_type=None, include_corrections=True) -> npt.NDArray: return np.array(trials[(trials.response.notnull()) & (trials.cue_onset.notnull())].response_time.values) -def summary(path: str) -> dict: - """Summarise session from a JSON or trials.csv file +def summary(path: str | pd.DataFrame, lever_durations: str | None = None) -> dict: + """Summarise session from a JSON or trials.csv file, or from a trials dataframe built by `get_trials` Args: - path (str): Path to (preprocessed) trials.csv or raw JSON + path (str | pd.DataFrame): Path to (preprocessed) trials.csv or raw JSON, or a trials dataframe + lever_durations (str, optional): Path to a lever push durations CSV, for a raw JSON only. Defaults to the + sibling file found by `find_lever_durations`. Returns: - dict: Summary dictionary + dict: Summary dictionary. Lever duration fields are only included if the trials have a `lever_duration` column. """ - if path.endswith(".json"): + if isinstance(path, str) and path.endswith(".json"): metadata = get_metadata(path) - df = get_trials(path) + df = get_trials(path, lever_durations=lever_durations) else: - df = pd.read_csv(path) + df = path if isinstance(path, pd.DataFrame) else pd.read_csv(path) metadata = { - "animal_id": df["animal_id"][0], - "session_date": df["session_date"][0], - "protocol": df["protocol"][0], - "environment": df["environment"][0], - "experiment": df["experiment"][0], + "animal_id": df["animal_id"].iloc[0], + "session_date": df["session_date"].iloc[0], + "protocol": df["protocol"].iloc[0], + "environment": df["environment"].iloc[0], + "experiment": df["experiment"].iloc[0], } # Trial counts @@ -471,6 +586,7 @@ def summary(path: str) -> dict: "rt_iqr_wc": rt_iqr_wc, "rt_idr": rt_idr, "rt_idr_wc": rt_idr_wc, + **(_lever_duration_summary(df) if LEVER_DURATION in df.columns else {}), } @@ -565,12 +681,27 @@ def _read_aligned_timestamps(path: str, session_start_time: datetime.datetime) - return np.asarray(dataset[()], dtype=np.float64) -def generate_report(path: str, output_dir: str = ".") -> str: +def generate_report( + path: str, output_dir: str = ".", lever_durations: str | None = None, trials: pd.DataFrame | None = None +) -> str: + """Render the session HTML report. + + Args: + path (str): Path to Visiomode JSON file. + output_dir (str, optional): Output directory for the report. Defaults to ".". + lever_durations (str, optional): Path to a lever push durations CSV. Ignored if `trials` is given. Defaults to + the sibling file found by `find_lever_durations`. + trials (pd.DataFrame, optional): Trials dataframe from `get_trials`, to avoid parsing the session again. + + Returns: + str: Path to the generated report. + """ template = env.get_template(SESSION_REPORT_TEMPLATE) metadata = get_metadata(path) - trials = get_trials(path) - session_summary = summary(path=path) + if trials is None: + trials = get_trials(path, lever_durations=lever_durations) + session_summary = summary(trials) template_identifiers = { "subject_id": metadata.get("animal_id"), @@ -605,31 +736,31 @@ def generate_report(path: str, output_dir: str = ".") -> str: session_summary.get("cued_wc"), session_summary.get("precued"), as_html=True ), "fig_rt_median": plots.plot_rt_median( - get_rts(path=path, include_corrections=False), + _select_rts(trials, include_corrections=False), stimulus_duration=metadata.get("stimulus_duration", 4000) / 1000, as_html=True, ) if is_targetonly(metadata.get("protocol")) else plots.plot_rt_medians_from_dict( { - "all": get_rts(path=path, include_corrections=False), - "hits": get_rts(path=path, sdt_type="hit", include_corrections=False), - "false_alarms": get_rts(path=path, sdt_type="false_alarm", include_corrections=False), + "all": _select_rts(trials, include_corrections=False), + "hits": _select_rts(trials, sdt_type="hit", include_corrections=False), + "false_alarms": _select_rts(trials, sdt_type="false_alarm", include_corrections=False), }, stimulus_duration=metadata.get("stimulus_duration", 4000) / 1000, as_html=True, ), "fig_rt_median_wc": plots.plot_rt_median( - get_rts(path=path, include_corrections=True), + _select_rts(trials, include_corrections=True), stimulus_duration=metadata.get("stimulus_duration", 4000) / 1000, as_html=True, ) if is_targetonly(metadata.get("protocol")) else plots.plot_rt_medians_from_dict( { - "all": get_rts(path=path, include_corrections=True), - "hits": get_rts(path=path, sdt_type="hit", include_corrections=True), - "false_alarms": get_rts(path=path, sdt_type="false_alarm", include_corrections=True), + "all": _select_rts(trials, include_corrections=True), + "hits": _select_rts(trials, sdt_type="hit", include_corrections=True), + "false_alarms": _select_rts(trials, sdt_type="false_alarm", include_corrections=True), }, stimulus_duration=metadata.get("stimulus_duration", 4000) / 1000, as_html=True, @@ -676,6 +807,22 @@ def generate_report(path: str, output_dir: str = ".") -> str: "fig_sdt_timeseries": plots.plot_trial_timeseries(trials=trials, use_sdt=True, as_html=True), } + if LEVER_DURATION in trials.columns and trials[LEVER_DURATION].notna().any(): + for fig_key, include_corrections in ( + ("fig_lever_duration_median", False), + ("fig_lever_duration_median_wc", True), + ): + subsets = _lever_duration_subsets(trials, include_corrections=include_corrections) + template_identifiers[fig_key] = plots.plot_lever_duration_medians_from_dict( + { + "all": subsets[""], + "hits": subsets["hits"], + **({} if is_targetonly(metadata.get("protocol")) else {"false_alarms": subsets["false_alarms"]}), + "precued": subsets["precued"], + }, + as_html=True, + ) + out_path = Path( f"{output_dir}{os.sep}sub-{metadata.get('animal_id')}_exp-{metadata.get('experiment')}_ses-{str(metadata.get('session_date')).replace('-', '')}_behaviour-{metadata.get('protocol')}_report-session.html" ) diff --git a/src/visiomode_analysis/session/plots.py b/src/visiomode_analysis/session/plots.py index 131c435..9421208 100644 --- a/src/visiomode_analysis/session/plots.py +++ b/src/visiomode_analysis/session/plots.py @@ -102,12 +102,12 @@ def plot_correction_pie(num_random, num_correction, as_html=False) -> str | go.F return fig -def _rt_stats(rts) -> tuple[float, float, float]: - """Median and 25th/75th percentiles of a reaction-time array, NaN for an empty one.""" - rts = np.asarray(rts, dtype=float) - if rts.size == 0: +def _median_iqr(values) -> tuple[float, float, float]: + """Median and 25th/75th percentiles of an array, NaN for an empty one.""" + values = np.asarray(values, dtype=float) + if values.size == 0: return np.nan, np.nan, np.nan - return float(np.median(rts)), float(np.percentile(rts, 25)), float(np.percentile(rts, 75)) + return float(np.median(values)), float(np.percentile(values, 25)), float(np.percentile(values, 75)) def _rt_yaxis_range(stimulus_duration) -> list[float] | None: @@ -118,15 +118,15 @@ def _rt_yaxis_range(stimulus_duration) -> list[float] | None: def plot_rt_median(rts, stimulus_duration=4, as_html=False) -> str | go.Figure: - median, q25, q75 = _rt_stats(rts) + median, q25, q75 = _median_iqr(rts) fig = go.Figure( go.Scatter( y=[median], error_y=dict( type="data", symmetric=False, - array=[q75], - arrayminus=[q25], + array=[q75 - median], + arrayminus=[median - q25], ), ), layout=go.Layout( @@ -142,23 +142,31 @@ def plot_rt_median(rts, stimulus_duration=4, as_html=False) -> str | go.Figure: def plot_rt_medians_from_dict(rt_dict: dict, stimulus_duration: int = 4, as_html=True) -> str | go.Figure: - stats = [_rt_stats(rt) for rt in rt_dict.values()] + return _plot_medians_from_dict(rt_dict, yaxis_range=_rt_yaxis_range(stimulus_duration), as_html=as_html) + + +def _plot_medians_from_dict( + data: dict, yaxis_range: list[float] | None = None, yaxis_title: str | None = None, as_html=True +) -> str | go.Figure: + """One marker per key at the median of its values, with error bars spanning the IQR.""" + stats = [_median_iqr(values) for values in data.values()] fig = go.Figure( go.Scatter( y=[median for median, _, _ in stats], - x=[key for key in rt_dict.keys()], + x=[key for key in data.keys()], error_y=dict( type="data", symmetric=False, - array=[q75 for _, _, q75 in stats], - arrayminus=[q25 for _, q25, _ in stats], + array=[q75 - median for median, _, q75 in stats], + arrayminus=[median - q25 for median, q25, _ in stats], ), mode="markers", ), layout=go.Layout( margin={"l": 20, "r": 20, "t": 20, "b": 20}, + yaxis_title=yaxis_title, ), - layout_yaxis_range=_rt_yaxis_range(stimulus_duration), + layout_yaxis_range=yaxis_range, ) if as_html: @@ -169,6 +177,11 @@ def plot_rt_medians_from_dict(rt_dict: dict, stimulus_duration: int = 4, as_html def plot_rt_distribution(rts, as_html=False) -> str | go.Figure: ... +def plot_lever_duration_medians_from_dict(duration_dict: dict, as_html=True) -> str | go.Figure: + """Median lever push duration (ms) per trial subset, with error bars spanning the IQR.""" + return _plot_medians_from_dict(duration_dict, yaxis_title="Duration (ms)", as_html=as_html) + + def plot_single_yvalue(value, ymin=0.0, ymax=1.0, as_html=False): fig = go.Figure( go.Scatter( diff --git a/tests/conftest.py b/tests/conftest.py index a42967b..31d593b 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -78,13 +78,42 @@ def _write(path, timestamps=None, session_start_time=None, behaviour_session=Non return _write +@pytest.fixture +def gonogo_session_with_lever_durations(tmp_path, gonogo_session_json_path): + """A BIDS-named copy of `gonogo_session_json_path` with a sibling `_lever-durations.csv` holding one + duration per lever push trial (precued, hits and false alarms), numbered 100, 101, ... in trial order. + Returns the paths of the JSON and of the durations file.""" + json_path = tmp_path / "sub-MM229_exp-expX_ses-20220101_behaviour-hf.json" + json_path.write_text(pathlib.Path(gonogo_session_json_path).read_text()) + + trials = session.get_trials(gonogo_session_json_path) + num_pushes = int(((trials["outcome"] == "precued") | trials["sdt_type"].isin(["hit", "false_alarm"])).sum()) + durations_path = tmp_path / "sub-MM229_exp-expX_ses-20220101_lever-durations.csv" + pd.DataFrame({"push_id": range(num_pushes), "duration": 100.0 + np.arange(num_pushes)}).to_csv( + durations_path, index=False + ) + return str(json_path), str(durations_path) + + @pytest.fixture def write_trials_csv(): """Factory fixture that writes a minimal preprocessed trials.csv with just enough columns - for `session.summary` (and therefore `subject.collate_sessions`) to run on.""" - - def _write(directory, filename, animal_id, session_date, protocol, experiment, environment="unknown"): - rows = [ + for `session.summary` (and therefore `subject.collate_sessions`) to run on. `lever_durations=True` + adds a `lever_duration` column to the default trials; pass `rows` to replace them (include a + `lever_duration` key in each row for durations).""" + + def _write( + directory, + filename, + animal_id, + session_date, + protocol, + experiment, + environment="unknown", + lever_durations=False, + rows=None, + ): + default_rows = [ dict(outcome="correct", correction=False, sdt_type="hit", response_time=0.5, response="leverpush"), dict( outcome="incorrect", correction=False, sdt_type="false_alarm", response_time=0.3, response="leverpush" @@ -95,12 +124,17 @@ def _write(directory, filename, animal_id, session_date, protocol, experiment, e dict(outcome="incorrect", correction=True, sdt_type=None, response_time=np.nan, response=None), dict(outcome="precued", correction=False, sdt_type=None, response_time=np.nan, response=None), ] - df = pd.DataFrame(rows) + df = pd.DataFrame(default_rows if rows is None else rows) df["animal_id"] = animal_id df["session_date"] = session_date df["protocol"] = protocol df["experiment"] = experiment df["environment"] = environment + if lever_durations: + if rows is not None: + raise ValueError("lever_durations=True only applies to the default rows") + # Lever pushes are the hit, false alarm and precued trials. + df["lever_duration"] = [80.0, 120.0, np.nan, np.nan, 60.0] path = os.path.join(directory, filename) df.to_csv(path, index=False) diff --git a/tests/test_cli.py b/tests/test_cli.py index 64402de..4f6d2df 100644 --- a/tests/test_cli.py +++ b/tests/test_cli.py @@ -38,6 +38,26 @@ def test_session_cmd_handles_an_all_miss_legacy_singletarget_session(runner, wri assert any(name.endswith("_behaviour-singletarget_report-session.html") for name in written) +def test_session_cmd_writes_lever_durations_from_explicit_path( + runner, gonogo_session_with_lever_durations, tmp_path +): + json_path, durations_path = gonogo_session_with_lever_durations + out_dir = tmp_path / "out" + out_dir.mkdir() + moved_durations_path = tmp_path / "elsewhere.csv" + os.rename(durations_path, moved_durations_path) + + result = runner.invoke( + session.session_cmd, + [json_path, "-o", str(out_dir), "--no-report", "--lever-durations", str(moved_durations_path)], + ) + + assert result.exit_code == 0, result.output + (trials_csv,) = [name for name in os.listdir(out_dir) if name.endswith("_trials.csv")] + trials = pd.read_csv(out_dir / trials_csv) + assert trials["lever_duration"].notna().any() + + def test_session_cmd_no_report_skips_report_generation(runner, gonogo_session_json_path, tmp_path): result = runner.invoke(session.session_cmd, [gonogo_session_json_path, "-o", str(tmp_path), "--no-report"]) @@ -103,6 +123,23 @@ def test_regressors_cmd_writes_regressors_file( assert any(name.endswith("_regressors.npz") for name in os.listdir(tmp_path)) +def test_regressors_cmd_is_not_broken_by_a_malformed_lever_durations_file( + runner, gonogo_session_with_lever_durations, gonogo_regressor_timestamps_csv_path, tmp_path +): + json_path, durations_path = gonogo_session_with_lever_durations + pd.read_csv(durations_path).rename(columns={"push_id": "id"}).to_csv(durations_path, index=False) + out_dir = tmp_path / "out" + out_dir.mkdir() + + result = runner.invoke( + session.regressors_cmd, + [json_path, "-o", str(out_dir), "--regressor-timestamps", gonogo_regressor_timestamps_csv_path], + ) + + assert result.exit_code == 0, result.output + assert any(name.endswith("_regressors.npz") for name in os.listdir(out_dir)) + + def test_regressors_cmd_accepts_h5_timestamps(runner, gonogo_session_json_path, write_aligned_h5, tmp_path): h5_path = write_aligned_h5(tmp_path / "aligned.h5") diff --git a/tests/test_plots.py b/tests/test_plots.py index 793b867..dc04789 100644 --- a/tests/test_plots.py +++ b/tests/test_plots.py @@ -56,9 +56,11 @@ def test_plot_rt_median_reports_median_and_iqr_of_rts(): fig = plots.plot_rt_median(rts, stimulus_duration=5) scatter = fig.data[0] - assert scatter.y[0] == pytest.approx(np.median(rts)) - assert scatter.error_y.array[0] == pytest.approx(np.percentile(rts, 75)) - assert scatter.error_y.arrayminus[0] == pytest.approx(np.percentile(rts, 25)) + median = np.median(rts) + assert scatter.y[0] == pytest.approx(median) + # Error bars are offsets from the median, so together they span the IQR. + assert scatter.error_y.array[0] == pytest.approx(np.percentile(rts, 75) - median) + assert scatter.error_y.arrayminus[0] == pytest.approx(median - np.percentile(rts, 25)) assert fig.layout.yaxis.range == pytest.approx((0, 5)) @@ -71,6 +73,10 @@ def test_plot_rt_medians_from_dict_preserves_key_order_and_per_key_stats(): assert list(scatter.x) == ["all", "hits"] assert scatter.y[0] == pytest.approx(np.median(rt_dict["all"])) assert scatter.y[1] == pytest.approx(np.median(rt_dict["hits"])) + for i, rts in enumerate(rt_dict.values()): + median = np.median(rts) + assert scatter.error_y.array[i] == pytest.approx(np.percentile(rts, 75) - median) + assert scatter.error_y.arrayminus[i] == pytest.approx(median - np.percentile(rts, 25)) assert fig.layout.yaxis.range == pytest.approx((0, 4)) @@ -223,6 +229,22 @@ def test_plot_trial_timeseries_falls_back_to_violet_for_unknown_category(): assert fig.data[0].marker.color == "violet" +def test_plot_lever_duration_medians_from_dict_reports_median_and_iqr_offsets(): + duration_dict = {"all": [60.0, 80.0, 100.0, 120.0, 140.0], "false_alarms": []} + + fig = plots.plot_lever_duration_medians_from_dict(duration_dict, as_html=False) + scatter = fig.data[0] + + assert list(scatter.x) == ["all", "false_alarms"] + median = np.median(duration_dict["all"]) + assert scatter.y[0] == pytest.approx(median) + # Error bars are offsets from the median, so together they span the IQR. + assert scatter.error_y.array[0] == pytest.approx(np.percentile(duration_dict["all"], 75) - median) + assert scatter.error_y.arrayminus[0] == pytest.approx(median - np.percentile(duration_dict["all"], 25)) + assert np.isnan(scatter.y[1]) + assert fig.layout.yaxis.range is None + + # -- Every plotting function also supports as_html=True, returning an embeddable HTML div. -- @@ -235,6 +257,7 @@ def test_plot_trial_timeseries_falls_back_to_violet_for_unknown_category(): lambda: plots.plot_sdt_pie(1, 2, 3, 4, as_html=True), lambda: plots.plot_rt_median([0.1, 0.2, 0.3], as_html=True), lambda: plots.plot_rt_medians_from_dict({"all": [0.1, 0.2]}, as_html=True), + lambda: plots.plot_lever_duration_medians_from_dict({"all": [80.0, 90.0]}, as_html=True), lambda: plots.plot_single_yvalue(0.5, as_html=True), lambda: plots.plot_roc(0.5, 0.5, as_html=True), lambda: plots.plot_dprime(1.0, as_html=True), diff --git a/tests/test_session.py b/tests/test_session.py index 4a6b6f8..3ff9f76 100644 --- a/tests/test_session.py +++ b/tests/test_session.py @@ -1,5 +1,6 @@ import datetime import json +import warnings import pathlib import shutil @@ -261,3 +262,253 @@ def test_generate_regressors_rejects_unsupported_protocol(tmp_path): with pytest.raises(NotImplementedError): session.generate_regressors(trials_df, metadata, str(timestamps_path), output_dir=str(tmp_path)) + + +# -- Lever push durations (from a sibling `_lever-durations.csv`) -- + + +def _lever_push_mask(trials): + return (trials["outcome"] == "precued") | trials["sdt_type"].isin(["hit", "false_alarm"]) + + +def test_find_lever_durations_pairs_bids_named_json_with_sibling_csv(gonogo_session_with_lever_durations): + json_path, durations_path = gonogo_session_with_lever_durations + + assert session.find_lever_durations(json_path) == durations_path + + +def test_find_lever_durations_returns_none_without_a_sibling_csv(gonogo_session_json_path): + assert session.find_lever_durations(gonogo_session_json_path) is None + + +def test_get_trials_assigns_lever_durations_to_push_trials_in_order(gonogo_session_with_lever_durations): + json_path, _ = gonogo_session_with_lever_durations + + trials = session.get_trials(json_path) + + is_push = _lever_push_mask(trials) + assert list(trials.loc[is_push, "lever_duration"]) == [100.0 + i for i in range(is_push.sum())] + assert trials.loc[~is_push, "lever_duration"].isna().all() + + +def test_get_trials_omits_lever_duration_column_without_a_durations_file(gonogo_session_json_path): + trials = session.get_trials(gonogo_session_json_path) + + assert "lever_duration" not in trials.columns + + +def test_get_trials_uses_explicit_lever_durations_path_over_sibling(gonogo_session_with_lever_durations, tmp_path): + json_path, durations_path = gonogo_session_with_lever_durations + durations = pd.read_csv(durations_path) + override_path = tmp_path / "override.csv" + durations.assign(duration=durations["duration"] + 1000).to_csv(override_path, index=False) + + trials = session.get_trials(json_path, lever_durations=str(override_path)) + + assert trials["lever_duration"].min() == 1100.0 + + +def test_get_trials_matches_lever_durations_by_push_id_not_file_order(gonogo_session_with_lever_durations): + json_path, durations_path = gonogo_session_with_lever_durations + durations = pd.read_csv(durations_path) + durations.iloc[::-1].to_csv(durations_path, index=False) + + trials = session.get_trials(json_path) + + assert list(trials.loc[_lever_push_mask(trials), "lever_duration"]) == list(durations["duration"]) + + +def test_get_trials_warns_and_leaves_lever_durations_empty_on_count_mismatch(gonogo_session_with_lever_durations): + json_path, durations_path = gonogo_session_with_lever_durations + pd.read_csv(durations_path).iloc[:-1].to_csv(durations_path, index=False) + + with pytest.warns(UserWarning, match="lever pushes"): + trials = session.get_trials(json_path) + + assert "lever_duration" in trials.columns + assert trials["lever_duration"].isna().all() + + +def _write_lever_duration_trials_csv(tmp_path, write_trials_csv): + rows = [ + dict(outcome="correct", correction=False, sdt_type="hit", lever_duration=80.0), + dict(outcome="correct", correction=False, sdt_type="hit", lever_duration=100.0), + dict(outcome="correct", correction=True, sdt_type="hit", lever_duration=200.0), + dict(outcome="incorrect", correction=False, sdt_type="false_alarm", lever_duration=120.0), + dict(outcome="correct", correction=False, sdt_type="correct_rejection", lever_duration=np.nan), + dict(outcome="precued", correction=False, sdt_type=None, lever_duration=60.0), + ] + for row in rows: + row["response"] = "leverpush" if not np.isnan(row["lever_duration"]) else None + row["response_time"] = 0.5 if row["sdt_type"] in ("hit", "false_alarm") else np.nan + return write_trials_csv(tmp_path, "trials.csv", "A1", "2022-01-01", "gonogo", "expX", rows=rows) + + +def _iqr(values): + return np.percentile(values, 75) - np.percentile(values, 25) + + +@pytest.mark.parametrize( + "subset, durations, durations_wc", + [ + ("", [80.0, 100.0, 120.0, 60.0], [80.0, 100.0, 200.0, 120.0, 60.0]), + ("_cued", [80.0, 100.0, 120.0], [80.0, 100.0, 200.0, 120.0]), + ("_hits", [80.0, 100.0], [80.0, 100.0, 200.0]), + ("_false_alarms", [120.0], [120.0]), + ("_precued", [60.0], [60.0]), + ], +) +def test_summary_reports_lever_duration_stats_per_subset(tmp_path, write_trials_csv, subset, durations, durations_wc): + result = session.summary(_write_lever_duration_trials_csv(tmp_path, write_trials_csv)) + + for suffix, values in (("", durations), ("_wc", durations_wc)): + assert result[f"lever_duration_mean{subset}{suffix}"] == pytest.approx(np.mean(values)) + assert result[f"lever_duration_median{subset}{suffix}"] == pytest.approx(np.median(values)) + assert result[f"lever_duration_iqr{subset}{suffix}"] == pytest.approx(_iqr(values)) + if len(values) > 1: + assert result[f"lever_duration_sd{subset}{suffix}"] == pytest.approx(np.std(values, ddof=1)) + else: + assert np.isnan(result[f"lever_duration_sd{subset}{suffix}"]) + + +def test_summary_reports_nan_lever_duration_stats_for_empty_subsets(tmp_path, write_trials_csv): + csv_path = _write_lever_duration_trials_csv(tmp_path, write_trials_csv) + df = pd.read_csv(csv_path) + df[df["sdt_type"] != "false_alarm"].to_csv(csv_path, index=False) + + result = session.summary(csv_path) + + for stat in session.LEVER_DURATION_STATS: + assert np.isnan(result[f"lever_duration_{stat}_false_alarms"]) + assert np.isnan(result[f"lever_duration_{stat}_false_alarms_wc"]) + + +def test_summary_has_forty_lever_duration_fields_when_durations_are_present(gonogo_session_with_lever_durations): + json_path, _ = gonogo_session_with_lever_durations + + result = session.summary(json_path) + + # 4 stats x 5 subsets x with/without corrections. + assert len([key for key in result if key.startswith("lever_duration_")]) == 40 + + +def test_summary_omits_lever_duration_fields_without_durations(gonogo_session_json_path): + result = session.summary(gonogo_session_json_path) + + assert not any(key.startswith("lever_duration") for key in result) + + +def test_generate_report_includes_lever_duration_plot_only_when_durations_are_present( + gonogo_session_with_lever_durations, gonogo_session_json_path, tmp_path +): + json_path, _ = gonogo_session_with_lever_durations + + with_durations = pathlib.Path(session.generate_report(json_path, output_dir=str(tmp_path))) + without_durations = pathlib.Path(session.generate_report(gonogo_session_json_path, output_dir=str(tmp_path))) + + assert "Median lever push duration" in with_durations.read_text(encoding="utf-8") + assert "Median lever push duration" not in without_durations.read_text(encoding="utf-8") + + +def test_get_trials_warns_and_leaves_lever_durations_empty_without_expected_columns( + gonogo_session_with_lever_durations, +): + json_path, durations_path = gonogo_session_with_lever_durations + pd.read_csv(durations_path).rename(columns={"push_id": "id"}).to_csv(durations_path, index=False) + + with pytest.warns(UserWarning, match="missing column.*push_id"): + trials = session.get_trials(json_path) + + assert trials["lever_duration"].isna().all() + + +@pytest.mark.parametrize("push_ids", ["duplicated", "gapped"]) +def test_get_trials_warns_and_leaves_lever_durations_empty_on_bad_push_ids( + gonogo_session_with_lever_durations, push_ids +): + json_path, durations_path = gonogo_session_with_lever_durations + durations = pd.read_csv(durations_path) + # Same number of pushes as lever push trials, but the IDs can't be trusted to line up with them. + durations.loc[1, "push_id"] = 0 if push_ids == "duplicated" else len(durations) + durations.to_csv(durations_path, index=False) + + with pytest.warns(UserWarning, match="push_id"): + trials = session.get_trials(json_path) + + assert trials["lever_duration"].isna().all() + + +def test_summary_accepts_a_trials_dataframe(gonogo_session_with_lever_durations): + json_path, _ = gonogo_session_with_lever_durations + + from_dataframe = session.summary(session.get_trials(json_path)) + from_json = session.summary(json_path) + + assert from_dataframe.keys() == from_json.keys() + for key, value in from_json.items(): + # NaN != NaN, so missing values have to be compared separately. + assert from_dataframe[key] == value or (pd.isna(from_dataframe[key]) and pd.isna(value)), key + + +def test_preprocess_session_uses_lever_durations_override_throughout(gonogo_session_with_lever_durations, tmp_path): + json_path, durations_path = gonogo_session_with_lever_durations + override_path = tmp_path / "override.csv" + pd.read_csv(durations_path).to_csv(override_path, index=False) + # A stale sibling that would warn if anything fell back to auto-detection. + pd.read_csv(durations_path).iloc[:-1].to_csv(durations_path, index=False) + out_dir = tmp_path / "out" + out_dir.mkdir() + + with warnings.catch_warnings(): + warnings.simplefilter("error", UserWarning) + session.preprocess_session(json_path, output_dir=str(out_dir), lever_durations=str(override_path)) + + (report_path,) = out_dir.glob("*_report-session.html") + assert "Median lever push duration" in report_path.read_text(encoding="utf-8") + + +def test_generate_report_omits_lever_duration_plot_when_durations_could_not_be_matched( + gonogo_session_with_lever_durations, tmp_path +): + json_path, durations_path = gonogo_session_with_lever_durations + pd.read_csv(durations_path).iloc[:-1].to_csv(durations_path, index=False) + + with pytest.warns(UserWarning): + report_path = pathlib.Path(session.generate_report(json_path, output_dir=str(tmp_path))) + + assert "Median lever push duration" not in report_path.read_text(encoding="utf-8") + + +def test_get_trials_warns_and_leaves_lever_durations_empty_on_empty_file(gonogo_session_with_lever_durations): + json_path, durations_path = gonogo_session_with_lever_durations + pathlib.Path(durations_path).write_text("") + + with pytest.warns(UserWarning, match="Can't match"): + trials = session.get_trials(json_path) + + assert trials["lever_duration"].isna().all() + + +def test_get_trials_warns_and_leaves_lever_durations_empty_on_non_numeric_duration( + gonogo_session_with_lever_durations, +): + json_path, durations_path = gonogo_session_with_lever_durations + durations = pd.read_csv(durations_path) + # Right number of pushes with valid IDs, so the non-numeric value is the only problem. + durations["duration"] = durations["duration"].astype(object) + durations.loc[1, "duration"] = "not-a-number" + durations.to_csv(durations_path, index=False) + + with pytest.warns(UserWarning, match="not-a-number"): + trials = session.get_trials(json_path) + + assert trials["lever_duration"].isna().all() + + +def test_summary_reads_metadata_from_a_filtered_trials_dataframe(gonogo_session_json_path): + trials = session.get_trials(gonogo_session_json_path) + # Dropping the first trial leaves an index that no longer starts at 0. + result = session.summary(trials.iloc[1:]) + + assert result["animal_id"] == "MM229" + assert result["protocol"] == "gonogo" diff --git a/tests/test_subject.py b/tests/test_subject.py index d1f7c3f..13665f9 100644 --- a/tests/test_subject.py +++ b/tests/test_subject.py @@ -145,3 +145,27 @@ def test_collate_sessions_raises_when_all_sessions_are_excluded_and_not_ignoring with pytest.raises(FileNotFoundError, match="not marked as excluded"): subject.collate_sessions(directory=str(tmp_path), output_dir=None, ignore=False) + + +# -- Lever push durations -- + + +def test_collate_sessions_includes_lever_duration_metrics_when_present(tmp_path, write_trials_csv): + write_trials_csv(tmp_path, "a_trials.csv", "A1", "2022-01-01", "gonogo", "expX", lever_durations=True) + write_trials_csv(tmp_path, "b_trials.csv", "A1", "2022-01-02", "gonogo", "expX") + + subject_df = subject.collate_sessions(directory=str(tmp_path), output_dir=None) + + assert subject_df["lever_duration_median"].iloc[0] == pytest.approx(80.0) + assert subject_df["lever_duration_median_hits"].iloc[0] == pytest.approx(80.0) + # Sessions without durations are left empty. + assert pd.isna(subject_df["lever_duration_median"].iloc[1]) + + +def test_collate_sessions_omits_lever_duration_metrics_when_no_session_has_them(tmp_path, write_trials_csv): + write_trials_csv(tmp_path, "a_trials.csv", "A1", "2022-01-01", "gonogo", "expX") + write_trials_csv(tmp_path, "b_trials.csv", "A1", "2022-01-02", "gonogo", "expX") + + subject_df = subject.collate_sessions(directory=str(tmp_path), output_dir=None) + + assert not any(column.startswith("lever_duration") for column in subject_df.columns)