Python ICA.find_sources_epochs方法代码示例

# 需要导入模块: from mne.preprocessing import ICA [as 别名]
# 或者: from mne.preprocessing.ICA import find_sources_epochs [as 别名]

#.........这里部分代码省略.........
        assert_true(isinstance(ica_read.info, Info))

        assert_raises(RuntimeError, ica_read.decompose_raw, raw)
        sources = ica.get_sources_raw(raw)
        sources2 = ica_read.get_sources_raw(raw)
        assert_array_almost_equal(sources, sources2)

        _raw1 = ica.pick_sources_raw(raw, exclude=[1])
        _raw2 = ica_read.pick_sources_raw(raw, exclude=[1])
        assert_array_almost_equal(_raw1[:, :][0], _raw2[:, :][0])

    os.remove(test_ica_fname)
    # check scrore funcs
    for name, func in score_funcs.items():
        if name in score_funcs_unsuited:
            continue
        scores = ica.find_sources_raw(raw, target='EOG 061', score_func=func,
                                      start=0, stop=10)
        assert_true(ica.n_components_ == len(scores))

    # check univariate stats
    scores = ica.find_sources_raw(raw, score_func=stats.skew)
    # check exception handling
    assert_raises(ValueError, ica.find_sources_raw, raw,
                  target=np.arange(1))

    params = []
    params += [(None, -1, slice(2), [0, 1])]  # varicance, kurtosis idx params
    params += [(None, 'MEG 1531')]  # ECG / EOG channel params
    for idx, ch_name in product(*params):
        ica.detect_artifacts(raw, start_find=0, stop_find=50, ecg_ch=ch_name,
                             eog_ch=ch_name, skew_criterion=idx,
                             var_criterion=idx, kurt_criterion=idx)
    ## score funcs epochs ##

    # check score funcs
    for name, func in score_funcs.items():
        if name in score_funcs_unsuited:
            continue
        scores = ica.find_sources_epochs(epochs_eog, target='EOG 061',
                                         score_func=func)
        assert_true(ica.n_components_ == len(scores))

    # check univariate stats
    scores = ica.find_sources_epochs(epochs, score_func=stats.skew)

    # check exception handling
    assert_raises(ValueError, ica.find_sources_epochs, epochs,
                  target=np.arange(1))

    # ecg functionality
    ecg_scores = ica.find_sources_raw(raw, target='MEG 1531',
                                      score_func='pearsonr')

    with warnings.catch_warnings(record=True):  # filter attenuation warning
        ecg_events = ica_find_ecg_events(raw,
                                         sources[np.abs(ecg_scores).argmax()])

    assert_true(ecg_events.ndim == 2)

    # eog functionality
    eog_scores = ica.find_sources_raw(raw, target='EOG 061',
                                      score_func='pearsonr')
    with warnings.catch_warnings(record=True):  # filter attenuation warning
        eog_events = ica_find_eog_events(raw,
                                         sources[np.abs(eog_scores).argmax()])

    assert_true(eog_events.ndim == 2)

    # Test ica fiff export
    ica_raw = ica.sources_as_raw(raw, start=0, stop=100)
    assert_true(ica_raw.last_samp - ica_raw.first_samp == 100)
    assert_true(len(ica_raw._filenames) == 0)  # API consistency
    ica_chans = [ch for ch in ica_raw.ch_names if 'ICA' in ch]
    assert_true(ica.n_components_ == len(ica_chans))
    test_ica_fname = op.join(op.abspath(op.curdir), 'test-ica_raw.fif')
    ica.n_components = np.int32(ica.n_components)
    ica_raw.save(test_ica_fname, overwrite=True)
    ica_raw2 = io.Raw(test_ica_fname, preload=True)
    assert_allclose(ica_raw._data, ica_raw2._data, rtol=1e-5, atol=1e-4)
    ica_raw2.close()
    os.remove(test_ica_fname)

    # Test ica epochs export
    ica_epochs = ica.sources_as_epochs(epochs)
    assert_true(ica_epochs.events.shape == epochs.events.shape)
    sources_epochs = ica.get_sources_epochs(epochs)
    assert_array_equal(ica_epochs.get_data(), sources_epochs)
    ica_chans = [ch for ch in ica_epochs.ch_names if 'ICA' in ch]
    assert_true(ica.n_components_ == len(ica_chans))
    assert_true(ica.n_components_ == ica_epochs.get_data().shape[1])
    assert_true(ica_epochs.raw is None)
    assert_true(ica_epochs.preload is True)

    # test float n pca components
    ica.pca_explained_variance_ = np.array([0.2] * 5)
    ica.n_components_ = 0
    for ncomps, expected in [[0.3, 1], [0.9, 4], [1, 1]]:
        ncomps_ = _check_n_pca_components(ica, ncomps)
        assert_true(ncomps_ == expected)

