deconv

deconv ¶

Matrix-based deconvolution with regularization.

This module provides deconvolution operations using convolution matrix construction and regularized pseudo-inverse. Supports both TSVD (Truncated SVD) and Tikhonov regularization methods.

Also includes batch SVD deconvolution functions for efficient multi-voxel processing.

References

- dcmri (https://github.com/dcmri/dcmri) utils.deconv()
- Wu et al. (2003). Tracer arrival timing-insensitive technique.
  Magn Reson Med. 50:164-174.
- Hansen PC (1998). Rank-Deficient and Discrete Ill-Posed Problems.

deconv ¶

deconv(
    c,
    aif,
    t,
    *,
    method="tsvd",
    tol=0.1,
    lambda_reg=None,
    circulant=False,
)

Deconvolve a tissue concentration curve by an arterial input function.

Computes the residue function R(t) such that: C(t) = AIF(t) * R(t)

where * denotes convolution. This is the inverse problem of pharmacokinetic modeling.

PARAMETER	DESCRIPTION
`c`	Tissue concentration curve. Shape: (n_times,). TYPE: `ndarray`
`aif`	Arterial input function. Shape: (n_times,). TYPE: `ndarray`
`t`	Time points in seconds. Shape: (n_times,). Must be monotonically increasing. TYPE: `ndarray`
`method`	Regularization method: - "tsvd": Truncated SVD. Robust and widely used in DSC-MRI. - "tikhonov": Tikhonov regularization. Smoother but may over-smooth. TYPE: `('tsvd', 'tikhonov')` DEFAULT: `"tsvd"`
`tol`	Regularization tolerance. For TSVD, singular values below tol * max(s) are truncated. For Tikhonov, controls the regularization strength if lambda_reg is not provided. TYPE: `float` DEFAULT: `0.1`
`lambda_reg`	Tikhonov regularization parameter. Only used if method="tikhonov". If not provided, computed as tol * max(singular_values). TYPE: `float` DEFAULT: `None`
`circulant`	If True, use block-circulant matrix for delay-insensitive deconvolution (cSVD/oSVD approach). This makes the result insensitive to timing differences between AIF and tissue curves. TYPE: `bool` DEFAULT: `False`

RETURNS	DESCRIPTION
`ndarray`	Residue function R(t). Shape: (n_times,). The residue function satisfies R(0) = 1 for a single-compartment model with no delay.

Notes

The deconvolution problem is ill-posed due to noise amplification. Regularization is essential for stable solutions.

TSVD (Truncated SVD): - Sets small singular values to zero - Simple and robust - May introduce ringing artifacts (Gibbs phenomenon)

Tikhonov regularization: - Damps small singular values smoothly - Produces smoother solutions - May over-smooth genuine features

Circulant matrix (cSVD/oSVD): - Makes deconvolution insensitive to AIF timing - Standard approach in DSC-MRI for CBF estimation - Requires uniform time sampling

Examples:

>>> import numpy as np
>>> from osipy.common.convolution import deconv
>>> t = np.linspace(0, 60, 61)
>>> aif = np.exp(-t / 10) * (1 - np.exp(-t / 2))  # Gamma variate AIF
>>> # Create synthetic tissue curve with known residue function
>>> from osipy.common.convolution import conv
>>> true_irf = np.exp(-t / 20)  # Exponential residue
>>> c = conv(aif, true_irf, t)
>>> # Deconvolve to recover residue function
>>> recovered_irf = deconv(c, aif, t, method="tsvd", tol=0.1)

References

.. [1] dcmri library: https://github.com/dcmri/dcmri .. [2] Wu et al. (2003). Magn Reson Med. 50:164-174. .. [3] Ostergaard et al. (1996). Magn Reson Med. 36:715-725.

deconv_osvd ¶

deconv_osvd(c, aif, t, *, osc_threshold=0.035)

Oscillation-index SVD deconvolution.

Implements the oSVD algorithm which selects the optimal truncation threshold based on minimizing oscillations in the residue function.

PARAMETER	DESCRIPTION
`c`	Tissue concentration curve. Shape: (n_times,). TYPE: `ndarray`
`aif`	Arterial input function. Shape: (n_times,). TYPE: `ndarray`
`t`	Time points in seconds. Shape: (n_times,). TYPE: `ndarray`
`osc_threshold`	Target oscillation index threshold. The algorithm finds the smallest truncation that keeps oscillations below this. TYPE: `float` DEFAULT: `0.035`

RETURNS	DESCRIPTION
`ndarray`	Residue function with minimized oscillations. Shape: (n_times,).

Notes

The oscillation index (OI) is defined as: OI = (1/N) * sum(|R[i] - R[i-1]|) / max(R)

oSVD iteratively increases the truncation threshold until OI < osc_threshold.

References

.. [1] Wu et al. (2003). Magn Reson Med. 50:164-174.

deconvolve_svd ¶

deconvolve_svd(ct, aif, dt=1.0, threshold=0.1)

Deconvolve tissue curve with AIF using SVD.

Solves the inverse problem to recover the impulse response function from the measured tissue curve and AIF.

PARAMETER	DESCRIPTION
`ct`	Tissue concentration curve, shape (n_timepoints,) or (n_timepoints, n_voxels). TYPE: `NDArray`
`aif`	Arterial input function, shape (n_timepoints,). TYPE: `NDArray`
`dt`	Time step in seconds. Default is 1.0. TYPE: `float` DEFAULT: `1.0`
`threshold`	SVD truncation threshold as fraction of maximum singular value. Default is 0.1 (10%). TYPE: `float` DEFAULT: `0.1`

RETURNS	DESCRIPTION
`irf`	Recovered impulse response function. TYPE: `NDArray`
`residue`	Residue function (cumulative integral of IRF). TYPE: `NDArray`

deconvolve_svd_batch ¶

deconvolve_svd_batch(ct, aif, dt=1.0, threshold=0.1)

Deconvolve multiple tissue curves with a single AIF.

Optimized batch deconvolution: the convolution matrix SVD is computed only once and applied to all voxels.

PARAMETER	DESCRIPTION
`ct`	Tissue concentration curves, shape (n_timepoints, n_voxels). TYPE: `NDArray`
`aif`	Arterial input function, shape (n_timepoints,). TYPE: `NDArray`
`dt`	Time step in seconds. Default is 1.0. TYPE: `float` DEFAULT: `1.0`
`threshold`	SVD truncation threshold. Default is 0.1. TYPE: `float` DEFAULT: `0.1`

RETURNS	DESCRIPTION
`irf`	Recovered impulse response functions, shape (n_timepoints, n_voxels). TYPE: `NDArray`
`residue`	Residue functions, shape (n_timepoints, n_voxels). TYPE: `NDArray`