MNE-Python compatibility (AMICAICA)

AMICAICA fits AMICA directly from an MNE-Python Raw/Epochs and hands the result back through the standard MNE ICA surface. It is additive: the scikit-learn-style AMICA interface and the EEGLAB output are unchanged; this is a second entry point for MNE users, not a replacement.

MNE is an optional dependency, so import pamica never requires it. Install the extra and import the wrapper explicitly:

uv add "pamica[mne]"        # or, from a source checkout: uv sync --extra mne

import mne
from pamica.mne_compat import AMICAICA

raw = mne.io.read_raw_eeglab("subject.set", preload=True)

ica = AMICAICA(n_mix=3, random_state=42).fit(raw, picks="eeg", max_iter=100)

sources = ica.get_sources(raw)     # an mne.io.RawArray of component activations
maps = ica.get_components()        # scalp maps, shape (n_channels, n_components)
ica.plot_components()              # native mne.viz topographies

# Reconstruct with some components removed:
clean = ica.apply(raw.copy(), exclude=[0, 3])

fit accepts a Raw or Epochs (epochs are concatenated along time, as MNE's own ICA does), any MNE picks selector, and forwards remaining keywords (max_iter, lrate, do_newton, ...) to AMICA.fit, so a setting left out takes the backend's default, as with AMICA (for example lrate=0.1). It rejects non-finite input, supports principal component analysis (PCA) reduction (pcakeep/pcadb) and rank-deficient data (see Rank-reduced fits and the PCA residual), and a degenerate (diverged) fit is refused by the consumer methods rather than emitting NaNs.

Choosing the backend

AMICAICA fits through AMICA, so it takes the same backend parameter: the PyTorch backend by default, or the Apple-GPU MLX backend with backend="mlx" (issue #313). Everything on this page works on both.

# raw: 32 average-referenced EEG channels, so rank 31
ica = AMICAICA(backend="mlx", random_state=42).fit(raw, picks="eeg", pcakeep=31)
clean = ica.apply(raw.copy(), exclude=[0])

As with AMICA, device (and a dtype fit keyword) apply to the PyTorch backend only and raise ValueError with backend="mlx", and backend="mlx" without MLX installed raises ImportError.

The export reads the fitted mean, sphere and per-model centers through the backend-agnostic float64 accessors AMICA.get_mean(), get_sphere() and get_model_center(), so both backends take the same code path. Precision follows the backend:

• A PyTorch fit (float64 by default) exports at float64 parity.
• An MLX fit computes in float32, so its unmixing matrix, mean and centers carry float32 rounding and the export is float32-consistent rather than float64-parity. get_sources agrees with ica.amica_.transform (which runs in float32) within float32 tolerance, and excluding a component changes the data by that component's back-projection to the same tolerance.
• Reconstruction does not depend on the backend's precision. MNE's mixing is the float64 pseudo-inverse of the exported unmixing and the PCA basis is orthonormal, so apply with nothing excluded returns the input to float64 round-off on either backend, residual included. Precision on the MLX backend gives the measured figures.

Interoperating with mne.preprocessing.ICA

to_mne_ica() returns a fully-populated mne.preprocessing.ICA, so the entire MNE ICA ecosystem (plotting, find_bads_eog/_ecg, exclusion workflows) works on an AMICA decomposition:

mne_ica = ica.to_mne_ica()
eog_idx, scores = mne_ica.find_bads_eog(raw)
mne_ica.plot_scores(scores)

The wrapper's get_sources, apply, get_components, plot_components and plot_sources delegate to this object, so they reproduce AMICA.transform exactly: MNE computes sources as unmixing_matrix_ @ pca_components_[:n_components_] @ (X / pre_whitener_ - pca_mean_), and the export maps pamica's mean, symmetric-ZCA sphere and unmixing into those matrices (writing the sphere as V diag(1/√e) Vᵀ with V orthonormal so MNE's scalp maps come out in channel space). The equivalence to_mne_ica().get_sources(raw) == AMICA.transform(X) is pinned by the test suite on real sample EEG.

Rank-reduced fits and the PCA residual

A fit can be rank-reduced by pcakeep/pcadb, or by automatic rank detection on Maxwell-filtered, average-referenced or interpolated data. AMICA then models only the retained PCA subspace, and n_components_ is that rank. The export still carries the full PCA basis, so apply restores the residual the reduction discarded, as MNE's own ICA does (issue #322). With nothing excluded, apply returns the input; excluding a component removes that component and nothing else.

This goes beyond the Fortran reference, whose back-projection reconstructs only the retained subspace. MNE's own n_pca_components gives that behavior back:

ica = AMICAICA(random_state=42).fit(raw, pcakeep=20, max_iter=100)

# Default: the residual is restored.
clean = ica.apply(raw.copy(), exclude=[0])

# Reference behavior: rank-reduced reconstruction.
reduced = ica.apply(raw.copy(), exclude=[0], n_pca_components=ica.n_components_)

Sources, component maps and the log-likelihood are the same either way; only reconstruction differs. The fitted basis is ica.pca_components_ (retained rows first, then the residual) with ica.pca_explained_variance_, and pamica vs. AMICA has the full rationale.

Multi-model fits

AMICA can learn a mixture of ICA models (n_models > 1). MNE's ICA represents only one unmixing matrix, so each model is exported as its own single-model mne.preprocessing.ICA, and the per-sample model dominance (which model best explains each timepoint) is exposed directly, since MNE has no concept for it:

ica = AMICAICA(n_models=2, random_state=42).fit(raw, max_iter=100)

# Per-model: model_idx selects the model on every consumer method.
sources_m1 = ica.get_sources(raw, model_idx=1)
ica.plot_components(model_idx=1)
model1 = ica.to_mne_ica(model_idx=1)   # a standard ICA for model 1

# Model dominance over time (P(model | sample), columns sum to 1):
prob = ica.get_model_probability(raw)  # (n_models, n_samples)
ica.plot_model_probability(raw)        # per-model probability + best-model LL

get_model_probability/plot_model_probability build on the public AMICA.model_loglik/model_probability accessors, which score arbitrary data through the fitted sphere and mean. Each per-model export folds that model's data-space center into pca_mean_, so to_mne_ica(model_idx=h).get_sources(raw) reproduces AMICA.transform(X, model_idx=h) (with X the picked channel array) for every model, not just the first.

Inspecting pamica-specific metadata

An mne.preprocessing.ICA has no field for AMICA's adaptive source densities or component sharing, so rather than drop them, the wrapper exposes them directly:

from pamica.mne_compat import AMICAICA, PDFTYPE_NAMES

ica = AMICAICA(n_models=2, random_state=42).fit(raw, max_iter=100)

families = ica.get_pdftype(model_idx=0)        # (n_components,) codes 0-4
names = [PDFTYPE_NAMES[c] for c in families]    # e.g. "generalized_gaussian"
rho = ica.get_rho(model_idx=0)                  # (n_mix, n_components) GG shape
shared = ica.shared_components()                # [(model, comp), ...] groups

get_pdftype returns each component's density family (0 generalized Gaussian, 1 super-Gaussian cosh, 2 Gaussian, 3 logistic, 4 sub-Gaussian cosh; they differ per component only under the adaptive switcher pdftype=1). get_rho is the generalized-Gaussian shape (meaningful for pdftype=0). shared_components lists components merged across models by share_comps (empty otherwise). The same three accessors exist on the scikit-learn-style AMICA.

Separation-quality metrics

The MIR/PMI metrics are available directly on an MNE object, so MNE-side users get the same separation-quality numbers as EEGLAB-side users:

mir_nats, variance = ica.mir(raw, model_idx=0)   # mutual information reduced
mi_matrix = ica.pmi(raw, model_idx=0)            # pairwise MI between sources

Both extract the fitted channels from the Raw/Epochs and delegate to AMICA.mir/pmi, reproducing the array API exactly.

pamica.mne_compat.AMICAICA

Fit AMICA from MNE objects and interoperate with mne.preprocessing.ICA.

The wrapper fits pamica's natural-gradient AMICA backend on the data of an MNE :class:~mne.io.Raw or :class:~mne.Epochs and lets MNE consume the result: :meth:get_sources, :meth:apply, :meth:get_components, :meth:plot_components and :meth:plot_sources all delegate to a real :class:mne.preprocessing.ICA built by :meth:to_mne_ica.

For a multi-model fit (n_models > 1) each model is exported as its own single-model MNE ICA (to_mne_ica(model_idx=...) / the model_idx argument on the consumer methods), and the per-sample model dominance -- which MNE's ICA cannot represent -- is exposed directly by :meth:get_model_probability / :meth:plot_model_probability.

Separation-quality metrics (issue #133) are available directly on an MNE object: :meth:mir (Mutual Information Reduction) and :meth:pmi (pairwise mutual information between sources). The pamica-specific fitted metadata MNE cannot hold -- source-density family, GG shape, component sharing -- is inspectable via :meth:get_pdftype / :meth:get_rho / :meth:shared_components.

The fit runs on the PyTorch backend by default and on the Apple-GPU MLX backend with backend="mlx" (issue #313); see the precision note in the Notes below.

Parameters:

Name	Type	Description	Default
n_models	int	Number of ICA models to learn (AMICA n_models).	1
n_mix	int	Number of mixture components per source (AMICA n_mix).	3
random_state	int or None	Seed for the AMICA fit (passed through as the backend seed) and stored on the exported :class:~mne.preprocessing.ICA.	None
device	str or device	Torch device for the fit (None = auto; the float64 parity backend falls back to CPU when auto-selection lands on MPS). See :class:AMICA. PyTorch backend only: with backend="mlx" it must stay None.	None
verbose	bool	Whether the underlying :class:AMICA prints fit progress.	True
backend	(torch, mlx)	Which :class:AMICA backend fits the data, forwarded to :class:AMICA. "mlx" is the Apple-GPU backend, float32 only; it requires MLX (ImportError here otherwise) and takes neither device nor a dtype fit keyword.	"torch"

Raises:

Type	Description
ValueError	If backend is not "torch" or "mlx", or device is set with backend="mlx".
ImportError	If backend="mlx" and MLX is not installed.

Attributes:

Name	Type	Description
amica_	AMICA	The fitted pamica model (holds all n_models models).
info_	Info	The picked measurement info the fit was run on (channel subset only).
ch_names_	list of str	Names of the fitted channels, in order.
n_components_	int	Number of ICA components. Equals the number of fitted channels unless the data are rank-deficient, in which case the model is sized to the detected numerical rank (issue #223).
pre_whitener_	np.ndarray of shape (n_channels, 1)	Per-channel-type scaling applied before fitting, following MNE's own ICA convention (one std per channel type, applied as X / pre_whitener_).
pca_components_	np.ndarray of shape (n_channels, n_channels) or None	The full orthonormal PCA basis of the fit, in pre-whitened channel space, computed once at fit time and exported as the MNE ICA's pca_components_. The first n_components_ rows span the subspace AMICA modeled. For a rank-reduced fit the remaining rows span the PCA residual that reduction discarded, which MNE's apply restores by default (issue #322); a full-rank fit has no such rows. Shared by every model of a multi-model fit (one sphere per fit). None after a degenerate fit, which exports nothing.
pca_explained_variance_	np.ndarray of shape (n_channels,) or None	Variance of the pre-whitened fit data along each row of pca_components_, in descending order: the retained rows carry the largest covariance eigenvalues and the residual rows the ones the reduction discarded (clipped at 0 against round-off). None after a degenerate fit.
reject_by_annotation_	bool	Whether the last Raw fit dropped bad_*-annotated samples (issue #251). Always False for an Epochs fit.
good_sample_mask_	np.ndarray of bool or None	For a Raw fit, which samples of the recording's timeline the fit consumed: True exactly for the fitted columns, False for samples covered by bad_* annotations (when rejection is on) and for samples outside a start/stop range. None for an Epochs fit.
converged_	bool	Whether the last fit ended usable (not degenerate). A degenerate fit is kept for inspection but refused by the consumer methods (issue #50).
stop_reason_	str or None	Why the backend fit stopped (e.g. "max_iter", "nan_ll").

Notes

Model h's AMICA transform is S = W_fort @ (sphere @ (X - mean) - c_h), where c_h is that model's data-space center (identically zero for a single model, since the c update is gated to n_models > 1). MNE computes sources as S = unmixing_matrix_ @ pca_components_[:n_components_] @ (X / pre_whitener_ - pca_mean_). X is scaled by channel type before fitting, exactly as MNE's own ICA does, so the two pipelines agree; AMICA's sphering absorbs a global rescale, so this changes nothing for single-channel-type data. Writing the symmetric-ZCA sphere as V @ diag(1/sqrt(e)) @ V.T with V orthonormal, the exported ICA for model h uses pca_components_ = V.T, unmixing_matrix_ = W_fort @ sphere @ V and pca_mean_ = mean + inv(sphere) @ c_h (which reduces to mean when c_h is zero). Keeping pca_components_ orthonormal is what makes MNE's get_components (scalp maps inv(sphere) @ inv(W_fort)) come out right, since MNE assumes orthonormal PCA rows. The mapping is pinned by a round-trip test (to_mne_ica(model_idx=h).get_sources(raw) equals amica_.transform(X, model_idx=h)).

Rank-deficient input -- Maxwell-filtered MEG, average referencing, channel interpolation, or an explicit pcakeep/pcadb -- is supported. The sphere is then (n_kept, n_channels) and has no eigendecomposition, so the retained basis is its right singular vectors instead. n_components_ reports the retained rank.

The discarded PCA subspace (the residual) is never part of the decomposition, and the export keeps it: pca_components_ holds the full (n_channels, n_channels) basis, retained rows first, and MNE's ICA.apply keeps rows past n_components_ as residual PCA components. So :meth:apply restores the residual by default, as MNE's own ICA does, and excluding a component removes only that component (issue #322). This goes beyond the Fortran reference, whose output has no representation of the residual (the written sphere's rows past numeigs are zero), so any back-projection from it drops the residual. Sources, component maps and the log-likelihood are unchanged; only reconstruction differs. For the reference's rank-reduced reconstruction pass MNE's own n_pca_components=ica.n_components_ to :meth:apply. See docs/guides/amica-differences.md and ADR 0005.

Precision: the export always reads the fit through the same float64 accessors (AMICA.get_sphere/get_mean/get_model_center), so both backends go through one code path. A PyTorch fit (float64 by default) exports at float64 parity. An MLX fit computes in float32, so its unmixing matrix, mean and centers carry float32 rounding (its sphere is the float64 one it was built from), and the export is float32-consistent rather than float64-parity: :meth:get_sources agrees with amica_.transform (which runs in float32) to float32 tolerance, and excluding a component changes the data by that component's back-projection to the same tolerance. Reconstruction does not depend on that precision: MNE's mixing is the float64 pseudo-inverse of the exported unmixing and the PCA basis is orthonormal to float64 round-off, so :meth:apply with nothing excluded returns the input to float64 round-off on either backend, PCA residual included.

Source code in pamica/mne_compat/core.py

class AMICAICA:
    """Fit AMICA from MNE objects and interoperate with ``mne.preprocessing.ICA``.

    The wrapper fits pamica's natural-gradient AMICA backend on the data of an
    MNE :class:`~mne.io.Raw` or :class:`~mne.Epochs` and lets MNE consume the
    result: :meth:`get_sources`, :meth:`apply`, :meth:`get_components`,
    :meth:`plot_components` and :meth:`plot_sources` all delegate to a real
    :class:`mne.preprocessing.ICA` built by :meth:`to_mne_ica`.

    For a multi-model fit (``n_models > 1``) each model is exported as its own
    single-model MNE ICA (``to_mne_ica(model_idx=...)`` / the ``model_idx``
    argument on the consumer methods), and the per-sample model dominance --
    which MNE's ``ICA`` cannot represent -- is exposed directly by
    :meth:`get_model_probability` / :meth:`plot_model_probability`.

    Separation-quality metrics (issue #133) are available directly on an MNE
    object: :meth:`mir` (Mutual Information Reduction) and :meth:`pmi` (pairwise
    mutual information between sources). The pamica-specific fitted metadata MNE
    cannot hold -- source-density family, GG shape, component sharing -- is
    inspectable via :meth:`get_pdftype` / :meth:`get_rho` / :meth:`shared_components`.

    The fit runs on the PyTorch backend by default and on the Apple-GPU MLX
    backend with ``backend="mlx"`` (issue #313); see the precision note in the
    Notes below.

    Parameters
    ----------
    n_models : int, default=1
        Number of ICA models to learn (AMICA ``n_models``).
    n_mix : int, default=3
        Number of mixture components per source (AMICA ``n_mix``).
    random_state : int or None, default=None
        Seed for the AMICA fit (passed through as the backend ``seed``) and
        stored on the exported :class:`~mne.preprocessing.ICA`.
    device : str or torch.device, optional
        Torch device for the fit (``None`` = auto; the float64 parity backend
        falls back to CPU when auto-selection lands on MPS). See :class:`AMICA`.
        PyTorch backend only: with ``backend="mlx"`` it must stay ``None``.
    verbose : bool, default=True
        Whether the underlying :class:`AMICA` prints fit progress.
    backend : {"torch", "mlx"}, default="torch"
        Which :class:`AMICA` backend fits the data, forwarded to
        :class:`AMICA`. ``"mlx"`` is the Apple-GPU backend, float32 only; it
        requires MLX (``ImportError`` here otherwise) and takes neither
        ``device`` nor a ``dtype`` fit keyword.

    Raises
    ------
    ValueError
        If ``backend`` is not ``"torch"`` or ``"mlx"``, or ``device`` is set
        with ``backend="mlx"``.
    ImportError
        If ``backend="mlx"`` and MLX is not installed.

    Attributes
    ----------
    amica_ : AMICA
        The fitted pamica model (holds all ``n_models`` models).
    info_ : mne.Info
        The picked measurement info the fit was run on (channel subset only).
    ch_names_ : list of str
        Names of the fitted channels, in order.
    n_components_ : int
        Number of ICA components. Equals the number of fitted channels unless the
        data are rank-deficient, in which case the model is sized to the detected
        numerical rank (issue #223).
    pre_whitener_ : np.ndarray of shape (n_channels, 1)
        Per-channel-type scaling applied before fitting, following MNE's own ICA
        convention (one ``std`` per channel type, applied as ``X / pre_whitener_``).
    pca_components_ : np.ndarray of shape (n_channels, n_channels) or None
        The full orthonormal PCA basis of the fit, in pre-whitened channel
        space, computed once at fit time and exported as the MNE ICA's
        ``pca_components_``. The first ``n_components_`` rows span the
        subspace AMICA modeled. For a rank-reduced fit the remaining rows span
        the PCA residual that reduction discarded, which MNE's ``apply``
        restores by default (issue #322); a full-rank fit has no such rows.
        Shared by every model of a multi-model fit (one sphere per fit).
        ``None`` after a degenerate fit, which exports nothing.
    pca_explained_variance_ : np.ndarray of shape (n_channels,) or None
        Variance of the pre-whitened fit data along each row of
        ``pca_components_``, in descending order: the retained rows carry the
        largest covariance eigenvalues and the residual rows the ones the
        reduction discarded (clipped at 0 against round-off). ``None`` after
        a degenerate fit.
    reject_by_annotation_ : bool
        Whether the last ``Raw`` fit dropped ``bad_*``-annotated samples
        (issue #251). Always ``False`` for an ``Epochs`` fit.
    good_sample_mask_ : np.ndarray of bool or None
        For a ``Raw`` fit, which samples of the recording's timeline the fit
        consumed: ``True`` exactly for the fitted columns, ``False`` for
        samples covered by ``bad_*`` annotations (when rejection is on) and
        for samples outside a ``start``/``stop`` range. ``None`` for an
        ``Epochs`` fit.
    converged_ : bool
        Whether the last fit ended usable (not degenerate). A degenerate fit is
        kept for inspection but refused by the consumer methods (issue #50).
    stop_reason_ : str or None
        Why the backend fit stopped (e.g. ``"max_iter"``, ``"nan_ll"``).

    Notes
    -----
    Model ``h``'s AMICA transform is ``S = W_fort @ (sphere @ (X - mean) - c_h)``,
    where ``c_h`` is that model's data-space center (identically zero for a
    single model, since the ``c`` update is gated to ``n_models > 1``). MNE
    computes sources as ``S = unmixing_matrix_ @ pca_components_[:n_components_]
    @ (X / pre_whitener_ - pca_mean_)``.
    ``X`` is scaled by channel type before fitting, exactly as MNE's own ICA does,
    so the two pipelines agree; AMICA's sphering absorbs a global rescale, so this
    changes nothing for single-channel-type data. Writing the symmetric-ZCA ``sphere`` as
    ``V @ diag(1/sqrt(e)) @ V.T`` with ``V`` orthonormal, the exported ICA for
    model ``h`` uses ``pca_components_ = V.T``,
    ``unmixing_matrix_ = W_fort @ sphere @ V`` and
    ``pca_mean_ = mean + inv(sphere) @ c_h`` (which reduces to ``mean`` when
    ``c_h`` is zero). Keeping ``pca_components_`` orthonormal is what makes MNE's
    ``get_components`` (scalp maps ``inv(sphere) @ inv(W_fort)``) come out right,
    since MNE assumes orthonormal PCA rows. The mapping is pinned by a round-trip
    test (``to_mne_ica(model_idx=h).get_sources(raw)`` equals
    ``amica_.transform(X, model_idx=h)``).

    Rank-deficient input -- Maxwell-filtered MEG, average referencing, channel
    interpolation, or an explicit ``pcakeep``/``pcadb`` -- is supported. The sphere
    is then ``(n_kept, n_channels)`` and has no eigendecomposition, so the
    retained basis is its right singular vectors instead. ``n_components_``
    reports the retained rank.

    The discarded PCA subspace (the residual) is never part of the
    decomposition, and the export keeps it: ``pca_components_`` holds the full
    ``(n_channels, n_channels)`` basis, retained rows first, and MNE's
    ``ICA.apply`` keeps rows past ``n_components_`` as residual PCA components.
    So :meth:`apply` restores the residual by default, as MNE's own ICA does,
    and excluding a component removes only that component (issue #322). **This
    goes beyond the Fortran reference**, whose output has no representation of
    the residual (the written sphere's rows past ``numeigs`` are zero), so any
    back-projection from it drops the residual. Sources, component maps and the
    log-likelihood are unchanged; only reconstruction differs. For the
    reference's rank-reduced reconstruction pass MNE's own
    ``n_pca_components=ica.n_components_`` to :meth:`apply`. See
    ``docs/guides/amica-differences.md`` and ADR 0005.

    Precision: the export always reads the fit through the same float64
    accessors (``AMICA.get_sphere``/``get_mean``/``get_model_center``), so
    both backends go through one code path. A PyTorch fit (float64 by
    default) exports at float64 parity. An MLX fit computes in float32, so its
    unmixing matrix, mean and centers carry float32 rounding (its sphere is
    the float64 one it was built from), and the export is float32-consistent
    rather than float64-parity: :meth:`get_sources` agrees with
    ``amica_.transform`` (which runs in float32) to float32 tolerance, and
    excluding a component changes the data by that component's
    back-projection to the same tolerance. Reconstruction does not depend on
    that precision: MNE's mixing is the float64 pseudo-inverse of the exported
    unmixing and the PCA basis is orthonormal to float64 round-off, so
    :meth:`apply` with nothing excluded returns the input to float64
    round-off on either backend, PCA residual included.
    """

    def __init__(
        self,
        n_models: int = 1,
        n_mix: int = 3,
        random_state: Optional[int] = None,
        device: Optional[Union[str, torch.device]] = None,
        verbose: bool = True,
        backend: str = "torch",
    ):
        # The same checks, errors and messages as AMICA's own constructor.
        _check_backend(backend, device)
        self.n_models = n_models
        self.n_mix = n_mix
        self.random_state = random_state
        self.device = device
        self.verbose = verbose
        self.backend = backend

        self.amica_: Optional[AMICA] = None
        self.info_ = None
        self.ch_names_: Optional[list] = None
        self.n_components_: Optional[int] = None
        self.pre_whitener_: Optional[np.ndarray] = None
        self.pca_components_: Optional[np.ndarray] = None
        self.pca_explained_variance_: Optional[np.ndarray] = None
        self.reject_by_annotation_: bool = False
        self.good_sample_mask_: Optional[np.ndarray] = None
        self.converged_: bool = False
        self.stop_reason_: Optional[str] = None
        self._n_samples: Optional[int] = None
        self._fit_kind: Optional[str] = None
        # Per-model export cache (one mne.preprocessing.ICA per model_idx).
        self._mne_ica_cache: dict = {}

    # ------------------------------------------------------------------
    # Fit
    # ------------------------------------------------------------------
    def fit(
        self,
        inst,
        picks=None,
        start: Optional[int] = None,
        stop: Optional[int] = None,
        reject_by_annotation: bool = True,
        **fit_kwargs,
    ) -> "AMICAICA":
        """Fit AMICA to the data of an MNE ``Raw`` or ``Epochs``.

        Parameters
        ----------
        inst : mne.io.BaseRaw | mne.BaseEpochs
            The data to decompose. ``Epochs`` are concatenated along time
            (``np.hstack``), matching how MNE's own ICA fits epoched data.
        picks : str | list | slice | None, default=None
            Channels to fit, in any form MNE accepts (e.g. ``"eeg"``,
            ``"data"``, a name list). ``None`` selects all good data channels
            (bads excluded), matching MNE's ICA default.
        start, stop : int | None
            Sample range for ``Raw`` input, passed to ``get_data``; ``None`` uses
            the full recording. Not supported for ``Epochs`` (raises).
        reject_by_annotation : bool, default=True
            If ``True``, omit samples covered by annotations whose description
            starts with ``"bad"`` when fitting from :class:`~mne.io.Raw`,
            matching MNE's own ``ICA.fit`` default (issue #251). The samples
            actually fitted are recorded in ``good_sample_mask_``. Has no
            effect for ``Epochs`` (MNE drops bad-annotated epochs at epoching).
        **fit_kwargs
            Forwarded to :meth:`AMICA.fit` (e.g. ``max_iter``, ``lrate``,
            ``do_newton``) and the backend constructor (e.g. ``block_size``).
            ``pcakeep``/``pcadb`` (PCA reduction) are supported: the export
            builds the retained basis from the reduced sphere (issue #225),
            and ``fit`` also computes the basis of the discarded residual
            from the fit data's covariance, so :meth:`apply` restores it by
            default (issue #322; ``n_pca_components=ica.n_components_``
            gives the reference's rank-reduced reconstruction). See
            ``pca_components_``.

        Returns
        -------
        self : AMICAICA

        Raises
        ------
        TypeError
            If ``inst`` is not an MNE ``Raw``/``Epochs``, or (from
            :meth:`AMICA.fit`) a keyword in ``fit_kwargs`` is neither an
            :meth:`AMICA.fit` parameter nor a constructor keyword of the
            selected backend.
        ValueError
            If ``start``/``stop`` are given for ``Epochs``, ``stop`` exceeds the
            recording length, the selected data is non-finite, or no samples
            remain to fit (an empty ``start``/``stop`` range, or ``bad``
            annotations covering the entire selected range); or (from
            :meth:`AMICA.fit`) ``fit_kwargs`` carries ``dtype`` or
            ``device`` with ``backend="mlx"``.
        """
        if not isinstance(inst, (mne.io.BaseRaw, mne.BaseEpochs)):
            raise TypeError(
                "AMICAICA.fit expects an mne.io.Raw or mne.Epochs, got "
                f"{type(inst).__name__}."
            )
        is_raw = isinstance(inst, mne.io.BaseRaw)
        if not is_raw and (start is not None or stop is not None):
            raise ValueError("start/stop are only supported for Raw input, not Epochs.")

        picked = inst.copy().pick("data" if picks is None else picks, exclude="bads")
        ch_names = list(picked.ch_names)

        if is_raw:
            if stop is not None and stop > inst.n_times:
                raise ValueError(
                    f"stop={stop} exceeds the recording length ({inst.n_times} "
                    "samples)."
                )
            range_kw = {}
            if start is not None:
                range_kw["start"] = start
            if stop is not None:
                range_kw["stop"] = stop
            X, good_sample_mask = _raw_data_and_mask(
                picked, reject_by_annotation, **range_kw
            )
            if X.shape[1] == 0:
                # Distinguish an empty start/stop range from annotation
                # rejection, so the error points at the actual cause.
                range_n = (inst.n_times if stop is None else stop) - (
                    0 if start is None else start
                )
                if range_n <= 0:
                    raise ValueError(
                        f"AMICAICA.fit: the start/stop range selects no samples "
                        f"(start={start!r}, stop={stop!r})."
                    )
                raise ValueError(
                    "AMICAICA.fit: no samples left to fit ('bad' annotations "
                    "cover the entire selected start/stop range)."
                )
            fit_kind = "raw"
        else:
            # Concatenate epochs along time, as MNE's ICA does for fitting.
            good_sample_mask = None
            X = np.hstack(picked.get_data())
            fit_kind = "epochs"

        X = np.ascontiguousarray(X, dtype=np.float64)
        if not np.isfinite(X).all():
            bad = [ch_names[i] for i in np.flatnonzero(~np.isfinite(X).all(axis=1))]
            raise ValueError(
                "AMICAICA.fit: input contains non-finite (NaN/Inf) samples in "
                f"channel(s) {bad}; clean bad segments/interpolation before "
                "fitting."
            )

        # Channel-type scaling, matching MNE's own ICA._compute_pre_whitener:
        # one std per channel type, applied as data / pre_whitener. Required for
        # mixed magnetometer/gradiometer data, whose physical units differ by
        # orders of magnitude (issue #225). AMICA's sphering absorbs a *global*
        # rescale exactly, so for single-channel-type data (ordinary EEG) this
        # leaves the sources unchanged; it bites only across types, where it
        # decides which directions survive rank reduction.
        pre_whitener = _compute_pre_whitener(X, picked.info)
        X = X / pre_whitener

        if isinstance(self.random_state, (int, np.integer)):
            fit_kwargs.setdefault("seed", int(self.random_state))

        amica = AMICA(
            n_models=self.n_models,
            n_mix=self.n_mix,
            device=self.device,
            verbose=self.verbose,
            backend=self.backend,
        )
        amica.fit(X, **fit_kwargs)
        backend = amica.model_
        if backend is None:
            raise RuntimeError("AMICAICA.fit: AMICA.fit returned without a model.")

        # The full PCA basis (retained rows, then the residual a rank-reduced
        # fit discarded), computed once here from the fitted sphere and the fit
        # data so nothing is recomputed from data at apply time (issue #322).
        # The sphere does not depend on the model, so one basis serves every
        # per-model export. Read through the backend-agnostic float64 accessor
        # (issue #313), which refuses a degenerate fit: such a fit is kept for
        # inspection but never exported (_check_fitted refuses it), so it gets
        # no basis.
        pca_components: Optional[np.ndarray] = None
        pca_explained_variance: Optional[np.ndarray] = None
        if amica.converged_:
            pca_components, pca_explained_variance = _pca_basis(amica.get_sphere(), X)

        # Publish to self only after every fallible step above succeeds, so a
        # failed (re)fit leaves the previously fitted state intact rather than a
        # mix of the old model and the new attempt's metadata (mirrors the
        # local-first pattern in AMICA.fit).
        self.info_ = picked.info
        self.ch_names_ = ch_names
        self.pre_whitener_ = pre_whitener
        self.pca_components_ = pca_components
        self.pca_explained_variance_ = pca_explained_variance
        self.good_sample_mask_ = good_sample_mask
        self.reject_by_annotation_ = bool(reject_by_annotation) if is_raw else False
        # The fitted model dimension, which is the input channel count unless
        # rank reduction shrank it (issue #223). MNE represents this natively:
        # the first n_components_ rows of pca_components_ are the ICA subspace
        # and any further rows are residual PCA components (issue #322).
        self.n_components_ = backend.n_channels
        self._n_samples = X.shape[1]
        self._fit_kind = fit_kind
        self.amica_ = amica
        self.converged_ = amica.converged_
        self.stop_reason_ = amica.stop_reason_
        self._mne_ica_cache = {}  # invalidate any cached exports
        return self

    # ------------------------------------------------------------------
    # Export
    # ------------------------------------------------------------------
    def to_mne_ica(self, model_idx: int = 0) -> _MNEICA:
        """Build (and cache) a fully-populated :class:`mne.preprocessing.ICA`.

        The returned object is a genuine MNE ICA: ``get_sources``, ``apply``,
        ``get_components``, ``save`` and the ``mne.viz`` plotters operate on it
        natively. See the class :class:`Notes <AMICAICA>` for the
        mean/sphere/unmixing to ``pca_mean_``/``pca_components_``/
        ``unmixing_matrix_`` mapping.

        For a multi-model fit each model is exported as its own single-model
        MNE ICA (MNE has no multi-model concept); pass ``model_idx`` to pick
        one. The per-model exports are cached and returned by reference:
        mutating one (for example setting ``.exclude``) persists across
        subsequent :meth:`apply`/:meth:`get_sources` calls for that model until
        the next :meth:`fit`.

        The export carries the full ``(n_channels, n_channels)`` PCA basis
        computed at fit time (see ``pca_components_``), with ``n_components_``
        the retained rank and ``n_pca_components`` left at ``None``. For a
        rank-reduced fit the rows past ``n_components_`` span the discarded
        PCA residual, which MNE keeps as residual PCA components, so the
        exported ICA's ``apply`` restores it by default (issue #322). This goes
        beyond the Fortran reference, whose output has no representation of
        the residual; ``apply(..., n_pca_components=ica.n_components_)`` gives
        the reference's rank-reduced reconstruction. When a residual is
        exported, its dimension and that opt-out are logged at INFO.

        Parameters
        ----------
        model_idx : int, default=0
            Which AMICA model to export (``0..n_models-1``).

        Returns
        -------
        ica : mne.preprocessing.ICA
        """
        self._check_fitted("build the MNE ICA")
        model_idx = self._check_model_idx(model_idx)
        if model_idx in self._mne_ica_cache:
            return self._mne_ica_cache[model_idx]

        amica = self.amica_
        if (
            amica is None
            or amica.model_ is None
            or self.ch_names_ is None
            or self.pca_components_ is None
            or self.pca_explained_variance_ is None
        ):
            raise RuntimeError(
                "AMICAICA: internal state is inconsistent; refit before to_mne_ica()."
            )
        # The backend-agnostic float64 accessors (issue #313), identical for
        # the PyTorch and MLX backends; see the class Notes on precision.
        mean = amica.get_mean()
        sphere = amica.get_sphere()
        w_fort = amica.get_unmixing_matrix(model_idx=model_idx)
        c = amica.get_model_center(model_idx=model_idx)  # sphered space
        n_ch = sphere.shape[0]

        # The full PCA basis was computed once at fit time (see _pca_basis): the
        # first n_ch rows span the ICA subspace, and any rows after them are the
        # residual a rank-reduced fit discarded. MNE keeps those as residual PCA
        # components in apply (n_pca_components=None keeps every row), which is
        # what restores the residual (issue #322).
        pca_components = self.pca_components_.copy()
        cov_evals = self.pca_explained_variance_.copy()
        n_residual = pca_components.shape[0] - n_ch
        v = self.pca_components_[:n_ch].T  # retained basis, as columns
        unmixing = w_fort @ sphere @ v
        # Fold the per-model center c (in sphered space) into pca_mean via the
        # data-space offset pinv(sphere) @ c, so MNE's (X - pca_mean) reproduces
        # AMICA's W(sphere(X - mean) - c). c is identically zero for a single
        # model, leaving pca_mean == mean bit-for-bit. pinv rather than solve: the
        # sphere is non-square under rank reduction, and for a square sphere the
        # two agree to round-off. pinv(sphere) @ c lies in the sphere's row
        # space, so the fold never reaches the residual rows.
        pca_mean = mean + np.linalg.pinv(sphere) @ c if np.any(c) else mean

        ica = _MNEICA(
            n_components=n_ch,
            method="infomax",
            max_iter="auto",
            random_state=self.random_state,
        )
        # `method` is inert here: the fitted attributes are hand-populated and
        # ICA.fit (the only place `method` branches) is never called, but
        # ICA.__init__ still requires a valid method name.
        ica.info = self.info_
        ica.ch_names = list(self.ch_names_)
        ica.n_components_ = n_ch
        ica.pca_mean_ = pca_mean
        ica.pca_components_ = pca_components
        ica.pca_explained_variance_ = cov_evals
        ica.unmixing_matrix_ = unmixing
        ica.pre_whitener_ = self.pre_whitener_
        # The number of iterations the fit ran: ``iteration`` is the 0-based
        # index of the last one whose E-step ran (issue #339), which is also the
        # reference's own ``iter`` count at a stop.
        ica.n_iter_ = int(getattr(amica.model_, "iteration", 0)) + 1
        # MNE's own fit sets these; read_ica_eeglab (the precedent for building
        # an ICA from an external decomposition) sets reject_=None. Without them
        # ICA.save()/plot_properties raise AttributeError.
        ica.reject_ = None
        ica.n_samples_ = int(self._n_samples) if self._n_samples is not None else 0
        ica._update_mixing_matrix()
        ica._update_ica_names()
        ica.current_fit = self._fit_kind

        if n_residual:
            logger.info(
                "AMICAICA.to_mne_ica: exporting %d ICA components plus a "
                "%d-dimensional PCA residual (the subspace rank reduction "
                "discarded before fitting). ICA.apply restores the residual by "
                "default, which goes beyond the Fortran reference; pass "
                "n_pca_components=%d (ica.n_components_) to apply for the "
                "reference's rank-reduced reconstruction (issue #322).",
                n_ch,
                n_residual,
                n_ch,
            )

        self._mne_ica_cache[model_idx] = ica
        return ica

    # ------------------------------------------------------------------
    # MNE consumer surface (delegates to the exported ICA)
    # ------------------------------------------------------------------
    def get_sources(self, inst, *args, model_idx: int = 0, **kwargs):
        """Sources for ``inst`` from model ``model_idx`` (see ``ICA.get_sources``).

        ``model_idx`` is keyword-only so positional arguments pass straight
        through to MNE's ``ICA.get_sources`` (e.g. ``add_channels``,
        ``start``/``stop``).
        """
        return self.to_mne_ica(model_idx).get_sources(inst, *args, **kwargs)

    def apply(self, inst, *args, model_idx: int = 0, **kwargs):
        """Remove selected components of model ``model_idx`` and back-project.

        ``model_idx`` is keyword-only so positional arguments pass straight
        through to MNE's ``ICA.apply``. Pass ``exclude=[...]`` (or set it on the
        exported ICA) to drop components; with no exclusions this reconstructs
        the input.

        For a rank-reduced fit (``pcakeep``/``pcadb`` or automatic rank
        detection) the PCA residual, the subspace the reduction discarded and
        AMICA never modeled, is **restored by default**: excluding a component
        removes that component's back-projection and nothing else, as with
        MNE's own ICA (issue #322). This differs from the Fortran reference,
        whose back-projection reconstructs only the retained subspace. To get
        that rank-reduced reconstruction, pass
        ``n_pca_components=ica.n_components_``; a float ``n_pca_components``
        keeps residual rows by cumulative explained variance, following MNE.
        A full-rank fit has no residual, so the option changes nothing there.
        """
        return self.to_mne_ica(model_idx).apply(inst, *args, **kwargs)

    def get_components(self, *, model_idx: int = 0) -> np.ndarray:
        """Scalp maps of model ``model_idx``, ``(n_channels, n_components)``."""
        return self.to_mne_ica(model_idx).get_components()

    def plot_components(self, *args, model_idx: int = 0, **kwargs):
        """Plot model ``model_idx`` component topographies (see ``ICA.plot_components``)."""
        return self.to_mne_ica(model_idx).plot_components(*args, **kwargs)

    def plot_sources(self, inst, *args, model_idx: int = 0, **kwargs):
        """Plot model ``model_idx`` component time courses (see ``ICA.plot_sources``)."""
        return self.to_mne_ica(model_idx).plot_sources(inst, *args, **kwargs)

    # ------------------------------------------------------------------
    # Multi-model dominance (issue #141)
    # ------------------------------------------------------------------
    def get_model_probability(
        self, inst, *, reject_by_annotation: bool = True
    ) -> np.ndarray:
        """Per-sample posterior probability of each model on ``inst`` (dominance).

        Returns ``P(model | sample)`` as ``(n_models, n_samples)`` via
        :meth:`AMICA.model_probability` on ``inst``'s data (``Epochs`` are
        concatenated along time). Each column sums to 1; all ones for a single
        model. MNE's own ``ICA`` has no multi-model concept, so this is exposed
        here rather than through the exported per-model ICA objects.

        Parameters
        ----------
        inst : mne.io.BaseRaw | mne.BaseEpochs
            Data to score.
        reject_by_annotation : bool, default=True
            For ``Raw`` input, evaluate only samples not covered by ``inst``'s
            ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
            default (issue #251). The output keeps the full ``n_times``
            timeline -- rejected columns are ``NaN`` -- so it stays aligned
            with ``inst`` without any manual re-indexing. ``False`` scores
            every sample (the model evaluates fine on artifact segments; their
            probabilities are simply dominated by whatever model claims the
            artifact).
        """
        self._check_fitted("compute the model probability")
        if self.amica_ is None:
            raise RuntimeError(
                "AMICAICA: internal state is inconsistent; refit before scoring."
            )
        X, good = self._data_for(inst, reject_by_annotation=reject_by_annotation)
        return _expand_to_timeline(self.amica_.model_probability(X), good)

    def plot_model_probability(
        self,
        inst,
        *,
        srate: Optional[float] = None,
        reject_by_annotation: bool = True,
        **kwargs,
    ):
        """Plot per-model probability + best-model log-likelihood over ``inst``.

        Delegates to :func:`pamica.viz.plot_model_probability` with the live
        per-model log-likelihood (:meth:`AMICA.model_loglik`) on ``inst``'s
        data. ``srate`` defaults to the fitted recording's sampling rate, so the
        x-axis is in seconds; extra keywords (``smooth_sec``, ``window_sec``,
        ``axes``) pass through.

        With ``reject_by_annotation`` (default, ``Raw`` only), samples covered
        by ``inst``'s ``bad_*`` annotations plot as gaps: the log-likelihood is
        evaluated on the good samples only and rejected columns are ``NaN``, so
        the time axis stays aligned with ``inst`` (issue #251). Note
        ``smooth_sec`` widens the gaps by the smoothing window, since the
        Hanning smoothing propagates ``NaN`` across its support.
        """
        from ..viz import plot_model_probability as _plot_model_probability

        self._check_fitted("plot the model probability")
        if self.amica_ is None or self.info_ is None:
            raise RuntimeError(
                "AMICAICA: internal state is inconsistent; refit before plotting."
            )
        X, good = self._data_for(inst, reject_by_annotation=reject_by_annotation)
        lht = _expand_to_timeline(self.amica_.model_loglik(X), good)
        if srate is None:
            srate = float(self.info_["sfreq"])
        return _plot_model_probability(lht=lht, srate=srate, **kwargs)

    # ------------------------------------------------------------------
    # Separation-quality metrics (issue #143, on top of #133)
    # ------------------------------------------------------------------
    def mir(
        self,
        inst,
        *,
        model_idx: int = 0,
        nbins: Optional[int] = None,
        reject_by_annotation: bool = True,
    ) -> tuple:
        """Mutual Information Reduction of model ``model_idx`` on ``inst``.

        How much mutual information the fitted unmixing removes from the data,
        in nats (issue #133). Delegates to :meth:`AMICA.mir` on ``inst``'s
        fitted-channel data; MIR is shift-invariant, so mean/``c`` centering is
        irrelevant.

        Parameters
        ----------
        inst : mne.io.BaseRaw | mne.BaseEpochs
            Data to score (``Epochs`` concatenated along time).
        model_idx : int, default=0
            Which model's unmixing to use.
        nbins : int, optional
            Histogram bin count; see :func:`pamica.metrics.mir`.
        reject_by_annotation : bool, default=True
            For ``Raw`` input, score only samples not covered by ``inst``'s
            ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
            default (issue #251).

        Returns
        -------
        mir_nats : float
        variance : float
        """
        self._check_fitted("compute MIR")
        model_idx = self._check_model_idx(model_idx)
        assert self.amica_ is not None
        X, _ = self._data_for(inst, reject_by_annotation=reject_by_annotation)
        return self.amica_.mir(X, model_idx=model_idx, nbins=nbins)

    def pmi(
        self,
        inst,
        *,
        model_idx: int = 0,
        nbins: Optional[int] = None,
        reject_by_annotation: bool = True,
    ) -> np.ndarray:
        """Pairwise Mutual Information between model ``model_idx``'s sources on ``inst``.

        The residual pairwise dependence between fitted sources, in nats
        (issue #133). Delegates to :meth:`AMICA.pmi` on ``inst``'s fitted-channel
        data.

        Parameters
        ----------
        inst : mne.io.BaseRaw | mne.BaseEpochs
            Data to score (``Epochs`` concatenated along time).
        model_idx : int, default=0
            Which model's sources to use.
        nbins : int, optional
            Histogram bin count; see :func:`pamica.metrics.pairwise_mi`.
        reject_by_annotation : bool, default=True
            For ``Raw`` input, score only samples not covered by ``inst``'s
            ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
            default (issue #251).

        Returns
        -------
        mi_matrix : np.ndarray of shape (n_components, n_components)
            Symmetric; the diagonal is each source's own entropy.
        """
        self._check_fitted("compute PMI")
        model_idx = self._check_model_idx(model_idx)
        assert self.amica_ is not None
        X, _ = self._data_for(inst, reject_by_annotation=reject_by_annotation)
        return self.amica_.pmi(X, model_idx=model_idx, nbins=nbins)

    # ------------------------------------------------------------------
    # pamica-specific metadata (issue #142)
    #
    # MNE's ``ICA`` carries no source-density family, GG shape, or
    # component-sharing state, so these are exposed here rather than silently
    # dropped by the ``mne.preprocessing.ICA`` export.
    # ------------------------------------------------------------------
    def get_pdftype(self, *, model_idx: int = 0) -> np.ndarray:
        """Per-component source-density family code for model ``model_idx``.

        One integer per ICA component (0-4); map to names with
        :data:`pamica.mne_compat.PDFTYPE_NAMES`. All components share one family
        unless the adaptive switcher (``pdftype=1``) moved them (issue #26).
        """
        self._check_fitted("get the density family")
        model_idx = self._check_model_idx(model_idx)
        assert self.amica_ is not None
        return self.amica_.get_pdftype(model_idx=model_idx)

    def get_rho(self, *, model_idx: int = 0) -> np.ndarray:
        """Generalized-Gaussian shape ``rho`` for model ``model_idx``.

        Shape ``(n_mix, n_components)``; ``rho == 2`` is Gaussian-shaped,
        ``rho == 1`` Laplacian, ``rho < 1`` heavier-tailed. Meaningful only for
        the generalized-Gaussian family (``pdftype=0``).
        """
        self._check_fitted("get rho")
        model_idx = self._check_model_idx(model_idx)
        assert self.amica_ is not None
        return self.amica_.get_rho(model_idx=model_idx)

    def shared_components(self) -> list:
        """Components shared across models by ``share_comps`` (issue #60).

        One group of ``(model_idx, component_idx)`` pairs per shared component;
        empty when nothing is shared (always so for a single model or a default
        multi-model fit with ``share_comps`` off).
        """
        self._check_fitted("get the shared components")
        assert self.amica_ is not None
        return self.amica_.shared_components()

    # ------------------------------------------------------------------
    # Helpers
    # ------------------------------------------------------------------
    def _data_for(self, inst, *, reject_by_annotation: bool = False):
        """The fitted-channel data of ``inst`` as ``(n_channels, n_samples)``.

        Selects the exact channels the fit used (by name) so the array aligns
        with the stored sphere/unmixing; ``Epochs`` are concatenated along time.
        The channel-type pre-whitener is applied, because the backend was fitted
        on scaled data and would otherwise be handed a differently-scaled array
        (issue #225).

        With ``reject_by_annotation`` (``Raw`` only), samples covered by
        ``inst``'s own ``bad_*`` annotations are omitted -- the evaluation-time
        analogue of MNE's ``ICA.score_sources`` behavior (issue #251). The
        rejection follows the *passed* instance's annotations, not the fit-time
        mask, so it stays correct for any Raw handed in.

        Returns
        -------
        X : np.ndarray of shape (n_channels, n_kept)
        good_mask : np.ndarray of bool, shape (inst.n_times,), or None
            ``True`` for the timeline positions the columns of ``X`` occupy;
            ``None`` when nothing was rejected (``Epochs``, or rejection off).
        """
        if not isinstance(inst, (mne.io.BaseRaw, mne.BaseEpochs)):
            raise TypeError(
                f"expected an mne.io.Raw or mne.Epochs, got {type(inst).__name__}."
            )
        picked = inst.copy().pick(self.ch_names_)
        good_mask = None
        if isinstance(inst, mne.io.BaseRaw):
            if reject_by_annotation:
                X, good_mask = _raw_data_and_mask(picked, True)
                if X.shape[1] == 0:
                    # One consistent hard error for all four scoring methods
                    # (get_model_probability would otherwise return silent
                    # all-NaN and mir/pmi an opaque zero-size reduction error).
                    raise ValueError(
                        "AMICAICA: every sample of the provided Raw is covered "
                        "by 'bad' annotations, so there is nothing to score. "
                        "Pass reject_by_annotation=False to score annotated "
                        "segments anyway."
                    )
            else:
                X = picked.get_data()
        else:
            X = np.hstack(picked.get_data())
        X = np.ascontiguousarray(X, dtype=np.float64)
        if self.pre_whitener_ is not None:
            X = X / self.pre_whitener_
        return X, good_mask

    def _check_model_idx(self, model_idx: int) -> int:
        if not isinstance(model_idx, (int, np.integer)):
            raise TypeError(
                f"model_idx must be an int, got {type(model_idx).__name__}."
            )
        # Bound against the fitted backend's model count (the source of truth),
        # not the mutable constructor hyperparameter self.n_models.
        n = (
            self.amica_.model_.n_models
            if self.amica_ is not None and self.amica_.model_ is not None
            else self.n_models
        )
        if not (0 <= model_idx < n):
            raise ValueError(
                f"model_idx={model_idx} out of range for a {n}-model fit "
                f"(valid: 0..{n - 1})."
            )
        return int(model_idx)

    def _check_fitted(self, action: str = "this call") -> None:
        """Raise if no usable model is available.

        ``ValueError`` when never fitted (matching :class:`AMICA`'s convention);
        ``RuntimeError`` when the fit ended degenerate (non-finite parameters,
        issue #50), so the failure surfaces here rather than as opaque NaNs
        downstream.
        """
        if self.amica_ is None:
            raise ValueError(f"AMICAICA must be fitted before {action}; run fit().")
        if not self.converged_:
            raise RuntimeError(
                f"Refusing to {action}: the AMICA fit ended degenerate "
                f"(stop_reason={self.stop_reason_!r}), so it holds non-finite "
                "parameters and would produce NaN output. Lower lrate, disable "
                "Newton, or check data conditioning, then refit."
            )

    def __repr__(self) -> str:
        config = (
            f"backend={self.backend!r}, n_models={self.n_models}, n_mix={self.n_mix}"
        )
        if self.amica_ is None:
            return f"<AMICAICA (unfitted, {config})>"
        if not self.converged_:
            return (
                f"<AMICAICA (degenerate fit, stop_reason={self.stop_reason_!r}, "
                f"{config}, {self._fit_kind})>"
            )
        return (
            f"<AMICAICA (fitted: {self.n_components_} components, {config}, "
            f"{self._fit_kind})>"
        )

fit(inst, picks=None, start=None, stop=None, reject_by_annotation=True, **fit_kwargs)

Fit AMICA to the data of an MNE Raw or Epochs.

Parameters:

Name	Type	Description	Default
inst	BaseRaw | BaseEpochs	The data to decompose. Epochs are concatenated along time (np.hstack), matching how MNE's own ICA fits epoched data.	required
picks	str | list | slice | None	Channels to fit, in any form MNE accepts (e.g. "eeg", "data", a name list). None selects all good data channels (bads excluded), matching MNE's ICA default.	None
start	int | None	Sample range for Raw input, passed to get_data; None uses the full recording. Not supported for Epochs (raises).	None
stop	int | None	Sample range for Raw input, passed to get_data; None uses the full recording. Not supported for Epochs (raises).	None
reject_by_annotation	bool	If True, omit samples covered by annotations whose description starts with "bad" when fitting from :class:~mne.io.Raw, matching MNE's own ICA.fit default (issue #251). The samples actually fitted are recorded in good_sample_mask_. Has no effect for Epochs (MNE drops bad-annotated epochs at epoching).	True
**fit_kwargs		Forwarded to :meth:AMICA.fit (e.g. max_iter, lrate, do_newton) and the backend constructor (e.g. block_size). pcakeep/pcadb (PCA reduction) are supported: the export builds the retained basis from the reduced sphere (issue #225), and fit also computes the basis of the discarded residual from the fit data's covariance, so :meth:apply restores it by default (issue #322; n_pca_components=ica.n_components_ gives the reference's rank-reduced reconstruction). See pca_components_.	{}

Returns:

Name	Type	Description
self	AMICAICA

Raises:

Type	Description
TypeError	If inst is not an MNE Raw/Epochs, or (from :meth:AMICA.fit) a keyword in fit_kwargs is neither an :meth:AMICA.fit parameter nor a constructor keyword of the selected backend.
ValueError	If start/stop are given for Epochs, stop exceeds the recording length, the selected data is non-finite, or no samples remain to fit (an empty start/stop range, or bad annotations covering the entire selected range); or (from :meth:AMICA.fit) fit_kwargs carries dtype or device with backend="mlx".

Source code in pamica/mne_compat/core.py

def fit(
    self,
    inst,
    picks=None,
    start: Optional[int] = None,
    stop: Optional[int] = None,
    reject_by_annotation: bool = True,
    **fit_kwargs,
) -> "AMICAICA":
    """Fit AMICA to the data of an MNE ``Raw`` or ``Epochs``.

    Parameters
    ----------
    inst : mne.io.BaseRaw | mne.BaseEpochs
        The data to decompose. ``Epochs`` are concatenated along time
        (``np.hstack``), matching how MNE's own ICA fits epoched data.
    picks : str | list | slice | None, default=None
        Channels to fit, in any form MNE accepts (e.g. ``"eeg"``,
        ``"data"``, a name list). ``None`` selects all good data channels
        (bads excluded), matching MNE's ICA default.
    start, stop : int | None
        Sample range for ``Raw`` input, passed to ``get_data``; ``None`` uses
        the full recording. Not supported for ``Epochs`` (raises).
    reject_by_annotation : bool, default=True
        If ``True``, omit samples covered by annotations whose description
        starts with ``"bad"`` when fitting from :class:`~mne.io.Raw`,
        matching MNE's own ``ICA.fit`` default (issue #251). The samples
        actually fitted are recorded in ``good_sample_mask_``. Has no
        effect for ``Epochs`` (MNE drops bad-annotated epochs at epoching).
    **fit_kwargs
        Forwarded to :meth:`AMICA.fit` (e.g. ``max_iter``, ``lrate``,
        ``do_newton``) and the backend constructor (e.g. ``block_size``).
        ``pcakeep``/``pcadb`` (PCA reduction) are supported: the export
        builds the retained basis from the reduced sphere (issue #225),
        and ``fit`` also computes the basis of the discarded residual
        from the fit data's covariance, so :meth:`apply` restores it by
        default (issue #322; ``n_pca_components=ica.n_components_``
        gives the reference's rank-reduced reconstruction). See
        ``pca_components_``.

    Returns
    -------
    self : AMICAICA

    Raises
    ------
    TypeError
        If ``inst`` is not an MNE ``Raw``/``Epochs``, or (from
        :meth:`AMICA.fit`) a keyword in ``fit_kwargs`` is neither an
        :meth:`AMICA.fit` parameter nor a constructor keyword of the
        selected backend.
    ValueError
        If ``start``/``stop`` are given for ``Epochs``, ``stop`` exceeds the
        recording length, the selected data is non-finite, or no samples
        remain to fit (an empty ``start``/``stop`` range, or ``bad``
        annotations covering the entire selected range); or (from
        :meth:`AMICA.fit`) ``fit_kwargs`` carries ``dtype`` or
        ``device`` with ``backend="mlx"``.
    """
    if not isinstance(inst, (mne.io.BaseRaw, mne.BaseEpochs)):
        raise TypeError(
            "AMICAICA.fit expects an mne.io.Raw or mne.Epochs, got "
            f"{type(inst).__name__}."
        )
    is_raw = isinstance(inst, mne.io.BaseRaw)
    if not is_raw and (start is not None or stop is not None):
        raise ValueError("start/stop are only supported for Raw input, not Epochs.")

    picked = inst.copy().pick("data" if picks is None else picks, exclude="bads")
    ch_names = list(picked.ch_names)

    if is_raw:
        if stop is not None and stop > inst.n_times:
            raise ValueError(
                f"stop={stop} exceeds the recording length ({inst.n_times} "
                "samples)."
            )
        range_kw = {}
        if start is not None:
            range_kw["start"] = start
        if stop is not None:
            range_kw["stop"] = stop
        X, good_sample_mask = _raw_data_and_mask(
            picked, reject_by_annotation, **range_kw
        )
        if X.shape[1] == 0:
            # Distinguish an empty start/stop range from annotation
            # rejection, so the error points at the actual cause.
            range_n = (inst.n_times if stop is None else stop) - (
                0 if start is None else start
            )
            if range_n <= 0:
                raise ValueError(
                    f"AMICAICA.fit: the start/stop range selects no samples "
                    f"(start={start!r}, stop={stop!r})."
                )
            raise ValueError(
                "AMICAICA.fit: no samples left to fit ('bad' annotations "
                "cover the entire selected start/stop range)."
            )
        fit_kind = "raw"
    else:
        # Concatenate epochs along time, as MNE's ICA does for fitting.
        good_sample_mask = None
        X = np.hstack(picked.get_data())
        fit_kind = "epochs"

    X = np.ascontiguousarray(X, dtype=np.float64)
    if not np.isfinite(X).all():
        bad = [ch_names[i] for i in np.flatnonzero(~np.isfinite(X).all(axis=1))]
        raise ValueError(
            "AMICAICA.fit: input contains non-finite (NaN/Inf) samples in "
            f"channel(s) {bad}; clean bad segments/interpolation before "
            "fitting."
        )

    # Channel-type scaling, matching MNE's own ICA._compute_pre_whitener:
    # one std per channel type, applied as data / pre_whitener. Required for
    # mixed magnetometer/gradiometer data, whose physical units differ by
    # orders of magnitude (issue #225). AMICA's sphering absorbs a *global*
    # rescale exactly, so for single-channel-type data (ordinary EEG) this
    # leaves the sources unchanged; it bites only across types, where it
    # decides which directions survive rank reduction.
    pre_whitener = _compute_pre_whitener(X, picked.info)
    X = X / pre_whitener

    if isinstance(self.random_state, (int, np.integer)):
        fit_kwargs.setdefault("seed", int(self.random_state))

    amica = AMICA(
        n_models=self.n_models,
        n_mix=self.n_mix,
        device=self.device,
        verbose=self.verbose,
        backend=self.backend,
    )
    amica.fit(X, **fit_kwargs)
    backend = amica.model_
    if backend is None:
        raise RuntimeError("AMICAICA.fit: AMICA.fit returned without a model.")

    # The full PCA basis (retained rows, then the residual a rank-reduced
    # fit discarded), computed once here from the fitted sphere and the fit
    # data so nothing is recomputed from data at apply time (issue #322).
    # The sphere does not depend on the model, so one basis serves every
    # per-model export. Read through the backend-agnostic float64 accessor
    # (issue #313), which refuses a degenerate fit: such a fit is kept for
    # inspection but never exported (_check_fitted refuses it), so it gets
    # no basis.
    pca_components: Optional[np.ndarray] = None
    pca_explained_variance: Optional[np.ndarray] = None
    if amica.converged_:
        pca_components, pca_explained_variance = _pca_basis(amica.get_sphere(), X)

    # Publish to self only after every fallible step above succeeds, so a
    # failed (re)fit leaves the previously fitted state intact rather than a
    # mix of the old model and the new attempt's metadata (mirrors the
    # local-first pattern in AMICA.fit).
    self.info_ = picked.info
    self.ch_names_ = ch_names
    self.pre_whitener_ = pre_whitener
    self.pca_components_ = pca_components
    self.pca_explained_variance_ = pca_explained_variance
    self.good_sample_mask_ = good_sample_mask
    self.reject_by_annotation_ = bool(reject_by_annotation) if is_raw else False
    # The fitted model dimension, which is the input channel count unless
    # rank reduction shrank it (issue #223). MNE represents this natively:
    # the first n_components_ rows of pca_components_ are the ICA subspace
    # and any further rows are residual PCA components (issue #322).
    self.n_components_ = backend.n_channels
    self._n_samples = X.shape[1]
    self._fit_kind = fit_kind
    self.amica_ = amica
    self.converged_ = amica.converged_
    self.stop_reason_ = amica.stop_reason_
    self._mne_ica_cache = {}  # invalidate any cached exports
    return self

to_mne_ica(model_idx=0)

Build (and cache) a fully-populated :class:mne.preprocessing.ICA.

The returned object is a genuine MNE ICA: get_sources, apply, get_components, save and the mne.viz plotters operate on it natively. See the class :class:Notes <AMICAICA> for the mean/sphere/unmixing to pca_mean_/pca_components_/ unmixing_matrix_ mapping.

For a multi-model fit each model is exported as its own single-model MNE ICA (MNE has no multi-model concept); pass model_idx to pick one. The per-model exports are cached and returned by reference: mutating one (for example setting .exclude) persists across subsequent :meth:apply/:meth:get_sources calls for that model until the next :meth:fit.

The export carries the full (n_channels, n_channels) PCA basis computed at fit time (see pca_components_), with n_components_ the retained rank and n_pca_components left at None. For a rank-reduced fit the rows past n_components_ span the discarded PCA residual, which MNE keeps as residual PCA components, so the exported ICA's apply restores it by default (issue #322). This goes beyond the Fortran reference, whose output has no representation of the residual; apply(..., n_pca_components=ica.n_components_) gives the reference's rank-reduced reconstruction. When a residual is exported, its dimension and that opt-out are logged at INFO.

Parameters:

Name	Type	Description	Default
model_idx	int	Which AMICA model to export (0..n_models-1).	0

Returns:

Name	Type	Description
ica	ICA

Source code in pamica/mne_compat/core.py

def to_mne_ica(self, model_idx: int = 0) -> _MNEICA:
    """Build (and cache) a fully-populated :class:`mne.preprocessing.ICA`.

    The returned object is a genuine MNE ICA: ``get_sources``, ``apply``,
    ``get_components``, ``save`` and the ``mne.viz`` plotters operate on it
    natively. See the class :class:`Notes <AMICAICA>` for the
    mean/sphere/unmixing to ``pca_mean_``/``pca_components_``/
    ``unmixing_matrix_`` mapping.

    For a multi-model fit each model is exported as its own single-model
    MNE ICA (MNE has no multi-model concept); pass ``model_idx`` to pick
    one. The per-model exports are cached and returned by reference:
    mutating one (for example setting ``.exclude``) persists across
    subsequent :meth:`apply`/:meth:`get_sources` calls for that model until
    the next :meth:`fit`.

    The export carries the full ``(n_channels, n_channels)`` PCA basis
    computed at fit time (see ``pca_components_``), with ``n_components_``
    the retained rank and ``n_pca_components`` left at ``None``. For a
    rank-reduced fit the rows past ``n_components_`` span the discarded
    PCA residual, which MNE keeps as residual PCA components, so the
    exported ICA's ``apply`` restores it by default (issue #322). This goes
    beyond the Fortran reference, whose output has no representation of
    the residual; ``apply(..., n_pca_components=ica.n_components_)`` gives
    the reference's rank-reduced reconstruction. When a residual is
    exported, its dimension and that opt-out are logged at INFO.

    Parameters
    ----------
    model_idx : int, default=0
        Which AMICA model to export (``0..n_models-1``).

    Returns
    -------
    ica : mne.preprocessing.ICA
    """
    self._check_fitted("build the MNE ICA")
    model_idx = self._check_model_idx(model_idx)
    if model_idx in self._mne_ica_cache:
        return self._mne_ica_cache[model_idx]

    amica = self.amica_
    if (
        amica is None
        or amica.model_ is None
        or self.ch_names_ is None
        or self.pca_components_ is None
        or self.pca_explained_variance_ is None
    ):
        raise RuntimeError(
            "AMICAICA: internal state is inconsistent; refit before to_mne_ica()."
        )
    # The backend-agnostic float64 accessors (issue #313), identical for
    # the PyTorch and MLX backends; see the class Notes on precision.
    mean = amica.get_mean()
    sphere = amica.get_sphere()
    w_fort = amica.get_unmixing_matrix(model_idx=model_idx)
    c = amica.get_model_center(model_idx=model_idx)  # sphered space
    n_ch = sphere.shape[0]

    # The full PCA basis was computed once at fit time (see _pca_basis): the
    # first n_ch rows span the ICA subspace, and any rows after them are the
    # residual a rank-reduced fit discarded. MNE keeps those as residual PCA
    # components in apply (n_pca_components=None keeps every row), which is
    # what restores the residual (issue #322).
    pca_components = self.pca_components_.copy()
    cov_evals = self.pca_explained_variance_.copy()
    n_residual = pca_components.shape[0] - n_ch
    v = self.pca_components_[:n_ch].T  # retained basis, as columns
    unmixing = w_fort @ sphere @ v
    # Fold the per-model center c (in sphered space) into pca_mean via the
    # data-space offset pinv(sphere) @ c, so MNE's (X - pca_mean) reproduces
    # AMICA's W(sphere(X - mean) - c). c is identically zero for a single
    # model, leaving pca_mean == mean bit-for-bit. pinv rather than solve: the
    # sphere is non-square under rank reduction, and for a square sphere the
    # two agree to round-off. pinv(sphere) @ c lies in the sphere's row
    # space, so the fold never reaches the residual rows.
    pca_mean = mean + np.linalg.pinv(sphere) @ c if np.any(c) else mean

    ica = _MNEICA(
        n_components=n_ch,
        method="infomax",
        max_iter="auto",
        random_state=self.random_state,
    )
    # `method` is inert here: the fitted attributes are hand-populated and
    # ICA.fit (the only place `method` branches) is never called, but
    # ICA.__init__ still requires a valid method name.
    ica.info = self.info_
    ica.ch_names = list(self.ch_names_)
    ica.n_components_ = n_ch
    ica.pca_mean_ = pca_mean
    ica.pca_components_ = pca_components
    ica.pca_explained_variance_ = cov_evals
    ica.unmixing_matrix_ = unmixing
    ica.pre_whitener_ = self.pre_whitener_
    # The number of iterations the fit ran: ``iteration`` is the 0-based
    # index of the last one whose E-step ran (issue #339), which is also the
    # reference's own ``iter`` count at a stop.
    ica.n_iter_ = int(getattr(amica.model_, "iteration", 0)) + 1
    # MNE's own fit sets these; read_ica_eeglab (the precedent for building
    # an ICA from an external decomposition) sets reject_=None. Without them
    # ICA.save()/plot_properties raise AttributeError.
    ica.reject_ = None
    ica.n_samples_ = int(self._n_samples) if self._n_samples is not None else 0
    ica._update_mixing_matrix()
    ica._update_ica_names()
    ica.current_fit = self._fit_kind

    if n_residual:
        logger.info(
            "AMICAICA.to_mne_ica: exporting %d ICA components plus a "
            "%d-dimensional PCA residual (the subspace rank reduction "
            "discarded before fitting). ICA.apply restores the residual by "
            "default, which goes beyond the Fortran reference; pass "
            "n_pca_components=%d (ica.n_components_) to apply for the "
            "reference's rank-reduced reconstruction (issue #322).",
            n_ch,
            n_residual,
            n_ch,
        )

    self._mne_ica_cache[model_idx] = ica
    return ica

get_sources(inst, *args, model_idx=0, **kwargs)

Sources for inst from model model_idx (see ICA.get_sources).

model_idx is keyword-only so positional arguments pass straight through to MNE's ICA.get_sources (e.g. add_channels, start/stop).

Source code in pamica/mne_compat/core.py

def get_sources(self, inst, *args, model_idx: int = 0, **kwargs):
    """Sources for ``inst`` from model ``model_idx`` (see ``ICA.get_sources``).

    ``model_idx`` is keyword-only so positional arguments pass straight
    through to MNE's ``ICA.get_sources`` (e.g. ``add_channels``,
    ``start``/``stop``).
    """
    return self.to_mne_ica(model_idx).get_sources(inst, *args, **kwargs)

apply(inst, *args, model_idx=0, **kwargs)

Remove selected components of model model_idx and back-project.

model_idx is keyword-only so positional arguments pass straight through to MNE's ICA.apply. Pass exclude=[...] (or set it on the exported ICA) to drop components; with no exclusions this reconstructs the input.

For a rank-reduced fit (pcakeep/pcadb or automatic rank detection) the PCA residual, the subspace the reduction discarded and AMICA never modeled, is restored by default: excluding a component removes that component's back-projection and nothing else, as with MNE's own ICA (issue #322). This differs from the Fortran reference, whose back-projection reconstructs only the retained subspace. To get that rank-reduced reconstruction, pass n_pca_components=ica.n_components_; a float n_pca_components keeps residual rows by cumulative explained variance, following MNE. A full-rank fit has no residual, so the option changes nothing there.

Source code in pamica/mne_compat/core.py

def apply(self, inst, *args, model_idx: int = 0, **kwargs):
    """Remove selected components of model ``model_idx`` and back-project.

    ``model_idx`` is keyword-only so positional arguments pass straight
    through to MNE's ``ICA.apply``. Pass ``exclude=[...]`` (or set it on the
    exported ICA) to drop components; with no exclusions this reconstructs
    the input.

    For a rank-reduced fit (``pcakeep``/``pcadb`` or automatic rank
    detection) the PCA residual, the subspace the reduction discarded and
    AMICA never modeled, is **restored by default**: excluding a component
    removes that component's back-projection and nothing else, as with
    MNE's own ICA (issue #322). This differs from the Fortran reference,
    whose back-projection reconstructs only the retained subspace. To get
    that rank-reduced reconstruction, pass
    ``n_pca_components=ica.n_components_``; a float ``n_pca_components``
    keeps residual rows by cumulative explained variance, following MNE.
    A full-rank fit has no residual, so the option changes nothing there.
    """
    return self.to_mne_ica(model_idx).apply(inst, *args, **kwargs)

get_components(*, model_idx=0)

Scalp maps of model model_idx, (n_channels, n_components).

Source code in pamica/mne_compat/core.py

def get_components(self, *, model_idx: int = 0) -> np.ndarray:
    """Scalp maps of model ``model_idx``, ``(n_channels, n_components)``."""
    return self.to_mne_ica(model_idx).get_components()

plot_components(*args, model_idx=0, **kwargs)

Plot model model_idx component topographies (see ICA.plot_components).

Source code in pamica/mne_compat/core.py

def plot_components(self, *args, model_idx: int = 0, **kwargs):
    """Plot model ``model_idx`` component topographies (see ``ICA.plot_components``)."""
    return self.to_mne_ica(model_idx).plot_components(*args, **kwargs)

plot_sources(inst, *args, model_idx=0, **kwargs)

Plot model model_idx component time courses (see ICA.plot_sources).

Source code in pamica/mne_compat/core.py

def plot_sources(self, inst, *args, model_idx: int = 0, **kwargs):
    """Plot model ``model_idx`` component time courses (see ``ICA.plot_sources``)."""
    return self.to_mne_ica(model_idx).plot_sources(inst, *args, **kwargs)

get_model_probability(inst, *, reject_by_annotation=True)

Per-sample posterior probability of each model on inst (dominance).

Returns P(model | sample) as (n_models, n_samples) via :meth:AMICA.model_probability on inst's data (Epochs are concatenated along time). Each column sums to 1; all ones for a single model. MNE's own ICA has no multi-model concept, so this is exposed here rather than through the exported per-model ICA objects.

Parameters:

Name	Type	Description	Default
inst	BaseRaw | BaseEpochs	Data to score.	required
reject_by_annotation	bool	For Raw input, evaluate only samples not covered by inst's bad_* annotations, mirroring MNE's ICA.score_sources default (issue #251). The output keeps the full n_times timeline -- rejected columns are NaN -- so it stays aligned with inst without any manual re-indexing. False scores every sample (the model evaluates fine on artifact segments; their probabilities are simply dominated by whatever model claims the artifact).	True

Source code in pamica/mne_compat/core.py

def get_model_probability(
    self, inst, *, reject_by_annotation: bool = True
) -> np.ndarray:
    """Per-sample posterior probability of each model on ``inst`` (dominance).

    Returns ``P(model | sample)`` as ``(n_models, n_samples)`` via
    :meth:`AMICA.model_probability` on ``inst``'s data (``Epochs`` are
    concatenated along time). Each column sums to 1; all ones for a single
    model. MNE's own ``ICA`` has no multi-model concept, so this is exposed
    here rather than through the exported per-model ICA objects.

    Parameters
    ----------
    inst : mne.io.BaseRaw | mne.BaseEpochs
        Data to score.
    reject_by_annotation : bool, default=True
        For ``Raw`` input, evaluate only samples not covered by ``inst``'s
        ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
        default (issue #251). The output keeps the full ``n_times``
        timeline -- rejected columns are ``NaN`` -- so it stays aligned
        with ``inst`` without any manual re-indexing. ``False`` scores
        every sample (the model evaluates fine on artifact segments; their
        probabilities are simply dominated by whatever model claims the
        artifact).
    """
    self._check_fitted("compute the model probability")
    if self.amica_ is None:
        raise RuntimeError(
            "AMICAICA: internal state is inconsistent; refit before scoring."
        )
    X, good = self._data_for(inst, reject_by_annotation=reject_by_annotation)
    return _expand_to_timeline(self.amica_.model_probability(X), good)

plot_model_probability(inst, *, srate=None, reject_by_annotation=True, **kwargs)

Plot per-model probability + best-model log-likelihood over inst.

Delegates to :func:pamica.viz.plot_model_probability with the live per-model log-likelihood (:meth:AMICA.model_loglik) on inst's data. srate defaults to the fitted recording's sampling rate, so the x-axis is in seconds; extra keywords (smooth_sec, window_sec, axes) pass through.

With reject_by_annotation (default, Raw only), samples covered by inst's bad_* annotations plot as gaps: the log-likelihood is evaluated on the good samples only and rejected columns are NaN, so the time axis stays aligned with inst (issue #251). Note smooth_sec widens the gaps by the smoothing window, since the Hanning smoothing propagates NaN across its support.

Source code in pamica/mne_compat/core.py

def plot_model_probability(
    self,
    inst,
    *,
    srate: Optional[float] = None,
    reject_by_annotation: bool = True,
    **kwargs,
):
    """Plot per-model probability + best-model log-likelihood over ``inst``.

    Delegates to :func:`pamica.viz.plot_model_probability` with the live
    per-model log-likelihood (:meth:`AMICA.model_loglik`) on ``inst``'s
    data. ``srate`` defaults to the fitted recording's sampling rate, so the
    x-axis is in seconds; extra keywords (``smooth_sec``, ``window_sec``,
    ``axes``) pass through.

    With ``reject_by_annotation`` (default, ``Raw`` only), samples covered
    by ``inst``'s ``bad_*`` annotations plot as gaps: the log-likelihood is
    evaluated on the good samples only and rejected columns are ``NaN``, so
    the time axis stays aligned with ``inst`` (issue #251). Note
    ``smooth_sec`` widens the gaps by the smoothing window, since the
    Hanning smoothing propagates ``NaN`` across its support.
    """
    from ..viz import plot_model_probability as _plot_model_probability

    self._check_fitted("plot the model probability")
    if self.amica_ is None or self.info_ is None:
        raise RuntimeError(
            "AMICAICA: internal state is inconsistent; refit before plotting."
        )
    X, good = self._data_for(inst, reject_by_annotation=reject_by_annotation)
    lht = _expand_to_timeline(self.amica_.model_loglik(X), good)
    if srate is None:
        srate = float(self.info_["sfreq"])
    return _plot_model_probability(lht=lht, srate=srate, **kwargs)

mir(inst, *, model_idx=0, nbins=None, reject_by_annotation=True)

Mutual Information Reduction of model model_idx on inst.

How much mutual information the fitted unmixing removes from the data, in nats (issue #133). Delegates to :meth:AMICA.mir on inst's fitted-channel data; MIR is shift-invariant, so mean/c centering is irrelevant.

Parameters:

Name	Type	Description	Default
inst	BaseRaw | BaseEpochs	Data to score (Epochs concatenated along time).	required
model_idx	int	Which model's unmixing to use.	0
nbins	int	Histogram bin count; see :func:pamica.metrics.mir.	None
reject_by_annotation	bool	For Raw input, score only samples not covered by inst's bad_* annotations, mirroring MNE's ICA.score_sources default (issue #251).	True

Returns:

Name	Type	Description
mir_nats	float
variance	float

Source code in pamica/mne_compat/core.py

def mir(
    self,
    inst,
    *,
    model_idx: int = 0,
    nbins: Optional[int] = None,
    reject_by_annotation: bool = True,
) -> tuple:
    """Mutual Information Reduction of model ``model_idx`` on ``inst``.

    How much mutual information the fitted unmixing removes from the data,
    in nats (issue #133). Delegates to :meth:`AMICA.mir` on ``inst``'s
    fitted-channel data; MIR is shift-invariant, so mean/``c`` centering is
    irrelevant.

    Parameters
    ----------
    inst : mne.io.BaseRaw | mne.BaseEpochs
        Data to score (``Epochs`` concatenated along time).
    model_idx : int, default=0
        Which model's unmixing to use.
    nbins : int, optional
        Histogram bin count; see :func:`pamica.metrics.mir`.
    reject_by_annotation : bool, default=True
        For ``Raw`` input, score only samples not covered by ``inst``'s
        ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
        default (issue #251).

    Returns
    -------
    mir_nats : float
    variance : float
    """
    self._check_fitted("compute MIR")
    model_idx = self._check_model_idx(model_idx)
    assert self.amica_ is not None
    X, _ = self._data_for(inst, reject_by_annotation=reject_by_annotation)
    return self.amica_.mir(X, model_idx=model_idx, nbins=nbins)

pmi(inst, *, model_idx=0, nbins=None, reject_by_annotation=True)

Pairwise Mutual Information between model model_idx's sources on inst.

The residual pairwise dependence between fitted sources, in nats (issue #133). Delegates to :meth:AMICA.pmi on inst's fitted-channel data.

Parameters:

Name	Type	Description	Default
inst	BaseRaw | BaseEpochs	Data to score (Epochs concatenated along time).	required
model_idx	int	Which model's sources to use.	0
nbins	int	Histogram bin count; see :func:pamica.metrics.pairwise_mi.	None
reject_by_annotation	bool	For Raw input, score only samples not covered by inst's bad_* annotations, mirroring MNE's ICA.score_sources default (issue #251).	True

Returns:

Name	Type	Description
mi_matrix	np.ndarray of shape (n_components, n_components)	Symmetric; the diagonal is each source's own entropy.

Source code in pamica/mne_compat/core.py

def pmi(
    self,
    inst,
    *,
    model_idx: int = 0,
    nbins: Optional[int] = None,
    reject_by_annotation: bool = True,
) -> np.ndarray:
    """Pairwise Mutual Information between model ``model_idx``'s sources on ``inst``.

    The residual pairwise dependence between fitted sources, in nats
    (issue #133). Delegates to :meth:`AMICA.pmi` on ``inst``'s fitted-channel
    data.

    Parameters
    ----------
    inst : mne.io.BaseRaw | mne.BaseEpochs
        Data to score (``Epochs`` concatenated along time).
    model_idx : int, default=0
        Which model's sources to use.
    nbins : int, optional
        Histogram bin count; see :func:`pamica.metrics.pairwise_mi`.
    reject_by_annotation : bool, default=True
        For ``Raw`` input, score only samples not covered by ``inst``'s
        ``bad_*`` annotations, mirroring MNE's ``ICA.score_sources``
        default (issue #251).

    Returns
    -------
    mi_matrix : np.ndarray of shape (n_components, n_components)
        Symmetric; the diagonal is each source's own entropy.
    """
    self._check_fitted("compute PMI")
    model_idx = self._check_model_idx(model_idx)
    assert self.amica_ is not None
    X, _ = self._data_for(inst, reject_by_annotation=reject_by_annotation)
    return self.amica_.pmi(X, model_idx=model_idx, nbins=nbins)

get_pdftype(*, model_idx=0)

Per-component source-density family code for model model_idx.

One integer per ICA component (0-4); map to names with :data:pamica.mne_compat.PDFTYPE_NAMES. All components share one family unless the adaptive switcher (pdftype=1) moved them (issue #26).

Source code in pamica/mne_compat/core.py

def get_pdftype(self, *, model_idx: int = 0) -> np.ndarray:
    """Per-component source-density family code for model ``model_idx``.

    One integer per ICA component (0-4); map to names with
    :data:`pamica.mne_compat.PDFTYPE_NAMES`. All components share one family
    unless the adaptive switcher (``pdftype=1``) moved them (issue #26).
    """
    self._check_fitted("get the density family")
    model_idx = self._check_model_idx(model_idx)
    assert self.amica_ is not None
    return self.amica_.get_pdftype(model_idx=model_idx)

get_rho(*, model_idx=0)

Generalized-Gaussian shape rho for model model_idx.

Shape (n_mix, n_components); rho == 2 is Gaussian-shaped, rho == 1 Laplacian, rho < 1 heavier-tailed. Meaningful only for the generalized-Gaussian family (pdftype=0).

Source code in pamica/mne_compat/core.py

def get_rho(self, *, model_idx: int = 0) -> np.ndarray:
    """Generalized-Gaussian shape ``rho`` for model ``model_idx``.

    Shape ``(n_mix, n_components)``; ``rho == 2`` is Gaussian-shaped,
    ``rho == 1`` Laplacian, ``rho < 1`` heavier-tailed. Meaningful only for
    the generalized-Gaussian family (``pdftype=0``).
    """
    self._check_fitted("get rho")
    model_idx = self._check_model_idx(model_idx)
    assert self.amica_ is not None
    return self.amica_.get_rho(model_idx=model_idx)

shared_components()

Components shared across models by share_comps (issue #60).

One group of (model_idx, component_idx) pairs per shared component; empty when nothing is shared (always so for a single model or a default multi-model fit with share_comps off).

Source code in pamica/mne_compat/core.py

def shared_components(self) -> list:
    """Components shared across models by ``share_comps`` (issue #60).

    One group of ``(model_idx, component_idx)`` pairs per shared component;
    empty when nothing is shared (always so for a single model or a default
    multi-model fit with ``share_comps`` off).
    """
    self._check_fitted("get the shared components")
    assert self.amica_ is not None
    return self.amica_.shared_components()