Pareto-Based Diagnostics and Selection for Mechanics–Synergy Trade-Offs in Unmeasured Muscle Activation Reconstruction

TL;DR

This paper introduces a method to reconstruct unmeasured muscle activations by balancing mechanical accuracy and synergy patterns, improving diagnostics in sparse measurement scenarios.

Contribution

A synergy-informed, physics-constrained framework with Pareto-based diagnostics for muscle activation reconstruction under sparse measurements.

Findings

Using J1-dominant selections from the Pareto set reduces joint-moment error from 0.154 to ≈0.138 in synthetic Arm26 cases.

Pareto diagnostics reveal identifiability and selection sensitivity in scenarios with sparse measurements and no ground truth.

Abstract

Background: Reconstructing full muscle activation trajectories from sparse measurements is underdetermined: many activation patterns can explain similar joint moments, and purely mechanical inverse formulations can yield non-physiological solutions. Methods: We propose a synergy-informed, physics-constrained framework to reconstruct unmeasured muscle activations when only a subset of muscles is observed. A synergy reconstruction prior (SynRc) is obtained by identifying a synergy basis from proxy activations via non-negative matrix factorization (NMF) and estimating time-varying synergy excitations from measured channels. Unmeasured activations are then solved via bound-constrained multi-objective optimization that jointly minimizes (i) normalized joint-moment error between OpenSim forward-computed moments and inverse-dynamics moments and (ii) deviation from the SynRc prior, with an optional smoothness refinement stage. Results: Verification on synthetic OpenSim Arm26 (2-DOF) cases with known ground truth shows that J1-dominant selections from the stage-I Pareto set reduce normalized joint-moment error from 0.154 (SynRc-only) to ≈0.138, at the cost of larger deviation from the synergy prior. These Pareto diagnostics expose identifiability and selection sensitivity under sparse measurements when ground truth is unavailable. Conclusions: The proposed framework makes mechanics–synergy trade-offs explicit and provides structured diagnostics and selection guidance for sparse-measurement scenarios.

Figures

Full text

1. Introduction

Estimating muscle activation trajectories remains difficult when experimental sensing is sparse (e.g., sparse EMG). Reconstructing a full set of muscle activations from a small number of measured channels is fundamentally underdetermined: different activation patterns can yield similar net joint moments due to muscular redundancy and the nonlinear mapping from activation to musculotendon force [1,2]. As a result, purely mechanics-based inverse formulations can admit many feasible solutions, and additional structure (physiological priors, measurement constraints, or regularization) is required to obtain interpretable estimates.

Mechanics-based approaches such as static optimization (SO) and musculoskeletal optimal control enforce joint-moment consistency through a forward model [3,4,5,6]. These methods work well when model fidelity is high, but can be sensitive to modeling and measurement errors, and when only a subset of muscles is instrumented, additional assumptions are still necessary to infer unmeasured activations. EMG-driven and hybrid formulations incorporate measured excitations to improve physiological plausibility [7,8,9], but they still require design choices about how strongly to trust EMG versus physics when measurements are sparse or mismatched.

Synergy-based extrapolation leverages low-dimensional structure to infer unmeasured channels and can be effective when synergy structure is transferable [10,11,12,13,14,15,16,17,18]. However, synergy reconstruction alone does not guarantee joint-moment consistency in a musculoskeletal forward model, and fixed scalarization weights can obscure how reconstructed activations vary with modeling assumptions.

We propose a synergy-informed, physics-constrained approach that treats joint-moment consistency and neural consistency as competing objectives. First, a synergy reconstruction prior (SynRc) is computed by identifying a synergy basis from proxy activations via non-negative matrix factorization (NMF) and estimating time-varying synergy excitations from measured channels via a bounded quadratic program with temporal regularization. Second, unmeasured activations are estimated by solving a bound-constrained multi-objective inverse problem that trades off (i) normalized joint-moment error $[eqn]$ and (ii) deviation from the SynRc prior $[eqn]$ , followed by an optional local refinement stage with a smoothness penalty. The resulting stage-I Pareto set makes the mechanics–synergy trade-off explicit and enables diagnostics and task-dependent solution selection when ground truth is unavailable.

Scope and claims: This paper does not claim universally accurate activation reconstruction across subjects and tasks, nor does it provide experimental validation or real-time deployment. Instead, it contributes a diagnostic and selection framework that exposes mechanics–synergy trade-offs via Pareto geometry and provides practical selection guidance under sparse measurements.

Intended use: The following framework targets offline analyses where interpretability and joint-moment fidelity matter under sparse observations (e.g., retrospective simulation studies or research workflows with sparse EMG) and where users need to assess identifiability and selection sensitivity rather than obtain a single black-box estimate.

The main contributions are as follows:

• A two-stage reconstruction framework combining a synergy reconstruction prior (SynRc) and a mechanics-constrained multi-objective inverse problem with optional smoothness refinement;
• A Pareto-based diagnostic and selection workflow (utopia-closest + neighborhood sensitivity) and a practical reporting template for applications without ground truth;
• Synthetic verification and stress tests illustrating how trade-off geometry changes across observability regimes (i.e., how sensitive joint moments are to unmeasured activations) and under noise/mismatch;
• A minimal reduced-dimension parameterization demo that motivates scalable implementations.

We emphasize that the present study is a verification on synthetic data with known ground truth [19]. The synthetic setting enables controlled analysis of identifiability and trade-offs but can underestimate the impact of experimental noise, modeling mismatch, and subject-specific variability; these limitations are discussed in Section 4. The remainder of this paper is organized as follows: Section 2 presents the proposed formulation and solver, Section 3 reports results and diagnostics, Section 4 discusses implications and limitations, and Section 5 concludes.

2. Methods

2.1. Nomenclature

