Python data.Data类代码示例

def test_add_conditional_manually():
    """Enforce the conditioning on additional variables."""
    settings = {'cmi_estimator': 'JidtKraskovCMI',
                'max_lag': 5,
                'n_perm_max_stat': 21,
                'n_perm_min_stat': 21,
                'n_perm_mi': 21}
    data = Data()
    data.generate_mute_data(10, 3)
    ais = ActiveInformationStorage()

    # Add a conditional with a lag bigger than the max_lag requested above
    settings['add_conditionals'] = (8, 0)
    with pytest.raises(IndexError):
        ais.analyse_single_process(settings=settings, data=data, process=0)

    # Add valid conditionals and test if they were added
    settings['add_conditionals'] = [(0, 1), (1, 3)]
    ais._initialise(settings, data, 0)
    # Get list of conditionals after intialisation and convert absolute samples
    # back to lags for comparison.
    cond_list = ais._idx_to_lag(ais.selected_vars_full)
    assert settings['add_conditionals'][0] in cond_list, (
        'First enforced conditional is missing from results.')
    assert settings['add_conditionals'][1] in cond_list, (
        'Second enforced conditional is missing from results.')

def test_ActiveInformationStorage_init():
    """Test instance creation for ActiveInformationStorage class."""
    # Test error on missing estimator
    settings = {'max_lag': 5}
    data = Data()
    data.generate_mute_data(10, 3)
    ais = ActiveInformationStorage()
    with pytest.raises(RuntimeError):
        ais.analyse_single_process(settings, data, process=0)

    # Test tau larger than maximum lag
    settings['cmi_estimator'] = 'JidtKraskovCMI'
    settings['tau'] = 10
    with pytest.raises(RuntimeError):
        ais.analyse_single_process(settings, data, process=0)
    # Test negative tau and maximum lag
    settings['tau'] = -10
    with pytest.raises(RuntimeError):
        ais.analyse_single_process(settings, data, process=0)
    settings['tau'] = 1
    settings['max_lag'] = -5
    with pytest.raises(RuntimeError):
        ais.analyse_single_process(settings, data, process=0)

    # Invalid: process is not an int
    settings['max_lag'] = 5
    with pytest.raises(RuntimeError):  # no int
        ais.analyse_single_process(settings, data, process=1.5)
    with pytest.raises(RuntimeError):  # negative
        ais.analyse_single_process(settings, data, process=-1)
    with pytest.raises(RuntimeError):  # not in data
        ais.analyse_single_process(settings, data, process=10)
    with pytest.raises(RuntimeError):  # wrong type
        ais.analyse_single_process(settings, data, process={})

def test_bivariate_te_one_realisation_per_replication():
    """Test boundary case of one realisation per replication."""
    # Create a data set where one pattern fits into the time series exactly
    # once, this way, we get one realisation per replication for each variable.
    # This is easyer to assert/verify later. We also test data.get_realisations
    # this way.
    settings = {
        'cmi_estimator': 'JidtKraskovCMI',
        'n_perm_max_stat': 21,
        'max_lag_target': 5,
        'max_lag_sources': 5,
        'min_lag_sources': 4}
    target = 0
    data = Data(normalise=False)
    n_repl = 10
    n_procs = 2
    n_points = n_procs * (settings['max_lag_sources'] + 1) * n_repl
    data.set_data(np.arange(n_points).reshape(
                                        n_procs,
                                        settings['max_lag_sources'] + 1,
                                        n_repl), 'psr')
    nw = BivariateTE()
    nw._initialise(settings, data, 'all', target)
    assert (not nw.selected_vars_full)
    assert (not nw.selected_vars_sources)
    assert (not nw.selected_vars_target)
    assert ((nw._replication_index == np.arange(n_repl)).all())
    assert (nw._current_value == (target, max(
           settings['max_lag_sources'], settings['max_lag_target'])))
    assert (nw._current_value_realisations[:, 0] ==
            data.data[target, -1, :]).all()

def test_multivariate_te_mute():
    """Test multivariate TE estimation on the MUTE example network.

    Test data comes from a network that is used as an example in the paper on
    the MuTE toolbox (Montalto, PLOS ONE, 2014, eq. 14). The network has the
    following (non-linear) couplings:

    0 -> 1, u = 2
    0 -> 2, u = 3
    0 -> 3, u = 2 (non-linear)
    3 -> 4, u = 1
    4 -> 3, u = 1

    The maximum order of any single AR process is never higher than 2.
    """
    data = Data()
    data.generate_mute_data(n_samples=1000, n_replications=10)
    settings = {
        'cmi_estimator':  'JidtKraskovCMI',
        'max_lag_sources': 3,
        'min_lag_sources': 1,
        'max_lag_target': 3,
        'n_perm_max_stat': 21,
        'n_perm_min_stat': 21,
        'n_perm_omnibus': 21,
        'n_perm_max_seq': 21}  # this should be equal to the min stats b/c we
                               # reuse the surrogate table from the min stats

    network_analysis = MultivariateTE()
    network_analysis.analyse_network(settings, data, targets=[1, 2])

