Quantification of the Mechanical Properties in the Human–Exoskeleton Upper Arm Interface During Overhead Work Postures in Healthy Young Adults

TL;DR

This study measures how stiff the interface between the human arm and an exoskeleton is during overhead work, providing new data for improving exoskeleton design.

Contribution

A novel methodology is introduced to quantify the non-linear stiffness of the human–exoskeleton upper arm interface.

Findings

Non-linear pHEI stiffness is best modeled using exponential functions.

Translational stiffness ranges from 2.1 N/mm to 4.5 N/mm depending on the direction.

Rotational stiffness around the arm's longitudinal axis is significantly lower at 0.03 N·m/°.

Abstract

Exoskeletons transfer loads to the human body via physical human–exoskeleton interfaces (pHEI). However, the human–exoskeleton interaction remains poorly understood, and the mechanical properties of the pHEI are not well characterized. Therefore, we present a novel methodology to precisely characterize pHEI interaction stiffnesses under various loading conditions. Forces and torques were applied in three orthogonal axes to the upper arm pHEI of 21 subjects using an electromechanical apparatus. Interaction loads and displacements were measured, and stiffness data were derived as well as mathematically described using linear and non-linear regression models, yielding all the diagonal elements of the stiffness tensor. We find that the non-linear nature of pHEI stiffness is best described using exponential functions, though we also provide linear approximations for simplified modeling. We identify statistically significant differences between loading conditions and report median translational stiffnesses between 2.1 N/mm along and 4.5 N/mm perpendicular to the arm axis, as well as rotational stiffnesses of 0.2 N·m/° perpendicular to the arm, while rotations around the longitudinal axis are almost an order of magnitude smaller (0.03 N·m/°). The resulting stiffness models are suitable for use in digital human–exoskeleton models, potentially leading to more accurate estimations of biomechanical efficacy and discomfort of exoskeletons.

Figures

Full text

1. Introduction

Exoskeletons are wearable devices that assist, enhance, or substitute the motor abilities of users [1]. Their applications range from supporting physically impaired individuals or serving as rehabilitation tools to alleviating the physical burden on workers or soldiers during strenuous activities [2,3]. Occupational exoskeletons can reduce muscle and joint strain, aiming to prevent overexertion and long-term musculoskeletal disorders [4,5]. To achieve this, they must apply forces through physical human–exoskeleton interfaces (pHEI) [6]. In addition to supportive forces, undesirable loads may arise at the pHEI, such as parasitic shear forces or torques [7]. These can result from misalignment of human and technical joint axes [7,8], movement restrictions, or design issues [9]. In fact, interaction forces and torques can occur in all directions during movement (see [9,10,11]), with their magnitude largely determined by the exoskeleton’s mechanical architecture [9]. Consequently, efficacy and efficiency decrease when components of the acting force or torque cannot be utilized for joint torque reduction [8] or when energy is lost through tissue deformation (with soft tissues acting as energy sinks) [12,13]. Additionally, these loads can lead to discomfort, pain, or injuries [7,8]. Also, the compliance of soft tissue allows relative movements of body attachments to the human body [14,15,16], which can exacerbate joint axis misalignment and related issues [7,14]. If body attachments slip, shift, or rotate, the entire exoskeleton moves relative to the body, altering lever arms and force directions, thus affecting interaction loads at the pHEI and stresses on joints and muscles. Therefore, the mechanical properties of the pHEI directly influence the efficacy and efficiency of the support, as well as the perception of discomfort [17].

To evaluate the biomechanical efficacy of exoskeletons, experimental methods are employed, such as the investigation of muscle activity or subjective strain [18,19,20]. In addition, musculoskeletal human models (MHM), e.g., the Anybody Modeling System (AMS, AnyBody Technology, Aalborg, Denmark) or OpenSim (National Institute of Health, Bethesda, MD, USA), are utilized, allowing for biomechanical analyses to be conducted earlier during the development process [17]. Furthermore, these models enable the evaluation of parameters that cannot easily be captured experimentally (e.g., joint kinetics) [17]. The models are used to analyze the biomechanical efficacy of exoskeletons [9,21,22,23,24,25,26] or to compare exoskeleton concepts and optimize their design parameters [9,25,27]. As previously described, the mechanical properties of the pHEI directly influence the efficacy and discomfort of the exoskeleton. Consequently, contact modeling in multibody systems significantly impacts the predictive capabilities of simulation models [17].