$[eqn]$ muscle activation vector at time t (bounded in $[eqn]$ ). $[eqn]$ measured and unmeasured activations. $[eqn]$ measured-channel SynRc reconstruction error. $[eqn]$ synergy excitation at time t (bounded in $[eqn]$ ). $[eqn]$ joint-moment-consistency, neural-consistency, and smoothness objectives.Observability regime: The degree to which joint moments are informative about unmeasured activations under fixed kinematics. High observability means joint moments are sensitive to unmeasured channels; low observability means many unmeasured-activation patterns produce similar joint moments.M number of muscles. $[eqn]$ number of joint degrees of freedom (dimension of the joint-moment vector).r synergy dimension (rank) in NMF.T number of time samples. $[eqn]$ synergy weight matrix. $[eqn]$ Euclidean norm. $[eqn]$ weights for $[eqn]$ and $[eqn]$ in local refinement. $[eqn]$ inverse-dynamics and forward-computed joint moments.

Throughout, we use “activation” to denote the muscle-level control signal (the OpenSim activation state) and reserve “synergy excitation” for the low-dimensional coefficients $[eqn]$ . In the synergy/EMG literature, activation and excitation are sometimes used interchangeably after preprocessing and normalization; we adopt activation-centric notation for consistency with musculoskeletal simulation conventions. We use “joint moment” as the primary term for $[eqn]$ and refer to the mechanics objective as “joint-moment-consistency”.

2.2. Notation and Problem Setup

Let M denote the number of muscles and T the number of time samples. Let $[eqn]$ be the muscle activation vector at time t, and define the activation matrix as $[eqn]$ . Muscles are partitioned into a measured set $[eqn]$ and an unmeasured set $[eqn]$ , with sizes $[eqn]$ and $[eqn]$ , respectively. Measured activations $[eqn]$ are treated as fixed inputs; the unmeasured activations $[eqn]$ are decision variables.

Given joint kinematics $[eqn]$ and inverse-dynamics moments $[eqn]$ , the forward musculoskeletal model defines a (generally nonlinear) mapping

$$[eqn]$$

where $[eqn]$ denotes the net muscle-induced joint moments computed by OpenSim under fixed kinematics (Section 2.6). Both $[eqn]$ and $[eqn]$ are vectors in $[eqn]$ . The inverse problem is used to estimate $[eqn]$ , such that $[eqn]$ matches $[eqn]$ while remaining physiologically plausible.

2.3. Algorithm Overview

For clarity, we summarize the proposed reconstruction procedure:

1. NMF: Identify a synergy basis $[eqn]$ from proxy activations using NMF (Section 2.4).
2. SynRc: Estimate synergy excitations $[eqn]$ from the measured channels via a bounded quadratic program (QP) to obtain a synergy reconstruction prior (SynRc) for unmeasured activations (Section 2.5).
3. Stage I: Approximate a stage-I Pareto set for $[eqn]$ over the unmeasured activation time series using derivative-free multi-objective search (Section 2.8).
4. Select trade-off point: Select a stage-I point using the recommended heuristics and use-case guidance (Section 3.6).
5. Stage II: Perform local refinement using a scalarized objective with an additional smoothness penalty (Section 2.8).

Figure 1 provides a schematic overview of the two-stage workflow.

2.4. Synergy Basis Identification from Proxy Activations

We adopt a standard synergy model where activations are approximated by a low-rank non-negative factorization:

$$[eqn]$$

where $[eqn]$ is a proxy activation matrix (e.g., from a generic controller or simulated library), $[eqn]$ contains synergy weights, and $[eqn]$ contains time-varying synergy excitations. The synergy dimension r is chosen based on reconstruction quality (e.g., VAF) and interpretability. We compute $[eqn]$ using NMF [20] with non-negativity constraints.

Practical Sources of Proxy Activations

In experimental settings, $[eqn]$ is not a “ground truth” signal but a practical device to construct a task-relevant synergy basis $[eqn]$ . This practice is motivated by prior observations that synergy structure can be partially shared across healthy subjects and remain relatively stable across repetitions for similar tasks [13,14,15]. Common sources include the following: (i) static optimization (SO) activations computed from the same kinematics and inverse-dynamics moments; (ii) activations from computed muscle control or dynamic optimization on similar tasks; (iii) EMG-derived muscle excitations after preprocessing and calibration (mapped to the activation scale as needed); (iv) a synergy library pooled from prior subjects and tasks in a publicly available dataset or a lab-specific database. The present study focuses on synthetic verification and does not evaluate all proxy sources above; we list them to contextualize how $[eqn]$ may be obtained in experimental workflows. The proxy library can be subject-specific (preferred when available) or generic; in either case, $[eqn]$ should be treated as a prior with uncertainty rather than a fixed truth. In practice, SynRc provides an internal sanity check: if the measured channels cannot be reconstructed well from $[eqn]$ (high measured-channel error $[eqn]$ , Equation (5)), then the proxy synergy basis is likely mismatched to the target task/subject (as illustrated in the proxy-mismatch example discussed below), and subsequent stage-I selection should rely more heavily on joint-moment consistency (lower $[eqn]$ ) or the synergy basis should be updated. In implementation, we apply row-wise normalization of $[eqn]$ (each muscle divided by its maximum over time) to reduce amplitude bias across muscles, perform NMF in the normalized space, and then rescale reconstructions to the original magnitude. To select r, we report variance accounted for (VAF) both globally and per muscle:

$$[eqn]$$

where $[eqn]$ is the mean activation of muscle m over time; local VAF is computed analogously for each muscle (row). Unless otherwise stated, VAF is computed in the normalized space (consistent with the NMF fit). Unless otherwise stated, NMF is solved using the alternating least squares algorithm with multiple random initializations, and the best factorization (minimum reconstruction error) is retained.

2.5. Synergy Reconstruction Prior from Sparse Observations (SynRc)

Given $[eqn]$ and measured activations $[eqn]$ , we estimate synergy excitations $[eqn]$ using a bounded quadratic program. Let $[eqn]$ be the rows of $[eqn]$ corresponding to $[eqn]$ . For each time step t, we solve

$$[eqn]$$

