SurfaceEMG#

class myogen.simulator.SurfaceEMG(muscle_model, sampling_frequency__Hz=2048.0, mean_conduction_velocity__mm_s=4.1, sampling_points_in_t_and_z_domains=256, sampling_points_in_theta_domain=180, MUs_to_simulate=None, mean_fiber_length__mm=32, var_fiber_length__mm=3, radius_bone__mm=0, fat_thickness__mm=1, skin_thickness__mm=1, muscle_conductivity_radial__S_m=0.09, muscle_conductivity_longitudinal__S_m=0.4, fat_conductivity__S_m=0.041, skin_conductivity__S_m=1, electrode_grid_inclination_angle__deg=0, electrode_grid_dimensions__rows_cols=(8, 8), inter_electrode_distance__mm=4, electrode_radius__mm=0.75, electrode_grid_center_positions=None)[source]#

Bases: object

Surface Electromyography (sEMG) Simulation.

This class provides a biophysical simulation framework for generating surface electromyography signals from the muscle. It implements a multi-layered cylindrical volume conductor model [1] that accounts for anatomical accuracy, physiological variability, and measurement noise characteristics typical of clinical and research EMG recordings.

The simulation is built on established biophysical principles and validated anatomical data, making it suitable for algorithm development, method validation, and educational purposes in EMG signal processing and motor control research.