# 需要导入模块: from mne.preprocessing import ICA [as 别名]
# 或者: from mne.preprocessing.ICA import find_sources_epochs [as 别名]
mne.viz.plot_events(events, raw.info["sfreq"])

event_ids = {"faces": 1, "scrambled": 2}

tmin, tmax = -0.2, 0.6
baseline = None  # no baseline as high-pass is applied
reject = dict(mag=1.5e-12)

epochs = mne.Epochs(raw, events, event_ids, tmin, tmax, picks=picks, baseline=baseline, preload=True, reject=reject)

# Fit ICA, find and remove major artifacts

ica = ICA(None, 50).decompose_epochs(epochs, decim=2)

for ch_name in ["MRT51-2908", "MLF14-2908"]:  # ECG, EOG contaminated chs
    scores = ica.find_sources_epochs(epochs, ch_name, "pearsonr")
    ica.exclude += list(np.argsort(np.abs(scores))[-2:])

ica.plot_topomap(np.unique(ica.exclude))  # plot components found


# select ICA sources and reconstruct MEG signals, compute clean ERFs

epochs = ica.pick_sources_epochs(epochs)

evoked = [epochs[k].average() for k in event_ids]

contrast = evoked[1] - evoked[0]

evoked.append(contrast)

# 需要导入模块: from mne.preprocessing import ICA [as 别名]
# 或者: from mne.preprocessing.ICA import find_sources_epochs [as 别名]