with $[eqn]$ . The first term fits measured activations, the second provides Tikhonov regularization, and the third promotes temporal smoothness in synergy excitation trajectories. Because Equation (4) is a convex bound-constrained quadratic program at each time step, it can be solved efficiently with standard QP solvers (MATLAB R2024b quadprog).

In practice, $[eqn]$ and $[eqn]$ can be selected by a small grid search that minimizes the measured-channel reconstruction error, while $[eqn]$ controls the smoothness of $[eqn]$ and can be set based on expected measurement noise. We report the measured-channel error as a diagnostic of proxy-synergy transferability:

$$[eqn]$$

The reconstructed full activation prior is

$$[eqn]$$

and we extract the unmeasured prior $[eqn]$ for use in the multi-objective problem.

2.6. Forward Computation of Muscle-Induced Joint Moments

For a given activation trajectory and fixed kinematics, OpenSim computes musculotendon equilibrium and the associated tendon forces. For each coordinate k, the net muscle-induced joint moment can be expressed as

$$[eqn]$$

where $[eqn]$ is the musculotendon force and $[eqn]$ is the moment arm of muscle m about coordinate k. Equation (7) includes only muscle-induced contributions; if reserve actuators or other generalized forces are present, their moment contributions should be included in $[eqn]$ when computing $[eqn]$ . This computation is implemented through the OpenSim Java API by setting $[eqn]$ , assigning activations, equilibrating muscles, and accumulating moment contributions. The underlying musculotendon formulation follows standard Hill-type muscle models as commonly implemented in OpenSim [21]. In the Arm26 benchmark used here, we compute moments for shoulder elevation and elbow flexion ( $[eqn]$ ). If muscle equilibration fails at a time step (e.g., due to numerical issues), the forward moment evaluation is treated as invalid and the corresponding objective evaluation is penalized, which is important for robust derivative-free optimization.

2.7. Multi-Objective Inverse Problem for Unmeasured Activations

We parameterize the unmeasured activations over time as decision variables,

$$[eqn]$$

with bound constraints $[eqn]$ . Measured activations are held fixed to their observed values, and $[eqn]$ fills only the unmeasured channels.

Why a Second-Stage Optimization?

SynRc yields a physiologically structured prior that is consistent with measured channels and a low-dimensional synergy model, but it does not enforce that the resulting full activations reproduce inverse-dynamics joint moments through the nonlinear forward map $[eqn]$ in Equation (1). The second stage explicitly enforces joint-moment consistency and quantifies the trade-off between mechanical fidelity and neural consistency. When inverse dynamics are unreliable or only a fast prior-only estimate is needed, the SynRc-only baseline provides a lightweight reconstruction that does not enforce joint-moment consistency by design. Mechanically sensitive studies benefit from the full multi-objective stage (torque-comparison results).

We define the following objectives:

1. Joint-moment-consistency objective:

$$[eqn]$$

Normalization makes $[eqn]$ scale-invariant across tasks and facilitates comparison across Pareto solutions.Interpretation note:  $[eqn]$ compares the inverse-dynamics target $[eqn]$ , which represents net joint moments, against the muscle-induced forward moments $[eqn]$ produced by OpenSim under fixed kinematics. Therefore, even the synthetic computed muscle control (CMC) reference is not guaranteed to minimize $[eqn]$ , and local refinement can achieve a lower $[eqn]$ by fitting residual mismatch between $[eqn]$ and $[eqn]$ . We accordingly interpret $[eqn]$ as a joint-moment-consistency diagnostic rather than a guarantee of physiological correctness. 2. Neural-consistency objective (deviation from SynRc prior):

$$[eqn]$$

3. Smoothness penalty (used in local refinement):

$$[eqn]$$

2.8. Two-Stage Solver: Pareto Front Search and Local Refinement

The primary trade-off is between $[eqn]$ and $[eqn]$ , so we first approximate the Pareto front for $[eqn]$ under bound constraints using a derivative-free direct search method [22], and report the resulting non-dominated output as a stage-I Pareto set. Initial points include the SynRc prior $[eqn]$ and additional random feasible points to encourage exploration. As a reproducible default when no strong preference is available, we select the stage-I point using the utopia-closest criterion on min–max normalized objectives (Equation (12)); task-dependent alternatives are discussed in Section 3.6. Specifically, we select

$$[eqn]$$

After selecting a point from the stage-I Pareto set based on task requirements, we perform local refinement with a scalarized objective that includes smoothness:

$$[eqn]$$

where $[eqn]$ is derived from the local slope of the stage-I Pareto set in the $[eqn]$ plane and $[eqn]$ is a user-specified smoothness weight. In practice, given a selected stage-I index $[eqn]$ , we approximate

$$[eqn]$$

using neighboring stage-I points (with appropriate boundary handling). The Pareto stage and local refinement are implemented with MATLAB paretosearch and patternsearch, respectively. In our implementation, we target 30 non-dominated solutions and use a function-evaluation budget of 90,000; objective evaluations are parallelized when available. Additional solver settings are available from the corresponding author upon reasonable request (within the scope necessary to reproduce the reported results). Here, ParetoSetSize specifies a target number of non-dominated solutions to return; the resulting set is an approximation under finite evaluation budgets and mesh tolerances. Accordingly, we refer to the output as a stage-I Pareto set (an approximate trade-off set) rather than the true Pareto front, and we report selection sensitivity by inspecting neighboring points and, when refinement is used, a small $[eqn]$ sweep.

2.9. Synthetic Benchmark Design and Implementation Details