Parameters:
  • muscle_model (Muscle) – Pre-computed muscle model containing motor unit territories, fiber distributions, and recruitment patterns. Must be created using myogen.simulator.Muscle class. The muscle model defines the spatial organization of motor units and their associated muscle fibers within the FDI muscle volume.

  • sampling_frequency__Hz (float, default 2048.0) –

    Temporal sampling frequency for MUAP and EMG signals in Hertz. Higher frequencies provide better temporal resolution but increase computational cost and memory usage.

    Typical ranges:
    • Clinical EMG: 1024-2048 Hz.

    • Research applications: 2048-4096 Hz.

    • High-resolution studies: >4096 Hz.

  • mean_conduction_velocity__mm_s (float, default 4.1) – Mean conduction velocity of action potentials along muscle fibers in mm/s. Based on experimental measurements from Botelho et al. (2019). This parameter significantly affects MUAP duration and propagation patterns. Typical range for FDI muscle: 3.5-4.5 mm/s.

  • sampling_points_in_t_and_z_domains (int, default 256) –

    Spatial and temporal discretization resolution for numerical integration. Controls the accuracy of the volume conductor calculations but significantly impacts computational cost (scales quadratically). Higher values provide better numerical accuracy at the expense of simulation time.

    • Fast simulations: 128-256 points.

    • Standard accuracy: 256-512 points.

    • High accuracy: 512-1024 points.

  • sampling_points_in_theta_domain (int, default 180) – Angular discretization for cylindrical coordinate system in degrees. 180 points provides 2° angular resolution, which is sufficient for most surface EMG applications. Higher values provide better spatial resolution for circumferential electrode placement studies.

  • MUs_to_simulate (list[int] or None, default None) – Indices of specific motor units to simulate. If None, simulates all motor units present in the muscle model. For computational efficiency, consider simulating subsets for initial analysis or algorithm development. Motor unit indices correspond to recruitment order (0 = first recruited).

  • mean_fiber_length__mm (float, default 32) – Mean muscle fiber length in millimeters based on anatomical measurements (Jacobson et al., 1992). This parameter affects MUAP duration and amplitude characteristics. Longer fibers produce longer-duration MUAPs with potentially higher amplitudes. Typical range for FDI: 28-36 mm.

  • var_fiber_length__mm (float, default 3) – Standard deviation of fiber length distribution in millimeters. Introduces realistic biological variability in MUAP shapes within motor units. Higher values increase MUAP shape variability but maintain physiological realism. Typical range: 2-5 mm.

  • radius_bone__mm (float, default 0) – Inner bone radius for cylindrical volume conductor model in millimeters. For surface EMG applications, bone effects are typically minimal, so default value of 0 is appropriate. Non-zero values may be used for intramuscular EMG or detailed anatomical studies.

  • fat_thickness__mm (float, default 1) – Subcutaneous fat layer thickness in millimeters. This parameter significantly affects MUAP amplitude and spatial distribution due to fat’s low electrical conductivity acting as an insulator. Based on ultrasound measurements (Störchle et al., 2018). Typical range: 0.5-3 mm for hand muscles.

  • skin_thickness__mm (float, default 1) – Skin layer thickness in millimeters. Generally less critical than fat thickness but affects high-frequency components of MUAPs and electrode-tissue coupling. Based on histological measurements (Brodar, 1960). Typical range: 0.5-2 mm.

  • muscle_conductivity_radial__S_m (float, default 0.09) – Muscle tissue electrical conductivity in radial direction (perpendicular to fibers) in Siemens per meter. Muscle tissue exhibits anisotropic conductivity due to fiber orientation, with lower conductivity perpendicular to fibers. Value based on Botelho et al. (2019).

  • muscle_conductivity_longitudinal__S_m (float, default 0.4) – Muscle tissue electrical conductivity in longitudinal direction (parallel to fibers) in Siemens per meter. Approximately 4-5x higher than radial conductivity due to fiber structure and orientation. This anisotropy significantly affects MUAP spatial patterns and propagation characteristics.

  • fat_conductivity__S_m (float, default 0.041) – Subcutaneous fat tissue electrical conductivity in Siemens per meter. Low conductivity makes fat act as an electrical insulator, spatially smoothing MUAP patterns and reducing amplitudes. Based on experimental measurements from tissue impedance studies.

  • skin_conductivity__S_m (float, default 1) – Skin tissue electrical conductivity in Siemens per meter. Relatively high conductivity provides good electrical coupling between underlying tissues and surface electrodes. Value based on Roeleveld et al. (1997).

  • electrode_grid_inclination_angle__deg (float, default 0) –

    Rotation angle of electrode grid relative to muscle fiber direction in degrees.
    • 0° = electrodes aligned parallel to muscle fibers

    • 90° = electrodes aligned perpendicular to muscle fibers

    This parameter affects MUAP propagation patterns and is important for studies of fiber direction estimation and conduction velocity analysis.

  • electrode_grid_dimensions__rows_cols (tuple[int, int], default (8, 8)) – Number of electrode rows and columns in the recording grid. Total number of electrodes = rows × columns. Higher densities provide better spatial resolution but increase data volume.

  • inter_electrode_distance__mm (float, default 4) – Distance between adjacent electrodes in millimeters. Affects spatial resolution and signal crosstalk between channels. Smaller distances provide higher spatial resolution but may increase correlated noise between adjacent channels.

  • electrode_radius__mm (float, default 0.75) – Physical radius of individual circular electrodes in millimeters. Affects spatial averaging (larger electrodes average over more tissue volume) and signal-to-noise ratio (larger electrodes typically have lower impedance). Standard clinical surface electrodes: 0.5-2 mm radius.

  • electrode_grid_center_positions (list[tuple[float | int, float | int]] or None, default None) –

    List of electrode grid center positions in cylindrical coordinates.

    Each position is specified as (z_position_mm, angle_radians) where:
    • z_position_mm: position along muscle fiber direction (+ = distal, - = proximal)

    • angle_radians: circumferential position around muscle (0 = anterior)

    If None, uses single position at anatomical muscle belly center computed from the center of mass of all motor unit territories.

    Multiple positions enable spatial analysis of EMG characteristics across different muscle regions.

    Example positions:
    • [(0, 0)]: Single center position

    • [(0, 0), (-10, 0), (10, 0)]: Proximal-center-distal sequence

    • [(0, 0), (0, π/2), (0, π)]: Circumferential positions

resulting_muaps_array#

Generated MUAP templates after calling simulate_muaps(). Shape: (n_positions, n_motor_units, n_electrode_rows, n_electrode_cols, n_time_samples)

Type:

MUAP_SHAPE__TENSOR

surface_emg__tensor#

Generated surface EMG signals after calling simulate_surface_emg(). Shape: (n_positions, n_pools, n_electrode_rows, n_electrode_cols, n_time_samples)

Type:

SURFACE_EMG__TENSOR

noisy_surface_emg#

Surface EMG with added noise after calling add_noise(). Same shape as surface_emg__tensor.

Type:

SURFACE_EMG__TENSOR

Raises:
  • TypeError – If any input parameter has incorrect type (enforced by beartype decorator).

  • ValueError – If muscle_model is empty or if electrode parameters are physically unrealistic (e.g., negative dimensions, impossible conductivity values).