Often, in such models, loads are applied directly to the rigid body segments (bones) or through rigidly connected structures, offering the least accurate representation of reality [17]. Sometimes, joint degrees of freedom and relative movements are allowed, or in newer approaches, contact elements with adjustable force transmission properties are implemented [9,17,23,26,28]. However, the mechanical properties of the pHEI cannot be defined, which is why other approaches utilize spring-damper models [17]. Model complexity ranges from two-dimensional human models with one-dimensional contact modeling [29,30] to three-dimensional human and contact models [22,31,32]. However, mostly general interaction models have been implemented, without values based on biomechanical measurements (except for Shafiei et al. [32], discussed later), possibly due to a lack of available data.

In general, the physical interaction between humans and exoskeletons is not well researched, with testing methods varying greatly and lacking harmonization, reproducibility, and realistic testing conditions [6]. So far, supportive normal forces and interaction torques have been measured in pHEI, for example, to recognize user intention and to control active exoskeletons [6,33,34,35,36]. In some cases, shear forces and interaction pressures have been assessed due to their relevance for discomfort and safety [6,10,11,37,38,39,40,41]. However, the mechanical properties of the pHEI have rarely been characterized, even though more accurate MHM could yield more precise estimations of efficacy and discomfort of exoskeleton concepts. Few studies have considered not only forces but also the movements of pHEI [13,32,42,43,44,45]. Mostly, multiple body interfaces were simultaneously and/or multiaxially loaded [42,43,44,45], complicating the characterization of individual pHEI and the application of models in MHM. Fewer studies have considered interaction loads and displacements in isolated loading directions. Langlois et al. determined forces and displacements in an upper-arm pHEI in the arm’s longitudinal direction and were able to derive interaction stiffness and damping parameters [13]. Yet, mechanics of other types of loads (torques) or loading directions (e.g., perpendicular to the arm) were not characterized in any of the studies [13,42,43,44,45]. Shafiei et al. applied forces and torques in different directions to the calf and thigh pHEI and measured loads and displacements with a load cell and motion capture (mocap) to derive a 6 × 6 stiffness tensor [32]. This approach allowed for full-dimensional characterization, although the stiffnesses were approximated linearly [32]. We assume that the anatomy and properties of different soft tissues primarily determine interaction stiffness. Typically, those tissues exhibit non-linear properties, following exponential, polynomial, or power functions [46,47,48]. Such trends were also noted in pHEI (e.g., [42,43]). Since Shafiei et al. applied loads randomly and manually [32], the measurement range for the linear approximation of each stiffness remains unclear. Furthermore, with only three subjects studied, case studies are feasible, but the statistical analysis of mechanical properties is not.

To our knowledge, no pHEI has been mechanically characterized in all direct loading directions with statistical validity. The upper arm interface—vital for shoulder exoskeletons—is particularly under-researched. Also, while there are indications that load–displacement curves progressively increase [42,43], accurate mathematical descriptions are still lacking. Hence, we present a novel methodology that enables the investigation of mechanical pHEI properties under different loading conditions (forces and torques in various axes), allowing controlled and repeatable application of loads to multiple subjects. The loads and displacements of the pHEI are measured simultaneously with high precision to derive mechanical characteristics. Also, relevant factors that could potentially influence the interaction (such as boundary conditions like attachment pressure, muscle activity, interaction speed, and body posture) can be defined within the pHEI, ensuring the acquisition of application-relevant data. We utilize this methodology to examine the upper arm pHEI and derive parameters for use in MHM. Since the mechanical properties of pHEI can be defined in three orthogonal translational and rotational axes in MHM (compare [9,22]), we establish relationships between respective loads and relative displacement in three axes. Hereinafter, we refer to these relationships as “interaction stiffness” for simplicity, even though minor damping and inertia components are included. We characterize non-linear stiffness profiles with mathematical functions that can be easily applied in MHM. We also investigate the magnitudes of interaction stiffness under different loading conditions as well as correlations with anthropometric parameters (e.g., body mass index, arm circumference). Finally, the results are discussed regarding plausibility and applicability to MHM.

2. Materials and Methods

2.1. Measurement Setup