We verify the method on a synthetic OpenSim arm benchmark (Arm26) [23] with two degrees of freedom (shoulder elevation and elbow flexion) and six muscles. Ground-truth kinematics and full muscle activations are generated using OpenSim computed muscle control (CMC) [24,25] and treated as the reference solution. Inverse-dynamics moments are computed from the same kinematics using OpenSim’s inverse dynamics tool [1,2]. To emulate typical preprocessing of experimental inverse dynamics, we filter the ID moment trajectories using a zero-phase 4th-order Butterworth low-pass filter with a 5 Hz cutoff (MATLAB butter/filtfilt). The filtered ID moments are then interpolated to the activation time grid for objective evaluation. We emulate sparse sensing by treating four muscles as measured—triceps long head (TRIlong), triceps medial head (TRImed), biceps long head (BIClong), and brachialis (BRA)—and two muscles as unmeasured—triceps lateral head (TRIlat) and biceps short head (BICshort). The measured activation channels are taken from the CMC reference and treated as the “observed” signals in this synthetic verification.

Synergies are identified via NMF with $[eqn]$ from proxy activations. For SynRc, we set $[eqn]$ for temporal regularization and selected $[eqn]$ and $[eqn]$ using a small grid search that minimizes measured-channel reconstruction error. In this verification study, the proxy activations are derived from the same synthetic pipeline, enabling controlled analysis with known ground truth while also representing an optimistic scenario in which synergy structure is well matched to the target task. The full dataset contains $[eqn]$ time samples spanning $[eqn]$ s. All computations are performed in MATLAB (R2024b) with the OpenSim Java API (v4.x). The OpenSim model files used in this study are publicly available from the OpenSim official website. The scripts and derived numerical results required to reproduce the benchmark are available from the corresponding author upon reasonable request (within the scope necessary for academic use and reproduction of the reported results). Solver settings and hyperparameters are summarized in Table 1.

2.9.1. Case 2: Alternative Measured Set(Setup)

We repeat the Arm26 benchmark under an alternative measured set while retaining the same model and motion. Specifically, we treat TRIlong, TRIlat, BIClong, and BICshort as measured and reconstruct TRImed and BRA as unmeasured. Results are reported in Section 3.2.

2.9.2. Case 3: Different Model and Task(Setup) (gait10dof18musc Walking)

To demonstrate portability beyond Arm26, we consider a different musculoskeletal model and motion. We use OpenSim’s distribution test asset gait10dof18musc_subject01 with a walking inverse-kinematics trajectory. Because this minimal asset does not include a full CMC/ID pipeline, we construct a fully synthetic benchmark with known ground truth. We generate smooth ground-truth activations for a right-leg muscle subset (hamstrings_r, rect_fem_r, vasti_r, gastroc_r, soleus_r, tib_ant_r). We treat rect_fem_r, vasti_r, soleus_r, and tib_ant_r as measured, and reconstruct hamstrings_r and gastroc_r as unmeasured. To form the mechanical objective, we precompute OpenSim moment arms for these muscles about hip flexion, knee angle, and ankle angle along the walking kinematics, and use a lightweight forward moment map $[eqn]$ (for each coordinate k), with synthetic targets $[eqn]$ generated from the ground truth. This Case 3 forward model is intentionally a linearized approximation: it uses OpenSim geometry (time-varying moment arms $[eqn]$ ) but replaces the full activation-to-force dynamics with a static proportional model using $[eqn]$ , thereby ignoring force–length/velocity effects, tendon compliance, passive forces, and activation dynamics. Accordingly, Case 3 is intended as a minimal cross-model demonstration that the Pareto/selection machinery remains meaningful under different musculoskeletal geometry, rather than as a quantitative gait prediction study. To avoid a degenerate trade-off in this small muscle subset, we additionally scale the maximum isometric force of the unmeasured muscles by a factor of 2 (a stress-test analogous to the observability scaling in Section 3.5, not a calibrated subject-specific parameter). Because $[eqn]$ is generated from the same forward map, this case does not probe ID–forward mismatch; ID–forward mismatch diagnostics are instead characterized by the Arm26 robustness tests in Section 3.4. Synergies are identified from separate proxy activations, and the SynRc prior is computed from the measured channels as in the Arm26 benchmark. Results are reported in Section 3.3.

2.10. Evaluation Metrics and Baselines

Because the benchmark is synthetic, ground-truth activations are known for all muscles. We evaluate reconstruction accuracy for each unmeasured muscle using the time-domain root-mean-square error (RMSE) between the estimate and the ground truth. To characterize mechanical consistency, we also report the normalized joint-moment error $[eqn]$ (Section 2.7) and visualize the stage-I Pareto set over $[eqn]$ , where $[eqn]$ measures deviation from the SynRc prior. We compare the following:

• SO baseline: OpenSim static optimization activations [1,3].
• SynRc-only baseline: The SynRc prior for the unmeasured activations (no joint-moment-consistency optimization).
• Joint-moment-only + smoothness baseline: A joint-moment-tracking reconstruction without a synergy prior, minimizing $[eqn]$ under a reduced-dimension parameterization.
• Proposed method: Two-stage multi-objective optimization with different stage-I selections and smoothness weights $[eqn]$ .

These baselines are not intended as a comprehensive comparison with state-of-the-art activation estimation pipelines. Rather, they isolate the effects of a synergy-based prior (SynRc-only) and mechanics-based fitting (SO and joint-moment-only) so that the paper’s main claim—making mechanics–synergy trade-offs explicit and diagnosable via Pareto geometry and selection sensitivity—can be evaluated without relying on a claim of superior reconstruction accuracy.

3. Results

3.1. Synergy Identification and SynRc Prior Quality

For $[eqn]$ , NMF yields a high global reconstruction quality (total VAF $[eqn]$ ) while preserving interpretability; local VAF remains high for most muscles but is lower for BRA (54.96%), indicating that this muscle is less well represented by the chosen synergy dimension. Figure 2 illustrates NMF reconstruction versus proxy activations. SynRc reconstructs TRIlat with very low RMSE (≈ $[eqn]$ ) but shows larger errors for BICshort (RMSE $[eqn]$ ), consistent with limited observability from the measured subset (Figure 3).

Robustness to Proxy-Synergy Mismatch