Parameters:
  • muscle_model (Muscle)

  • sampling_frequency__Hz (float)

  • mean_conduction_velocity__mm_s (float)

  • sampling_points_in_t_and_z_domains (int)

  • sampling_points_in_theta_domain (int)

  • MUs_to_simulate (list[int] | None)

  • mean_fiber_length__mm (float)

  • var_fiber_length__mm (float)

  • radius_bone__mm (float)

  • fat_thickness__mm (float)

  • skin_thickness__mm (float)

  • muscle_conductivity_radial__S_m (float)

  • muscle_conductivity_longitudinal__S_m (float)

  • fat_conductivity__S_m (float)

  • skin_conductivity__S_m (float)

  • electrode_grid_inclination_angle__deg (float)

  • electrode_grid_dimensions__rows_cols (tuple[int, int])

  • inter_electrode_distance__mm (float)

  • electrode_radius__mm (float)

  • electrode_grid_center_positions (list[tuple[float | int, float | int] | None] | None)

Examples

Standard Configuration:

>>> muscle = myogen.simulator.Muscle(...)
>>> surface_emg = SurfaceEMG(
...     muscle_model=muscle,
...     sampling_frequency__Hz=2048,
...     electrode_grid_dimensions__rows_cols=(8, 8),
...     MUs_to_simulate=list(range(10))  # First 10 motor units
... )

High-Resolution Configuration:

>>> surface_emg = SurfaceEMG(
...     muscle_model=muscle,
...     sampling_frequency__Hz=4096,  # Higher temporal resolution
...     sampling_points_in_t_and_z_domains=512,  # Higher spatial resolution
...     electrode_grid_dimensions__rows_cols=(16, 8),  # More electrodes
...     inter_electrode_distance__mm=2,  # Higher spatial sampling
... )

Multi-Position Spatial Analysis:

>>> positions = [
...     (0, 0),                    # Center
...     (-8, np.deg2rad(0)),      # 8mm proximal
...     (8, np.deg2rad(0)),       # 8mm distal
...     (0, np.deg2rad(90)),      # 90° circumferential
... ]
>>> surface_emg = SurfaceEMG(
...     muscle_model=muscle,
...     electrode_grid_center_positions=positions
... )

Thin Subcutaneous Layer Configuration:

>>> surface_emg = SurfaceEMG(
...     muscle_model=muscle,
...     fat_thickness__mm=0.5,    # Thin fat layer
...     skin_thickness__mm=0.8,   # Thin skin layer
... )  # Results in higher amplitude, more spatially localized MUAPs

Notes

The simulation framework is based on:
  • Volume Conductor Theory: Multi-layered cylindrical model representing muscle, fat,

and skin tissues with anisotropic conductivity properties.

  • Motor Unit Physiology: Realistic motor unit territories, fiber distributions, and

action potential propagation based on anatomical and physiological measurements.

  • Electrode Array Modeling: High-density surface electrode grids with configurable

spatial arrangements and electrode properties.

  • Signal Processing Pipeline: Complete workflow from single-fiber action potentials

to composite surface EMG signals with noise characteristics.

Key capabilities:
  • Multi-Position Electrode Arrays: Simulate EMG at multiple spatial locations

  • Motor Unit Action Potentials (MUAPs): Generate individual MUAP templates

  • Composite Surface EMG: Convolve MUAPs with spike trains for realistic EMG

  • Noise Modeling: Add measurement noise with specified signal-to-noise ratios

  • Computational Optimization: Efficient matrix reuse for large-scale simulations

  • Visualization Support: Built-in plotting capabilities for model validation

Typical workflow:
  1. Create muscle model using myogen.simulator.Muscle

  2. Initialize SurfaceEMG with desired parameters

  3. Generate MUAP templates using simulate_muaps()

  4. Create motor neuron pool using myogen.simulator.MotorNeuronPool

  5. Generate surface EMG using simulate_surface_emg()

  6. Optionally add noise using add_noise()

Performance considerations:
  • Memory usage scales as: n_positions × n_motor_units × grid_area × time_samples

  • Computational cost increases with: fiber count, electrode density, temporal resolution

  • For large simulations (>50 MUs, >64 electrodes), consider parameter optimization

  • GPU acceleration used automatically for convolution operations in surface EMG generation

Validation and quality assurance:
  • Parameters based on peer-reviewed anatomical and physiological studies

  • Built-in validation checks for parameter consistency

  • Extensive documentation with literature references

  • Type checking enforced via beartype decorators

See also

myogen.simulator.Muscle

For creating anatomically accurate muscle models

myogen.simulator.MotorNeuronPool

For generating motor neuron spike trains

myogen.utils.plotting.plot_surface_emg

For visualizing simulation results

References

