Phase 3: preserve pamica-specific metadata through the MNE wrapper

docs/api/mne-compat.md

@@ -86,4 +86,27 @@ 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:
+
+```python
+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`](amica.md).
+
 ::: pamica.mne_compat.AMICAICA

docs/changelog.md

@@ -42,6 +42,13 @@ MNE-Python compatibility layer (epic #139), additive: the scikit-learn-style
   folds each model's data-space center `c` into `pca_mean_`, so the round trip
   holds for the multi-model case too. `pamica.viz.plot_model_probability` now also
   accepts a live `lht` array, not only a written `AmicaOutput` (phase 2, #141).
+- pamica-specific fitted metadata is inspectable through the MNE wrapper rather
+  than silently dropped by the `mne.preprocessing.ICA` export: `get_pdftype(model_idx=...)`
+  returns each component's source-density family code (0-4, named by
+  `pamica.mne_compat.PDFTYPE_NAMES`), `get_rho(model_idx=...)` the
+  generalized-Gaussian shape parameters, and `shared_components()` the components
+  merged across models by `share_comps`. The same accessors are added to
+  `AMICA`/`AMICATorchNG` (phase 3, #142).
 
 ## 0.2.2
 

pamica/amica.py

@@ -477,6 +477,47 @@ def model_probability(self, X: np.ndarray) -> np.ndarray:
 
         return self.model_.model_probability(X)
 
+    def get_pdftype(self, model_idx: int = 0) -> np.ndarray:
+        """Per-source density-family code for model ``model_idx`` (issue #142).
+
+        Delegates to :meth:`AMICATorchNG.get_pdftype`. One integer per source
+        component (0-4; see :data:`pamica.torch_impl.PDFTYPE_NAMES`).
+
+        Returns
+        -------
+        np.ndarray of int, shape (n_sources,)
+        """
+        self._check_usable("get the density family")
+        assert self.model_ is not None
+
+        return self.model_.get_pdftype(model_idx=model_idx)
+
+    def get_rho(self, model_idx: int = 0) -> np.ndarray:
+        """Generalized-Gaussian shape ``rho`` for model ``model_idx`` (issue #142).
+
+        Delegates to :meth:`AMICATorchNG.get_rho`.
+
+        Returns
+        -------
+        np.ndarray of float, shape (n_mix, n_sources)
+        """
+        self._check_usable("get rho")
+        assert self.model_ is not None
+
+        return self.model_.get_rho(model_idx=model_idx)
+
+    def shared_components(self) -> list:
+        """Components shared across models by ``share_comps`` (issue #142).
+
+        Delegates to :meth:`AMICATorchNG.shared_components`: one group of
+        ``(model_idx, source_idx)`` pairs per shared column; empty when nothing
+        is shared.
+        """
+        self._check_usable("get the shared components")
+        assert self.model_ is not None
+
+        return self.model_.shared_components()
+
     def variance_order(
         self, model_idx: int = 0, return_svar: bool = False
     ) -> Union[np.ndarray, tuple]:

pamica/mne_compat/__init__.py

@@ -10,6 +10,6 @@
 pattern as the optional ``mlx`` backend).
 """
 
-from .core import AMICAICA
+from .core import AMICAICA, PDFTYPE_NAMES
 
-__all__ = ["AMICAICA"]
+__all__ = ["AMICAICA", "PDFTYPE_NAMES"]

pamica/mne_compat/core.py

@@ -26,6 +26,9 @@
 from mne.preprocessing import ICA as _MNEICA  # ty: ignore[unresolved-import]
 
 from ..amica import AMICA
+from ..torch_impl import PDFTYPE_NAMES
+
+__all__ = ["AMICAICA", "PDFTYPE_NAMES"]
 
 
 class AMICAICA:
@@ -415,6 +418,48 @@ def plot_model_probability(self, inst, *, srate: Optional[float] = None, **kwarg
             srate = float(self.info_["sfreq"])
         return _plot_model_probability(lht=lht, srate=srate, **kwargs)
 
+    # ------------------------------------------------------------------
+    # 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 column;
+        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
     # ------------------------------------------------------------------

pamica/tests/mne_tests/test_mne_wrapper.py

@@ -18,7 +18,7 @@
 
 mne = pytest.importorskip("mne")
 
-from pamica.mne_compat import AMICAICA  # noqa: E402  (after importorskip)
+from pamica.mne_compat import AMICAICA, PDFTYPE_NAMES  # noqa: E402  (after importorskip)
 
 mne.set_log_level("ERROR")
 
@@ -434,3 +434,34 @@ def test_degenerate_fit_refused_on_dominance_methods(raw):
         ica.get_model_probability(raw)
     with pytest.raises(RuntimeError, match="degenerate"):
         ica.plot_model_probability(raw)
+
+
+# --- pamica-specific metadata (issue #142) ----------------------------------
+def test_get_pdftype_and_rho_per_model(fitted_2m):
+    for h in range(2):
+        pdf = fitted_2m.get_pdftype(model_idx=h)
+        assert pdf.shape == (fitted_2m.n_components_,)
+        assert np.array_equal(np.unique(pdf), [0])  # default generalized Gaussian
+        assert PDFTYPE_NAMES[int(pdf[0])] == "generalized_gaussian"
+        rho = fitted_2m.get_rho(model_idx=h)
+        assert rho.shape == (fitted_2m.n_mix, fitted_2m.n_components_)
+        assert np.all(np.isfinite(rho))
+
+
+def test_shared_components_empty_by_default(fitted_2m):
+    assert fitted_2m.shared_components() == []
+
+
+def test_metadata_respects_model_idx_bounds(fitted_2m):
+    with pytest.raises(ValueError, match="out of range"):
+        fitted_2m.get_pdftype(model_idx=2)
+    with pytest.raises(ValueError, match="out of range"):
+        fitted_2m.get_rho(model_idx=5)
+
+
+def test_metadata_requires_fit():
+    ica = AMICAICA()
+    with pytest.raises(ValueError, match="must be fitted"):
+        ica.get_pdftype()
+    with pytest.raises(ValueError, match="must be fitted"):
+        ica.shared_components()

pamica/tests/torch_tests/test_amica_ng_wrapper.py

@@ -643,6 +643,12 @@ def test_unfitted_output_raises_not_fitted():
         model.mir(np.zeros((NW, 16)))
     with pytest.raises(ValueError, match="fitted"):
         model.pmi(np.zeros((NW, 16)))
+    with pytest.raises(ValueError, match="fitted"):
+        model.get_pdftype()
+    with pytest.raises(ValueError, match="fitted"):
+        model.get_rho()
+    with pytest.raises(ValueError, match="fitted"):
+        model.shared_components()
 
 
 def test_degenerate_fit_refuses_output(real_data, tmp_path, caplog):
@@ -675,6 +681,9 @@ def test_degenerate_fit_refuses_output(real_data, tmp_path, caplog):
         lambda: model.variance_order(),
         lambda: model.mir(real_data[:, :512]),
         lambda: model.pmi(real_data[:, :512]),
+        lambda: model.get_pdftype(),
+        lambda: model.get_rho(),
+        lambda: model.shared_components(),
     ):
         with pytest.raises(RuntimeError, match="degenerate.*nan_ll"):
             action()

pamica/tests/torch_tests/test_ng_metadata.py

@@ -0,0 +1,145 @@
+"""Fitted-parameter metadata accessors (issue #142).
+
+``get_pdftype``/``get_rho``/``shared_components`` expose the source-density
+family, GG shape, and component-sharing state that the EEGLAB/MNE exports do not
+carry. Real sample EEG only (no synthetic data).
+"""
+
+from pathlib import Path
+
+import numpy as np
+import pytest
+
+from pamica import AMICATorchNG
+from pamica.torch_impl import PDFTYPE_NAMES
+from pamica.torch_impl.utils import load_eeglab_data
+
+SAMPLE_DIR = Path(__file__).resolve().parents[2] / "sample_data"
+DATA_FILE = SAMPLE_DIR / "eeglab_data.fdt"
+NW = 32
+FIELD = 30504
+NMIX = 3
+SEED = 42
+
+pytestmark = pytest.mark.skipif(not DATA_FILE.exists(), reason="sample data missing")
+
+
+@pytest.fixture(scope="module")
+def real_data():
+    return load_eeglab_data(str(DATA_FILE), data_dim=NW, field_dim=FIELD).astype(
+        np.float64
+    )
+
+
+def _fit(X, n_models=1, max_iter=8, **kw):
+    m = AMICATorchNG(
+        n_channels=NW, n_models=n_models, n_mix=NMIX, seed=SEED, device="cpu", **kw
+    )
+    m.fit(X, max_iter=max_iter, verbose=False)
+    return m
+
+
+def test_pdftype_names_cover_all_codes():
+    assert set(PDFTYPE_NAMES) == {0, 1, 2, 3, 4}
+
+
+def test_get_pdftype_default_is_generalized_gaussian(real_data):
+    m = _fit(real_data)
+    pdf = m.get_pdftype()
+    assert pdf.shape == (NW,)
+    assert np.array_equal(np.unique(pdf), [0])  # pdftype=0 default
+    assert PDFTYPE_NAMES[int(pdf[0])] == "generalized_gaussian"
+
+
+def test_get_pdftype_reflects_fixed_family(real_data):
+    m = _fit(real_data, pdftype=2)  # Gaussian
+    assert np.array_equal(np.unique(m.get_pdftype()), [2])
+
+
+def test_get_rho_shape_and_bounds(real_data):
+    m = _fit(real_data)
+    rho = m.get_rho()
+    assert rho.shape == (NMIX, NW)
+    assert np.all(rho >= m.minrho) and np.all(rho <= m.maxrho)
+
+
+def test_shared_components_empty_without_sharing(real_data):
+    assert _fit(real_data, n_models=1).shared_components() == []
+    assert _fit(real_data, n_models=2).shared_components() == []  # share off default
+
+
+def test_shared_components_reports_partial_merges(real_data):
+    """A PARTIAL merge (some columns single-model) exercises both the grouping
+    and the exclude-filter: shared groups appear, non-shared columns do not."""
+    m = _fit(
+        real_data[:, :4096],
+        n_models=2,
+        max_iter=25,
+        block_size=1024,
+        share_comps=True,
+        share_start=8,
+        share_iter=10,
+        do_newton=True,
+    )
+    groups = m.shared_components()
+    # Partial: not zero (something merged), and fewer than n_sources (some
+    # columns stay single-model and are correctly excluded from the groups).
+    assert 0 < len(groups) < NW
+    for group in groups:
+        models = {h for h, _ in group}
+        assert len(models) >= 2  # a shared column spans >= 2 models
+        assert all(0 <= h < 2 and 0 <= i < NW for h, i in group)
+
+
+def test_get_rho_differs_per_model(real_data):
+    """get_rho indexes comp_list per model, so a 2-model fit's rho differs."""
+    m = _fit(real_data, n_models=2, max_iter=10)
+    assert not np.allclose(m.get_rho(model_idx=0), m.get_rho(model_idx=1))
+
+
+def test_get_pdftype_adaptive_switcher_is_mixed(real_data):
+    """The adaptive switcher (pdftype=1) moves sources into families {1, 4};
+    get_pdftype must surface those per-source codes, not a uniform array."""
+    m = AMICATorchNG(
+        n_channels=NW,
+        n_models=1,
+        n_mix=1,  # adaptive mode is single-component
+        seed=SEED,
+        device="cpu",
+        pdftype=1,
+        kurt_start=3,
+        num_kurt=5,
+        kurt_int=1,
+    )
+    m.fit(real_data, max_iter=20, verbose=False)
+    codes = set(np.unique(m.get_pdftype()).tolist())
+    assert codes.issubset({1, 4}) and codes  # only the two cosh families
+
+
+def test_metadata_accessors_reject_bad_model_idx(real_data):
+    m = _fit(real_data, n_models=2, max_iter=5)
+    for call in (m.get_pdftype, m.get_rho):
+        with pytest.raises(ValueError, match="out of range"):
+            call(model_idx=2)
+        with pytest.raises(ValueError, match="out of range"):
+            call(model_idx=-1)  # negative would silently wrap without the guard
+
+
+def test_metadata_accessors_require_fit():
+    m = AMICATorchNG(n_channels=NW, n_models=1, n_mix=NMIX, seed=SEED, device="cpu")
+    for call in (m.get_pdftype, m.get_rho, m.shared_components):
+        with pytest.raises(RuntimeError, match="fitted"):
+            call()
+
+
+def test_amica_wrapper_delegates_metadata(real_data):
+    """The scikit-learn-style AMICA wrapper exposes the same accessors."""
+    from pamica import AMICA
+
+    a = AMICA(n_models=2, n_mix=NMIX, device="cpu", verbose=False)
+    a.fit(real_data, max_iter=8, seed=SEED)
+    assert a.get_pdftype(model_idx=1).shape == (NW,)
+    assert a.get_rho(model_idx=1).shape == (NMIX, NW)
+    assert a.shared_components() == []  # share off
+    with pytest.raises(ValueError, match="out of range"):
+        a.get_pdftype(model_idx=-1)  # guard propagates through the wrapper

pamica/torch_impl/__init__.py

@@ -5,11 +5,12 @@
 with support for CUDA, ROCm, and Apple Silicon (MPS) backends.
 """
 
-from .core import AMICATorchNG
+from .core import AMICATorchNG, PDFTYPE_NAMES
 from .utils import setup_device, check_numerical_stability
 
 __all__ = [
     "AMICATorchNG",
+    "PDFTYPE_NAMES",
     "setup_device",
     "check_numerical_stability",
 ]