The preceding analysis uses proxy synergies derived from the same synthetic pipeline as the target motion (an optimistic setting). To emulate a mismatch between the proxy synergy library and the target task, we identify the synergy basis from static-optimization (SO) activations and then recompute the SynRc prior using the same measured channels from the CMC reference. Figure 4 shows that a mismatch in the proxy basis can substantially degrade prior accuracy for the unmeasured muscles (most strongly for BICshort), even though the measured channels remain well reconstructed. Quantitatively, the measured-channel RMSE increases from ≈ $[eqn]$ (matched proxy) to ≈ $[eqn]$ (SO proxy), while unmeasured-muscle RMSE increases from ≈ $[eqn]$ to ≈ $[eqn]$ for TRIlat and from ≈ $[eqn]$ to ≈ $[eqn]$ for BICshort. This suggests a practical diagnostic in applications without ground truth: the measured-channel reconstruction error of SynRc provides an internal check on whether the proxy synergy basis is transferable, and can inform how strongly $[eqn]$ should be trusted relative to $[eqn]$ during stage-I selection.

3.2. Case 2: Alternative Measured Set

We next evaluate the Case 2 configuration with an alternative measured set (Section 2.9.1). For this configuration, the SynRc diagnostic remains low ( $[eqn]$ ), and the stage-I trade-off set collapses to a narrow band, indicating that the two objectives are largely aligned for this motion and measured set. Figure 5 and Table 2 show that the utopia-closest selection coincides with the SynRc prior and that stage-II refinement produces only negligible changes, illustrating how the framework can diagnose when additional joint-moment-driven optimization provides little benefit beyond the synergy-based reconstruction.

3.3. Case 3: Different Model and Task (gait10dof18musc Walking)

Case 3 evaluates the method on a different model and task (Section 2.9.2). Figure 6 shows a nontrivial $[eqn]$ geometry under this different model and task, and Table 3 summarizes representative solutions and unmeasured-muscle RMSE against the known ground truth.

3.4. Synthetic Robustness Tests Under Noise and Mismatch

The preceding synthetic benchmark is optimistic because the measured activations, proxy synergies, and inverse-dynamics moments are generated within a consistent OpenSim pipeline. To better characterize when the proposed selection rules are reliable (and when they can fail), we perform four synthetic robustness tests by introducing controlled perturbations to (i) inverse-dynamics (ID) moments, (ii) measured activation channels, (iii) proxy synergy mismatch, and (iv) forward-model parameter mismatch. Because re-running the full Pareto search for every perturbation is computationally expensive, we adopt a post-hoc protocol: we keep the baseline stage-I Pareto decision set $[eqn]$ fixed (returned by paretosearch, $[eqn]$ ) and re-evaluate objective values under each perturbation to assess the stability of the $[eqn]$ trade-off geometry and the resulting point selection. For a deterministic representative choice, we apply the same utopia-closest rule (Equation (12)) on the perturbed objectives and report the selected index $[eqn]$ . This protocol evaluates selection sensitivity within a fixed candidate set and should be interpreted as a robustness diagnostic rather than a re-computation of the true Pareto front under perturbation.

Supplementary Figure S1 summarizes the results. ID perturbations include additive Gaussian noise (5% of the root mean square (RMS) of each degree-of-freedom moment) and a small phase bias (5 ms) applied to the ID target. Measured-channel perturbations include per-channel gain error ( $[eqn]$ ) and additive noise (standard deviation 0.02 in activation units, followed by clipping to $[eqn]$ ); this directly increases the measured-channel SynRc reconstruction error $[eqn]$ and can degrade the proxy-based prior. Finally, proxy-synergy mismatch is graded by mixing the synergy basis with a random basis (mismatch level $[eqn]$ ), which monotonically increases both $[eqn]$ and the unmeasured-muscle prior RMSE. To emulate a common source of musculoskeletal model uncertainty, we also perturb forward-model muscle strength by scaling each muscle’s maximum isometric force by a random factor in $[eqn]$ and recompute forward moments under the perturbed model. Perturbation magnitudes are chosen to be moderate on the normalized scales of this benchmark and are intended as stress-tests rather than as calibrated experimental noise models. Supplementary Table S1 reports the corresponding objective values and unmeasured-muscle RMSE at the selected point. To account for stochastic variability, we repeat the ID-noise and proxy-mismatch perturbations over 20 random seeds, and we repeat the measured-channel perturbation and the model-parameter mismatch over 5 random seeds. Supplementary Table S2 reports mean ± SD statistics at the selected point and the empirical range of the selected index $[eqn]$ . Overall, these tests highlight that (a) even modest ID bias can substantially alter $[eqn]$ (and thus, the selected point), and (b) increases in $[eqn]$ track degradation of the synergy prior under measured-signal noise and proxy mismatch, motivating a shift toward $[eqn]$ -dominant selections when $[eqn]$ is high.

3.5. Multi-Case Synthetic Study: Observability Regimes

Method behavior depends on the degree to which the mechanical objective is informative about the unmeasured channels. To provide a small set of representative synthetic cases (beyond a single benchmark), we construct three “observability regimes” by scaling the maximum isometric force of the unmeasured muscles (TRIlat and BICshort) while keeping the measured channels and kinematics fixed. Intuitively, when unmeasured muscles are weak, net joint moments become less sensitive to their activations, producing a flatter trade-off geometry; when unmeasured muscles are strong, the sensitivity increases and the Pareto set becomes more structured. For each regime, we keep the baseline stage-I decision set fixed and re-evaluate $[eqn]$ under the modified forward model, while $[eqn]$ (deviation from the SynRc prior) is unchanged by construction. Figure 7 shows that the resulting $[eqn]$ geometry and the utopia-closest selection can change substantially across regimes, highlighting that a narrow $[eqn]$ spread is not a solver failure but a diagnostic of limited mechanical observability under sparse measurements. Table 4 reports the selected index and the empirical $[eqn]$ range within the candidate set for each regime.

3.6. Pareto Front and Task-Dependent Solution Selection