A surface EMG generation model with multilayer cylindrical description of the volume conductor. IEEE Transactions on Biomedical Engineering 51, 415–426. https://doi.org/10.1109/TBME.2003.820998

Methods

__init__

add_noise

Add noise to the surface EMG signal to achieve a specified Signal-to-Noise Ratio (SNR).

simulate_muaps

Simulate MUAPs for the muscle model.

simulate_surface_emg

Generate realistic surface EMG signals by convolving MUAP templates with motor neuron spike trains.

add_noise(snr_db, noise_type='gaussian')[source]#

Add noise to the surface EMG signal to achieve a specified Signal-to-Noise Ratio (SNR).

This function calculates the appropriate noise level based on the signal power and desired SNR, then adds noise to the simulated surface EMG signals. This is useful for studying the effects of measurement noise on EMG analysis algorithms.

Parameters:
  • snr_db (float) –

    Target Signal-to-Noise Ratio in decibels (dB). Higher values mean less noise.

    Typical values for surface EMG:
    • High quality: 20-30 dB

    • Medium quality: 10-20 dB

    • Low quality: 0-10 dB

    • Very noisy: < 0 dB

  • noise_type (str, default "gaussian") –

    Type of noise to add.

    Currently supports:
    • ”gaussian”: Additive white Gaussian noise (most common for EMG).

Returns:

noisy_surface_emg – Surface EMG tensor with added noise, same shape as self.surface_emg__tensor. The original signal is preserved in self.surface_emg__tensor.

Return type:

SURFACE_EMG__TENSOR

Raises:

ValueError – If surface EMG has not been simulated yet (call simulate_surface_emg first) If noise_type is not supported

Notes

SNR Calculation: SNR_dB = 10 * log10(P_signal / P_noise)

Where:
  • P_signal = mean power of the signal across all channels

  • P_noise = power of the additive noise

The noise power is calculated as: P_noise = P_signal / (10^(SNR_dB/10))

Noise Characteristics:
  • Gaussian noise: Zero-mean with variance calculated to achieve target SNR

  • Noise is added independently to each electrode channel

  • Noise power is calculated based on the RMS power of the signal

Examples

>>> # Simulate surface EMG first
>>> surface_emg = emg_simulator.simulate_surface_emg(motor_neuron_pool)
>>>
>>> # Add moderate noise (15 dB SNR)
>>> noisy_emg = emg_simulator.add_noise(snr_db=15)
>>>
>>> # Add high noise (5 dB SNR)
>>> very_noisy_emg = emg_simulator.add_noise(snr_db=5)
>>>
>>> print(f"Original EMG shape: {surface_emg.shape}")
>>> print(f"Noisy EMG shape: {noisy_emg.shape}")
simulate_muaps(show_plots=False, verbose=True)[source]#

Simulate MUAPs for the muscle model.

Parameters:
  • show_plots (bool, optional) – Whether to show plots of the MUAPs. Default is False.

  • verbose (bool, optional) – Whether to print verbose output. Default is True.

Returns:

resulting_muaps_array – Generated MUAP templates after calling simulate_muaps(). Shape: (n_positions, n_motor_units, n_electrode_rows, n_electrode_cols, n_time_samples)

Return type:

MUAP_SHAPE__TENSOR

simulate_surface_emg(motor_neuron_pool)[source]#

Generate realistic surface EMG signals by convolving MUAP templates with motor neuron spike trains.

This method implements the final step of surface EMG simulation by combining the previously computed Motor Unit Action Potential (MUAP) templates with physiologically realistic motor neuron firing patterns. The process involves temporal resampling to match simulation timesteps, convolution of spike trains with MUAP shapes, and summation across all active motor units for each electrode channel.

The simulation can optionally use GPU acceleration (via CuPy) for efficient parallel processing of convolution operations across multiple electrode positions, motor units, and time points. If CuPy is not available, the simulation automatically falls back to CPU computation with NumPy.

Parameters:

motor_neuron_pool (MotorNeuronPool) –

Motor neuron pool containing spike trains and simulation parameters.

Must contain:
  • spike_trains: Binary spike train matrix [n_pools, n_motor_units, n_time_samples]

  • timestep__ms: Simulation timestep in milliseconds

  • active_neuron_indices: List of active motor unit indices for each pool

The motor neuron pool defines the temporal firing patterns that will be convolved with the MUAP templates to generate surface EMG signals.

Returns:

surface_emg__tensor – Surface EMG tensor with shape (n_positions, n_pools, n_electrodes, n_time_samples)

