diff --git a/README.md b/README.md index cb6e64b..7a3b862 100644 --- a/README.md +++ b/README.md @@ -35,7 +35,7 @@ parameter inference techniques based on [sbi](https://sbi-dev.github.io/sbi)) wi inverse modelling pipeline. Designed for modularity, it adapts to your needs without constraints. -# Installation +# Installation on Linux `ncpi` requires Python 3.10 or higher. To install `ncpi`, you can use pip. The package is available on PyPI, so you can install it directly from there. @@ -73,6 +73,31 @@ pip uninstall scikit-learn numpy -y pip install scikit-learn==1.5.0 numpy ``` +# Installation on Windows + +To be able to install all dependencies of `ncpi` in Windows, first you have to install [Windows Subsystem for Linux (WSL)](https://documentation.ubuntu.com/wsl/stable/howto/install-ubuntu-wsl2/). + +After the WSL installation, we strongly recommend to install the latest updates by running the following commands within the Ubuntu terminal: + +```bash +$ sudo apt update +$ sudo apt upgrade -y +``` + +Once Ubuntu is up and running in WSL, [conda can be installed there](https://docs.conda.io/projects/conda/en/stable/user-guide/install/linux.html). Now you can follow the rest of the instructions of the [installation of `ncpi` in Linux](#Installation-on-Linux). + +If you encounter the error *Failed building wheel for pycatch22* when installing `pip install ncpi`, install: + +```bash +$ conda install -c conda-forge pycatch22 +``` + +If the error still persists, install the following dependencies: + +```bash +$ sudo apt install -y build-essential python3-dev +``` + # Folder Structure - `ncpi/`: Contains the source code for the library, organized into modules and classes. diff --git a/examples/EEG_AD/EEG_AD.py b/examples/EEG_AD/EEG_AD.py index 5bcabc6..a0e113f 100644 --- a/examples/EEG_AD/EEG_AD.py +++ b/examples/EEG_AD/EEG_AD.py @@ -1,245 +1,245 @@ -import os -import pickle -import time -import pandas as pd -import scipy -import numpy as np -import shutil -import ncpi -from ncpi import tools -from ncpi.tools import timer - -# Select either raw EEG data or source-reconstructed EEG data. This study used the raw EEG data for all analyses. -raw = True -if raw: - data_path = os.path.join('DATA', 'empirical_datasets', 'POCTEP_data', 'CLEAN', 'SENSORS') # Specify you data path here -else: - data_path = os.path.join('DATA', 'empirical_datasets', 'POCTEP_data', 'CLEAN', 'SOURCES', 'dSPM', 'DK') # Specify you data path here - -# Choose to either download data from Zenodo (True) or load it from a local path (False). -# Important: the zenodo downloads will take a while, so if you have already downloaded the data, set this to False and -# configure the zenodo_dir variables to point to the local paths where the data is stored. -zenodo_dw_sim = True # simulation data - -# Zenodo URL that contains the simulation data and ML models (used if zenodo_dw_sim is True) -zenodo_URL_sim = "https://zenodo.org/api/records/15351118" - -# Paths to zenodo simulation files -zenodo_dir_sim = "zenodo_sim_files" - -# Methods used to compute the features -all_methods = ['catch22','power_spectrum_parameterization_1'] - -# ML model used to compute the predictions -ML_model = 'MLPRegressor' - - -def create_POCTEP_dataframe(data_path): - ''' - Load the POCTEP dataset and create a DataFrame with the data. - - Parameters - ---------- - data_path: str - Path to the directory containing the POCTEP data files. - - Returns - ------- - df : pandas.DataFrame - DataFrame containing the POCTEP data. - ''' - - - # List files in the directory - ldir = os.listdir(data_path) - - ID = [] - group = [] - epoch = [] - sensor = [] - EEG = [] - - for pt,file in enumerate(ldir): - print(f'\rProcessing {file} - {pt + 1 }/{len(ldir)}', end="", flush=True) - - # load data - data = scipy.io.loadmat(data_path + '/' + file)['data'] - signal = data['signal'][0, 0] - - # get sampling frequency - fs = data['cfg'][0, 0]['fs'][0, 0][0, 0] - - # Electrodes (raw data)/ regions (if source data) - regions = np.arange(signal.shape[1]) - - # get channels - ch_names = data['cfg'][0, 0]['channels'][0, 0][0] - ch_names = [ch_names[ll][0] for ll in range(len(ch_names))] - - # 5-second epochs - epochs = np.arange(0, signal.shape[0], int(fs * 5)) - - for i in range(len(epochs) - 1): - ep = signal[epochs[i]:epochs[i + 1], :] - # z-score normalization - ep = (ep - np.mean(ep, axis=0)) / np.std(ep, axis=0) - - # Append data - for rg in regions: - ID.append(pt) - group.append(file.split('_')[0]) - epoch.append(i) - sensor.append(ch_names[rg]) - EEG.append(ep[:, rg]) - - # Create the Pandas DataFrame - df = pd.DataFrame({'ID': ID, - 'Group': group, - 'Epoch': epoch, - 'Sensor': sensor, - 'Data': EEG}) - df['Recording'] = 'EEG' - df['fs'] = fs - - return df - -@timer("Loading simulation data.") -def load_model_features(method, zenodo_dir_sim): - try: - with open(os.path.join(zenodo_dir_sim, 'data', method, 'sim_theta'), 'rb') as file: - theta = pickle.load(file) - with open(os.path.join(zenodo_dir_sim, 'data', method, 'sim_X'), 'rb') as file: - X = pickle.load(file) - except Exception as e: - print(f"Error loading simulation data: {e}") - - # Print info - print('theta:') - for key, value in theta.items(): - if isinstance(value, np.ndarray): - print(f'--Shape of {key}: {value.shape}') - else: - print(f'--{key}: {value}') - print(f'Shape of X: {X.shape}') - - return X, theta - -@timer("Feature extraction") -def feature_extraction(method, df): - # Parameters of the feature extraction method - if method == 'catch22': - params = None - elif method == 'power_spectrum_parameterization_1': - fooof_setup_emp = {'peak_threshold': 1., - 'min_peak_height': 0., - 'max_n_peaks': 5, - 'peak_width_limits': (10., 50.)} - params={'fs': df['fs'][0], - 'fmin': 5., - 'fmax': 45., - 'fooof_setup': fooof_setup_emp, - 'r_squared_th':0.9} - - features = ncpi.Features(method=method if method == 'catch22' else 'power_spectrum_parameterization', - params=params) - emp_data = features.compute_features(df) - - # Keep only the aperiodic exponent (1/f slope) - if method == 'power_spectrum_parameterization_1': - emp_data['Features'] = emp_data['Features'].apply(lambda x: x[1]) - - return emp_data - - -@timer('Computing predictions...') -def compute_predictions_neural_circuit(emp_data, method, ML_model, X, theta, zenodo_dir_sim): - # Add "Predictions" column to later store the parameters infered - emp_data['Predictions'] = np.nan - - # List of sensors - sensor_list = [ - 'Fp1','Fp2','F3','F4','C3','C4','P3','P4','O1', - 'O2','F7','F8','T3','T4','T5','T6','Fz','Cz','Pz' - ] - - # Create folder to save results - if not os.path.exists('data'): - os.makedirs('data') - if not os.path.exists(os.path.join('data', method)): - os.makedirs(os.path.join('data', method)) - - # Create inference object - inference = ncpi.Inference(model=ML_model) - # Not sure if this is really needed - inference.add_simulation_data(X, theta['data']) - - for s, sensor in enumerate(sensor_list): - print(f'--- Sensor: {sensor}') - - shutil.copy( - os.path.join(zenodo_dir_sim, 'ML_models/EEG', sensor, method, 'model'), - os.path.join('data', 'model.pkl') - ) - - shutil.copy( - os.path.join(zenodo_dir_sim, 'ML_models/EEG', sensor, method, 'scaler'), - os.path.join('data', 'scaler.pkl') - ) - - sensor_df = emp_data[emp_data['Sensor'].isin([sensor, s])] - - predictions = inference.predict(np.array(sensor_df['Features'].to_list())) - - sensor_df['Predictions'] = [list(pred) for pred in predictions] - emp_data.update(sensor_df['Predictions']) - - return emp_data - -def save_data(emp_data, method): - # Save the data including predictions of all parameters - emp_data.to_pickle(os.path.join('data', method, 'emp_data_all.pkl')) - - # Replace parameters of recurrent synaptic conductances with the ratio (E/I)_net - E_I_net = emp_data['Predictions'].apply(lambda x: (x[0]/x[2]) / (x[1]/x[3])) - others = emp_data['Predictions'].apply(lambda x: x[4:]) - emp_data['Predictions'] = (np.concatenate((E_I_net.values.reshape(-1,1), - np.array(others.tolist())), axis=1)).tolist() - - # Save the data including predictions of (E/I)_net - emp_data.to_pickle(os.path.join('data', method, 'emp_data_reduced.pkl')) - -if __name__ == "__main__": - # Check scikit-learn version - if not tools.ensure_module('scikit-learn', 'scikit-learn==1.3.2'): - print("Failed to install required scikit-learn version 1.3.2. Please install it manually.") - - # Download simulation data and ML models - if zenodo_dw_sim: - tools.timer("Downloading simulation data and ML models from Zenodo")( - tools.download_zenodo_record - )(zenodo_URL_sim, download_dir=zenodo_dir_sim) - - # Go through all methods - for method in all_methods: - print(f'\n\n--- Method: {method}') - - # Load parameters of the model (theta) and features (X) from simulation data - X, theta = load_model_features(method, zenodo_dir_sim) - - # Load empirical data and create DataFrame - df = create_POCTEP_dataframe(data_path=data_path) - - ########################## - # FEATURE EXTRACTION # - ########################## - - emp_data = feature_extraction(method, df) - - ####################################################### - # PREDICTIONS OF PARAMETERS OF THE NEURAL CIRCUIT # - ####################################################### - - emp_data = compute_predictions_neural_circuit(emp_data, method, ML_model, X, theta, zenodo_dir_sim) - - # Save the data including predictions of all parameters - save_data(emp_data, method) +import os +import pickle +import time +import pandas as pd +import scipy +import numpy as np +import shutil +import ncpi +from ncpi import tools +from ncpi.tools import timer + +# Select either raw EEG data or source-reconstructed EEG data. This study used the raw EEG data for all analyses. +raw = True +if raw: + data_path = os.path.join(os.sep, 'DATA', 'empirical_datasets', 'POCTEP_data', 'CLEAN', 'SENSORS') # Specify you data path here +else: + data_path = os.path.join(os.sep, 'DATA', 'empirical_datasets', 'POCTEP_data', 'CLEAN', 'SOURCES', 'dSPM', 'DK') # Specify you data path here + +# Choose to either download data from Zenodo (True) or load it from a local path (False). +# Important: the zenodo downloads will take a while, so if you have already downloaded the data, set this to False and +# configure the zenodo_dir variables to point to the local paths where the data is stored. +zenodo_dw_sim = True # simulation data + +# Zenodo URL that contains the simulation data and ML models (used if zenodo_dw_sim is True) +zenodo_URL_sim = "https://zenodo.org/api/records/15351118" + +# Paths to zenodo simulation files +zenodo_dir_sim = "zenodo_sim_files" + +# Methods used to compute the features +all_methods = ['catch22','power_spectrum_parameterization_1'] + +# ML model used to compute the predictions +ML_model = 'MLPRegressor' + + +def create_POCTEP_dataframe(data_path): + ''' + Load the POCTEP dataset and create a DataFrame with the data. + + Parameters + ---------- + data_path: str + Path to the directory containing the POCTEP data files. + + Returns + ------- + df : pandas.DataFrame + DataFrame containing the POCTEP data. + ''' + + + # List files in the directory + ldir = os.listdir(data_path) + + ID = [] + group = [] + epoch = [] + sensor = [] + EEG = [] + + for pt,file in enumerate(ldir): + print(f'\rProcessing {file} - {pt + 1 }/{len(ldir)}', end="", flush=True) + + # load data + data = scipy.io.loadmat(data_path + '/' + file)['data'] + signal = data['signal'][0, 0] + + # get sampling frequency + fs = data['cfg'][0, 0]['fs'][0, 0][0, 0] + + # Electrodes (raw data)/ regions (if source data) + regions = np.arange(signal.shape[1]) + + # get channels + ch_names = data['cfg'][0, 0]['channels'][0, 0][0] + ch_names = [ch_names[ll][0] for ll in range(len(ch_names))] + + # 5-second epochs + epochs = np.arange(0, signal.shape[0], int(fs * 5)) + + for i in range(len(epochs) - 1): + ep = signal[epochs[i]:epochs[i + 1], :] + # z-score normalization + ep = (ep - np.mean(ep, axis=0)) / np.std(ep, axis=0) + + # Append data + for rg in regions: + ID.append(pt) + group.append(file.split('_')[0]) + epoch.append(i) + sensor.append(ch_names[rg]) + EEG.append(ep[:, rg]) + + # Create the Pandas DataFrame + df = pd.DataFrame({'ID': ID, + 'Group': group, + 'Epoch': epoch, + 'Sensor': sensor, + 'Data': EEG}) + df['Recording'] = 'EEG' + df['fs'] = fs + + return df + +@timer("Loading simulation data.") +def load_model_features(method, zenodo_dir_sim): + try: + with open(os.path.join(zenodo_dir_sim, 'data', method, 'sim_theta'), 'rb') as file: + theta = pickle.load(file) + with open(os.path.join(zenodo_dir_sim, 'data', method, 'sim_X'), 'rb') as file: + X = pickle.load(file) + except Exception as e: + print(f"Error loading simulation data: {e}") + + # Print info + print('theta:') + for key, value in theta.items(): + if isinstance(value, np.ndarray): + print(f'--Shape of {key}: {value.shape}') + else: + print(f'--{key}: {value}') + print(f'Shape of X: {X.shape}') + + return X, theta + +@timer("Feature extraction") +def feature_extraction(method, df): + # Parameters of the feature extraction method + if method == 'catch22': + params = None + elif method == 'power_spectrum_parameterization_1': + fooof_setup_emp = {'peak_threshold': 1., + 'min_peak_height': 0., + 'max_n_peaks': 5, + 'peak_width_limits': (10., 50.)} + params={'fs': df['fs'][0], + 'fmin': 5., + 'fmax': 45., + 'fooof_setup': fooof_setup_emp, + 'r_squared_th':0.9} + + features = ncpi.Features(method=method if method == 'catch22' else 'power_spectrum_parameterization', + params=params) + emp_data = features.compute_features(df) + + # Keep only the aperiodic exponent (1/f slope) + if method == 'power_spectrum_parameterization_1': + emp_data['Features'] = emp_data['Features'].apply(lambda x: x[1]) + + return emp_data + + +@timer('Computing predictions...') +def compute_predictions_neural_circuit(emp_data, method, ML_model, X, theta, zenodo_dir_sim): + # Add "Predictions" column to later store the parameters infered + emp_data['Predictions'] = np.nan + + # List of sensors + sensor_list = [ + 'Fp1','Fp2','F3','F4','C3','C4','P3','P4','O1', + 'O2','F7','F8','T3','T4','T5','T6','Fz','Cz','Pz' + ] + + # Create folder to save results + if not os.path.exists('data'): + os.makedirs('data') + if not os.path.exists(os.path.join('data', method)): + os.makedirs(os.path.join('data', method)) + + # Create inference object + inference = ncpi.Inference(model=ML_model) + # Not sure if this is really needed + inference.add_simulation_data(X, theta['data']) + + for s, sensor in enumerate(sensor_list): + print(f'--- Sensor: {sensor}') + + shutil.copy( + os.path.join(zenodo_dir_sim, 'ML_models/EEG', sensor, method, 'model'), + os.path.join('data', 'model.pkl') + ) + + shutil.copy( + os.path.join(zenodo_dir_sim, 'ML_models/EEG', sensor, method, 'scaler'), + os.path.join('data', 'scaler.pkl') + ) + + sensor_df = emp_data[emp_data['Sensor'].isin([sensor, s])] + + predictions = inference.predict(np.array(sensor_df['Features'].to_list())) + + sensor_df['Predictions'] = [list(pred) for pred in predictions] + emp_data.update(sensor_df['Predictions']) + + return emp_data + +def save_data(emp_data, method): + # Save the data including predictions of all parameters + emp_data.to_pickle(os.path.join('data', method, 'emp_data_all.pkl')) + + # Replace parameters of recurrent synaptic conductances with the ratio (E/I)_net + E_I_net = emp_data['Predictions'].apply(lambda x: (x[0]/x[2]) / (x[1]/x[3])) + others = emp_data['Predictions'].apply(lambda x: x[4:]) + emp_data['Predictions'] = (np.concatenate((E_I_net.values.reshape(-1,1), + np.array(others.tolist())), axis=1)).tolist() + + # Save the data including predictions of (E/I)_net + emp_data.to_pickle(os.path.join('data', method, 'emp_data_reduced.pkl')) + +if __name__ == "__main__": + # Check scikit-learn version + if not tools.ensure_module('sklearn', 'scikit-learn==1.3.2'): + print("Failed to install required scikit-learn version 1.3.2. Please install it manually.") + + # Download simulation data and ML models + if zenodo_dw_sim: + tools.timer("Downloading simulation data and ML models from Zenodo")( + tools.download_zenodo_record + )(zenodo_URL_sim, download_dir=zenodo_dir_sim) + + # Go through all methods + for method in all_methods: + print(f'\n\n--- Method: {method}') + + # Load parameters of the model (theta) and features (X) from simulation data + X, theta = load_model_features(method, zenodo_dir_sim) + + # Load empirical data and create DataFrame + df = create_POCTEP_dataframe(data_path=data_path) + + ########################## + # FEATURE EXTRACTION # + ########################## + + emp_data = feature_extraction(method, df) + + ####################################################### + # PREDICTIONS OF PARAMETERS OF THE NEURAL CIRCUIT # + ####################################################### + + emp_data = compute_predictions_neural_circuit(emp_data, method, ML_model, X, theta, zenodo_dir_sim) + + # Save the data including predictions of all parameters + save_data(emp_data, method) diff --git a/examples/EEG_AD/figures/EEG_predictions.py b/examples/EEG_AD/figures/EEG_predictions.py index d1d5c44..25b2666 100644 --- a/examples/EEG_AD/figures/EEG_predictions.py +++ b/examples/EEG_AD/figures/EEG_predictions.py @@ -6,7 +6,7 @@ import ncpi # Set the path to the results folder -results_path = '../data' +results_path = os.path.join('..', 'data') # Select the statistical analysis method ('cohen', 'lmer') statistical_analysis = 'lmer' @@ -82,7 +82,7 @@ def append_lmer_results(lmer_results, group, elec, p_value_th, data_lmer): if row == 2 or row == 3: method = 'power_spectrum_parameterization_1' try: - data = pd.read_pickle(os.path.join('../data', method, 'emp_data_reduced.pkl')) + data = pd.read_pickle(os.path.join(results_path, method, 'emp_data_reduced.pkl')) except Exception as e: print(f'Error loading data for {method}: {e}') @@ -238,6 +238,6 @@ def append_lmer_results(lmer_results, group, elec, p_value_th, data_lmer): fig1.add_artist(tauexc_line_f) fig1.add_artist(tauinh_line_f) - fig1.savefig('EEG_predictions.png') + fig1.savefig(f'EEG_predictions_{statistical_analysis}.png') diff --git a/examples/LFP_developing_brain/LFP_developing_brain.py b/examples/LFP_developing_brain/LFP_developing_brain.py index e41e643..5dc2331 100644 --- a/examples/LFP_developing_brain/LFP_developing_brain.py +++ b/examples/LFP_developing_brain/LFP_developing_brain.py @@ -22,8 +22,8 @@ zenodo_URL_emp = "https://zenodo.org/api/records/15382047" # Paths to zenodo files -zenodo_dir_sim = "zenodo_sim_files" -zenodo_dir_emp= "zenodo_emp_files" +zenodo_dir_sim = os.path.join("zenodo_sim_files") +zenodo_dir_emp= os.path.join("zenodo_emp_files") # Methods used to compute the features all_methods = ['catch22','power_spectrum_parameterization_1'] @@ -79,18 +79,18 @@ def load_inference_data(method, X, theta, zenodo_dir_sim, ML_model): folder = 'SBI' shutil.copy( - os.path.join(zenodo_dir_sim, 'ML_models/4_param', folder, method, 'scaler'), + os.path.join(zenodo_dir_sim, 'ML_models', '4_param', folder, method, 'scaler'), os.path.join('data', 'scaler.pkl') ) shutil.copy( - os.path.join(zenodo_dir_sim, 'ML_models/4_param', folder, method, 'model'), + os.path.join(zenodo_dir_sim, 'ML_models', '4_param', folder, method, 'model'), os.path.join('data', 'model.pkl') ) if ML_model == 'NPE': shutil.copy( - os.path.join(zenodo_dir_sim, 'ML_models/4_param', folder, method, + os.path.join(zenodo_dir_sim, 'ML_models', '4_param', folder, method, 'density_estimator'), os.path.join('data', 'density_estimator.pkl') ) @@ -101,13 +101,13 @@ def load_inference_data(method, X, theta, zenodo_dir_sim, ML_model): @timer("Loading LFP data.") def load_empirical_data(zenodo_dir_emp): # Load empirical data - file_list = os.listdir(os.path.join(zenodo_dir_emp, 'development_EI_decorrelation/baseline/LFP')) + file_list = os.listdir(os.path.join(zenodo_dir_emp, 'development_EI_decorrelation', 'baseline', 'LFP')) emp_data = {'LFP': [], 'fs': [], 'age': []} for i,file_name in enumerate(file_list): print(f'\r Progress: {i+1} of {len(file_list)} files loaded', end='', flush=True) structure = scipy.io.loadmat(os.path.join(os.path.join(zenodo_dir_emp, - 'development_EI_decorrelation/baseline/LFP'), + 'development_EI_decorrelation', 'baseline', 'LFP'), file_name)) LFP = structure['LFP']['LFP'][0,0] sum_LFP = np.sum(LFP, axis=0) # sum LFP across channels diff --git a/examples/LFP_developing_brain/figures/LFP_predictions.py b/examples/LFP_developing_brain/figures/LFP_predictions.py index fadd0eb..24785df 100644 --- a/examples/LFP_developing_brain/figures/LFP_predictions.py +++ b/examples/LFP_developing_brain/figures/LFP_predictions.py @@ -7,17 +7,17 @@ import ncpi # Folder with parameters of LIF model simulations -sys.path.append(os.path.join(os.path.dirname(__file__), '../../simulation/Hagen_model/simulation/params')) +sys.path.append(os.path.join(os.path.dirname(__file__), '..', '..', 'simulation', 'Hagen_model', 'simulation', 'params')) # Path to the folder with prediction results -pred_results = '../data' +pred_results = os.path.join('..', 'data') # Calculate new firing rates (True) or load them from file if they already exist (False). If firing rates do not # exist, they will not be plotted. compute_firing_rate = False # Path to saved firing rates -fr_path = './data' +fr_path = os.path.join('.', 'data') # Number of samples to draw from the predictions for computing the firing rates n_samples = 50 @@ -118,23 +118,23 @@ LIF_params['J_ext'] = J_ext # Create a Simulation object - sim = ncpi.Simulation(param_folder='../../simulation/Hagen_model/simulation/params', - python_folder='../../simulation/Hagen_model/simulation/python', - output_folder='../../simulation/Hagen_model/simulation/output') + sim = ncpi.Simulation(param_folder = os.path.join('../../simulation/Hagen_model/simulation/params'), + python_folder = os.path.join('../../simulation/Hagen_model/simulation/python'), + output_folder = os.path.join('../../simulation/Hagen_model/simulation/output')) # Save parameters to a pickle file - with open(os.path.join('../../simulation/Hagen_model/simulation/output', 'network.pkl'), 'wb') as f: + with open(os.path.join('..', '..', 'simulation', 'Hagen_model', 'simulation', 'output', 'network.pkl'), 'wb') as f: pickle.dump(LIF_params, f) # Run the simulation sim.simulate('simulation.py', 'simulation_params.py') # Load spike times - with open(os.path.join('../../simulation/Hagen_model/simulation/output', 'times.pkl'), 'rb') as f: + with open(os.path.join('..', '..', 'simulation', 'Hagen_model', 'simulation', 'output', 'times.pkl'), 'rb') as f: times = pickle.load(f) # Load tstop - with open(os.path.join('../../simulation/Hagen_model/simulation/output', 'tstop.pkl'), 'rb') as f: + with open(os.path.join('..', '..', 'simulation', 'Hagen_model', 'simulation', 'output', 'tstop.pkl'), 'rb') as f: tstop = pickle.load(f) # Transient period @@ -155,15 +155,15 @@ if compute_firing_rate: if not os.path.exists('data'): os.makedirs('data') - with open('data/firing_rates_preds.pkl', 'wb') as f: + with open('data', 'firing_rates_preds.pkl', 'wb') as f: pickle.dump(firing_rates, f) - with open('data/IDs.pkl', 'wb') as f: + with open('data', 'IDs.pkl', 'wb') as f: pickle.dump(IDs, f) else: try: - with open(os.path.join(fr_path,'firing_rates_preds.pkl'), 'rb') as f: + with open(os.path.join(fr_path, 'firing_rates_preds.pkl'), 'rb') as f: firing_rates = pickle.load(f) - with open(os.path.join(fr_path,'IDs.pkl'), 'rb') as f: + with open(os.path.join(fr_path, 'IDs.pkl'), 'rb') as f: IDs = pickle.load(f) except FileNotFoundError: print('Firing rates not found.') @@ -383,5 +383,5 @@ ax.text(0.5, 0.49, y_labels[1], color = 'blue', alpha = 0.5, fontsize = 10, ha='center') # Save the figure -plt.savefig('LFP_predictions.png', bbox_inches='tight') +plt.savefig(f'LFP_predictions_{statistical_analysis}.png', bbox_inches='tight') # plt.show() diff --git a/examples/LFP_developing_brain/figures/emp_features.py b/examples/LFP_developing_brain/figures/emp_features.py index c004ada..21a3414 100644 --- a/examples/LFP_developing_brain/figures/emp_features.py +++ b/examples/LFP_developing_brain/figures/emp_features.py @@ -3,7 +3,7 @@ import matplotlib.pyplot as plt # Path to the folder with prediction results -pred_results = '../data' +pred_results = os.path.join('..', 'data') # Names of catch22 features try: diff --git a/examples/simulation/Hagen_model/figures/SBI_results.py b/examples/simulation/Hagen_model/figures/SBI_results.py index a3954cc..2238dae 100644 --- a/examples/simulation/Hagen_model/figures/SBI_results.py +++ b/examples/simulation/Hagen_model/figures/SBI_results.py @@ -17,7 +17,7 @@ compute_metrics = True # Path to the local directory where the metrics and posteriors will be saved -result_folder = 'SBI_results' +result_folder = os.path.join('SBI_results') # Choose whether to use a held-out dataset or folds from RepeatedKFold use_held_out_data = True @@ -34,7 +34,7 @@ zenodo_URL_sim = "https://zenodo.org/api/records/15351118" # Paths to zenodo files -zenodo_dir_sim = "/DATA/zenodo_sim_files" +zenodo_dir_sim = os.path.join(os.sep, 'DATA', 'zenodo_sim_files') # Download simulation data and ML models if zenodo_dw_sim: @@ -105,10 +105,10 @@ def create_white_to_color_cmap(color): # Path to ML models trained based on a held-out dataset approach if use_held_out_data: - ML_path = os.path.join(zenodo_dir_sim, 'ML_models/held_out_data_models') + ML_path = os.path.join(zenodo_dir_sim, 'ML_models', 'held_out_data_models') # Path to ML models trained based on a RepeatedKFold approach else: - ML_path = os.path.join(zenodo_dir_sim, 'ML_models/4_param') + ML_path = os.path.join(zenodo_dir_sim, 'ML_models', '4_param') # Limits of histograms lims = [[-15, 15], [-2, 5], [-2, 12], [0, 60]] diff --git a/examples/simulation/Hagen_model/figures/example_full_pipeline.py b/examples/simulation/Hagen_model/figures/example_full_pipeline.py index f103605..e69ea62 100644 --- a/examples/simulation/Hagen_model/figures/example_full_pipeline.py +++ b/examples/simulation/Hagen_model/figures/example_full_pipeline.py @@ -10,7 +10,7 @@ from ncpi import tools # Path to parameters of the LIF network model -sys.path.append(os.path.join(os.path.dirname(__file__), '../simulation/params')) +sys.path.append(os.path.join(os.path.dirname(__file__), '..', 'simulation', 'params')) # Choose to either download files and precomputed outputs used in simulations of the reference multicompartment neuron # network model (True) or load them from a local path (False) @@ -20,7 +20,7 @@ zenodo_URL_mult = "https://zenodo.org/api/records/15429373" # Zenodo directory where the data is stored (must be an absolute path to correctly load morphologies in NEURON) -zenodo_dir = '/DATA/multicompartment_neuron_network' +zenodo_dir = os.path.join(os.sep, 'DATA', 'multicompartment_neuron_network') # Set to True to run new simulations of the LIF network model, or False to load precomputed results from a pickle file # located in a 'data' folder. diff --git a/examples/simulation/Hagen_model/figures/sim_features.py b/examples/simulation/Hagen_model/figures/sim_features.py index 5bf9fd2..0721dbd 100644 --- a/examples/simulation/Hagen_model/figures/sim_features.py +++ b/examples/simulation/Hagen_model/figures/sim_features.py @@ -17,7 +17,7 @@ zenodo_URL_sim = "https://zenodo.org/api/records/15351118" # Paths to zenodo files -zenodo_dir_sim = "zenodo_sim_files" +zenodo_dir_sim = os.path.join("zenodo_sim_files") # Download simulation data and ML models if zenodo_dw_sim: diff --git a/examples/simulation/Hagen_model/figures/sim_predictions.py b/examples/simulation/Hagen_model/figures/sim_predictions.py index 59b2395..5f1e467 100644 --- a/examples/simulation/Hagen_model/figures/sim_predictions.py +++ b/examples/simulation/Hagen_model/figures/sim_predictions.py @@ -16,7 +16,7 @@ zenodo_URL_sim = "https://zenodo.org/api/records/15351118" # Paths to zenodo files -zenodo_dir_sim = "zenodo_sim_files" +zenodo_dir_sim = os.path.join("zenodo_sim_files") # ML model used to compute the predictions (MLPRegressor or Ridge) ML_model = 'MLPRegressor' @@ -69,10 +69,10 @@ def hellinger_distance(p, q): all_theta = {} for method in all_methods: try: - with open(os.path.join(zenodo_dir_sim,'ML_models/held_out_data_models', folder, method, + with open(os.path.join(zenodo_dir_sim,'ML_models', 'held_out_data_models', folder, method, 'predictions'), 'rb') as file: all_preds[method] = np.array(pickle.load(file)) - with open(os.path.join(zenodo_dir_sim,'ML_models/held_out_data_models', 'datasets', method, + with open(os.path.join(zenodo_dir_sim,'ML_models', 'held_out_data_models', 'datasets', method, 'held_out_dataset'), 'rb') as file: X_test, theta_test = pickle.load(file) all_theta[method] = np.array(theta_test) diff --git a/examples/simulation/Hagen_model/simulation/compute_features.py b/examples/simulation/Hagen_model/simulation/compute_features.py index 419323b..bcfb57b 100644 --- a/examples/simulation/Hagen_model/simulation/compute_features.py +++ b/examples/simulation/Hagen_model/simulation/compute_features.py @@ -6,10 +6,10 @@ import ncpi # Path to the folder containing the processed data -sim_file_path = '/DATOS/pablomc/data/Hagen_model_v1' +sim_file_path = os.path.join(os.sep, 'DATOS', 'pablomc', 'data', 'Hagen_model_v1') # Path to the folder where the features will be saved -features_path = '/DATOS/pablomc/data' +features_path = os.path.join(os.sep, 'DATOS', 'pablomc', 'data') # Set to True if features should be computed for the EEG data instead of the CDM data compute_EEG = False diff --git a/examples/simulation/Hagen_model/simulation/preprocessing.py b/examples/simulation/Hagen_model/simulation/preprocessing.py index b193268..f4d9ce0 100644 --- a/examples/simulation/Hagen_model/simulation/preprocessing.py +++ b/examples/simulation/Hagen_model/simulation/preprocessing.py @@ -5,10 +5,10 @@ from tqdm import tqdm # Path to the folders containing the simulation data -sim_file_path_pre = '/DATOS/pablomc/Hagen_model_v1' +sim_file_path_pre = os.path.join(os.sep, 'DATOS', 'pablomc', 'Hagen_model_v1') # Path to the folder containing the processed data -sim_file_path_post = '/DATOS/pablomc/data/Hagen_model_v1' +sim_file_path_post = os.path.join(os.sep, 'DATOS', 'pablomc', 'data', 'Hagen_model_v1') def process_batch(ldir): diff --git a/examples/simulation/Hagen_model/simulation/python/analysis.py b/examples/simulation/Hagen_model/simulation/python/analysis.py index 1096940..f146812 100644 --- a/examples/simulation/Hagen_model/simulation/python/analysis.py +++ b/examples/simulation/Hagen_model/simulation/python/analysis.py @@ -17,7 +17,7 @@ zenodo_URL_mult = "https://zenodo.org/api/records/15429373" # Zenodo directory where the data is stored (must be an absolute path to correctly load morphologies in neuron) -zenodo_dir = '/DATA/multicompartment_neuron_network' +zenodo_dir = os.path.join(os.sep, 'DATA', 'multicompartment_neuron_network') # Download data if zenodo_dw_mult: diff --git a/examples/simulation/Hagen_model/simulation/sim_statistics.py b/examples/simulation/Hagen_model/simulation/sim_statistics.py index 0f9b08b..458fa42 100644 --- a/examples/simulation/Hagen_model/simulation/sim_statistics.py +++ b/examples/simulation/Hagen_model/simulation/sim_statistics.py @@ -15,7 +15,7 @@ zenodo_URL_sim = "https://zenodo.org/api/records/15351118" # Paths to zenodo files -zenodo_dir_sim = "zenodo_sim_files" +zenodo_dir_sim = os.path.join("zenodo_sim_files") # Download simulation data and ML models if zenodo_dw_sim: diff --git a/examples/simulation/Hagen_model/train/RepeatedKFold.py b/examples/simulation/Hagen_model/train/RepeatedKFold.py index 27e6499..da4ecd0 100644 --- a/examples/simulation/Hagen_model/train/RepeatedKFold.py +++ b/examples/simulation/Hagen_model/train/RepeatedKFold.py @@ -6,7 +6,7 @@ import ncpi # Path to folder where simulation features are stored -sim_file_path = 'zenodo_sim_files/data' +sim_file_path = os.path.join('zenodo_sim_files', 'data') # Choose whether to use a held-out dataset or the full dataset held_out_dataset = False diff --git a/pyproject.toml b/pyproject.toml index eef7320..a0e74ec 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -52,12 +52,14 @@ utils = [ statistics = [ "rpy2", + "statsmodels" ] visualization = [ "seaborn", "mne", - "mpl_toolkits" + "mpl_toolkits", + "pdf2image" ] feats_hctsa = [