The multi-objective formulation yields a Pareto front over $[eqn]$ . How to read the Pareto plot: each point is a candidate reconstruction in the stage-I set. Moving toward lower $[eqn]$ improves joint-moment consistency (mechanically dominant), whereas moving toward lower $[eqn]$ increases neural consistency by staying closer to the synergy prior (neurally dominant). In applications, the stage-I point can be selected based on task-specific confidence in ID moments versus synergy priors. We report representative cases where the stage-I point and the smoothness weight $[eqn]$ are varied to illustrate the trade-offs. For this synthetic dataset, the stage-I Pareto set returned by paretosearch spans a narrow range of joint-moment errors ( $[eqn]$ to $[eqn]$ ) but a broader range of neural deviations ( $[eqn]$ to $[eqn]$ ). Within this explored stage-I set, the dominant trade-off is thus in closeness to the SynRc prior while producing only modest changes in net joint moments. This narrow $[eqn]$ spread reflects the redundancy of the benchmark: for this motion and measured set, the net joint moments are relatively insensitive to variations in the unmeasured channels within the explored region, so many distinct activation patterns yield similar moments. We interpret stage-I point 1 as $[eqn]$ -dominant (minimum $[eqn]$ on the stage-I set, largest $[eqn]$ ), stage-I point 30 as neural-dominant (minimum $[eqn]$ ), and stage-I point 18 as an intermediate stage-I choice used in the subsequent reconstruction tests. Because the Pareto stage is an approximation under finite budgets and mesh tolerances, and because stage-II refinement minimizes a different scalarized objective that also includes $[eqn]$ , the refined solution is not guaranteed to remain on the stage-I Pareto set in $[eqn]$ . In particular, $[eqn]$ -dominant selections can move away from the stage-I Pareto set while reducing $[eqn]$ beyond the minimum observed in stage I (Table 5).

3.6.1. Neighborhood Sensitivity

To quantify local sensitivity around the utopia-closest selection, Supplementary Table S4 reports objectives and unmeasured-muscle RMSE for $[eqn]$ and a small index neighborhood ( $[eqn]$ ) on the stage-I set. In this Arm26 benchmark, the neighborhood shows small $[eqn]$ variation but non-negligible changes in smoothness ( $[eqn]$ ) and BICshort RMSE, so reporting $[eqn]$ sensitivity provides a simple robustness check when ground truth is unavailable.

3.6.2. Scalarization Baseline

A common alternative to explicit Pareto analysis is to scalarize objectives with a fixed weighted sum. To provide a simple baseline, we apply a weighted-sum rule on the stage-I Pareto set: $[eqn]$ , where $[eqn]$ is a weight ratio and $[eqn]$ are min–max normalized objectives over the stage-I set (same normalization as Equation (12)). Figure 8 shows that, even after normalization, the selected index can be piecewise constant over wide ranges of $[eqn]$ and can jump abruptly as $[eqn]$ changes, effectively selecting only a few discrete solutions. This sensitivity makes weight selection difficult to justify a priori and motivates reporting the full trade-off set and neighborhood sensitivity rather than committing to a single fixed scalarization.

Figure 9 shows the stage-I $[eqn]$ set, the utopia-closest point (star; Equation (12)), and representative refined solutions (cross).

3.7. Unmeasured Activation Reconstruction Results

Figure 5 and Figure 6 illustrate two representative trade-off geometries (near-degenerate versus non-degenerate) to clarify when the proposed Pareto machinery is informative; Figure 9 reports the main-case trade-off and refined off-Pareto behavior.

Table 6 summarizes RMSE for unmeasured muscles under SO, SynRc-only (prior-only; no joint-moment-consistency optimization), a joint-moment-only baseline with smoothness regularization, and the proposed method with four configurations (stage-I point and $[eqn]$ ). Relative to SO, the synergy-based prior (SynRc-only) reduces BICshort RMSE by ≈ $[eqn]$ ; the multi-objective stage primarily improves joint-moment tracking and influences TRIlat depending on the stage-I point and $[eqn]$ . The joint-moment-only baseline improves joint-moment tracking compared with SynRc-only (Table 5) but yields TRIlat reconstruction that is substantially less accurate than the synergy-driven prior, illustrating the role of $[eqn]$ in preserving neural-consistent structure when joint-moment information alone is insufficient to uniquely determine unmeasured activations. In this benchmark, BICshort remains close to the SynRc prior across configurations, whereas TRIlat is more sensitive to the stage-I selection; Supplementary Figure S3 provides a brief sensitivity check for the refinement smoothness weight $[eqn]$ . To connect activation reconstruction to the multi-objective formulation, Table 5 reports the corresponding joint-moment error $[eqn]$ , prior deviation $[eqn]$ , and smoothness $[eqn]$ for the synthetic reference (CMC), the SynRc-only baseline (direct synergy reconstruction without joint-moment-consistency optimization), a joint-moment-only smoothness-regularized baseline, and the optimized solutions (Opt). $[eqn]$ -dominant selections (stage-I point 1) reduce $[eqn]$ from 0.154 (SynRc-only) and 0.141 (CMC) to ≈ $[eqn]$ , but increase deviation from the SynRc prior ( $[eqn]$ to $[eqn]$ ) and can degrade TRIlat accuracy. Conversely, the more neural-consistent selection (stage-I point 18) keeps $[eqn]$ close to zero and yields near-ground-truth TRIlat reconstruction, at the expense of a larger $[eqn]$ comparable to the synergy baseline. Notably, even the synthetic CMC reference has nonzero $[eqn]$ because inverse-dynamics moments represent net joint moments and are not, in general, identical to muscle-induced moments computed by the forward model under fixed kinematics. Figure 10 visualizes the resulting shoulder and elbow moments across representative solutions, highlighting how $[eqn]$ -dominant refinement primarily alters the elbow-flexion moment near the terminal phase to reduce mismatch with the inverse-dynamics target. Supplementary Figure S2 further compares reconstructed unmeasured activations across the four representative optimized configurations, illustrating how stage-I selection and $[eqn]$ jointly shape TRIlat and BICshort trajectories relative to SO and ground truth.