#.........这里部分代码省略.........
        ica_read = read_ica(test_ica_fname)
        assert_true(ica.n_pca_components ==
                    ica_read.n_pca_components)
        ica.n_pca_components = 4
        ica_read.n_pca_components = 4

        ica.exclude = []
        ica.save(test_ica_fname)
        ica_read = read_ica(test_ica_fname)

        assert_true(ica.ch_names == ica_read.ch_names)

        assert_true(np.allclose(ica.mixing_matrix_, ica_read.mixing_matrix_,
                                rtol=1e-16, atol=1e-32))
        assert_array_equal(ica.pca_components_,
                           ica_read.pca_components_)
        assert_array_equal(ica.pca_mean_, ica_read.pca_mean_)
        assert_array_equal(ica.pca_explained_variance_,
                           ica_read.pca_explained_variance_)
        assert_array_equal(ica._pre_whitener, ica_read._pre_whitener)

        # assert_raises(RuntimeError, ica_read.decompose_raw, raw)
        sources = ica.get_sources_raw(raw)
        sources2 = ica_read.get_sources_raw(raw)
        assert_array_almost_equal(sources, sources2)

        _raw1 = ica.pick_sources_raw(raw, exclude=[1])
        _raw2 = ica_read.pick_sources_raw(raw, exclude=[1])
        assert_array_almost_equal(_raw1[:, :][0], _raw2[:, :][0])

    os.remove(test_ica_fname)
    # score funcs raw, with catch since "ties preclude exact" warning
    # XXX this should be fixed by a future PR...
    with warnings.catch_warnings(True) as w:
        sfunc_test = [ica.find_sources_raw(raw, target='EOG 061',
                score_func=n, start=0, stop=10)
                for n, f in score_funcs.items()]
    # score funcs raw

    # check lenght of scores
    [assert_true(ica.n_components_ == len(scores)) for scores in sfunc_test]

    # check univariate stats
    scores = ica.find_sources_raw(raw, score_func=stats.skew)
    # check exception handling
    assert_raises(ValueError, ica.find_sources_raw, raw,
                  target=np.arange(1))

    ## score funcs epochs ##

    # check lenght of scores
    # XXX this needs to be fixed, some of the score funcs don't seem to be
    # suited for the testing data.
    with warnings.catch_warnings(True) as w:
        sfunc_test = [ica.find_sources_epochs(epochs_eog, target='EOG 061',
                score_func=n)
                for n, f in score_funcs.items()]

    # check lenght of scores
    [assert_true(ica.n_components_ == len(scores)) for scores in sfunc_test]

    # check univariat stats
    scores = ica.find_sources_epochs(epochs, score_func=stats.skew)

    # check exception handling
    assert_raises(ValueError, ica.find_sources_epochs, epochs,
                  target=np.arange(1))

    # ecg functionality
    ecg_scores = ica.find_sources_raw(raw, target='MEG 1531',
                                      score_func='pearsonr')

    ecg_events = ica_find_ecg_events(raw, sources[np.abs(ecg_scores).argmax()])

    assert_true(ecg_events.ndim == 2)

    # eog functionality
    eog_scores = ica.find_sources_raw(raw, target='EOG 061',
                                      score_func='pearsonr')
    eog_events = ica_find_eog_events(raw, sources[np.abs(eog_scores).argmax()])

    assert_true(eog_events.ndim == 2)

    # Test ica fiff export
    ica_raw = ica.sources_as_raw(raw, start=0, stop=100)
    assert_true(ica_raw.last_samp - ica_raw.first_samp == 100)
    ica_chans = [ch for ch in ica_raw.ch_names if 'ICA' in ch]
    assert_true(ica.n_components_ == len(ica_chans))
    test_ica_fname = op.join(op.abspath(op.curdir), 'test_ica.fif')
    ica_raw.save(test_ica_fname)
    ica_raw2 = fiff.Raw(test_ica_fname, preload=True)
    assert_array_almost_equal(ica_raw._data, ica_raw2._data)
    ica_raw2.close()
    os.remove(test_ica_fname)

    # regression test for plot method
    assert_raises(ValueError, ica.plot_sources_raw, raw,
                  order=np.arange(50))
    assert_raises(ValueError, ica.plot_sources_epochs, epochs,
                  order=np.arange(50))

# 需要导入模块: from mne.preprocessing import ICA [as 别名]
# 或者: from mne.preprocessing.ICA import find_sources_epochs [as 别名]

# In[6]:

mne.viz.plot_drop_log(epochs.drop_log)


#### ICA

ica = ICA(n_components=0.90, n_pca_components=64,
          max_pca_components=100,
          noise_cov=None)

ica.decompose_epochs(epochs)

eog_scores_1_normal = ica.find_sources_epochs(epochs, target="EOG001",
                                              score_func="pearsonr")
eog_scores_2_normal = ica.find_sources_epochs(epochs, target="EOG003",
                                              score_func="pearsonr")

# get maximum correlation index for EOG
eog_source_idx_1_normal = np.abs(eog_scores_1_normal).argmax()
eog_source_idx_2_normal = np.abs(eog_scores_2_normal).argmax()

source_idx = range(0, ica.n_components_)
ica.plot_topomap(source_idx, ch_type="mag")


# select ICA sources and reconstruct MEG signals, compute clean ERFs
# Add detected artefact sources to exclusion list
# We now add the eog artefacts to the ica.exclusion list
if eog_source_idx_1_normal == eog_source_idx_2_normal:

# 需要导入模块: from mne.preprocessing import ICA [as 别名]
# 或者: from mne.preprocessing.ICA import find_sources_epochs [as 别名]