The experimental setup consists of two main structures: a support frame and a cart equipped with a drive unit (Figure 1a). The support frame features various fixtures for positioning and securing the subject’s shoulder and arm, adjustable according to the subject’s anthropometry (Figure 1b). The drive unit can apply both forces and torques (load types) in two directions (positive and negative load directions) to the physical human–exoskeleton interface (pHEI) under investigation (e.g., upper arm attachment, see Figure 1c) via a linkage (Figure 1b). To prevent overloading the subject, several safety mechanisms have been implemented, including a safety controller, mechanical end stops, and emergency switches. Adjustments on the cart allow the application of loads from three orthogonal axes onto the pHEI (Figure 1c, Appendix A.1.2 Figure A1a). The first load axis aligns approximately with the arm’s longitudinal axis (distoproximal axis), while the other two are perpendicular to it: one pointing lateromedially towards the pHEI and the other caudocranially from below. Various pHEI configurations can be tested in this setup. In this instance, a pneumatic upper arm attachment was utilized to control the attachment pressure via padded blood pressure cuffs and a compressed air supply. To measure the forces and torques acting on the pHEI, a load cell (6-axis force-torque sensor K27 × 53, measuring range up to 500 N and 20 N·m, with a maximum measurement error of 1% according to the manufacturer, SENSIX S.A.S., Poitiers, France) is positioned in the linkage between the drive unit and the upper arm attachment (Figure 1b and Figure A1b). Movements of the arm and upper arm attachment are captured using an optical motion capture system with 13 cameras (seven Oqus 700+, four Oqus 400, and one Oqus 310+ for marker capture, and one Oqus 210c for video recording, Qualisys AB, Gothenburg, Sweden). Six mocap markers are attached to the upper arm attachment and the subjects at bony landmarks of the acromioclavicular joint and medial and lateral epicondyles of the humerus (Figure 1b,c). Further details regarding the experimental setup can be found in Appendix A.1.

2.2. Measurement Procedure

For the experiments, 21 subjects (10 female, 11 male) aged between 20 and 34 years were recruited from employees at the Fraunhofer IPA. Subjects with sensory disorders or those taking analgesics were excluded from the study. The subjects were informed about the experimental procedures, safety features, and potential risks. Subsequently, their height was recorded, and the length and circumference of the left upper arm (with 90° elbow and shoulder flexion) were measured using a measuring tape. Body mass index (BMI), body weight, body fat percentage, and muscle mass were recorded using a digital body analysis scale (Model HBF-511B-E, OMRON Medizintechnik Handelsgesellschaft mbH, Mannheim, Germany). Mocap markers were attached to the left elbow and shoulder, and arm sleeves (Unisex Armlinge, Brubeck Body Guard, Toronto, ON, Canada) were placed over the upper arm to standardize friction properties and simulate realistic working conditions, as workers typically wear clothing beneath the arm attachment. The subject was positioned in the measurement setup (shoulder: 90° flexion, 45° horizontal extension, 0° rotation; elbow: 90° flexion, 0° pronation/supination), and the setup was adjusted to the first configuration (first selected load axis) to match the subject’s anthropometry. The motion limits of the drives were recorded to prevent collisions between the upper arm attachment and fixtures during the experiments. Additionally, it was ensured (by monitoring a live load cell signal) that the subject did not exert significant pre-tension on the upper arm attachment while in the starting position with a relaxed arm, and this was regularly checked throughout the experiment. Finally, the subject’s arm was secured using the fixtures.

During the experiments, two types of trials were conducted: one experiment to determine discomfort thresholds, which is not presented here, and the subsequently presented trials aimed at assessing interaction stiffness. At the beginning of each trial, we closed the upper arm attachment, told the subject to relax the upper arm, applied a constant attachment pressure of 65 mbar, and initiated an automated measurement in a predefined load axis, type, and direction. The drives applied varying, defined interaction forces between 20 N and 140 N in seven approximately equal increments at a speed of 10 mm/s, or interaction torques between 2 N·m and 10 N·m in five approximately equal increments at a speed of 8°/s. Additionally, a previously determined subject-specific discomfort threshold (range: 0–150 N, 0–10 N·m) was applied, resulting in six or eight different stimuli being presented five times each (in total 30 or 40 stimuli per loading condition). Subjects were first presented with all loads in a defined order, followed by the same loads in four repeated measurements in randomized order. The drives moved the upper arm attachment until the defined load (controlled via predefined current limits) or a predefined motion limit was reached. The movement was held for two seconds, the arm attachment was returned to the starting position, and the next movement was initiated. After completing one set, we repeated measurements for the other load type and direction in the same loading axis. We then rearranged the experimental setup, readjusted it to the subject anthropometry, and repeated the trials in the two other loading axes. The order in which loading axes were selected and the order of loading types and directions within each loading axis were randomized. Short breaks were provided between each change of the loading condition (axis, type or direction), during which the pressure from the upper arm attachment was released, allowing the subject to remove their arm and stand up if needed.