def test_visualise_multivariate_te():
    """Visualise output of multivariate TE estimation."""
    data = Data()
    data.generate_mute_data(100, 5)
    settings = {
        'cmi_estimator':  'JidtKraskovCMI',
        'max_lag_sources': 5,
        'min_lag_sources': 4,
        'n_perm_max_stat': 25,
        'n_perm_min_stat': 25,
        'n_perm_omnibus': 50,
        'n_perm_max_seq': 50,
        }
    network_analysis = MultivariateTE()
    results = network_analysis.analyse_network(settings, data,
                                               targets=[0, 1, 2])
    # generate graph plots
    visualise_graph.plot_selected_vars(results, target=1, sign_sources=False)
    plt.show()
    visualise_graph.plot_network(results, fdr=False)
    plt.show()
    visualise_graph.plot_network(results, fdr=True)
    plt.show()
    visualise_graph.plot_selected_vars(results, target=1, sign_sources=True)
    plt.show()

def test_multivariate_te_mute():
    """Test multivariate TE estimation on the MUTE example network.

    Test data comes from a network that is used as an example in the paper on
    the MuTE toolbox (Montalto, PLOS ONE, 2014, eq. 14). The network has the
    following (non-linear) couplings:

    0 -> 1, u = 2
    0 -> 2, u = 3
    0 -> 3, u = 2 (non-linear)
    3 -> 4, u = 1
    4 -> 3, u = 1

    The maximum order of any single AR process is never higher than 2.
    """
    dat = Data()
    dat.generate_mute_data(n_samples=1000, n_replications=10)
    analysis_opts = {
        'cmi_calc_name': 'jidt_kraskov',
        'n_perm_max_stat': 21,
        'n_perm_min_stat': 21,
        'n_perm_omnibus': 21,
        'n_perm_max_seq': 21,  # this should be equal to the min stats b/c we
                               # reuse the surrogate table from the min stats
        }

    network_analysis = Multivariate_te(max_lag_sources=3, min_lag_sources=1,
                                       max_lag_target=3, options=analysis_opts)
    res = network_analysis.analyse_network(dat, targets=[1, 2])

def test_calculate_mean():
    """Test if mean over CMI estimates is calculated correctly."""
    data = Data()
    data.generate_mute_data(100, 5)
    res_0 = np.load(os.path.join(os.path.dirname(__file__),
                    'data/mute_results_0.p'))
    comp_settings = {
        'cmi_estimator': 'JidtKraskovCMI',
        'n_perm_max_stat': 50,
        'n_perm_min_stat': 50,
        'n_perm_omnibus': 200,
        'n_perm_max_seq': 50,
        'tail': 'two',
        'n_perm_comp': 6,
        'alpha_comp': 0.2,
        'stats_type': 'dependent'
        }
    comp = NetworkComparison()
    comp._initialise(comp_settings)
    comp._create_union(res_0)
    cmi = comp._calculate_cmi_all_links(data)
    cmi_mean = comp._calculate_mean([cmi, cmi])
    for t in comp.union.targets_analysed:
        assert (cmi_mean[t] == cmi[t]).all(), ('Error in mean of CMI for '
                                               'target {0}'.format(t))

def test_return_local_values():
    """Test estimation of local values."""
    max_lag = 5
    settings = {
        'cmi_estimator': 'JidtKraskovCMI',
        'local_values': True,  # request calculation of local values
        'n_perm_max_stat': 21,
        'n_perm_min_stat': 21,
        'n_perm_mi': 21,
        'max_lag': max_lag,
        'tau': 1}
    data = Data()
    data.generate_mute_data(100, 3)
    ais = ActiveInformationStorage()
    processes = [1, 2]
    results = ais.analyse_network(settings, data, processes)

    for p in processes:
        lais = results.get_single_process(p, fdr=False)['ais']
        if lais is np.nan:
            continue
        assert type(lais) is np.ndarray, (
            'LAIS estimation did not return an array of values: {0}'.format(
                lais))
        assert lais.shape[0] == data.n_replications, (
            'Wrong dim (no. replications) in LAIS estimate: {0}'.format(
                lais.shape))
        assert lais.shape[1] == data.n_realisations_samples((0, max_lag)), (
            'Wrong dim (no. samples) in LAIS estimate: {0}'.format(lais.shape))

