calab ===== .. py:module:: calab .. autoapi-nested-parse:: CaLab: calcium imaging analysis tools — deconvolution and data preparation. Attributes ---------- .. autoapisummary:: calab.DriftModel calab.SpikeModel Classes ------- .. autoapisummary:: calab.BiexpFitResult calab.CaDeconResult calab.CellGroundTruth calab.DeconConfig calab.DeconvolutionResult calab.HeadlessBrowser calab.KernelConfig calab.MarkovConfig calab.NoiseConfig calab.PhotobleachingConfig calab.PoissonConfig calab.RandomWalkDrift calab.SaturationConfig calab.SimulationConfig calab.SimulationResult calab.SinusoidalDrift calab.SolveTraceResult calab.presets Functions --------- .. autoapisummary:: calab.bandpass_filter calab.build_kernel calab.compute_lipschitz calab.compute_upsample_factor calab.decon calab.deconvolve_from_export calab.estimate_kernel calab.fit_biexponential calab.load_caiman calab.load_export_params calab.load_minian calab.load_tuning_data calab.run_deconvolution calab.run_deconvolution_full calab.save_for_tuning calab.simulate calab.solve_trace calab.tau_to_ar2 calab.tune Package Contents ---------------- .. py:class:: BiexpFitResult Bases: :py:obj:`NamedTuple` Result from bi-exponential kernel fitting. .. attribute:: tau_rise Slow-component rise time constant (seconds). :type: float .. attribute:: tau_decay Slow-component decay time constant (seconds). :type: float .. attribute:: beta Slow-component amplitude. :type: float .. attribute:: residual Fit residual (lower is better). :type: float .. attribute:: tau_rise_fast Fast-component rise time constant (seconds), 0 if single-component. :type: float .. attribute:: tau_decay_fast Fast-component decay time constant (seconds), 0 if single-component. :type: float .. attribute:: beta_fast Fast-component amplitude, 0 if single-component. :type: float .. py:attribute:: beta :type: float .. py:attribute:: beta_fast :type: float .. py:attribute:: residual :type: float .. py:attribute:: tau_decay :type: float .. py:attribute:: tau_decay_fast :type: float .. py:attribute:: tau_rise :type: float .. py:attribute:: tau_rise_fast :type: float .. py:class:: CaDeconResult Bases: :py:obj:`NamedTuple` Full result from CaDecon (automated deconvolution via InDeCa algorithm). .. attribute:: activity Deconvolved activity matrix, shape ``(n_cells, n_timepoints)``, float32. :type: np.ndarray .. attribute:: alphas Per-cell scaling factors, shape ``(n_cells,)``, float64. :type: np.ndarray .. attribute:: baselines Per-cell baseline estimates, shape ``(n_cells,)``, float64. :type: np.ndarray .. attribute:: pves Per-cell proportion of variance explained, shape ``(n_cells,)``, float64. :type: np.ndarray .. attribute:: kernel_slow Slow biexponential kernel waveform, float32. :type: np.ndarray .. attribute:: kernel_fast Fast biexponential kernel waveform, float32 (empty if single-component). :type: np.ndarray .. attribute:: fs Sampling rate in Hz. :type: float .. attribute:: metadata Extensible dict with biexp params, convergence info, h_free, etc. :type: dict .. py:attribute:: activity :type: numpy.ndarray .. py:attribute:: alphas :type: numpy.ndarray .. py:attribute:: baselines :type: numpy.ndarray .. py:attribute:: fs :type: float .. py:attribute:: kernel_fast :type: numpy.ndarray .. py:attribute:: kernel_slow :type: numpy.ndarray .. py:attribute:: metadata :type: dict .. py:attribute:: pves :type: numpy.ndarray .. py:class:: CellGroundTruth Bases: :py:obj:`pydantic.BaseModel` Ground truth for a single simulated cell. .. py:attribute:: alpha :type: float .. py:attribute:: clean_calcium :type: numpy.ndarray .. py:attribute:: model_config .. py:attribute:: snr :type: float .. py:attribute:: spikes :type: numpy.ndarray .. py:attribute:: tau_decay_s :type: float .. py:attribute:: tau_rise_s :type: float .. py:class:: DeconConfig Bases: :py:obj:`pydantic.BaseModel` Configuration for automated CaDecon deconvolution via the bridge. Fields map 1:1 to the BridgeConfig TypeScript interface in @calab/io. Only ``autorun`` is always serialized; all other fields use ``exclude_none=True`` so absent values fall through to browser defaults. .. py:attribute:: aspect_ratio :type: float | None .. py:attribute:: autorun :type: bool .. py:attribute:: convergence_tol :type: float | None .. py:attribute:: hp_filter_enabled :type: bool | None .. py:attribute:: lp_filter_enabled :type: bool | None .. py:attribute:: max_iterations :type: int | None .. py:attribute:: num_subsets :type: int | None .. py:attribute:: seed :type: int | None .. py:attribute:: target_coverage :type: float | None .. py:attribute:: upsample_target :type: int | None .. py:class:: DeconvolutionResult Bases: :py:obj:`NamedTuple` Full result from FISTA deconvolution. .. attribute:: activity Non-negative deconvolved activity estimates, same shape as input traces. :type: np.ndarray .. attribute:: baseline Estimated scalar baseline (per-trace if multi-trace input). :type: float | np.ndarray .. attribute:: reconvolution K*activity + baseline, the model fit to the trace. :type: np.ndarray .. attribute:: iterations Number of FISTA iterations run (per-trace if multi-trace input). :type: int | np.ndarray .. attribute:: converged Whether convergence criterion was met (per-trace if multi-trace input). :type: bool | np.ndarray .. py:attribute:: activity :type: numpy.ndarray .. py:attribute:: baseline :type: float | numpy.ndarray .. py:attribute:: converged :type: bool | numpy.ndarray .. py:attribute:: iterations :type: int | numpy.ndarray .. py:attribute:: reconvolution :type: numpy.ndarray .. py:class:: HeadlessBrowser(*, visible: bool = False) Manages a Playwright browser for headless CaDecon / CaTune runs. Can be used as a context manager for single runs, or kept alive across multiple ``decon()`` / ``tune()`` calls for batch processing. .. rubric:: Examples Single run:: with HeadlessBrowser() as hb: result = calab.decon(traces, fs, headless=hb, autorun=True) Batch (reuses browser across datasets):: with HeadlessBrowser() as hb: for traces in datasets: result = calab.decon(traces, fs, headless=hb, autorun=True) .. py:method:: close() -> None Shut down page, context, browser, and Playwright. Each teardown runs in its own try/except so a crashed browser (where ``context.close()`` raises) still lets ``browser.close()`` and ``pw.stop()`` run. Without this, a single teardown error would leak the remaining resources — most visibly an orphaned Chromium process that lingers until the Python process exits. .. py:method:: navigate(url: str) -> None Navigate the managed page to *url*. If the page is already on a different URL, it navigates there (effectively reusing the same tab). .. py:method:: start() -> None Launch the browser. Called automatically by ``__enter__``. If any initialization step raises after earlier steps succeeded, ``close()`` runs before re-raising. This prevents partial-init failures (e.g. Chromium launch fails after the Playwright driver started) from leaking driver processes and Chromium instances across retries. .. py:property:: is_alive :type: bool Whether the browser is still connected. .. py:property:: page :type: playwright.sync_api.Page The managed Playwright ``Page`` (raises if not started). .. py:class:: KernelConfig Bases: :py:obj:`pydantic.BaseModel` Double-exponential kernel: h(t) = exp(-t/tau_decay) - exp(-t/tau_rise). Attribution: standard calcium response model (CaImAn, Suite2p, CaLab). .. py:attribute:: tau_decay_cv :type: float .. py:attribute:: tau_decay_s :type: float .. py:attribute:: tau_rise_cv :type: float .. py:attribute:: tau_rise_s :type: float .. py:class:: MarkovConfig Bases: :py:obj:`pydantic.BaseModel` Two-state HMM spike generator (silent/active) with bursty firing. Attribution: CaLab web simulator Markov spike model. .. py:attribute:: model_type :type: Literal['markov'] :value: 'markov' .. py:attribute:: p_active_to_silent :type: float .. py:attribute:: p_silent_to_active :type: float .. py:attribute:: p_silent_to_active_cv :type: float .. py:attribute:: p_spike_when_active :type: float .. py:attribute:: p_spike_when_silent :type: float .. py:class:: NoiseConfig Bases: :py:obj:`pydantic.BaseModel` Noise model: Gaussian + optional Poisson (shot) noise. Attribution: Gaussian from CaLab web simulator. Shot noise from CASCADE (Rupprecht et al., 2021). .. py:attribute:: shot_noise_enabled :type: bool .. py:attribute:: shot_noise_fraction :type: float .. py:attribute:: snr :type: float .. py:attribute:: snr_spread :type: float .. py:class:: PhotobleachingConfig Bases: :py:obj:`pydantic.BaseModel` Exponential photobleaching: ``F(t) *= 1 - amp * (1 - exp(-t/tau))``. Attribution: NAOMi (Charles et al., 2019). .. py:attribute:: amplitude_cv :type: float .. py:attribute:: amplitude_fraction :type: float .. py:attribute:: decay_time_constant_s :type: float .. py:attribute:: enabled :type: bool .. py:class:: PoissonConfig Bases: :py:obj:`pydantic.BaseModel` Homogeneous Poisson spike generator. Attribution: standard model in OASIS (Friedrich et al., 2017) and CaImAn (Giovannucci et al., 2019). .. py:attribute:: model_type :type: Literal['poisson'] :value: 'poisson' .. py:attribute:: rate_hz :type: float .. py:class:: RandomWalkDrift Bases: :py:obj:`pydantic.BaseModel` Mean-reverting Gaussian random walk baseline drift (default). Models slow irregular baseline fluctuations from tissue movement, focus drift, and neuropil signal changes. From MLspike (Deneux et al., 2016, Nature Communications). .. py:attribute:: mean_reversion :type: float .. py:attribute:: model_type :type: Literal['random_walk'] :value: 'random_walk' .. py:attribute:: step_std_cv :type: float .. py:attribute:: step_std_fraction :type: float .. py:class:: SaturationConfig Bases: :py:obj:`pydantic.BaseModel` Hill equation indicator saturation: F_sat = F^n / (F^n + Kd^n). Attribution: MLspike (Deneux et al., 2016). .. py:attribute:: enabled :type: bool .. py:attribute:: hill_coefficient :type: float .. py:attribute:: k_d :type: float .. py:attribute:: k_d_cv :type: float .. py:class:: SimulationConfig Bases: :py:obj:`pydantic.BaseModel` Complete configuration for synthetic calcium trace generation. Per-cell variation (_cv fields) live on each config struct alongside the nominal value they modify. Alpha is here because it doesn't belong to any pipeline step. .. py:attribute:: alpha_cv :type: float .. py:attribute:: alpha_mean :type: float .. py:attribute:: drift :type: DriftModel .. py:attribute:: fs_hz :type: float .. py:attribute:: kernel :type: KernelConfig .. py:attribute:: noise :type: NoiseConfig .. py:attribute:: num_cells :type: int .. py:attribute:: num_timepoints :type: int .. py:attribute:: photobleaching :type: PhotobleachingConfig .. py:attribute:: saturation :type: SaturationConfig .. py:attribute:: seed :type: int .. py:attribute:: spike_model :type: SpikeModel .. py:attribute:: spike_sim_hz :type: float .. py:class:: SimulationResult Bases: :py:obj:`pydantic.BaseModel` Complete simulation result with observed traces and per-cell ground truth. .. py:attribute:: config :type: SimulationConfig .. py:attribute:: ground_truth :type: list[CellGroundTruth] .. py:attribute:: model_config .. py:attribute:: traces :type: numpy.ndarray .. py:class:: SinusoidalDrift Bases: :py:obj:`pydantic.BaseModel` Deterministic sinusoidal baseline drift. Useful as a simple test signal but not physically motivated. .. py:attribute:: amplitude_cv :type: float .. py:attribute:: amplitude_fraction :type: float .. py:attribute:: cycles_max :type: float .. py:attribute:: cycles_min :type: float .. py:attribute:: model_type :type: Literal['sinusoidal'] :value: 'sinusoidal' .. py:class:: SolveTraceResult Bases: :py:obj:`NamedTuple` Result from a single-trace InDeCa solve. .. attribute:: s_counts Spike counts at the original sampling rate, shape ``(n_timepoints,)``, float32. :type: np.ndarray .. attribute:: alpha Amplitude scaling factor. :type: float .. attribute:: baseline Estimated baseline. :type: float .. attribute:: threshold Spike threshold used. :type: float .. attribute:: pve Proportion of variance explained (0–1). :type: float .. attribute:: iterations Number of FISTA iterations run. :type: int .. attribute:: converged Whether the solver converged. :type: bool .. py:attribute:: alpha :type: float .. py:attribute:: baseline :type: float .. py:attribute:: converged :type: bool .. py:attribute:: iterations :type: int .. py:attribute:: pve :type: float .. py:attribute:: s_counts :type: numpy.ndarray .. py:attribute:: threshold :type: float .. py:class:: presets Built-in indicator presets. Each method returns a SimulationConfig. .. py:method:: clean(**overrides: object) -> SimulationConfig :staticmethod: Near-ideal traces: minimal noise, no drift. For algorithm debugging. .. py:method:: gcamp6f(**overrides: object) -> SimulationConfig :staticmethod: GCaMP6f. Time constants from Chen et al., 2013, Nature. .. py:method:: gcamp6m(**overrides: object) -> SimulationConfig :staticmethod: GCaMP6m. Moderate kinetics. Chen et al., 2013. .. py:method:: gcamp6s(**overrides: object) -> SimulationConfig :staticmethod: GCaMP6s. Slow kinetics, high SNR. Chen et al., 2013. .. py:method:: jgcamp8f(**overrides: object) -> SimulationConfig :staticmethod: jGCaMP8f. Fast indicator, noisier. Zhang et al., 2023. .. py:method:: ogb1(**overrides: object) -> SimulationConfig :staticmethod: OGB-1 synthetic dye. Stosiek et al., 2003. .. py:function:: bandpass_filter(trace: numpy.ndarray, tau_rise: float, tau_decay: float, fs: float) -> numpy.ndarray Apply FFT bandpass filter derived from kernel time constants. Delegates to Rust. .. py:function:: build_kernel(tau_rise: float, tau_decay: float, fs: float) -> numpy.ndarray Build double-exponential calcium kernel. Delegates to Rust. .. py:function:: compute_lipschitz(kernel: numpy.ndarray) -> float Compute Lipschitz constant. Delegates to Rust. .. py:function:: compute_upsample_factor(fs: float, target_fs: float) -> int Compute the upsample factor for a given sampling rate and target. Delegates to Rust. :param fs: Original sampling rate in Hz. :type fs: float :param target_fs: Target sampling rate in Hz. :type target_fs: float :returns: Upsampling multiplier (>= 1). :rtype: int .. py:function:: decon(traces: numpy.ndarray, fs: float = 30.0, timeout: float | None = None, port: int | None = None, app_url: str | None = None, open_browser: bool = True, headless: calab._bridge._headless.HeadlessBrowser | bool | None = None, *, autorun: bool = False, upsample_target: int | None = None, hp_filter_enabled: bool | None = None, lp_filter_enabled: bool | None = None, max_iterations: int | None = None, convergence_tol: float | None = None, num_subsets: int | None = None, target_coverage: float | None = None, aspect_ratio: float | None = None, seed: int | None = None) -> calab._compute.CaDeconResult | None Open CaDecon in the browser for automated deconvolution. Starts a localhost HTTP server serving the provided traces, opens CaDecon with a ``?bridge=`` parameter pointing to the server, and waits for the browser to export deconvolution results back. :param traces: Calcium traces, shape ``(n_cells, n_timepoints)`` or ``(n_timepoints,)``. :type traces: np.ndarray :param fs: Sampling rate in Hz. Default: 30.0. :type fs: float :param timeout: Seconds to wait for results. None = wait forever (until Ctrl-C). :type timeout: float, optional :param port: Port to bind to. None = auto-assign. :type port: int, optional :param app_url: Override CaDecon URL (for local dev). Default: GitHub Pages. :type app_url: str, optional :param open_browser: Whether to auto-open the browser. Default: True. :type open_browser: bool :param headless: ``None``/``False``: default (use ``webbrowser.open``). ``True``: create a temporary headless browser for this call. ``HeadlessBrowser``: reuse an existing browser instance (for batch). :type headless: HeadlessBrowser or bool or None :param autorun: If True, the solver starts automatically after loading. Default: False. :type autorun: bool :param upsample_target: Target sampling rate for upsampling. Must be > 0. :type upsample_target: int, optional :param hp_filter_enabled: Enable high-pass filter. :type hp_filter_enabled: bool, optional :param lp_filter_enabled: Enable low-pass filter. :type lp_filter_enabled: bool, optional :param max_iterations: Maximum solver iterations (1–200). :type max_iterations: int, optional :param convergence_tol: Convergence tolerance (0–1 exclusive). :type convergence_tol: float, optional :param num_subsets: Number of random subsets. Must be > 0. :type num_subsets: int, optional :param target_coverage: Target coverage fraction (0–1]. :type target_coverage: float, optional :param aspect_ratio: Subset aspect ratio. Must be > 0. :type aspect_ratio: float, optional :param seed: Random seed for subset placement. :type seed: int, optional :returns: Deconvolution results if received, None if timeout/cancelled. :rtype: CaDeconResult or None .. py:function:: deconvolve_from_export(traces: numpy.ndarray, params_path: str | pathlib.Path, return_full: bool = False) Run deconvolution using parameters from a CaTune export JSON. Loads parameters via :func:`load_export_params`, optionally applies the bandpass filter, and runs FISTA deconvolution. :param traces: Calcium traces, shape ``(n_timepoints,)`` or ``(n_cells, n_timepoints)``. :type traces: np.ndarray :param params_path: Path to the CaTune export JSON file. :type params_path: str or Path :param return_full: If True, return a :class:`~calab.DeconvolutionResult` namedtuple. Default False returns just the activity array. :type return_full: bool, optional :returns: Deconvolved activity (or full result if ``return_full=True``). :rtype: np.ndarray or DeconvolutionResult .. py:function:: estimate_kernel(traces_flat: numpy.ndarray, spikes_flat: numpy.ndarray, trace_lengths: numpy.ndarray, alphas: numpy.ndarray, baselines: numpy.ndarray, kernel_length: int, *, max_iters: int = 200, tol: float = 0.0001, warm_kernel: numpy.ndarray | None = None, smooth_lambda: float = 0.0) -> numpy.ndarray Estimate a free-form kernel from traces and spike trains. Delegates to Rust. :param traces_flat: Concatenated 1-D traces (all cells flattened). :type traces_flat: np.ndarray :param spikes_flat: Concatenated 1-D spike trains (matching traces_flat). :type spikes_flat: np.ndarray :param trace_lengths: Length of each individual trace in the flat arrays. :type trace_lengths: np.ndarray :param alphas: Per-trace amplitude scaling factors. :type alphas: np.ndarray :param baselines: Per-trace baseline estimates. :type baselines: np.ndarray :param kernel_length: Desired output kernel length in samples. :type kernel_length: int :param max_iters: Maximum FISTA iterations for kernel estimation. :type max_iters: int :param tol: Convergence tolerance. :type tol: float :param warm_kernel: Kernel from a previous iteration for warm-start. :type warm_kernel: np.ndarray, optional :param smooth_lambda: Total-variation smoothness penalty weight. :type smooth_lambda: float :returns: Estimated free-form kernel, shape ``(kernel_length,)``, float32. :rtype: np.ndarray .. py:function:: fit_biexponential(h_free: numpy.ndarray, fs: float, *, refine: bool = True, skip: int = 0, warm: BiexpFitResult | None = None) -> BiexpFitResult Fit a bi-exponential model to a free-form kernel. Delegates to Rust. :param h_free: Free-form kernel (1-D). :type h_free: np.ndarray :param fs: Sampling rate in Hz. :type fs: float :param refine: Whether to refine with a fast (second) component. :type refine: bool :param skip: Number of leading samples to skip in the fit. :type skip: int :param warm: Previous fit result for warm-start. :type warm: BiexpFitResult, optional :rtype: BiexpFitResult .. py:function:: load_caiman(path: str, trace_key: str = 'estimates/C', fs_key: str = 'params/data/fr', fs: float | None = None) -> tuple[numpy.ndarray, dict] Load traces from a CaImAn HDF5 results file. :param path: Path to the CaImAn HDF5 file (e.g., ``caiman_results.hdf5``). :type path: str :param trace_key: HDF5 key for the traces array. Default: ``"estimates/C"``. :type trace_key: str :param fs_key: HDF5 key for the sampling rate. Default: ``"params/data/fr"``. :type fs_key: str :param fs: Override sampling rate. If provided, ``fs_key`` is ignored. :type fs: float, optional :returns: * **traces** (*np.ndarray*) -- Traces array, shape ``(n_cells, n_timepoints)``, dtype float64. * **metadata** (*dict*) -- Metadata dict with keys: ``source``, ``sampling_rate_hz``, ``num_cells``, ``num_timepoints``. :raises ImportError: If h5py is not installed. :raises FileNotFoundError: If the HDF5 file does not exist. :raises KeyError: If ``trace_key`` is not found in the file. .. py:function:: load_export_params(path: str | pathlib.Path) -> dict Load deconvolution parameters from a CaTune export JSON. Parses the JSON export schema 1.1.0 (defined in src/lib/export.ts) and returns the parameters needed for deconvolution. :param path: Path to the CaTune export JSON file. :type path: str or Path :returns: Dictionary with keys: ``tau_rise``, ``tau_decay``, ``lambda_``, ``fs``, ``filter_enabled``. :rtype: dict :raises FileNotFoundError: If the JSON file does not exist. :raises KeyError: If required parameter fields are missing. .. py:function:: load_minian(path: str, trace_key: str = 'C', fs: float | None = None) -> tuple[numpy.ndarray, dict] Load traces from a Minian Zarr output directory. :param path: Path to the Minian Zarr directory (e.g., ``minian_output/``). :type path: str :param trace_key: Zarr key for the traces array. Default: ``"C"``. :type trace_key: str :param fs: Sampling rate in Hz. Minian does not store this, so it must be provided (or will default to None in metadata). :type fs: float, optional :returns: * **traces** (*np.ndarray*) -- Traces array, shape ``(n_cells, n_timepoints)``, dtype float64. * **metadata** (*dict*) -- Metadata dict with keys: ``source``, ``sampling_rate_hz``, ``num_cells``, ``num_timepoints``. :raises ImportError: If zarr is not installed. :raises FileNotFoundError: If the Zarr directory does not exist. :raises KeyError: If ``trace_key`` is not found in the store. .. py:function:: load_tuning_data(path: str | pathlib.Path) -> tuple[numpy.ndarray, dict] Load calcium traces and metadata saved by :func:`save_for_tuning`. :param path: Path stem (without extension), matching the ``path`` argument used in :func:`save_for_tuning`. :type path: str or Path :returns: * **traces** (*np.ndarray*) -- Loaded traces array, dtype float64. * **metadata** (*dict*) -- Metadata from the JSON sidecar. :raises FileNotFoundError: If either ``{path}.npy`` or ``{path}_metadata.json`` is missing. .. py:function:: run_deconvolution(traces: numpy.ndarray, fs: float, tau_r: float, tau_d: float, lam: float, max_iters: int = 2000, conv_mode: str = 'fft', constraint: str = 'nonneg') -> numpy.ndarray Run FISTA deconvolution on one or more calcium traces. Delegates to the Rust solver via calab._solver. :param traces: Input traces, shape ``(n_timepoints,)`` for a single trace or ``(n_cells, n_timepoints)`` for multiple traces. :type traces: np.ndarray :param fs: Sampling rate in Hz. :type fs: float :param tau_r: Rise time constant in seconds. :type tau_r: float :param tau_d: Decay time constant in seconds. :type tau_d: float :param lam: L1 penalty (sparsity regularization strength). :type lam: float :param max_iters: Maximum number of FISTA iterations, by default 2000. :type max_iters: int, optional :param conv_mode: Convolution mode: ``'fft'`` (default) or ``'banded'`` (O(T) AR2). :type conv_mode: str, optional :param constraint: Constraint type: ``'nonneg'`` (default, L1 + non-negative) or ``'box01'`` (box constraint [0, 1], no L1 penalty). :type constraint: str, optional :returns: Non-negative activity estimates, same shape as input ``traces``. :rtype: np.ndarray .. py:function:: run_deconvolution_full(traces: numpy.ndarray, fs: float, tau_r: float, tau_d: float, lam: float, max_iters: int = 2000, conv_mode: str = 'fft', constraint: str = 'nonneg') -> DeconvolutionResult Run FISTA deconvolution returning full results. :param traces: Input traces, shape ``(n_timepoints,)`` for a single trace or ``(n_cells, n_timepoints)`` for multiple traces. :type traces: np.ndarray :param fs: Sampling rate in Hz. :type fs: float :param tau_r: Rise time constant in seconds. :type tau_r: float :param tau_d: Decay time constant in seconds. :type tau_d: float :param lam: L1 penalty (sparsity regularization strength). :type lam: float :param max_iters: Maximum number of FISTA iterations, by default 2000. :type max_iters: int, optional :param conv_mode: Convolution mode: ``'fft'`` (default) or ``'banded'`` (O(T) AR2). :type conv_mode: str, optional :param constraint: Constraint type: ``'nonneg'`` (default, L1 + non-negative) or ``'box01'`` (box constraint [0, 1], no L1 penalty). :type constraint: str, optional :returns: Namedtuple with fields: ``activity``, ``baseline``, ``reconvolution``, ``iterations``, ``converged``. :rtype: DeconvolutionResult .. py:function:: save_for_tuning(traces: numpy.ndarray, fs: float, path: str | pathlib.Path, metadata: dict | None = None) -> None Save calcium traces in CaTune-compatible format. Creates two files: - ``{path}.npy`` -- Float64 array, shape ``(n_cells, n_timepoints)``, C-contiguous, little-endian (`` my_recording.npy + my_recording_metadata.json **Converting from Minian** (users have zarr/xarray):: import zarr import calab store = zarr.open("minian_output", mode="r") traces = store["C"][:] # shape: (n_cells, n_frames) fs = 30.0 # user must know their frame rate calab.save_for_tuning(traces, fs, "my_recording") :param traces: Calcium traces, shape ``(n_timepoints,)`` for a single trace or ``(n_cells, n_timepoints)`` for multiple traces. :type traces: np.ndarray :param fs: Sampling rate in Hz. :type fs: float :param path: Output path stem (without extension). Will create ``{path}.npy`` and ``{path}_metadata.json``. :type path: str or Path :param metadata: Additional metadata to include in the JSON sidecar. Keys are merged into the output; built-in keys take precedence. :type metadata: dict, optional :raises ValueError: If ``traces`` has more than 2 dimensions. .. py:function:: simulate(config: SimulationConfig | None = None, **kwargs: object) -> SimulationResult Generate synthetic calcium imaging traces with full ground truth. :param config: Full configuration. If None, a default config is created. :type config: SimulationConfig, optional :param \*\*kwargs: Override fields on the default/provided config (e.g., num_cells=50, seed=123). :returns: Contains traces array and per-cell ground truth. :rtype: SimulationResult .. py:function:: solve_trace(trace: numpy.ndarray, tau_rise: float, tau_decay: float, fs: float, *, upsample_factor: int = 1, max_iters: int = 500, tol: float = 0.0001, hp_enabled: bool = False, lp_enabled: bool = False, warm_counts: numpy.ndarray | None = None, lambda_: float = 0.0) -> SolveTraceResult Run the InDeCa pipeline on a single trace. Delegates to Rust. :param trace: 1-D calcium trace. :type trace: np.ndarray :param tau_rise: Time constants in seconds. :type tau_rise: float :param tau_decay: Time constants in seconds. :type tau_decay: float :param fs: Sampling rate in Hz. :type fs: float :param upsample_factor: Upsampling multiplier (1 = no upsampling). :type upsample_factor: int :param max_iters: Maximum FISTA iterations. :type max_iters: int :param tol: Convergence tolerance. :type tol: float :param hp_enabled: Enable high-pass / low-pass filtering. :type hp_enabled: bool :param lp_enabled: Enable high-pass / low-pass filtering. :type lp_enabled: bool :param warm_counts: Spike counts from a previous iteration (at original rate) for warm-start. :type warm_counts: np.ndarray, optional :param lambda_: L1 sparsity penalty. :type lambda_: float :rtype: SolveTraceResult .. py:function:: tau_to_ar2(tau_rise: float, tau_decay: float, fs: float) -> tuple[float, float, float, float] Derive AR(2) coefficients from tau parameters. Pure Python (trivial math, no solver needed). :returns: (g1, g2, d, r) where g1 = d + r, g2 = -(d * r), d = exp(-dt/tau_decay), r = exp(-dt/tau_rise). :rtype: tuple[float, float, float, float] .. py:function:: tune(traces: numpy.ndarray, fs: float = 30.0, timeout: float | None = None, port: int | None = None, app_url: str | None = None, open_browser: bool = True) -> dict | None Open CaTune in the browser for interactive parameter tuning. Starts a localhost HTTP server serving the provided traces, opens CaTune with a ``?bridge=`` parameter pointing to the server, and waits for the user to export parameters from the web app. :param traces: Calcium traces, shape ``(n_cells, n_timepoints)`` or ``(n_timepoints,)``. :type traces: np.ndarray :param fs: Sampling rate in Hz. Default: 30.0. :type fs: float :param timeout: Seconds to wait for params. None = wait forever (until Ctrl-C). :type timeout: float, optional :param port: Port to bind to. None = auto-assign. :type port: int, optional :param app_url: Override CaTune URL (for local dev). Default: GitHub Pages. :type app_url: str, optional :param open_browser: Whether to auto-open the browser. Default: True. :type open_browser: bool :returns: Exported parameters dict if received, None if timeout/cancelled. Keys: ``tau_rise``, ``tau_decay``, ``lambda_``, ``fs``, ``filter_enabled``. :rtype: dict or None .. py:data:: DriftModel .. py:data:: SpikeModel