diff --git a/.github/workflows/cron-scenario.yml b/.github/workflows/cron-scenario.yml index 4f4355a5..c1489828 100644 --- a/.github/workflows/cron-scenario.yml +++ b/.github/workflows/cron-scenario.yml @@ -167,6 +167,7 @@ jobs: interpolated_trajectories session_metrics latency_tests + np_sync summary_emails inputs_videos decision @@ -275,6 +276,10 @@ jobs: dataframe_count="$(docker exec vr4mice-mysql-restore mysql -N -uroot -psimple -e 'SELECT COUNT(*) FROM base_analysis.__data_frame;')" summary_plot_count="$(docker exec vr4mice-mysql-restore mysql -N -uroot -psimple -e 'SELECT COUNT(*) FROM base_analysis.__summary_plots;')" summary_email_count="$(docker exec vr4mice-mysql-restore mysql -N -uroot -psimple -e 'SELECT COUNT(*) FROM summary_emails.summary_plot_email WHERE send_error IS NULL;')" + np_sync_barcode_sync_count="" + if grep -qx 'np_sync' "$RUNNER_TEMP/mysql-dump/schemas.txt"; then + np_sync_barcode_sync_count="$(docker exec vr4mice-mysql-restore mysql -N -uroot -psimple -e 'SELECT COUNT(*) FROM np_sync.__barcode_sync;')" + fi if [ "${summary_email_count}" -lt 1 ] 2>/dev/null; then echo "summary_emails tables:" >&2 docker exec vr4mice-mysql-restore mysql -uroot -psimple -e 'SHOW TABLES FROM summary_emails;' >&2 || true @@ -290,6 +295,11 @@ jobs: echo "base_analysis.__data_frame rows: ${dataframe_count}" echo "base_analysis.__summary_plots rows: ${summary_plot_count}" echo "summary_emails.summary_plot_email successful rows: ${summary_email_count}" + if [ -n "${np_sync_barcode_sync_count}" ]; then + echo "np_sync.__barcode_sync rows: ${np_sync_barcode_sync_count}" + else + echo "np_sync schema not present in dump; skipping np_sync.__barcode_sync row-count check" + fi test "${dataset_count}" -ge 1 test "${batch_count}" -ge 2 @@ -298,6 +308,9 @@ jobs: test "${state_count}" -ge 1 test "${dataframe_count}" -ge 1 test "${summary_plot_count}" -ge 1 + if [ -n "${np_sync_barcode_sync_count}" ]; then + test "${np_sync_barcode_sync_count}" -ge 1 + fi - name: Verify decision tables (--aws populated) run: | diff --git a/dj_pipeline/.env.compose.example b/dj_pipeline/.env.compose.example index 381ca147..f6e73f51 100644 --- a/dj_pipeline/.env.compose.example +++ b/dj_pipeline/.env.compose.example @@ -11,6 +11,9 @@ DB_DATA_PATH=/mnt/database/vr4mice/vr4mice_database/database/ SHARED_PATH=/mnt/database/shared DATA_PATH=/mnt/database/vr4mice/vr4mice_database/data SCREEN_RECORDINGS_PATH=/mnt/neuropixel_data/vr4mice/raw_screen_recordings +# Optional: point to a real FreelyMoving_Mouse_NP checkout to enable np_sync. +# If omitted, docker-compose defaults to ./np_pipeline_stub and np_sync is skipped. +NP_PIPELINE_SRC_PATH=./np_pipeline_stub JUPYTER_PORT=8887 CLIENT_IMAGE=mmathislab/vr4mice_app:0.1.0 CLIENT_CONTAINER_NAME=vr4mice_${USER} diff --git a/dj_pipeline/.env.example b/dj_pipeline/.env.example index 6d4da1a1..ef899ec0 100644 --- a/dj_pipeline/.env.example +++ b/dj_pipeline/.env.example @@ -1,5 +1,6 @@ # Pipeline / DataJoint credentials (loaded into the client container via docker-compose env_file). # DJ_HOST port must match DB_PORT in .env.compose when using the local Docker database. +# For the "np_sync" mode: DJ_HOST must match np_pipeline's DJ_HOST (same MySQL server). DJ_HOST=your-db-hostname-or-ip:3306 DJ_USER=your-db-user DJ_PWD=your-db-password diff --git a/dj_pipeline/Makefile b/dj_pipeline/Makefile index fb261f1c..faf35fbe 100644 --- a/dj_pipeline/Makefile +++ b/dj_pipeline/Makefile @@ -37,11 +37,16 @@ BUILD_ARGS := \ .PHONY: build_all git_info git-whoami git-user-env git-user-local up_all base_install down_all client_build client_up \ client_stop client_rm client_down ipython bash mysql aws notebook sync_build \ sync_up add-cron cron-local cron-aws cron-now cron-aws-only cron-reboot import-dumps export-dumps \ - check-compose-project + check-compose-project format check-compose-project: @bash docker/check_compose_conflict.sh "$(COMPOSE_PROJECT)" +# Matches CI's formatter.yml: black==22.6 against dj_pipeline. +format: + pip install --quiet black==22.6 + black . + # Fail with logs if the client service is not running (common after compose project / name mismatch). define verify_client_running @if [ "$$($(DOCKER_COMPOSE) ps --status running -q client 2>/dev/null | wc -l | tr -d ' ')" = "0" ]; then \ diff --git a/dj_pipeline/README.md b/dj_pipeline/README.md index de92e8bd..a94d95a6 100644 --- a/dj_pipeline/README.md +++ b/dj_pipeline/README.md @@ -298,6 +298,7 @@ When using the local Docker database, `DJ_HOST` port must match `DB_PORT` in `.e - `COMPOSE_PROJECT` (default `vr4mice`) — must match between `make`, cron scripts, and manual `docker compose -p …` calls - `DB_BIND_IP`, `DB_PORT`, `MYSQL_ROOT_PASSWORD` - `DB_DATA_PATH`, `SHARED_PATH`, `DATA_PATH`, `SCREEN_RECORDINGS_PATH` +- `NP_PIPELINE_SRC_PATH` (optional): path to a real `FreelyMoving_Mouse_NP` checkout. If unset, compose mounts `./np_pipeline_stub`, and only `run.py np_sync` is skipped while the rest of the pipeline remains functional. - `CLIENT_IMAGE`, `CLIENT_CONTAINER_NAME`, `DB_CONTAINER_NAME`, `CLIENT_NETWORK_MODE`, `JUPYTER_PORT` Remote/AWS DB credentials for scheduled AWS runs live in `.env-aws` (copy from `.env-aws.example`); this file is **not** committed. diff --git a/dj_pipeline/cron_scenario.py b/dj_pipeline/cron_scenario.py index 634f9c2c..a46dde35 100755 --- a/dj_pipeline/cron_scenario.py +++ b/dj_pipeline/cron_scenario.py @@ -295,6 +295,33 @@ def import_decision_schema(): return decision + def import_np_sync_schema(): + from vr4mice.schema import np_sync + + return np_sync + + # np_sync depends on the separate np_pipeline package/repo and should not + # fail the whole cron scenario when unavailable in a deployment/CI setup. + try: + np_sync = import_np_sync_schema() + except Exception as err: + logger.warning( + "Skipping np_sync: np_pipeline is not available (%s). " + "Behavioral analysis is unaffected.", + err, + ) + else: + from vr4mice.utils.populate_helpers import populate_pending + + run_step( + "np_sync.BarcodeSync.populate", + lambda: populate_pending( + np_sync.BarcodeSync, + np_sync.BarcodeSync().key_source, + logger=logger, + ), + ) + decision = run_import("import decision schema", import_decision_schema) if decision: run_step( diff --git a/dj_pipeline/docker-compose.yml b/dj_pipeline/docker-compose.yml index 8323340f..5fe01ab6 100644 --- a/dj_pipeline/docker-compose.yml +++ b/dj_pipeline/docker-compose.yml @@ -1,3 +1,5 @@ +name: ${COMPOSE_PROJECT:-vr4mice} + services: db: image: ${DB_IMAGE:-mysql:8.0} @@ -23,6 +25,7 @@ services: - MPLCONFIGDIR=/app/.cache/matplotlib - XDG_CONFIG_HOME=/app/.cache - PATH=/app/.local/bin:/usr/local/bin:/usr/bin:/bin + - PYTHONPATH=/app/np_pipeline/src:/app/np_pipeline:${PYTHONPATH} - TAG=${TAG} network_mode: ${CLIENT_NETWORK_MODE:-host} # Run via bash so the bind-mounted script does not need +x on shared/NFS checkouts. @@ -36,7 +39,8 @@ services: - "${SCREEN_RECORDINGS_PATH:-/mnt/neuropixel_data/vr4mice/raw_screen_recordings}:/vr4mice_screen_recordings" - "./base/base_min_schemas:/base_schemas" - "./base/base_actions:/base_actions" - ports: + - "${NP_PIPELINE_SRC_PATH:-./np_pipeline_stub}:/app/np_pipeline" + ports: - "${JUPYTER_PORT:-8887}:8888" working_dir: /app env_file: diff --git a/dj_pipeline/run.py b/dj_pipeline/run.py index 58059732..ff6f5214 100644 --- a/dj_pipeline/run.py +++ b/dj_pipeline/run.py @@ -25,6 +25,9 @@ "latency": compute latencies based on photodiode signals "inputs_videos": process input videos and extract frames "decision": analyze decision-making metrics + "np_sync": align VR time to Neuropixels time via barcodes (requires the + np_pipeline package/repo; skips gracefully, without affecting any other + mode, if it isn't installed or no dataset has a linked NP recording yet) "maintenance": rebuild DataJoint lineage tables (one-time setup) """ @@ -77,9 +80,10 @@ def check_folder_existence(folder_path): "sync_days", "inputs_videos", "decision", + "np_sync", "maintenance", ], - help="Mode to execute: 'connect', 'populate', 'summary', 'dlc', 'fetch', 'sync_days', 'analysis', 'inputs_videos', 'decision', 'maintenance'", + help="Mode to execute: 'connect', 'populate', 'summary', 'dlc', 'fetch', 'sync_days', 'analysis', 'inputs_videos', 'decision', 'np_sync', 'maintenance'", ) args = parser.parse_args() @@ -281,6 +285,29 @@ def check_folder_existence(folder_path): decision.PredictionModel10Windows().populate() decision.DecisionPoints10Windows().populate() + elif args.mode == "np_sync": + # np_sync depends on the separate np_pipeline package/repo, which may be + # absent (ModuleNotFoundError) or unreachable (e.g. its own DB connection + # fails while its schema modules activate). Either way, isolate the + # failure to this mode so every other (purely behavioral) mode is + # unaffected. + try: + from vr4mice.schema import np_sync + except Exception as err: + logger.warning( + "Skipping np_sync: np_pipeline is not available (%s). " + "Behavioral analysis is unaffected.", + err, + ) + else: + from vr4mice.utils.populate_helpers import populate_pending + + populate_pending( + np_sync.BarcodeSync, + np_sync.BarcodeSync().key_source, + logger=logger, + ) + elif args.mode == "fetch": from vr4mice.actions.fetch_data import fetch_data diff --git a/dj_pipeline/vr4mice/analysis/np_sync.py b/dj_pipeline/vr4mice/analysis/np_sync.py new file mode 100644 index 00000000..ee783a9d --- /dev/null +++ b/dj_pipeline/vr4mice/analysis/np_sync.py @@ -0,0 +1,78 @@ +"""Fit a VR-time-to-NP-time alignment from barcode values shared by both streams.""" + +from __future__ import annotations + +from dataclasses import dataclass + +import numpy as np +import scipy.interpolate +import scipy.stats + +# Number of leading (chronologically earliest) VR barcode events excluded from the +# regression fit by default. DLC-live starts receiving data slightly after the +# Unity/game stream does, which is why downstream analysis already drops +# `trial == 1` as a DLC-live initialization trial (see vr4mice/analysis/analysis.py). +# Barcodes in that same early window can carry an unreliable onset_time_unity, +# biasing the fit if used as tie points, so we drop the earliest few by default. +DEFAULT_SKIP_FIRST_N_BARCODES = 10 + + +@dataclass(frozen=True) +class BarcodeAlignmentFit: + """Linear fit + interpolator mapping VR time to NP time.""" + + slope: float + intercept: float + r2: float + interpol_func: scipy.interpolate.interp1d + shared_barcodes: np.ndarray + + +def align_barcodes( + vr_times: np.ndarray, + vr_values: np.ndarray, + np_times: np.ndarray, + np_values: np.ndarray, + skip_first_n_barcodes: int = 0, +) -> BarcodeAlignmentFit: + """Fit VR time -> NP time from barcode values shared between both streams. + + Alignment approach (intersect1d + linregress + interp1d) ported from + cross_analysis_schemas/schemas/vr_np_sync.py::BarcodeSync.align_barcodes in + https://github.com/AdaptiveMotorControlLab/auxPipelines-DataJoint_Mathis, + adapted for this repo's VR (vr4mice) / NP (np_pipeline) schemas. + + Args: + vr_times: VR-side barcode onset times, ordered by event index (chronological). + vr_values: VR-side barcode integer payloads, same order as `vr_times`. + np_times: NP-side barcode onset times, ordered by event index. + np_values: NP-side barcode integer payloads, same order as `np_times`. + skip_first_n_barcodes: number of leading (earliest) VR events to exclude + before matching, to avoid the DLC-live startup-lag window. + """ + if skip_first_n_barcodes: + vr_times = vr_times[skip_first_n_barcodes:] + vr_values = vr_values[skip_first_n_barcodes:] + + shared_barcodes, vr_index, np_index = np.intersect1d( + vr_values, np_values, return_indices=True + ) + + vr_shared_times = np.asarray(vr_times)[vr_index] + np_shared_times = np.asarray(np_times)[np_index] + + linreg = scipy.stats.linregress(vr_shared_times, np_shared_times) + interpol_func = scipy.interpolate.interp1d( + vr_shared_times, + np_shared_times, + bounds_error=False, + fill_value="extrapolate", + ) + + return BarcodeAlignmentFit( + slope=linreg.slope, + intercept=linreg.intercept, + r2=linreg.rvalue**2, + interpol_func=interpol_func, + shared_barcodes=shared_barcodes, + ) diff --git a/dj_pipeline/vr4mice/schema/np_sync.py b/dj_pipeline/vr4mice/schema/np_sync.py new file mode 100644 index 00000000..62ddcecf --- /dev/null +++ b/dj_pipeline/vr4mice/schema/np_sync.py @@ -0,0 +1,285 @@ +"""Alignment of VR (behavior) time to Neuropixels (NP) native time via shared barcodes. + +Not every VR ``Dataset`` has a corresponding NP recording. ``BarcodeSync.key_source`` +contains only strict VR/NP matches (via ``base.Base`` + +``session_link.RecordingSessionLink`` + successful NP probe barcode extraction), +so VR-only datasets are excluded from this table and are not visited by +``populate()``. + +Cross-repo foreign keys: ``BarcodeSync`` references ``np_pipeline`` tables +directly via ``-> ``. This requires ``vr4mice`` and ``np_pipeline`` to share one +MySQL server (``DJ_HOST``) and one process-global ``dj.conn()`` — both already +satisfied by this repo's connection setup. `np_sync` is imported as its own +isolated step in ``run.py``/``cron_scenario.py``, so a missing ``np_pipeline`` +only skips NP sync, not the rest of the pipeline. +""" + +import pickle + +import datajoint as dj +import numpy as np + +try: + from np_pipeline.schemas import barcodes as np_barcodes, session_link +except ModuleNotFoundError: + import sys + + # Keep np_pipeline imports resilient when process CWD is not /app. + sys.path.insert(0, "/app/np_pipeline/src") + sys.path.insert(0, "/app/np_pipeline") + from np_pipeline.schemas import barcodes as np_barcodes, session_link + +from vr4mice.analysis.np_sync import DEFAULT_SKIP_FIRST_N_BARCODES, align_barcodes +from vr4mice.schema import base +from vr4mice.schema import barcodes as vr_barcodes +from vr4mice.schema import vr4mice +from vr4mice.utils.logger import Logger +from vr4mice.utils.schema_config import get_schema + +schema_name = "np_sync" +schema = get_schema(schema_name, locals()) + +logger = Logger.get_logger() + + +def _candidate_relation(): + """Rows where one VR dataset has at least one NP probe with successful barcodes.""" + vr = ( + base.Base + * (vr_barcodes.TeensyBarcodes & 'extraction_status = "success"').proj() + ) + npx = ( + session_link.RecordingSessionLink + * (np_barcodes.ProbeBarcodeExtraction & 'extraction_status = "success"').proj() + ) + + return vr * npx + + +@schema +class BarcodeSync(dj.Computed): + """Linear fit + interpolator mapping VR Unity/game time to NP probe time. + + One row per related VR/NP barcode pair with successful extraction on both sides. + The table is explicitly keyed by ``TeensyBarcodes`` (VR side) and + ``ProbeBarcodeExtraction`` (NP side). ``key_source`` enforces strict matching + through ``base.Base`` + ``RecordingSessionLink`` + ``recording_id`` so + ``populate()`` only visits related rows; sessions with no strict NP match are + invisible to this table. + + Downstream code should not fetch ``slope``/``intercept``/``interpol_func`` + directly; use ``align_timepoints``/``align_timepoints_lin`` instead. + """ + + definition = """ + -> vr_barcodes.TeensyBarcodes + -> np_barcodes.ProbeBarcodeExtraction + --- + skip_first_n_barcodes: int32 # leading VR barcode events excluded from the fit + slope: float64 # Slope of the linear fit mapping VR time to NP time + intercept: float64 # Intercept of the linear fit mapping VR time to NP time + r2: float64 # R-squared value of the linear fit + interpol_func: # pickled scipy.interpolate.interp1d, VR time -> NP time + barcode_overlap: float64 # Fraction of NP barcodes also found on the VR side + """ + + @property + def key_source(self): + source = _candidate_relation() + unexpected = sorted(set(source.heading.primary_key) - set(self.primary_key)) + if unexpected: + raise dj.DataJointError( + f"key_source primary key has attributes BarcodeSync lacks: {unexpected}" + ) + return source + + skip_first_n_barcodes = DEFAULT_SKIP_FIRST_N_BARCODES + min_shared_barcodes = 20 + min_barcode_overlap = 0.90 + + def make(self, key): + """Fit and insert one VR-time-to-NP-time alignment. + + Fetches decoded barcode events for one dataset/recording/probe key from both + `vr_barcodes.TeensyBarcodes.Event` (VR side) and + `np_barcodes.ProbeBarcodeExtraction.Event` (NP side), fits the alignment via + `vr4mice.analysis.np_sync.align_barcodes`, and inserts the resulting fit + parameters and pickled interpolator. On failure, records the error in + `vr4mice.FailedSession` and logs a warning instead of raising, matching the + error-handling convention used by `barcodes.TeensyTTL`/`TeensyBarcodes`. + """ + if vr4mice.FailedSession.should_skip(key, self.__class__.__name__, logger): + return + + try: + candidate_rows = (_candidate_relation() & key).to_dicts() + if not candidate_rows: + reason = "No NP candidate matched key fields" + vr4mice.FailedSession().add_entry( + f"{key['dataset']}", f"{self.__class__.__name__}", reason + ) + logger.warning( + "%s %s for dataset %s", + self.__class__.__name__, + reason, + key["dataset"], + ) + return + if len(candidate_rows) > 1: + reason = "Ambiguous NP candidates matched key fields" + vr4mice.FailedSession().add_entry( + f"{key['dataset']}", f"{self.__class__.__name__}", reason + ) + logger.warning( + "%s %s for dataset %s", + self.__class__.__name__, + reason, + key["dataset"], + ) + return + + vr_values, vr_times = (vr_barcodes.TeensyBarcodes.Event & key).to_arrays( + "barcode_value", "onset_time_unity", order_by="barcode_index" + ) + np_values, np_times = ( + np_barcodes.ProbeBarcodeExtraction.Event & key + ).to_arrays("barcode_value", "onset_time", order_by="barcode_index") + + if len(np_values) == 0: + reason = ( + "No NP barcode events found for key at populate time " + "(events may have been removed after key_source selection)" + ) + vr4mice.FailedSession().add_entry( + f"{key['dataset']}", f"{self.__class__.__name__}", reason + ) + logger.warning( + "%s %s. key: %s", + self.__class__.__name__, + reason, + key, + ) + return + + fit = align_barcodes( + vr_times, + vr_values, + np_times, + np_values, + skip_first_n_barcodes=self.skip_first_n_barcodes, + ) + + # Measure overlap on full streams (independent of fit-time skipping) + # so DEFAULT_SKIP_FIRST_N_BARCODES does not bias this quality metric. + full_shared_barcodes = np.intersect1d(vr_values, np_values) + np_unique_barcodes = np.unique(np_values) + barcode_overlap = len(full_shared_barcodes) / len(np_unique_barcodes) + + if len(fit.shared_barcodes) < self.min_shared_barcodes: + reason = ( + "Insufficient shared barcodes for reliable NP-VR alignment " + f"(shared={len(fit.shared_barcodes)}, " + f"min_required={self.min_shared_barcodes})" + ) + vr4mice.FailedSession().add_entry( + f"{key['dataset']}", f"{self.__class__.__name__}", reason + ) + logger.warning( + "%s %s for dataset %s", + self.__class__.__name__, + reason, + key["dataset"], + ) + return + + if barcode_overlap <= self.min_barcode_overlap: + reason = ( + "Insufficient NP-VR barcode overlap for reliable alignment " + f"(overlap={barcode_overlap:.4f}, " + f"min_required>{self.min_barcode_overlap:.4f})" + ) + vr4mice.FailedSession().add_entry( + f"{key['dataset']}", f"{self.__class__.__name__}", reason + ) + logger.warning( + "%s %s for dataset %s", + self.__class__.__name__, + reason, + key["dataset"], + ) + return + + self.insert1( + { + **key, + "skip_first_n_barcodes": self.skip_first_n_barcodes, + "slope": fit.slope, + "intercept": fit.intercept, + "r2": fit.r2, + "interpol_func": pickle.dumps(fit.interpol_func), + "barcode_overlap": barcode_overlap, + } + ) + logger.info( + "%s aligned %d shared barcodes for %s", + self.__class__.__name__, + len(fit.shared_barcodes), + key, + ) + + except Exception as err: + dataset = key["dataset"] + vr4mice.FailedSession().add_entry( + f"{dataset}", f"{self.__class__.__name__}", str(err) + ) + logger.warning( + "Can't populate %s, key: %s. Error: %s.", + self.__class__.__name__, + key, + err, + ) + return None + + @classmethod + def align_timepoints(cls, key, timepoints: list): + """Convert a list of VR times to NP times via the fitted interpolator. + + Preferred over `align_timepoints_lin` for most uses: the interpolator is fit + only within the range spanned by the shared barcode events, so it captures + any small clock drift between the VR and NP streams rather than assuming a + perfectly constant offset/rate. + + Args: + key: A restriction identifying exactly one `BarcodeSync` row (e.g. a + dataset/recording/DAQ key). + timepoints: VR-side times (`onset_time_unity`-style values) to convert. + `None` entries pass through as `None`; DataJoint-`NULL`/`NaN`-producing + extrapolation misses become `None` as well. + + Returns: + A list of NP-side times (or `None`), same length and order as `timepoints`. + """ + interpol_func = pickle.loads((cls & key).fetch1("interpol_func")) + timepoints = np.array( + [np.nan if tx is None else tx for tx in timepoints], dtype=np.float64 + ) + return [ + float(tx) if not np.isnan(tx) else None for tx in interpol_func(timepoints) + ] + + @classmethod + def align_timepoints_lin(cls, key, timepoints: list): + """Convert a list of VR times to NP times via the fitted line only (`y = slope*x + intercept`). + + Faster than `align_timepoints` and fine for a single global rate/offset, but + ignores any local clock drift the interpolator would otherwise correct for. + + Args: + key: A restriction identifying exactly one `BarcodeSync` row. + timepoints: VR-side times to convert; `None` entries pass through as `None`. + + Returns: + A list of NP-side times (or `None`), same length and order as `timepoints`. + """ + slope, intercept = (cls & key).fetch1("slope", "intercept") + return [tx * slope + intercept if tx is not None else None for tx in timepoints] diff --git a/dj_pipeline/vr4mice/utils/maintenance.py b/dj_pipeline/vr4mice/utils/maintenance.py index 67c1c7c5..654915f7 100644 --- a/dj_pipeline/vr4mice/utils/maintenance.py +++ b/dj_pipeline/vr4mice/utils/maintenance.py @@ -41,6 +41,7 @@ def _schema_pairs() -> Iterable[Tuple[str, object]]: inputs_videos, interpolated_trajectories, latency_tests, + np_sync, session_metrics, summary_emails, vr4mice, @@ -55,6 +56,7 @@ def _schema_pairs() -> Iterable[Tuple[str, object]]: ("summary_emails", summary_emails.schema), ("dlc", dlc.schema), ("barcodes", barcodes.schema), + ("np_sync", np_sync.schema), ("session_metrics", session_metrics.schema), ("interpolated_trajectories", interpolated_trajectories.schema), ("latency_tests", latency_tests.schema), diff --git a/docs/software/datajoint.md b/docs/software/datajoint.md index 9874a45e..2ae410d2 100644 --- a/docs/software/datajoint.md +++ b/docs/software/datajoint.md @@ -396,4 +396,40 @@ Extracts the binary sync trace from the sync ROI video. class AlignedVideoFrame(dj.Computed) ``` **Depends on:** `VideoSyncSignal`, `vr4mice.State` -Aligns game steps to video frames using photodiode when available; stores QA metrics. \ No newline at end of file +Aligns game steps to video frames using photodiode when available; stores QA metrics. + +### `barcodes.py` + +*Schema for barcodes recorded through the Teensy synchronization channel.* + +```python +class TeensyTTL(dj.Imported) +``` +**Depends on:** `vr4mice.DLC` +Raw Teensy barcode channel imported from a DLC PROC file. + +```python +class TeensyBarcodes(dj.Computed) +``` +**Depends on:** `TeensyTTL`, `dlc.DLCProcessor`, `base_analysis.DataFrame` +Barcode events decoded from the TTL signal sampled by the Teensy. + +### `np_sync.py` + +*Alignment of VR (behavior) time to Neuropixels (NP) native time via shared barcodes.* + +```python +class BarcodeSync(dj.Computed) +``` +**Depends on:** `barcodes.TeensyBarcodes`, `np_pipeline.schemas.barcodes.ProbeBarcodeExtraction` (table key), plus `base.Base` and `np_pipeline.schemas.session_link.RecordingSessionLink` in `key_source` filtering +Fits a linear regression + interpolator mapping VR Unity/game time to NP probe +time, from barcode values shared between `barcodes.TeensyBarcodes.Event` and +`np_pipeline`'s `ProbeBarcodeExtraction.Event`. `key_source` is restricted to +datasets with a linked, successfully barcode-decoded NP recording, so `populate()` +is a no-op — not an error — for VR-only sessions with no neural data. Use +`align_timepoints`/`align_timepoints_lin` to convert VR times to NP times. + +> **Cross-repo foreign keys.** `BarcodeSync` references `np_pipeline` tables +> directly via `-> `. Requires `vr4mice` and `np_pipeline` to share one MySQL +> server (`DJ_HOST`). `np_sync` is imported as its own isolated step in +> `run.py`/`cron_scenario.py`, so a missing `np_pipeline` only skips NP sync. \ No newline at end of file diff --git a/tests/unit/test_np_sync.py b/tests/unit/test_np_sync.py new file mode 100644 index 00000000..991c7f29 --- /dev/null +++ b/tests/unit/test_np_sync.py @@ -0,0 +1,218 @@ +"""Unit tests for the VR-to-NP barcode alignment fit, plus guards that +behavioral analysis stays unaffected when np_pipeline is unavailable. +""" + +import importlib +import sys +import ast +import textwrap +import pickle +from pathlib import Path + +import numpy as np +import pytest + +from np_sync import align_barcodes + +REPO_ROOT = Path(__file__).parent.parent.parent +RUN_PY = REPO_ROOT / "dj_pipeline" / "run.py" +CRON_SCENARIO_PY = REPO_ROOT / "dj_pipeline" / "cron_scenario.py" +SCHEMA_NP_SYNC_PY = REPO_ROOT / "dj_pipeline" / "vr4mice" / "schema" / "np_sync.py" + + +def _linear_barcode_streams(*, n=20, slope=2.0, intercept=100.0, skip_vr=0, skip_np=0): + """Two barcode streams sharing all-but-`skip_*` values, related by a known linear fit.""" + values = np.arange(n) + vr_times = np.arange(n, dtype=np.float64) + np_times = slope * vr_times + intercept + + vr_values = values[skip_vr:] if skip_vr else values + vr_times = vr_times[skip_vr:] if skip_vr else vr_times + np_values = values[skip_np:] if skip_np else values + np_times = np_times[skip_np:] if skip_np else np_times + + return vr_times, vr_values, np_times, np_values + + +def test_align_barcodes_recovers_known_linear_fit(): + vr_times, vr_values, np_times, np_values = _linear_barcode_streams( + slope=2.0, intercept=100.0 + ) + + fit = align_barcodes(vr_times, vr_values, np_times, np_values) + + assert fit.slope == pytest.approx(2.0) + assert fit.intercept == pytest.approx(100.0) + assert fit.r2 == pytest.approx(1.0) + assert len(fit.shared_barcodes) == 20 + + +def test_align_barcodes_uses_only_shared_barcode_values(): + # VR side is missing the first 3 barcode values, NP side the last 3. + vr_times, vr_values, np_times, np_values = _linear_barcode_streams( + n=20, skip_vr=3 + ) + np_values = np_values[:-3] + np_times = np_times[:-3] + + fit = align_barcodes(vr_times, vr_values, np_times, np_values) + + assert len(fit.shared_barcodes) == 14 # 20 - 3 leading - 3 trailing + assert fit.slope == pytest.approx(2.0) + + +def test_align_barcodes_skip_first_n_excludes_leading_vr_events(): + vr_times, vr_values, np_times, np_values = _linear_barcode_streams(n=20) + # Corrupt the first two VR onset times so an unskipped fit would be biased. + vr_times = vr_times.copy() + vr_times[:2] += 1000.0 + + biased_fit = align_barcodes(vr_times, vr_values, np_times, np_values) + assert biased_fit.slope != pytest.approx(2.0) + + corrected_fit = align_barcodes( + vr_times, vr_values, np_times, np_values, skip_first_n_barcodes=2 + ) + assert corrected_fit.slope == pytest.approx(2.0) + assert corrected_fit.intercept == pytest.approx(100.0) + assert len(corrected_fit.shared_barcodes) == 18 + + +def test_align_barcodes_interpol_func_maps_vr_time_to_np_time(): + vr_times, vr_values, np_times, np_values = _linear_barcode_streams( + slope=2.0, intercept=100.0 + ) + + fit = align_barcodes(vr_times, vr_values, np_times, np_values) + + assert fit.interpol_func(5.0) == pytest.approx(110.0) + + +def test_analysis_np_sync_has_no_datajoint_or_np_pipeline_dependency(monkeypatch): + """The pure alignment math must import and run with np_pipeline/datajoint absent.""" + monkeypatch.setitem(sys.modules, "datajoint", None) + monkeypatch.setitem(sys.modules, "np_pipeline", None) + + import np_sync as np_sync_analysis + + importlib.reload(np_sync_analysis) + + vr_times, vr_values, np_times, np_values = _linear_barcode_streams( + slope=2.0, intercept=100.0 + ) + fit = np_sync_analysis.align_barcodes(vr_times, vr_values, np_times, np_values) + assert fit.slope == pytest.approx(2.0) + + +def _function_source(file_text: str, def_line: str) -> str: + """Return the source of a top-level-indented function, up to the next `def `.""" + start = file_text.index(def_line) + rest = file_text[start + len(def_line) :] + end = rest.find("\n def ") + return def_line + (rest if end == -1 else rest[:end]) + + +def _load_schema_align_timepoints_method(): + """Load BarcodeSync.align_timepoints from source without importing DB deps.""" + text = SCHEMA_NP_SYNC_PY.read_text() + module = ast.parse(text) + + method_src = None + for node in module.body: + if isinstance(node, ast.ClassDef) and node.name == "BarcodeSync": + for item in node.body: + if isinstance(item, ast.FunctionDef) and item.name == "align_timepoints": + method_src = ast.get_source_segment(text, item) + break + break + + if method_src is None: + raise ValueError("BarcodeSync.align_timepoints not found") + + namespace = {"np": np, "pickle": pickle} + exec(textwrap.dedent(method_src), namespace) + return namespace["align_timepoints"] + + +def test_schema_align_timepoints_preserves_none_entries(): + align_timepoints = _load_schema_align_timepoints_method() + + class _FakeRelation: + def fetch1(self, field_name): + assert field_name == "interpol_func" + return pickle.dumps(np.square) + + class _FakeCls: + def __and__(self, key): + assert key == {"dataset": "dummy"} + return _FakeRelation() + + aligned = align_timepoints(_FakeCls(), {"dataset": "dummy"}, [1.0, None, 2.0]) + assert aligned == [1.0, None, 4.0] + + +def test_cron_scenario_core_schemas_import_excludes_np_sync(): + """import_core_schemas() must not import np_sync, so a missing np_pipeline + can't null out the whole behavioral core_schemas tuple (see import_np_sync_schema, + which is deliberately a separate run_import() call).""" + text = CRON_SCENARIO_PY.read_text() + core_schemas_src = _function_source(text, " def import_core_schemas():") + + assert "np_sync" not in core_schemas_src + assert "def import_np_sync_schema():" in text + + +def test_run_py_np_sync_mode_catches_broad_exception(): + """The np_sync CLI mode must not let an import/connection failure crash the + process; it should catch broadly and warn, not just ModuleNotFoundError.""" + text = RUN_PY.read_text() + start = text.index('elif args.mode == "np_sync":') + end = text.index('elif args.mode == "fetch":') + block = text[start:end] + + assert "from vr4mice.schema import np_sync" in block + assert "except Exception as err:" in block + assert "except ModuleNotFoundError" not in block + + +def test_schema_make_handles_empty_np_events_with_clear_reason(): + text = SCHEMA_NP_SYNC_PY.read_text() + make_src = _function_source(text, " def make(self, key):") + + assert "if len(np_values) == 0:" in make_src + assert "No NP barcode events found for key at populate time" in make_src + + +def test_schema_make_uses_key_only_without_identity_parsing(): + text = SCHEMA_NP_SYNC_PY.read_text() + make_src = _function_source(text, " def make(self, key):") + + assert "_dataset_identity" not in make_src + assert "_identity_matches_row" not in make_src + assert "_select_rows_by_np_event_count" not in make_src + assert "No NP candidate matched key fields" in make_src + + +def test_schema_make_fails_ambiguous_matches_directly_without_event_count_tie_break(): + text = SCHEMA_NP_SYNC_PY.read_text() + make_src = _function_source(text, " def make(self, key):") + + assert "Ambiguous NP candidates matched key fields" in make_src + assert "event_count matching" not in make_src + + +def test_schema_make_has_quality_gate_for_min_shared_and_overlap(): + text = SCHEMA_NP_SYNC_PY.read_text() + make_src = _function_source(text, " def make(self, key):") + + assert "min_shared_barcodes = 20" in text + assert "min_barcode_overlap = 0.90" in text + assert "Insufficient shared barcodes for reliable NP-VR alignment" in make_src + assert "Insufficient NP-VR barcode overlap for reliable alignment" in make_src + + +def test_schema_module_docstring_states_vr_only_sessions_excluded_from_key_source(): + text = SCHEMA_NP_SYNC_PY.read_text() + + assert "so VR-only datasets are excluded from this table" in text + assert "are not visited by" in text