def test_multivariate_te_init():
    analysis_opts = {'cmi_calc_name': 'jidt_kraskov'}
    max_lag_target = 5
    max_lag_sources = 7
    min_lag_sources = 4
    target = 0
    sources = [2, 3, 4]
    dat = Data()
    dat.generate_mute_data(100, 5)
    nw_0 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    nw_0.analyse_single_target(dat, target, sources)

    # This should just run: Test what happens if the target max lag is bigger
    # than the source max lag
    max_lag_sources = 5
    max_lag_target = 7
    nw_1 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    nw_1.analyse_single_target(dat, target, sources)

    # The following should crash: min lag bigger than max lag
    max_lag_sources = 5
    min_lag_sources = 7
    nw_2 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    with pytest.raises(AssertionError):
        nw_2.analyse_single_target(dat, target, sources)

def test_check_source_set():
    """Test the method _check_source_set.

    This method sets the list of source processes from which candidates are
    taken for multivariate TE estimation.
    """

    dat = Data()
    dat.generate_mute_data(100, 5)
    max_lag_sources = 7
    min_lag_sources = 5
    max_lag_target = 5
    analysis_opts = {'cmi_calc_name': 'jidt_kraskov'}
    nw_0 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    sources = [1, 2, 3]
    nw_0._check_source_set(sources, dat.n_processes)

    # Assert that initialisation fails if the target is also in the source list
    sources = [0, 1, 2, 3]
    nw_0 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    nw_0.target = 0
    with pytest.raises(RuntimeError):
        nw_0._check_source_set(sources, dat.n_processes)

    sources = 1
    nw_0 = Multivariate_te(max_lag_sources, min_lag_sources, max_lag_target,
                           analysis_opts)
    nw_0._check_source_set(sources, dat.n_processes)
    assert (type(nw_0.source_set) is list)

def test_add_conditional_manually():
    """Enforce the conditioning on additional variables."""
    settings = {'cmi_estimator': 'JidtKraskovCMI',
                'max_lag_sources': 5,
                'min_lag_sources': 3,
                'max_lag_target': 7}
    nw = BivariateTE()
    data = Data()
    data.generate_mute_data()

    # Add a conditional with a lag bigger than the max_lag requested above
    settings['add_conditionals'] = (8, 0)
    with pytest.raises(IndexError):
        nw.analyse_single_target(settings=settings, data=data, target=0)

    # Add valid conditionals and test if they were added
    settings['add_conditionals'] = [(0, 1), (1, 3)]
    nw._initialise(settings=settings, data=data, target=0, sources=[1, 2])
    # Get list of conditionals after intialisation and convert absolute samples
    # back to lags for comparison.
    cond_list = nw._idx_to_lag(nw.selected_vars_full)
    assert settings['add_conditionals'][0] in cond_list, (
        'First enforced conditional is missing from results.')
    assert settings['add_conditionals'][1] in cond_list, (
        'Second enforced conditional is missing from results.')

def ft2idtxlconverter(filename, FTstructname, fileversion):
    # TODO: This will need better error handling !
    if fileversion == "v7.3":
#        try:
        print('Creating Python dictionary from FT data structure: ' + FTstructname)
        NPData = _ft_trial_2_numpyarray(filename, FTstructname)
        label = _ft_label_2_list(filename, FTstructname)
        NPfsample = _ft_fsample_2_float(filename, FTstructname)
        NPtime = _ft_time_2_numpyarray(filename, FTstructname)
        # convert data into IDTxl's Data class
        d = Data()
        # fieldtrip had "channel x timesamples" data,
        # but numpy sees the data as stored internally in the hdf5 file as:
        # "timesamples x channel"
        # we collected the replications
        # in the tirhd diemsnion --> dimension are:
        # s(amples) x p(rocesses) x r(eplications) = 'spr'
        d.set_data(NPData, 'spr')
        TXLdata = {"dataset" : d , "label" : label,
                   "time" : NPtime, "fsample" : NPfsample}

#        except(OSError, RuntimeError):
#            print('incorrect file version, the given file was not a MATLAB'
#                  ' m-file version 7.3')
#            return
    else:
        print('At present only m-files in format 7.3 are aupported,'
              'please consider reopening and resaving your m-file in that'
              'version')
    return TXLdata

