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

Electromagnetic stimulation probes and modulates the neural systems that control movement. Key to understanding their effects is the muscle recruitment curve, which maps evoked potential size against stimulation intensity. Current methods to estimate curve parameters require large samples; however, obtaining these is often impractical due to experimental constraints. Here, we present a hierarchical Bayesian framework that accounts for small samples, handles outliers, simulates high-fidelity data, and returns a posterior distribution over curve parameters that quantify estimation uncertainty. It uses a rectified-logistic function that estimates motor threshold and outperforms conventionally used sigmoidal alternatives in predictive performance, as demonstrated through cross-validation. In simulations, our method outperforms non-hierarchical models by reducing threshold estimation error on sparse data and requires fewer participants to detect shifts in threshold compared to frequentist testing. We present two common use cases involving electrical and electromagnetic stimulation data and provide an open-source library for Python, called hbMEP, for diverse applications.