#.........这里部分代码省略.........
        assert_true(np.allclose(ica.mixing_matrix_, ica_read.mixing_matrix_,
                                rtol=1e-16, atol=1e-32))
        assert_array_equal(ica.pca_components_,
                           ica_read.pca_components_)
        assert_array_equal(ica.pca_mean_, ica_read.pca_mean_)
        assert_array_equal(ica.pca_explained_variance_,
                           ica_read.pca_explained_variance_)
        assert_array_equal(ica._pre_whitener, ica_read._pre_whitener)

        # assert_raises(RuntimeError, ica_read.decompose_raw, raw)
        sources = ica.get_sources_raw(raw)
        sources2 = ica_read.get_sources_raw(raw)
        assert_array_almost_equal(sources, sources2)

        _raw1 = ica.pick_sources_raw(raw, exclude=[1])
        _raw2 = ica_read.pick_sources_raw(raw, exclude=[1])
        assert_array_almost_equal(_raw1[:, :][0], _raw2[:, :][0])

    os.remove(test_ica_fname)
    # check scrore funcs
    for name, func in score_funcs.items():
        if name in score_funcs_unsuited:
            continue
        scores = ica.find_sources_raw(raw, target='EOG 061', score_func=func,
                                      start=0, stop=10)
        assert_true(ica.n_components_ == len(scores))

    # check univariate stats
    scores = ica.find_sources_raw(raw, score_func=stats.skew)
    # check exception handling
    assert_raises(ValueError, ica.find_sources_raw, raw,
                  target=np.arange(1))

    params = []
    params += [(None, -1, slice(2), [0, 1])]  # varicance, kurtosis idx params
    params += [(None, 'MEG 1531')]  # ECG / EOG channel params
    for idx, ch_name in product(*params):
        ica.detect_artifacts(raw, start_find=0, stop_find=50, ecg_ch=ch_name,
                             eog_ch=ch_name, skew_criterion=idx,
                             var_criterion=idx, kurt_criterion=idx)
    ## score funcs epochs ##

    # check score funcs
    for name, func in score_funcs.items():
        if name in score_funcs_unsuited:
            continue
        scores = ica.find_sources_epochs(epochs_eog, target='EOG 061',
                                         score_func=func)
        assert_true(ica.n_components_ == len(scores))

    # check univariate stats
    scores = ica.find_sources_epochs(epochs, score_func=stats.skew)

    # check exception handling
    assert_raises(ValueError, ica.find_sources_epochs, epochs,
                  target=np.arange(1))

    # ecg functionality
    ecg_scores = ica.find_sources_raw(raw, target='MEG 1531',
                                      score_func='pearsonr')

    ecg_events = ica_find_ecg_events(raw, sources[np.abs(ecg_scores).argmax()])

    assert_true(ecg_events.ndim == 2)

    # eog functionality
    eog_scores = ica.find_sources_raw(raw, target='EOG 061',
                                      score_func='pearsonr')
    eog_events = ica_find_eog_events(raw, sources[np.abs(eog_scores).argmax()])

    assert_true(eog_events.ndim == 2)

    # Test ica fiff export
    ica_raw = ica.sources_as_raw(raw, start=0, stop=100)
    assert_true(ica_raw.last_samp - ica_raw.first_samp == 100)
    ica_chans = [ch for ch in ica_raw.ch_names if 'ICA' in ch]
    assert_true(ica.n_components_ == len(ica_chans))
    test_ica_fname = op.join(op.abspath(op.curdir), 'test_ica.fif')
    ica_raw.save(test_ica_fname)
    ica_raw2 = fiff.Raw(test_ica_fname, preload=True)
    assert_array_almost_equal(ica_raw._data, ica_raw2._data)
    ica_raw2.close()
    os.remove(test_ica_fname)

    # Test ica epochs export
    ica_epochs = ica.sources_as_epochs(epochs)
    assert_true(ica_epochs.events.shape == epochs.events.shape)
    sources_epochs = ica.get_sources_epochs(epochs)
    assert_array_equal(ica_epochs.get_data(), sources_epochs)
    ica_chans = [ch for ch in ica_epochs.ch_names if 'ICA' in ch]
    assert_true(ica.n_components_ == len(ica_chans))
    assert_true(ica.n_components_ == ica_epochs.get_data().shape[1])
    assert_true(ica_epochs.raw is None)
    assert_true(ica_epochs.preload == True)

    # regression test for plot method
    assert_raises(ValueError, ica.plot_sources_raw, raw,
                  order=np.arange(50))
    assert_raises(ValueError, ica.plot_sources_epochs, epochs,
                  order=np.arange(50))