Return type:

SURFACE_EMG__TENSOR

Raises:
  • AttributeError – If simulate_muaps() has not been called first to generate MUAP templates.

  • ValueError – If motor_neuron_pool contains no active motor units or if there’s a mismatch between simulated motor units and available spike trains.

  • RuntimeError – If GPU memory is insufficient for large simulations when using CuPy (will automatically fall back to CPU-based computation).

Notes

Computational Complexity:
  • Time complexity: O(n_positions × n_pools × n_MUs × n_electrodes × n_time_samples)

  • Space complexity: O(n_positions × n_pools × n_electrodes × n_time_samples)

  • GPU acceleration (if available) provides 10-50x speedup for large simulations

Memory Requirements:
  • Input MUAPs: n_positions × n_MUs × n_electrodes × n_time_samples × 8 bytes

  • Spike trains: n_pools × n_MUs × n_time_samples × 8 bytes

  • Output EMG: n_positions × n_pools × n_electrodes × n_time_samples × 8 bytes

  • Example: 4 positions × 10 pools × 100 MUs × 64 electrodes × 10000 samples ≈ 2 GB

Signal Characteristics:
  • Amplitude range: typically 10-1000 μV for surface EMG

  • Frequency content: 10-500 Hz (depending on muscle and electrode configuration)

  • Temporal patterns reflect motor unit recruitment and firing rate modulation

  • Spatial patterns depend on motor unit territories and electrode positioning

Validation Recommendations:
  • Verify EMG amplitude scaling with muscle activation level

  • Check frequency content matches physiological expectations

  • Ensure spatial patterns are consistent with motor unit territories

  • Compare multi-channel correlations with experimental data

Examples

Basic Surface EMG Simulation:

>>> # First generate MUAP templates
>>> muaps = surface_emg.simulate_muaps()
>>>
>>> # Create motor neuron pool with spike trains
>>> motor_pool = myogen.simulator.MotorNeuronPool(
...     spike_trains=spike_data,
...     timestep__ms=1.0,
...     active_neuron_indices=[[0, 1, 2, 3, 4]]  # First 5 MUs active
... )
>>>
>>> # Generate surface EMG
>>> emg_signals = surface_emg.simulate_surface_emg(motor_pool)
>>> print(f"EMG shape: {emg_signals.shape}")  # (1, 1, 8, 8, 10000)

Multi-Condition Simulation:

>>> # Motor pool with multiple activation conditions
>>> motor_pool = myogen.simulator.MotorNeuronPool(
...     spike_trains=multi_condition_spikes,  # [n_conditions, n_MUs, n_samples]
...     timestep__ms=0.5,
...     active_neuron_indices=[
...         [0, 1, 2],           # Low activation: 3 MUs
...         [0, 1, 2, 3, 4, 5],  # Medium activation: 6 MUs
...         list(range(10))      # High activation: 10 MUs
...     ]
... )
>>>
>>> emg_multi = surface_emg.simulate_surface_emg(motor_pool)
>>> print(f"Multi-condition EMG: {emg_multi.shape}")  # (n_pos, 3, 8, 8, n_samples)

Analysis of Specific Electrodes:

>>> emg = surface_emg.simulate_surface_emg(motor_pool)
>>>
>>> # Extract signal from center electrode
>>> center_row, center_col = 4, 4  # Center of 8x8 grid
>>> center_emg = emg[0, 0, center_row, center_col, :]
>>>
>>> # Compute RMS amplitude across grid
>>> rms_amplitudes = np.sqrt(np.mean(emg[0, 0, :, :, :] ** 2, axis=2))
>>> print(f"RMS amplitude map shape: {rms_amplitudes.shape}")  # (8, 8)

Multi-Position Comparison:

>>> # Compare EMG characteristics across electrode positions
>>> emg = surface_emg.simulate_surface_emg(motor_pool)
>>>
>>> for pos_idx in range(emg.shape[0]):
...     pos_emg = emg[pos_idx, 0, :, :, :]
...     max_amplitude = np.max(np.abs(pos_emg))
...     print(f"Position {pos_idx}: Max amplitude = {max_amplitude:.2f} μV")

See also

simulate_muaps

Generate MUAP templates (must be called first)

add_noise

Add measurement noise to simulated EMG signals

myogen.simulator.MotorNeuronPool

Create motor neuron spike trains

myogen.utils.plotting.plot_surface_emg

Visualize simulated EMG signals