4. Discussion

The results highlight several practical implications of synergy-informed, joint-moment-consistent reconstruction under sparse measurements. In the Arm26 benchmark, optimized solutions reduce $[eqn]$ from 0.154 (SynRc-only) to ≈ $[eqn]$ . More neural-consistent selections achieve TRIlat RMSE $[eqn]$ (Table 5 and Table 6). First, synergy reconstruction can be highly accurate for some muscles (TRIlat in this benchmark). However, the synergy prior alone does not guarantee joint-moment consistency. Enforcing joint-moment matching through the forward model reconciles neural priors with physics and makes the trade-off explicit via the stage-I Pareto set. Second, solution quality is sensitive to stage-I point selection. A $[eqn]$ -dominant point can degrade agreement with ground-truth activations even when the synergy prior is accurate, whereas a more neural-consistent point can recover near-ground-truth activations. This supports task-dependent selection rules based on relative confidence in inverse dynamics moments, model fidelity, and the transferability of synergy structure. Third, the smoothness weight $[eqn]$ acts as a user-specified sensitivity parameter in local refinement. Increasing $[eqn]$ reduces temporal fluctuations but can oversmooth and degrade accuracy. We therefore recommend reporting a small $[eqn]$ sweep (Supplementary Figure S3) and choosing $[eqn]$ in relation to sampling density and expected measurement noise, rather than treating it as a single fixed choice. Finally, the observability-regime cases show that the stage-I geometry can diagnose identifiability under sparse measurements. A very narrow $[eqn]$ spread indicates that joint-moment information is weakly informative about unmeasured channels within the explored region. In such regimes, selection should rely more on prior confidence and measured-channel diagnostics than on small differences in $[eqn]$ .

4.1. Selection Workflow Without Ground Truth

For applications without ground truth, we suggest reporting the following checklist to document diagnostics and selection choices. This checklist is intended as a reporting aid rather than a validated default pipeline:

1. Report the stage-I Pareto set $[eqn]$ and the SynRc diagnostic $[eqn]$ (Equation (5)) as a proxy-transfer indicator;
2. Report ID-target quality diagnostics (e.g., high-frequency ratio, spike score, estimated lag under synchronization checks, and the minimum achievable $[eqn]$ within the stage-I set) to contextualize how strongly $[eqn]$ should be trusted;
3. State the selection rationale (e.g., endpoint toward low $[eqn]$ or low $[eqn]$ , or a compromise rule such as utopia-closest) and report the selected index $[eqn]$ and $[eqn]$ ;
4. Report sensitivity by evaluating a small neighborhood on the stage-I set and, when refinement is used, a small $[eqn]$ sweep.

To make this procedure reproducible and less subjective, we recommend reporting a small, standardized “protocol summary” in every experiment (Table 7). Table 7 is designed as a fill-in template. Report each diagnostic, state the selection rule, and document neighborhood/ $[eqn]$ sensitivity. An objective default for the compromise option in Step 3 is the utopia-closest stage-I point using min–max normalized objectives (Equation (12)). When the empirical range of one objective is very narrow (e.g., $[eqn]$ in highly redundant or weakly observable regimes), min–max normalization can amplify small numerical differences. We therefore recommend interpreting $[eqn]$ together with neighborhood sensitivity (Supplementary Table S4) rather than as a definitive optimum.

4.2. Relation to Prior Work

Synergy extrapolation and related methods infer unmeasured activations from measured channels by exploiting low-dimensional structure [16,17,18,26]. Our formulation is complementary to these approaches. We treat the synergy-derived estimate (SynRc) as an explicit prior. Rather than committing to a single synergy-driven reconstruction, we place this prior in competition with the joint-moment-consistency objective $[eqn]$ computed from the musculoskeletal forward model.

This separation makes the trade-off transparent (via $[eqn]$ ) and avoids selecting fixed scalarization weights a priori. Compared with EMG-driven and hybrid approaches [7,8,9], the present study focuses on sparse-activation reconstruction under fixed kinematics. We do not model excitation-to-activation dynamics or perform EMG-driven parameter calibration. Integrating Pareto-based selection with EMG-driven pipelines (where measured channels are EMG-derived excitations) is a natural extension.

From an experimental perspective, EMG-based validation itself involves a spatial-resolution trade-off: surface EMG is noninvasive and practical but can be limited for deep or closely spaced muscles, whereas needle or fine-wire EMG provides higher spatial specificity at the cost of participant discomfort and constraints on natural movement. This trade-off further motivates sparse-measurement scenarios in which only a subset of muscles can be reliably instrumented.

4.3. Limitations and Future Directions

This study is a verification on synthetic data with known ground truth [19]. Ground-truth activations, kinematics, and inverse dynamics are all generated within a consistent OpenSim pipeline, so the benchmark can underestimate experimental uncertainty—particularly noise, soft-tissue artifacts, force-plate errors, and musculoskeletal parameter uncertainty. We partially address these concerns through synthetic stress-tests that perturb ID targets, measured-channel noise, proxy-synergy structure, and model parameters (Section 3.4), though these perturbations are not intended to fully replicate experimental conditions.

A further caveat concerns the inverse-dynamics targets themselves. The moments used as $[eqn]$ targets are net joint moments, which are not purely muscle-induced; even the CMC reference exhibits residual mismatch. $[eqn]$ -dominant optimization can exploit this mismatch by introducing sharp activation changes near motion endpoints. Mitigating such effects will require additional physiological constraints (e.g., activation dynamics and bounds on activation rate) and validation under calibrated experimental noise and modeling errors.

The Arm26 model is intentionally small (2 degrees of freedom, 6 muscles) and should be viewed as a controlled proof-of-concept. Establishing practical utility will require evaluation on larger models with stronger redundancy (e.g., lower-limb gait models with dozens of muscles). Case 3 provides a cross-model demonstration. However, it uses a linearized forward-moment map and a deliberate unmeasured $[eqn]$ scaling stress-test, and it generates $[eqn]$ from the same simplified map. Therefore, Case 3 should be interpreted as evidence that the Pareto/selection framework transfers to different musculoskeletal geometry, not as a validation of full gait dynamics or experimental inverse dynamics. We also provide a sensitivity analysis using SO-derived synergies. Broader mismatch across tasks/subjects remains important future work.