In all trials, movements of the pHEI and arm were recorded using the described camera system (see Section 2.1) and the software Qualisys Track Manager 2023.1 (Qualisys AB, Gothenburg, Sweden). Load cell data were also recorded via Qualisys Track Manager. Additionally, the drives’ positions, velocities, and accelerations were monitored using Elmo Application Studio II (Version Number 2.8.0.22, Elmo Motion Control Ltd., Petah Tikva, Israel). We calibrated the mocap system daily and reset the load cell without the arm present regularly (at least before each new loading condition) to avoid sensor drift and ensure that only interaction forces and torques were measured, excluding the weight of components (upper arm attachment, linkage).

2.3. Data Processing and Analysis

We excluded 14 of the 252 measurements (12 loading conditions, 21 subjects) because of missing data (e.g., recording stopped early) or erroneous load cell readings where the measured loads were shifted to unrealistically high values (e.g., baseline a few hundred Newtons above zero, even though no force was applied). The remaining data were imported and processed in an automated workflow using MATLAB (R2023b, The MathWorks, Inc., Natick, MA, USA). Load cell and mocap data were extracted in the respective main loading directions, and relative displacements were calculated as the difference between absolute motions of the arm and arm attachment (Figure 1d). The data were smoothed and subjected to several manual checks, including correcting obviously erroneous data in individual cases. In one instance, force data were offset by ~50 N due to preload during load cell resetting. In five cases, automated starting position detection from mocap data failed, requiring manual definition. Afterwards, the algorithm derived load–displacement curves for each stimulus (one single motion) and carried out an automated offset correction to ensure that force-displacement curves pass the zero point. For each loading condition of each subject, the ten (out of 30 or 40) individual stimuli with the longest load–displacement curves were selected for subsequent analyses of interaction stiffnesses. A detailed outline of the data processing can be found in Appendix A.2.

We fitted various functions to the datasets using the MATLAB function fitnlm() to analyze and mathematically describe the nature of the load–displacement curves. The functions included linear ( $[eqn]$ ), quadratic ( $[eqn]$ ), and exponential functions ( $[eqn]$ ), as well as power laws ( $[eqn]$ ), all defined to pass through the origin of load and displacement (force/torque: f(x), displacement: x, coefficients: a and b). We fitted functions to individual stimuli (individual stimulus regression, ISR), to all 10 selected repetitions per loading condition of a subject (subject-specific regression, SSR), and to all selected load–displacement curves of all subjects per loading condition (subject group regression, SGR). We then analyzed the adjusted coefficient of determination R^2^ for each case to facilitate qualitative comparisons of fitting accuracy. The slopes of the fitted first-order polynomials were utilized to compare the magnitudes of interaction stiffnesses across different loading conditions.

3. Results

3.1. Load–Displacement Profiles

Hereinafter, we examine the shape of the load–displacement curves, followed by an analysis of stiffness magnitudes under varying loading conditions. We fitted regression models to the load–displacement data for shape analysis, visually inspected the results, and compared the adjusted coefficients of determination (R^2^) across fittings. Figure 2 displays exemplary load–displacement curves for one subject under a specific loading condition (positive rotation about the lateromedial axis) and for all subjects under the same condition, each fitted with selected regression models. The data reveal that in this loading condition, interaction stiffnesses are more accurately described by exponential functions than by linear models, which is true for both individual subject fits (ISR, SSR) and the fittings to data from all subjects (SGR, see Figure 2).

Table 1 quantifies the loading conditions through the coefficients of determination of the regression models; Figure 2 illustrates one loading condition. The first section of the table lists the mean coefficients of determination for fittings to individual stimuli (ISR). Exponential functions best fit across all directions (mean: R^2^ = 0.99). Power laws and quadratic functions also yield high coefficients (mean: R^2^ = 0.98), while linear models show average R^2^ values of 0.86, indicating a significantly poorer fit than the other models. Notably, straight lines fit the primary data with slightly higher accuracy for rotations around the distoproximal axis and for all translational movements compared to rotations perpendicular to the distoproximal axis. This can be attributed to the steeper ascending load–displacement curves of the rotational movements perpendicular to the distoproximal axis.

