.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_tutorials/10_orientation_and_uncertainty_representations.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_tutorials_10_orientation_and_uncertainty_representations.py: 10. Orientation and Uncertainty Representations =============================================== This tutorial compares the uncertainty representations used in CaliBrain for fixed orientation, reduced free-orientation MEG, and free-orientation EEG. The main point is practical: the calibration method depends not only on the solver output, but also on the geometric object used to represent uncertainty. .. GENERATED FROM PYTHON SOURCE LINES 16-35 Scientific question ------------------- CaliBrain supports three main source configurations: - fixed orientation, where each source has one coefficient; - reduced free-orientation MEG, where each source has two coefficients; - free-orientation EEG, where each source has three coefficients. The uncertainty object changes with that representation: - fixed orientation uses scalar marginal intervals; - reduced free-orientation MEG can use either marginal intervals or full-covariance ellipses; - free-orientation EEG can use either marginal intervals or full-covariance ellipsoids. In the current workflow, calibration is usually evaluated after temporal aggregation, so this tutorial uses the aggregated coverage routines. .. GENERATED FROM PYTHON SOURCE LINES 35-46 .. code-block:: Python import matplotlib.pyplot as plt import numpy as np from calibrain import UncertaintyEstimator rng = np.random.default_rng(23) nominal_coverages = np.linspace(0.0, 1.0, 11) uncertainty = UncertaintyEstimator(nominal_coverages=nominal_coverages) .. GENERATED FROM PYTHON SOURCE LINES 47-52 Fixed orientation: scalar intervals ----------------------------------- For fixed orientation, each source has one value per time point. The uncertainty summary is therefore just a scalar variance per source. .. GENERATED FROM PYTHON SOURCE LINES 52-77 .. code-block:: Python n_sources_fixed = 36 n_times = 90 time = np.linspace(-0.1, 0.7, n_times) x_true_fixed = np.zeros((n_sources_fixed, n_times)) fixed_active = rng.choice(n_sources_fixed, size=4, replace=False) fixed_waveform = np.exp(-0.5 * ((time - 0.18) / 0.05) ** 2) x_true_fixed[fixed_active] = rng.normal(1.0, 0.12, size=(4, 1)) * fixed_waveform posterior_var_fixed = np.full(n_sources_fixed, 0.055**2) x_hat_fixed = x_true_fixed + rng.normal( scale=np.sqrt(posterior_var_fixed)[:, None], size=x_true_fixed.shape, ) curve_fixed = uncertainty.calibration_curve_intervals_aggregated( x_true=x_true_fixed, x_hat=x_hat_fixed, posterior_var=posterior_var_fixed, ) print("fixed interval type:", curve_fixed["interval_type"]) print("fixed empirical coverages:", np.round(curve_fixed["empirical_coverages"], 3)) .. rst-class:: sphx-glr-script-out .. code-block:: none fixed interval type: marginal fixed empirical coverages: [0. 0.056 0.139 0.222 0.417 0.444 0.583 0.694 0.833 0.917 1. ] .. GENERATED FROM PYTHON SOURCE LINES 78-88 Reduced free-orientation MEG: marginal intervals vs full-covariance ellipses ----------------------------------------------------------------------------- For reduced free-orientation MEG, each source has two coefficients. CaliBrain can evaluate calibration in two ways: - ``marginal``: component-wise intervals using only the diagonal variances; - ``full_cov``: two-dimensional ellipses using the full 2x2 covariance block. The second representation preserves the within-source covariance geometry. .. GENERATED FROM PYTHON SOURCE LINES 88-129 .. code-block:: Python n_sources_meg = 24 x_true_meg = np.zeros((n_sources_meg, 2, n_times)) meg_active = rng.choice(n_sources_meg, size=4, replace=False) for source_idx in meg_active: amp1, amp2 = rng.normal(loc=[1.0, 0.8], scale=[0.08, 0.08]) waveform = np.exp(-0.5 * ((time - rng.uniform(0.12, 0.22)) / 0.05) ** 2) x_true_meg[source_idx, 0] = amp1 * waveform x_true_meg[source_idx, 1] = amp2 * waveform x_hat_meg = np.zeros_like(x_true_meg) posterior_cov_meg = np.zeros((n_sources_meg, 2, 2)) V_tan = np.zeros((n_sources_meg, 3, 2)) V_tan[:, 0, 0] = 1.0 V_tan[:, 1, 1] = 1.0 x_true_meg_3d = np.zeros((n_sources_meg, 3, n_times)) x_true_meg_3d[:, :2, :] = x_true_meg for source_idx in range(n_sources_meg): cov_block = np.array([[0.030**2, 0.00045], [0.00045, 0.025**2]]) mean_error = rng.multivariate_normal(mean=[0.0, 0.0], cov=cov_block) x_hat_meg[source_idx, 0] = x_true_meg[source_idx, 0] + mean_error[0] x_hat_meg[source_idx, 1] = x_true_meg[source_idx, 1] + mean_error[1] posterior_cov_meg[source_idx] = cov_block curve_meg_marginal = uncertainty.calibration_curve_componentwise_meg_free_aggregated( x_true_2d=x_true_meg, x_hat_2d=x_hat_meg, posterior_uncert_2d=posterior_cov_meg, ) curve_meg_full = uncertainty.calibration_curve_ellipse_meg_free_aggregated( x_true_3d=x_true_meg_3d, x_hat_2d=x_hat_meg, posterior_cov_2d=posterior_cov_meg, V_tan=V_tan, ) print("MEG marginal empirical coverages:", np.round(curve_meg_marginal["empirical_coverages"], 3)) print("MEG full_cov empirical coverages:", np.round(curve_meg_full["empirical_coverages"], 3)) .. rst-class:: sphx-glr-script-out .. code-block:: none MEG marginal empirical coverages: [0. 0. 0. 0.021 0.042 0.062 0.062 0.062 0.062 0.104 1. ] MEG full_cov empirical coverages: [0. 0. 0. 0. 0. 0. 0. 0. 0. 0.042 1. ] .. GENERATED FROM PYTHON SOURCE LINES 130-136 Free-orientation EEG: marginal intervals vs full-covariance ellipsoids ---------------------------------------------------------------------- For free-orientation EEG, each source has three coefficients. The same logic applies, but the full-covariance representation is now a three-dimensional ellipsoid instead of a two-dimensional ellipse. .. GENERATED FROM PYTHON SOURCE LINES 136-179 .. code-block:: Python n_sources_eeg = 20 x_true_eeg = np.zeros((n_sources_eeg, 3, n_times)) eeg_active = rng.choice(n_sources_eeg, size=3, replace=False) for source_idx in eeg_active: amps = rng.normal(loc=[1.0, 0.75, 0.55], scale=[0.08, 0.08, 0.08]) waveform = np.exp(-0.5 * ((time - rng.uniform(0.10, 0.2)) / 0.045) ** 2) x_true_eeg[source_idx, 0] = amps[0] * waveform x_true_eeg[source_idx, 1] = amps[1] * waveform x_true_eeg[source_idx, 2] = amps[2] * waveform x_hat_eeg = np.zeros_like(x_true_eeg) posterior_cov_eeg = np.zeros((n_sources_eeg, 3, 3)) for source_idx in range(n_sources_eeg): cov_block = np.array( [ [0.032**2, 0.00035, 0.00025], [0.00035, 0.028**2, 0.00022], [0.00025, 0.00022, 0.024**2], ] ) mean_error = rng.multivariate_normal(mean=[0.0, 0.0, 0.0], cov=cov_block) x_hat_eeg[source_idx, 0] = x_true_eeg[source_idx, 0] + mean_error[0] x_hat_eeg[source_idx, 1] = x_true_eeg[source_idx, 1] + mean_error[1] x_hat_eeg[source_idx, 2] = x_true_eeg[source_idx, 2] + mean_error[2] posterior_cov_eeg[source_idx] = cov_block curve_eeg_marginal = uncertainty.calibration_curve_componentwise_eeg_free_aggregated( x_true=x_true_eeg, x_hat=x_hat_eeg, posterior_uncert=posterior_cov_eeg, ) curve_eeg_full = uncertainty.calibration_curve_ellipsoid_eeg_free_aggregated( x_true=x_true_eeg, x_hat=x_hat_eeg, posterior_cov=posterior_cov_eeg, ) print("EEG marginal empirical coverages:", np.round(curve_eeg_marginal["empirical_coverages"], 3)) print("EEG full_cov empirical coverages:", np.round(curve_eeg_full["empirical_coverages"], 3)) .. rst-class:: sphx-glr-script-out .. code-block:: none EEG marginal empirical coverages: [0. 0.017 0.017 0.033 0.05 0.05 0.05 0.083 0.117 0.15 1. ] EEG full_cov empirical coverages: [0. 0. 0. 0. 0. 0. 0. 0. 0. 0. 1.] .. GENERATED FROM PYTHON SOURCE LINES 180-190 Compare the calibration curves ------------------------------ This plot summarizes the main distinction: - fixed orientation has one natural marginal interval representation; - free orientation allows either a marginal representation or a full local covariance representation; - ``marginal`` and ``full_cov`` are therefore two different calibration diagnostics, not two labels for the same object. .. GENERATED FROM PYTHON SOURCE LINES 190-251 .. code-block:: Python fig, axes = plt.subplots(1, 3, figsize=(13.0, 4.0), sharex=True, sharey=True) axes[0].plot([0, 1], [0, 1], "--", color="0.5", label="perfect calibration") axes[0].plot( curve_fixed["nominal_coverages"], curve_fixed["empirical_coverages"], "o-", color="#4c72b0", label="fixed marginal", ) axes[0].set_title("Fixed orientation") axes[0].set_ylabel("Empirical coverage") axes[0].legend(loc="lower right") axes[1].plot([0, 1], [0, 1], "--", color="0.5") axes[1].plot( curve_meg_marginal["nominal_coverages"], curve_meg_marginal["empirical_coverages"], "o-", color="#55a868", label="marginal", ) axes[1].plot( curve_meg_full["nominal_coverages"], curve_meg_full["empirical_coverages"], "s-", color="#c44e52", label="full_cov", ) axes[1].set_title("Reduced free-orientation MEG") axes[1].legend(loc="lower right") axes[2].plot([0, 1], [0, 1], "--", color="0.5") axes[2].plot( curve_eeg_marginal["nominal_coverages"], curve_eeg_marginal["empirical_coverages"], "o-", color="#55a868", label="marginal", ) axes[2].plot( curve_eeg_full["nominal_coverages"], curve_eeg_full["empirical_coverages"], "s-", color="#c44e52", label="full_cov", ) axes[2].set_title("Free-orientation EEG") axes[2].legend(loc="lower right") for ax in axes: ax.set( xlabel="Nominal coverage", xlim=(0, 1), ylim=(0, 1), ) ax.grid(True, linestyle="--", alpha=0.35) fig.tight_layout() .. image-sg:: /auto_tutorials/images/sphx_glr_10_orientation_and_uncertainty_representations_001.png :alt: Fixed orientation, Reduced free-orientation MEG, Free-orientation EEG :srcset: /auto_tutorials/images/sphx_glr_10_orientation_and_uncertainty_representations_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 252-271 Practical interpretation ------------------------ The calibration workflow should always be read together with the uncertainty representation: - fixed orientation uses scalar source-wise intervals; - free orientation with ``marginal`` checks component-wise intervals; - free orientation with ``full_cov`` checks local ellipses or ellipsoids. In current CaliBrain workflows, all of these diagnostics are typically evaluated after averaging over time. The distinction between ``precal``, ``post_oracle``, ``post_pooled``, ``post_pooled_mismatch``, and ``post_fixed`` then acts on these coverage curves; it does not redefine the underlying uncertainty geometry. The next tutorial applies those calibration modes explicitly: - :doc:`Calibration Methods ` .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.274 seconds) .. _sphx_glr_download_auto_tutorials_10_orientation_and_uncertainty_representations.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: 10_orientation_and_uncertainty_representations.ipynb <10_orientation_and_uncertainty_representations.ipynb>` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: 10_orientation_and_uncertainty_representations.py <10_orientation_and_uncertainty_representations.py>` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: 10_orientation_and_uncertainty_representations.zip <10_orientation_and_uncertainty_representations.zip>`