DOI: 10.5281/zenodo.15875078
title: Dataset: Hierarchical Bayesian estimation of motor-evoked potential recruitment curves yields accurate and robust estimates

description: This repository contains data used in the manuscript:

Hierarchical Bayesian estimation of motor-evoked potential recruitment curves yields accurate and robust estimates

It includes three datasets consisting of recruitment curves:

rat_scs.csv - Rat Spinal Cord Stimulation (SCS)



Intensity column: pulse_amplitude

Feature columns: participant, compound_position

Response columns (AUC): LADM, LBiceps, LDeltoid, LECR, LFCR, LTriceps


human_tms.csv — Human Transcranial Magnetic Stimulation (TMS)



Intensity column: TMSInt

Feature columns: participant, participant_condition

Response columns (peak-to-peak): PKPK_ADM, PKPK_APB, PKPK_Biceps, PKPK_ECR, PKPK_FCR, PKPK_Triceps


human_scs.csv — Human Spinal Cord Stimulation (SCS)



Intensity column: sc_current

Feature columns: participant, sc_laterality

Response columns (AUC): ADM, APB, Biceps, Triceps


Visualization

A Python script is provided below to generate recruitment curve plots for all datasets. It uses pandas, matplotlib and seaborn.

Logic:

For each dataset, we group the data by the feature columns (e.g., participant and compound_position for rat_scs.csv). For each group, we plot recruitment curves by placing the intensity on the x-axis and the response on the y-axis. 

To run the script, simply update the DATA_DIR constant on line 11 to point to the directory containing the .csv files on your system. The generated plots are saved in the same DATA_DIR directory.

import os

import pandas as pd
import matplotlib.pyplot as plt
from matplotlib.figure import Figure
from matplotlib.backends.backend_pdf import PdfPages
import seaborn as sns

# Path to directory containing the CSV files
# This should be updated to the path on your system
DATA_DIR = "/home/vishu/hbmep-data"


def make_pdf(figures: list[Figure], output_path: str):
    """
    Save a list of matplotlib figures to a multi-page PDF.

    Args:
        figures (List[Figure]): List of figures to save.
        output_path (str): Path to the output PDF file.
    """
    with PdfPages(output_path) as pdf:
        for fig in figures:
            pdf.savefig(fig, bbox_inches='tight') 
            plt.close(fig)
    print(f"Saved to {output_path}")
    return


def plot(
    data_path: str,
    intensity: str,
    features: list[str],
    response: list[str],
    output_path: str
):
    """
    Plot response vs. intensity for each unique feature combination
    and save as PDF.

    Each page of the PDF contains scatter plots arranged in a grid:
    - Rows: unique combinations of features
    - Columns: response variables

    Args:
        data_path (str): Path to the input CSV file.
        intensity (str): Column name for stimulation intensity.
        features (list[str]): Feature columns to define unique groups.
        response (list[str]): List of response variable column names.
        output_path (str): Path to save the output PDF file.
    """
    num_rows_per_page = 10
    num_columns_per_page = len(response)
    colors = sns.color_palette("tab10", num_columns_per_page)

    df = pd.read_csv(data_path)
    df_features = df[features].apply(tuple, axis=1)
    combinations = sorted(df_features.unique().tolist())
    num_combinations = len(combinations)

    num_pages = num_combinations // num_rows_per_page
    if num_combinations % num_rows_per_page:
        num_pages += 1

    figures = []
    counter = 0
    print(f"Making {output_path} ...")
    for page_idx in range(num_pages):
        print(f"Page {page_idx + 1}/{num_pages} ...")
        num_current_page_rows = min(
            num_rows_per_page,
            num_combinations - page_idx * num_rows_per_page
        )
        figsize = (8, (11.69 / num_rows_per_page) * num_current_page_rows)
        fig, axes = plt.subplots(
            *(num_current_page_rows, num_columns_per_page),
            figsize=figsize,
            squeeze=False,
            constrained_layout=True
        )

        for i in range(num_current_page_rows):
            combination = combinations[counter]
            idx = df_features.isin([combination])
            ccdf = df[idx].reset_index(drop=True).copy()

            for j, response_name in enumerate(response):
                ax = axes[i, j]
                x = ccdf[intensity]
                y = ccdf[response_name]
                color = colors[j]
                sns.scatterplot(x=x, y=y, ax=ax, s=14, color=color)

                ax.set_xlabel("")
                ax.set_ylabel("")
                if not i: ax.set_title(response_name, fontsize=8)
                if not j: ax.set_ylabel("\n".join(combination), fontsize=8)
            
            counter += 1

        for i in range(num_current_page_rows):
            for j in range(num_columns_per_page):
                ax = axes[i, j]
                ax.spines[['top', 'right']].set_visible(False)
                ax.tick_params(axis='both', labelsize=6)
        
        fig.suptitle(f"Page {page_idx + 1}/{num_pages}")
        figures.append(fig)

    make_pdf(figures, output_path)
    return


def main():
    """
    Generate plots for all datasets and save them as PDFs.

    This function loads three datasets (Rat SCS, Human TMS, Human SCS),
    specifies their intensity, features, and response variables, and calls
    `plot` to generate and save the figures.
    """
    # Rat SCS
    data_path = os.path.join(DATA_DIR, "rat_scs.csv")
    intensity = 'pulse_amplitude'
    features = ['participant', 'compound_position']
    response = ['LADM', 'LBiceps', 'LDeltoid', 'LECR', 'LFCR', 'LTriceps']
    output_path = os.path.join(DATA_DIR, "rat_scs.pdf")
    plot(data_path, intensity, features, response, output_path)

    # Human TMS
    data_path = os.path.join(DATA_DIR, "human_tms.csv")
    intensity = 'TMSInt'
    features = ['participant', 'participant_condition']
    response = ['PKPK_ADM', 'PKPK_APB', 'PKPK_Biceps', 'PKPK_ECR', 'PKPK_FCR', 'PKPK_Triceps']
    output_path = os.path.join(DATA_DIR, "human_tms.pdf")
    plot(data_path, intensity, features, response, output_path)

    # Human SCS
    data_path = os.path.join(DATA_DIR, "human_scs.csv")
    intensity = 'sc_current'
    features = ['participant', 'sc_laterality']
    response = ['ADM', 'APB', 'Biceps', 'Triceps']
    output_path = os.path.join(DATA_DIR, "human_scs.pdf")
    plot(data_path, intensity, features, response, output_path)

    return


if __name__ == "__main__":
    main()