Finally, this paper does not include a head-to-head numerical comparison with alternative approaches (e.g., EMG-driven synergy extrapolation pipelines or full optimal-control formulations) because the scope is limited to synthetic verification. Such comparisons are an important next step once appropriate experimental datasets are available.

4.4. Computational Considerations

Computational cost is a major practical consideration because each objective evaluation requires OpenSim muscle equilibration across T time samples. Accordingly, the present implementation is intended for offline analyses where interpretability and mechanical fidelity justify longer runtimes (e.g., high-fidelity postoperative evaluation, retrospective simulation studies, or research workflows with sparse EMG). It is not intended for real-time control or closed-loop applications. For the Arm26 benchmark ( $[eqn]$ ), we use MATLAB paretosearch with ParetoSetSize=30 and a fixed function-evaluation budget. We use parallel objective evaluations when available. Local refinement is performed with patternsearch after Pareto-point selection. We report budgets for reproducibility because runtime depends strongly on hardware, OpenSim settings, and parallelization. Full solver settings are available from the corresponding author upon reasonable request (within the scope necessary to reproduce the reported results). As a concrete reference, on a workstation (Intel i7-10700K CPU (Intel, Santa Clara, CA, USA); MATLAB R2024b), the Pareto stage with a 90,000-evaluation budget took 201 min wall-clock using parallel objective evaluations (8 workers). This corresponds to an average OpenSim cost of ≈ $[eqn]$ s per objective evaluation (per worker), or ≈2 ms per time sample per evaluation in this Arm26 setting. Scaling beyond Arm26 is limited by both decision dimension and forward-model cost. Optimizing unmeasured activations at each sample yields $[eqn]$ decision variables. Wall-clock time grows roughly with $[eqn]$ (amortized across workers), where $[eqn]$ is the OpenSim forward cost per evaluation. If $[eqn]$ is approximately linear in both muscle count and time samples, per-evaluation cost scales roughly with $[eqn]$ . A lower-limb gait model with $[eqn]$ muscles and $[eqn]$ samples could therefore increase per-evaluation cost by ∼ $[eqn]$ . Naive 90,000-evaluation derivative-free searches would then be impractical. A practical path to scaling is therefore to reduce both (i) the decision dimension and (ii) the required number of evaluations. Practical options include warm-starting from SynRc/SO, optimizing directly in a reduced basis (synergy coefficients or spline knots), and using multi-fidelity/surrogate strategies to limit expensive OpenSim calls. As a small illustration, the joint-moment-only smoothness baseline (Table 5) is solved in a 20-dimensional knot parameterization with a 3000-evaluation budget in ≈ $[eqn]$ min on the same workstation, showing that reduced-dimension parameterizations can make time-coupled refinement problems substantially cheaper.

4.5. Minimal Extensibility Demo: Reduced-Dimension Parameterization

To illustrate that the framework is compatible with lower-dimensional decision variables (without switching to a larger musculoskeletal model), we construct a simple reduced-dimension parameterization for the unmeasured activations. Specifically, for each unmeasured muscle we represent $[eqn]$ by a small number of knot values and reconstruct the full time series by interpolation, reducing the decision dimension from $[eqn]$ to $[eqn]$ coefficients. Using the fixed stage-I solutions as reference, we quantify how this compression distorts the objectives by re-evaluating $[eqn]$ after compressing each candidate trajectory. Figure 11 shows that objective distortion decreases as the number of coefficients increases, demonstrating a concrete pathway to reduce decision dimension while preserving the trade-off geometry. This demo does not replace full re-optimization under a reduced parameterization. Instead, it provides a reproducible sanity check. It also motivates future work that solves the multi-objective problem directly in a lower-dimensional basis (i.e., re-running stage-I and/or stage-II over coefficients rather than over $[eqn]$ samples). While OpenSim evaluation cost per function call remains, reducing the decision dimension can improve the efficiency of derivative-free search and makes neighborhood and refinement sensitivity analyses more tractable in larger models.

5. Conclusions

We presented a synergy-informed, joint-moment-consistent framework for reconstructing unmeasured muscle activations under sparse observations. The method combines (i) a SynRc prior derived from NMF synergies and bounded quadratic programming and (ii) a bound-constrained multi-objective optimization that enforces agreement with inverse-dynamics joint moments while retaining neural consistency. The stage-I Pareto set provides an interpretable trade-off between objectives, enabling diagnostic, task-dependent solution selection. For sparse-measurement studies without ground truth, the Pareto geometry and the measured-only SynRc diagnostic $[eqn]$ provide an interpretable view of identifiability and prior transferability. The recommended reporting protocol (Table 7), together with neighborhood and parameter sensitivity checks (Supplementary Table S4 and Supplementary Figure S3), offers a reproducible way to document selection choices and their stability. Accordingly, the main practical output is not a single “best” reconstruction, but a set of diagnostics that makes assumptions and trade-offs transparent. Synthetic studies demonstrate robustness to noise and mismatch, reveal how the trade-off geometry changes across observability regimes, and illustrate a reduced-dimension parameterization pathway for scaling. In the Arm26 benchmark, $[eqn]$ -dominant selections reduced joint-moment error from 0.154 (SynRc-only) to ≈ $[eqn]$ , while more neural-consistent selections preserved the synergy prior and yielded near-ground-truth reconstruction for TRIlat (RMSE $[eqn]$ ). These findings should be interpreted in light of the study limitations: verification is synthetic, evaluations assume fixed kinematics, and inverse-dynamics targets represent net moments that may include non-muscle contributions. Future work will focus on experimental sparse-EMG applications with subject-specific models and on scalable implementations using reduced-dimension bases and surrogate strategies.