hladnik_grewe_heterogeneity/code/plots/figure2_analysis.py

import pandas as pd
import numpy as np
from IPython import embed

def get_feature_values(df, features):
    """Get the values of different 

    Parameters
    ----------
    df : [type]
        [description]
    features : [type]
        [description]

    Returns
    -------
    [type]
        [description]
    """
    datasets = df.dataset.unique()
    df = df.dropna()
    feature_dict = {}
    for f in features:
        feature_dict[f] = []

    for d in datasets:
        contrasts = df.contrast[(df.dataset == d) & (df.pop_size == 1)].unique()
        for c in contrasts:
            trials = df[(df.dataset == d) & (df.pop_size == 1) & (df.contrast == c)]
            for f in feature_dict.keys():
                feature_dict[f].append(np.mean(trials[f].values))
    for f in feature_dict.keys():
        feature_dict[f] = np.array(feature_dict[f])
    return feature_dict, datasets

def export_baseline_data(baseline_df, phase_df):
    selection = baseline_df[["FiringRate", "CvISI", "BurstFraction", "VectorStrengthGlobalEOD", "EodPeriod"]]
    selection = selection.rename(columns={"FiringRate": "firing_rate",
                                          "CvISI": "cv",
                                          "VectorStrengthGlobalEOD": "vector_strength", 
                                          "BurstFraction": "burst_fraction",
                                          "EodPeriod": "eod_period"}
                                )
    ids = [id.split("_")[0] for id in selection.index]
    selection.insert(0, "dataset_id", ids)
    joined = pd.merge(phase_df, selection)
    joined.to_csv("../../data/Figure2_baseline_properties.csv", sep=";")
    return joined


def export_driven_data(driven_df, mi_df):
    selection = driven_df[["datasetId", "responseModulation", "responseVariability", "normPrsLowerCutoffFreq","normPrsUpperCutoffFreq", "cutoffDiff", "normPrsPeakPower", "AutoRelRfPos"]]
    selection = selection.rename(columns={
                                "datasetId": "dataset_id",
                                "responseModulation": "response_modulation",
                                "responseVariability": "response_variability",
                                "normPrsLowerCutoffFreq": "gain_cutoff_lower",
                                "normPrsUpperCutoffFreq": "gain_cutoff_upper",
                                "BurstFraction": "burst_fraction",
                                "cutoffDiff": "gain_pass_band",
                                "normPrsPeakPower": "gain_maximum",
                                "AutoRelRfPos": "receptor_pos_relative"}
                                )
    features, mi_datasets = get_feature_values(mi_df, ["mi"])
    mutual_infos = []
    for d in driven_df.datasetId:
        index = np.argwhere(mi_datasets ==  d + "-invivo-1")[0]
        if len(index) == 1:
            mutual_infos.append(features["mi"][index][0])
        else:
            embed()
            mutual_infos.extend(features["mi"][index])
    mutual_infos = np.array(mutual_infos) 
    datasets = driven_df.datasetId.values 
    df = pd.DataFrame({"dataset_id": datasets, "mi": mutual_infos})
    joined = pd.merge(selection, df)
    joined.to_csv("../../data/Figure2_driven_properties.csv", sep=";")

    return joined


if __name__ == "__main__":
    data_frame = "../../data/Fig6_7.hdf5"
    phase_data_frame = "../../data/Figure2_phase_analysis.csv"

    driven_df = pd.read_hdf(data_frame, key='/dynamic_coding_5pc')
    driven_df['cutoffDiff'] = (driven_df.normPrsUpperCutoffFreq - driven_df.normPrsLowerCutoffFreq).values
    mi_df = pd.read_csv("../../data/PopulationCodingHomogeneous.csv", index_col=0)
    
    baseline_df = pd.read_hdf(data_frame, key='/baseline')
    baseline_df['EodPeriod'] = baseline_df.GlobalEODPeakTimes.apply(lambda x: np.mean(np.diff(x)))
    phase_df = pd.read_csv(phase_data_frame, sep=";", index_col=0)

    baseline_df = export_baseline_data(baseline_df, phase_df)

    export_driven_data(driven_df, mi_df)
    embed()