We observe similar trends in the analysis of direct fittings per subject and loading condition (SSR). Regression models with exponential functions show the highest coefficients of determination, while straight lines fit the data with less accuracy. Generally, rotational movements are modeled more accurately by exponential functions than translations due to the better fit of exponential functions to the more steeply curved trajectories in rotational movements. Overall, coefficients of determination are slightly lower than those in the ISR. However, the high overall R^2^ values indicate low variability and high repeatability of individual measurements per subject. Visual analysis of the plots further supports this finding, often showing negligible differences between the means of ISR and the corresponding SSR (see, for example, Figure 2a). Standard deviations are also low and often nearly indistinguishable from the mean values. This trend holds across all directions and almost all subjects. In a few exceptions, individual measurement trajectories deviate from the pattern, reducing SSR quality and decreasing R^2^ values.

Regression models fitted to subject-specific data (ISR, SSR) describe the qualitative nature of the interaction stiffness of the pHEI. However, fittings for all selected data from all subjects under one loading condition (SGR) may provide universally applicable stiffness models. Figure 2b and Figure A3 (Appendix B.1) present the load–displacement curves of all measurement data and the corresponding SGR. The coefficients of determination in Table 1 are significantly lower compared to ISR and SSR fittings, reflecting notable variability among subjects. On average, coefficients for quadratic, exponential, and power functions are nearly identical and only slightly higher than those for linear regression. R^2^ values for fittings to distoproximal translations are slightly lower due to greater scatter. Figure A3 in Appendix B.1 visually indicates increased variability for these directions as well. Perpendicular rotational loads (lateromedial and caudocranial rotations) exhibit lower variabilities and more pronounced curvature in primary data trajectories. Hence, progressive regression models (Q, E, and P) show slightly higher coefficients of determination for perpendicular rotations compared to other loading conditions (Table 1). The coefficients for the plotted regression models (see Figure A3) are presented in Table A1 (Appendix B.1) for linear and exponential regressions.

We found that SGR fittings often do not accurately represent the qualitative nature of the underlying data due to varying lengths of load–displacement curves across subjects. Smaller displacements yield more data points, dominating the fittings. An alternative could be curves of mean values from individual stimulus regression functions (ISR), presented in Figure A3 and corresponding model coefficients in Table A1 (Appendix B.1). These curves visually match the data better; however, valid R^2^ values are not available. Depending on subsequent modeling questions, one of the two functions may be more beneficial than the other.

3.2. Interaction Stiffness Magnitudes

To analyze the magnitudes of interaction stiffness and compare different loading conditions, we investigated the slopes of linear regression models, as a single parameter can represent stiffness. The slopes grouped by loading condition mainly did not exhibit a normal distribution (Anderson–Darling test for normality). However, they displayed similarly shaped distributions (tending to be slightly right-skewed). Thus, we performed Mann–Whitney tests to test the inequality of medians across different loading conditions. We found significant differences between almost all loading conditions, except for positive and negative caudocranial rotations relative to each other. Figure 3 presents means and box plots of all ISR linear regression model slopes.

Translational movements exhibited median interaction stiffnesses ranging from 2.1 N/mm to 4.5 N/mm. The lowest stiffnesses occurred in the distoproximal direction, followed by the negative lateromedial direction. The positive lateromedial direction showed significantly higher stiffness. Increased stiffnesses were also observed in the caudocranial direction, with the greatest interaction stiffness in the negative direction (Figure 3, Table A1). The interquartile range was approximately 27–39% of the respective median value, with the greatest deviations in the distoproximal and positive lateromedial directions. The lowest value was 1.1 N/mm, while the highest was 7.0 N/mm without outliers (the highest outlier: 10.9 N/mm).

Interaction stiffnesses for rotational movements can be categorized into two groups. Rotations around axes perpendicular to the arm (lateromedial, caudocranial) exhibited similar stiffnesses (medians between 0.18 N·m/° and 0.22 N·m/°, see Figure 3). Stiffnesses were slightly lower for rotations around the caudocranial axis. Interquartile ranges were between 21 and 34% of the respective median, with lower deviations in positive rotations. The lowest slope determined was 0.10 N·m/°, and the highest was 0.34 N·m/° without outliers (lowest outlier: 0.08 N·m/°, highest outlier: 0.38 N·m/°). The second group, rotational movements around the arm’s longitudinal axis (distoproximal), exhibited stiffnesses reduced by nearly an order of magnitude, with medians of 0.03 N·m/°. Interquartile ranges were 34–49% of the respective median, with the lowest value being 0.01 N·m/° and the highest 0.06 N·m/° without outliers (highest outlier: 0.07 N·m/°). The slopes of the linear regression models and the scatter of the underlying measurement data are also clearly depicted in Figure A3 and quantified in Table A1 (Appendix B.1).