def test_compare_jidt_open_cl_estimator():
    """Compare results from OpenCl and JIDT estimators for AIS calculation."""
    dat = Data()
    dat.generate_mute_data(100, 2)
    max_lag = 5
    analysis_opts = {
        'cmi_calc_name': 'opencl_kraskov',
        'n_perm_mi': 22,
        'alpha_mi': 0.05,
        'tail_mi': 'one',
        }
    processes = [2, 3]
    network_analysis = Single_process_storage(max_lag, analysis_opts, tau=1)
    res_opencl = network_analysis.analyse_network(dat, processes)
    analysis_opts['cmi_calc_name'] = 'jidt_kraskov'
    network_analysis = Single_process_storage(max_lag, analysis_opts, tau=1)
    res_jidt = network_analysis.analyse_network(dat, processes)
    # Note that I require equality up to three digits. Results become more exact for bigger
    # data sizes, but this takes too long for a unit test.
    np.testing.assert_approx_equal(res_opencl[2]['ais'], res_jidt[2]['ais'], significant=3, 
                                   err_msg='AIS results differ between OpenCl and JIDT estimator.')
    np.testing.assert_approx_equal(res_opencl[3]['ais'], res_jidt[3]['ais'], significant=3, 
                                   err_msg='AIS results differ between OpenCl and JIDT estimator.')
    print('AIS for MUTE data proc 2 - opencl: {0} and jidt: {1}'.format(res_opencl[2]['ais'], res_jidt[2]['ais']))
    print('AIS for MUTE data proc 3 - opencl: {0} and jidt: {1}'.format(res_opencl[3]['ais'], res_jidt[3]['ais']))

def test_max_statistic_sequential():
    data = Data()
    data.generate_mute_data(104, 10)
    settings = {
        'cmi_estimator': 'JidtKraskovCMI',
        'n_perm_max_stat': 21,
        'n_perm_min_stat': 21,
        'n_perm_omnibus': 21,
        'n_perm_max_seq': 21,
        'max_lag_sources': 5,
        'min_lag_sources': 1,
        'max_lag_target': 5
        }
    setup = MultivariateTE()
    setup._initialise(settings, data, sources=[0, 1], target=2)
    setup.current_value = (0, 4)
    setup.selected_vars_sources = [(1, 1), (1, 2)]
    setup.selected_vars_full = [(0, 1), (1, 1), (1, 2)]
    setup._selected_vars_realisations = np.random.rand(
                                    data.n_realisations(setup.current_value),
                                    len(setup.selected_vars_full))
    setup._current_value_realisations = np.random.rand(
                                    data.n_realisations(setup.current_value),
                                    1)
    [sign, p, te] = stats.max_statistic_sequential(analysis_setup=setup,
                                                   data=data)

def test_ais_fdr():
    settings = {'n_perm_max_seq': 1000, 'n_perm_mi': 1000}
    process_0 = {
        'selected_vars': [(0, 1), (0, 2), (0, 3)],
        'ais_pval': 0.0001,
        'ais_sign': True}
    process_1 = {
        'selected_vars': [(1, 0), (1, 1), (1, 2)],
        'ais_pval': 0.031,
        'ais_sign': True}
    process_2 = {
        'selected_vars': [],
        'ais_pval': 0.41,
        'ais_sign': False}
    res_1 = ResultsSingleProcessAnalysis(
        n_nodes=3, n_realisations=1000, normalised=True)
    res_1._add_single_result(process=0, settings=settings, results=process_0)
    res_1._add_single_result(process=1, settings=settings, results=process_1)
    res_2 = ResultsSingleProcessAnalysis(
        n_nodes=3, n_realisations=1000, normalised=True)
    res_2._add_single_result(process=2, settings=settings, results=process_2)

    settings = {
        'cmi_estimator': 'JidtKraskovCMI',
        'alpha_fdr': 0.05,
        'max_lag': 3}
    data = Data()
    data.generate_mute_data(n_samples=100, n_replications=3)
    analysis_setup = ActiveInformationStorage()
    analysis_setup._initialise(settings=settings, data=data, process=1)
    res_pruned = stats.ais_fdr(settings, res_1, res_2)
    assert (not res_pruned._single_process[2].selected_vars_sources), (
        'Process 2 has not been pruned from results.')

    alpha_fdr = res_pruned.settings.alpha_fdr
    for k in res_pruned.processes_analysed:
        if not res_pruned._single_process[k]['ais_sign']:
            assert (res_pruned._single_process[k]['ais_pval'] > alpha_fdr), (
                'P-value of non-sign. AIS is not 1.')
            assert (not res_pruned._single_process[k]['selected_vars']), (
                'List of significant past variables is not empty')
        else:
            assert (res_pruned._single_process[k]['ais_pval'] < 1), (
                'P-value of sign. AIS is not smaller 1.')
            assert (res_pruned._single_process[k]['selected_vars']), (
                'List of significant past variables is empty')

    # Test function call for single result
    res_pruned = stats.ais_fdr(settings, res_1)
    print('successful call on single result dict.')

    # Test None result for insufficient no. permutations, no FDR-corrected
    # results (the results class throws an error if no FDR-corrected results
    # exist).
    res_1.settings['n_perm_mi'] = 2
    res_2.settings['n_perm_mi'] = 2
    res_pruned = stats.ais_fdr(settings, res_1, res_2)
    with pytest.raises(RuntimeError):
        res_pruned.get_significant_processes(fdr=True)