4. Discussion

4.1. Load–Displacement Profiles

The results show that stiffness profiles in the pHEI progressively increase. Analyses of adjusted coefficients of determination for regression models on individual stimuli (ISR) and per subject (SSR) reveal that exponential functions most accurately fit the data (means of ISR: R^2^ = 0.99, SSR: R^2^ = 0.95), closely followed by quadratic functions and power laws. Visual inspections of graphs confirm this finding, with high similarity in form and trend evident in almost all cases. Only for small movements, especially rotations up to ~10°, do the models slightly underestimate the loads while closely aligning with the rest of the profiles (see Figure 2a). For most simulation-based assessments, these minor differences should not significantly impact. Isolated human soft tissues like skin, fat, and muscle also exhibit progressive stiffness profiles or stress-strain curves, which are flatter at small strains and steeper at higher strains [46,47,48,49]. Consequently, they are often characterized by power laws or quadratic, polynomial, or exponential functions [46,47,48]. This supports the assumption that pHEI stiffnesses are essentially defined by the soft tissue properties. Also, soft tissue stiffness profiles are often simplified using linear (or bi-linear) functions [46,48,49]. In our case, linear approximations did not accurately fit the curves, as expected; however, they still achieved high coefficients of determination, with a mean R^2^ value of 0.86. Fundamentally, linear functions tend to overestimate loads at small displacements and underestimate them at high displacements (Figure 2a). Nevertheless, they serve as a suitable means to intuitively compare interaction stiffness magnitudes, aid in rough calculations, or enable less computationally intensive simulation models. For precise personalized mechanical models of the pHEI, we recommend using exponential functions. However, universally applicable stiffness models would be more practical, as they eliminate the need for new data acquisition for each user.

With an average R^2^ of 0.71–0.76, fitted SGR matched the data less accurately than personalized fittings (ISR, SSR). All non-linear functions exhibited similar coefficients of determination, making them equally suitable for pHEI modeling. However, the models did not ideally follow the load–displacement curves but appeared somewhat too flat. This occurred because varying displacements of individual stimuli resulted in more data points at small displacements being included in the fittings than at high displacements, thus dominating the functions. The averaged ISR across all subjects appeared to follow the natural trend more closely, though it potentially overestimates forces/torques at higher displacements. Depending on the research question or investigation goal, one relationship may be more advantageous than the other, which is why we provide both (see Appendix B.1, Table A1). It should also be noted that in SGR, the differences between exponential (R^2^: 0.75) and linear fittings (R^2^: 0.71) were negligible. Consequently, linear stiffness models can be used within the examined measurement range without significant loss of accuracy. These models are easier to implement, intuitively understandable, and less computationally intensive. They are also suitable for estimating and comparing magnitudes (e.g., for different load directions, arm attachments, or boundary conditions).

4.2. Interaction Stiffness Magnitudes

Hereinafter, we interpret stiffness magnitudes before discussing the plausibility of the results. Distoproximal translational stiffnesses were significantly smaller than in perpendicular translations, as shear moduli, which are lower than E-moduli, dominate tissue mechanics (compare muscle tissue [50]). Positive distoproximal stiffness is slightly higher than negative because the upper arm thickens towards the shoulder. Negative lateromedial interaction stiffness is reduced compared to positive, since the relaxed biceps muscle tends to fall slightly medially and can be pushed past the humerus (shear moduli dominate mechanics) without being fully compressed against it, resulting in lower forces at the start of the movement (see Appendix B.1, Figure A3). In the opposite direction, the soft tissue compresses against the bone more quickly, as the humerus tends to lie laterally within the arm volume. Stiffness is lower when pressing against the underside of the arm compared to the top, due to the generally greater presence of adipose tissue in this area [51], which has lower E-moduli [52]. Rotational stiffness around the longitudinal arm axis is lower than for torques perpendicular to the arm, as distoproximal torques are dominated by shear moduli parallel to the muscle fibers. Additionally, lever arms and effective area decrease towards the bone, leading to increased tissue stress and deformation. In torques perpendicular to the arm, however, the tissue is compressed against an approximately constant bone area at the edges of the arm attachment. Positive lateromedial torques exhibit slightly higher stiffness than negative ones, because soft tissue layers are thinner proximally on top and distally on the bottom (see [53]), leading to earlier compression against the bone. Stiffness in caudocranial rotations is lower than in lateromedial rotations, as muscles primarily run above and below the humerus. Thus, in caudocranial rotations, the tissue is partially pushed past the humerus (shear moduli dominate mechanics), whereas in lateromedial rotations, it must be compressed more against the bone. Overall, all measured effects can be plausibly explained with anatomical and tissue knowledge. Additionally, the load cases were simulated in a highly simplified finite element (FE) model of an upper arm with arm attachment to roughly estimate interaction stiffness magnitudes (see Appendix B.2). Both the magnitude and ratios of stiffnesses (e.g., perpendicular vs. distoproximal rotations) were similar to the measured stiffnesses, with a maximum deviation of factor 2, confirming the plausibility of the measured magnitudes.

Comparisons with literature values further support this. Langlois et al. reported average stiffness in the longitudinal upper arm direction ranging from 1.2 to 2.4 N/mm at the same attachment pressure [13]. Despite differing experimental setups and arm attachments, these values align closely with ours (medians: 2.1–2.4). The slightly smaller values in Langlois et al. are likely due to shorter displacement distances (<20 mm vs. >30 mm), as the progressive stiffness curves cause a steeper slope of linear fittings in our study. No additional studies on the upper arm were identified, whereas thigh and calf interfaces have been examined more frequently [32,42,43]. Shafiei et al. determined distoproximal translational stiffnesses of 1.1–3.4 N/mm and perpendicular translational stiffnesses of 2.0–5.3 N/mm [32], closely aligning with our results. Rotational stiffnesses were slightly higher at 0.12–0.15 N·m/° (distoproximal) and 0.25–1.03 N·m/° (perpendicular) [32], likely due to anatomical differences. Since forces and torques were applied randomly and manually, deviations may also arise from varying measurement ranges, as the linear approximation of stiffnesses heavily depends on these ranges (see Section 3.1). Quinlivan et al. investigated stiffness at thigh interfaces, where parallel and perpendicular force components acted simultaneously [43]. At approximately 30 mm displacement, forces between 40 and 80 N were observed [43], yielding linearly approximated stiffnesses of about 1.3 to 2.7 N/mm. Asbeck et al. found similar values in a comparable setup with hip and thigh interfaces [42]. Both studies presented progressive increases in stiffness [42,43]. Although thigh and calf pHEI are not directly comparable to the upper arm, the similar magnitude of results supports the plausibility of the stiffnesses determined here. Furthermore, some approaches measure tissue stiffness rather than total pHEI stiffness directly. Yousaf et al. found tissue stiffness perpendicular to the forearm, ranging from 0.4 to 1.4 N/mm, using indenters [54]. The lower values compared to our study may result from smaller indentations (17 mm vs. >30 mm) and anatomical differences. Mathematical models, such as smooth orthogonal decomposition, have also estimated tissue stiffness from mocap marker trajectories on the calf and thigh, revealing similar translational stiffness magnitudes [55].

The results confirm the plausibility of both the shape and linearly approximated magnitudes of stiffnesses. Next, we discuss our methodology concerning the accuracy of measurement results and the applicability of the determined models.

4.3. Limitations and Applicability of Models

The results are subject to the following measurement inaccuracies. The calibrated load cell had a maximum error of 1% of the measurement range (approximately 5 N or 0.2 N·m). Regular load cell resetting ensured that only interaction loads were recorded. The mocap system was calibrated daily, capturing marker positions accurately to the millimeter. However, mocap markers cannot be directly attached to the bone, so skin movement relative to the bone under load may have distorted the measured position of the humerus and the resulting relative displacements by a few millimeters or degrees. Velocities were achieved within fractions of seconds, and the attachment pressure changed only marginally (a few mbar) under load, allowing boundary conditions to be considered constant. Data processing inaccuracies mainly arose from two sources.

First, minimal pre-tensions applied by the subject in the zero position could not be completely ruled out. However, it was assumed that no load can act at zero relative displacement, so the load–displacement curves were corrected to pass through the origin (minor corrections by a few millimeters or degrees). Second, in distoproximal rotations, absolute rotations were determined instead of relative rotations due to insufficient marker visibility. However, since the arm was very compliant under these loads (large rotations at low torques, see Appendix B.1, Figure A3) and no noticeable rotation of the humerus was observed during experiments, we assume that stiffness was only minimally distorted. Overall, inaccuracies from the measurement setup and data processing are considered minor, having only a marginal impact on the trends and magnitudes of the determined stiffnesses. The analysis of coefficients of determination (SSR) also shows a high repeatability of measurements and reproducibility of the determined stiffnesses. SSR demonstrated a high R^2^ (only slightly reduced compared to ISR), indicating low variation between individual stimuli. Slightly increased variations appeared in the positive translational directions (lower R^2^ of SSR, see Table 1). In these directions, subjects’ arms were less securely fixed (straps instead of rigid fixtures), requiring them to actively counteract motion or correct their posture occasionally. This may have led to brief increases in muscle activity or slight angle changes, affecting interaction stiffness. However, this impacted only a few subjects and often only single stimuli, rendering the impact on stiffness magnitudes and profiles negligible. Nevertheless, future efforts should focus on stabilizing subjects even better to optimize experiment repeatability.

When utilizing our magnitudes and models, note that we have only captured the main elements of a 6 × 6 stiffness tensor; thus, coupling deformations must be neglected when applying the models to rigid multi-body MHM. Furthermore, the “stiffness” results actually describe mechanical impedance, encompassing spring stiffness, damping, and inertia. However, under the investigated constant velocities (no relevant inertia), damping effects were negligibly small (estimated with damping coefficients from [13]). Furthermore, limitations exist in transferring the results to other applications due to the investigated measurement range, subject cohort, arm attachment, and chosen boundary conditions.

Relative displacements exceeding 30° or 30 mm were observed, and the models are expected to be valid within this range, but larger deviations may occur beyond it. The subjects represented body sizes typical of the German population (normally distributed, sizes close to the 5th to 95th percentile according to [56]), with roughly equal numbers of women and men. However, the age was relatively homogeneous, BMI was predominantly within the normal range, and the subjects’ physique did not necessarily match that of workers. In larger or differently composed groups, results may vary more widely or change in level and quality. The boundary conditions were chosen to reflect overhead work. Shoulder and elbow joint angles were aligned with those typical of overhead work [57], and attachment pressures were based on those found in commercial overhead exoskeletons in similar poses (59–72 mbar [58]). Although overhead work heavily strains the shoulder muscles, the upper arm is often minimally affected [21,59], so experiments were conducted with a relaxed upper arm. For other activities, such as heavy lifting or carrying, these boundary conditions may change, and different interaction stiffnesses are anticipated in such scenarios.

Finally, the results are fundamentally valid only for the specific arm attachment used in the study. To generate accurate models for other arm attachments, experiments must be repeated with those attachments, which is feasible with the current setup. However, since other studies using different upper arm pHEI [13] and even other body interfaces [32,42,43] have produced similar results, it can be concluded that the presented values might be applicable to other upper arm attachments with some loss of model accuracy. In general, the models may be less accurate in predicting outcomes for different users and applications. However, to our knowledge, published values for most loading conditions do not exist, so the functions described here can, for the first time, provide at least a valid magnitude of values for the upper arm interface. Future investigations should aim to quantify how altered boundary conditions and subject cohorts affect the results to derive more universally applicable models, possibly using correction factors or similar approaches. Variation between subjects may be greater than other influences combined. Our investigations have shown that interaction stiffnesses correlate with some easily measurable anthropometric characteristics (BMI, weight, upper arm circumference, gender, see Appendix B.3), which might allow for the future derivation of semi-personalized models that adjust general stiffness models (SGR) using correction functions for individual users.

5. Conclusions

In this study, we established a methodology for the comprehensive mechanical characterization of the pHEI under various boundary conditions. We found that the approach allows precise measurement of loads and displacements with high repeatability, yielding plausible interaction stiffnesses. To our knowledge, this is the first time that all diagonal elements of a 6 × 6 stiffness tensor for an upper arm pHEI have been characterized, enabling a full pHEI description in MHM (when neglecting coupling effects). In general, this may be the first instance of a pHEI being fully characterized with enough subjects to obtain statistically significant data. We identified significant stiffness differences between loading conditions. Forces and especially rotations along the arm’s longitudinal axis exhibited notably lower stiffness compared to other conditions. Beyond determining magnitudes, we derived precise mathematical models of pHEI stiffness profiles, best described by exponential functions to capture their non-linear nature. The results enable implementing general and highly accurate personalized interaction models of the pHEI in MHM, applicable to overhead work simulations. If some loss of accuracy is acceptable, our models can also be used for simulations with other arm attachments or application domains.

Future work could optimize the setup to achieve more rigid fixation of subjects. The methodology may be used for further studies to characterize commercial upper arm attachments or explore the impact of various boundary conditions (posture, attachment pressures, relative velocities, muscle activity), populations, or anthropometries on mechanics. Such investigations might foster the development of models with broader general applicability.