Numerical Identification of Deep Muscle Residual Tensions (Tones) Based on Multi-Directional Trunk Stiffness Data

Featured Application

This study presents a new method to estimate the residual forces (tones) of the trunk muscles at rest, including the deep ones that are difficult to reach with nonintrusive means. Compared to classical methods (intrusive means), this method is painless and safe for the subject, as it uses only trunk experimental stiffness in conjunction with numerical modeling. The present method can be used independently or to supplement classical means which provide direct measures for the tone of superficial muscles. With these measures the identification problem is reduced to finding the tones of only the deep muscles. The proposed method can then be applied, for the specific subject with further gain in accuracy and efficiency. The accuracy of results can be improved by considering more muscle classes in the identification methodology.

Abstract

This work proposes an identification methodology to estimate the residual tension values (tones) for the human trunk muscles, including the deep ones, using multidirectional trunk stiffness data in association with numerical modeling. The role of this residual muscle tension or contraction is to maintain posture and balance. Knowledge of the tone is useful for the diagnosis and treatment of several spinal diseases and is important for realistic modeling and numerical simulation of trunk behavior. Most of the existing techniques for the measurement and estimation of muscle tones, such as electromyography, are restricted to superficial muscles. Those designed for deep muscles are invasive and present risks of infection and pain. In contrast, the proposed identification approach is painless and safe for the subject. It proceeds by matching the experimental trunk stiffness with numerical upper and lower estimates of the stiffness, to construct an inclusive solution domain of possible tone values of superficial as well as deep trunk muscles. By dividing the trunk muscles into three classes, each assumed to share the same tone ratio, a reasonable solution domain is obtained that exhibits a significant overlap with ranges of muscle tones found in the literature.

1. Introduction

Muscle tone is one of the most common parameters governing human movement and posture [1,2]. It is described as residual muscle tension or contraction [3]. The knowledge of muscle tone is often needed for the diagnosis and treatment of movement disorders and spinal diseases [4,5]. However, data on tone levels are scarce in the literature because of the difficulties associated with tone measurement, especially for deep muscles. The human trunk is a complex musculoskeletal system and the muscles it comprises are diverse in size and type. For instance, muscles may be large or small, long or short, deep or superficial, global or local (multi-joint), and grouped or isolated. The muscle tone varies according to the muscle type and depends on other factors such as the body’s daily activity [6], gender, and age [7,8,9]). On the other hand, muscle tone levels influence the behavior of the trunk. Through a sensitivity analysis, [2,10] showed the model trunk flexural stiffness to be significantly influenced by the muscle tones, treated as model parameters. This finding suggests that accurate estimation of muscle tone is helpful for a better assessment of trunk behavior. Furthermore, the influence of muscle tones on trunk behavior provides a basis for their identification, which is the core idea of the present work.

In clinical practice and studies, different techniques are used to measure muscle tones [7,11,12,13]. For instance, electromyography (EMG) is a popular technique that uses electrodes to capture and record the electrical activity in muscles [14]. Then, the recorded signals are processed to quantify the tone levels. In muscle groups, the proximity of muscle elements amplifies the crosstalk phenomenon [15] that complicates the interpretation of the captured signals and reduces the accuracy of tone estimates [16]. Other techniques exist, such as myometry, which is based on oscillation frequency measurement [7,17], that combines the use of a portable device and mechanical measurements to quantify muscle tone. These methods are generally restricted to measuring the activity of superficial, or body member muscles, and are inappropriate for the identification of deep muscle tones.

In the present work, an identification methodology is proposed that delimits the ranges of possible muscle tone values. For the deep muscles, unlike traditional intrusive means of tone measurement, this method is painless and safe for the subject. The identification procedure is intended to capture the properties of the hardly accessible deep muscles, through the association of experimental trunk stiffness data with numerical simulation of the biomechanical trunk system behavior [10].

2. Materials and Methods

The proposed identification method makes use of trunk test data taken from the literature in conjunction with numerical modeling. The data consist of the experimental trunk stiffness measured along different directions for the sample population who participated in the original work [18]. A numerical biomechanical model [10] of the trunk is used to calculate the trunk stiffnesses depending on the muscle tones taken as input parameters. The identification procedure solves the boundaries of the domain of feasible tone values by matching the tone-dependent model stiffness with the experimental stiffness taken as a reference. In accordance with the common practice in biomechanical modeling, the muscle tone has been represented by its ratio to the maximum active force of the muscle, called tone ratio (RMTT) [2]. In the following, the RMTTs will be represented by a vector in a space denoted T.

Details on the experimental protocol relative to the multidirectional trunk test, and the relevant ethical and safety considerations, are given in [18]. A numerical musculoskeletal trunk model [10] is used to calculate the upper and lower bounds of the trunk flexural stiffnesses in five different directions. The domain of possible tone values, or equivalent RMTTs, is narrowed down by locating the experimental trunk stiffness [18] relative to those bounds.

The identification methodology is presented according to either of the following two simplifying hypotheses. In the first, the RMTT is assumed to be identical for all muscles. In the second, a more realistic hypothesis is adopted, where three classes of muscles are distinguished, each sharing a common RMTT value. One class comprises all the external muscles, whereas the internal muscles are divided into two classes: thin and thick muscles.

2.1. Experimental Trunk Stiffness

In most experimental studies on trunk stiffness, the investigations are restricted to flexion and extension in the sagittal plane [19,20]. In the work of [18], the trunk tests were conducted in several directions, including non-sagittal ones, providing abundant data for the given sample population. A group of fifteen healthy individuals, in sitting positions, were subjected to random position perturbations of 10 mm in eight different horizontal directions. The subjects’ trunk stiffnesses were evaluated for each direction based on a mathematical model of a second-order mass-spring-damper system fed with force and kinematics test data. For each of the eight perturbation directions, the group average trunk stiffness was calculated using the mean values from all subjects. The results show distinct stiffness values for all directions. The dissymmetry of the stiffness with respect to the frontal plane is understandable. However, it is less obvious with respect to the sagittal one, considering that the structure of the musculoskeletal system tends to be symmetric. In the present investigation, the trunk musculoskeletal system is assumed to be symmetric with respect to the sagittal plane. Therefore, the stiffness data for any pair of symmetric directions are merged together, reducing the analysis from eight to five distinct directions (Figure 1). The trunk stiffness values (

Table 1. Average experimental trunk stiffness [18] reduced to five directions.

Direction #	1	2	3	4	5
Trunk stiffness (N/mm)	1.30	1.46	2.11	1.90	1.58

) corresponding to the retained directions were deduced from graphic charts reported in [18] and served as reference values for the tone identification process.

2.2. Trunk Numerical Model

The trunk biomechanical model is a lumped type numerical model that enables calculation of the nonlinear static response of the torso, to an input set of applied loads, at a given posture of the body. The model (Figure 2) includes the spine, the rib cage, and 174 muscles [21,22]. The vertebrae, the intervertebral discs, the sternum, and the ribs are represented by beam elements and the ligaments are modeled as linear springs. Each muscle is represented by its tensile force and stiffness while taking into account the unilateral behavior of the muscle–tendon complex [10]. The maximum active force of a muscle is evaluated proportionally to its largest physiological cross-sectional area (PCSA).

The tone is taken proportional to the maximum active muscle force. The calculation of the model trunk stiffness, in the specified displacement direction, proceeds as follows [2,10]. An equilibrium problem (Equation (1)) is run for the trunk subject to a specified horizontal force, applied at the T10 level, with the muscle stiffnesses held constant. A trunk stiffness can be calculated (Equation (2)) as the ratio of the applied horizontal force to the induced displacement at T10. Often, the muscle activations are assumed to remain constant, at the tone level, prior to the induction of the paraspinal reflex. However, the actual muscle active stiffnesses may decrease during the displacement perturbation, due to excessive shortening or shortening velocity. Therefore, accurate calculation of trunk stiffness should take those changes in muscle stiffness into account.

2.3. Lower and Upper Stiffness Estimates

The trunk stiffness [2] measured through the trunk test is dependent on the changes in muscle activation and muscle–tendon stiffness, concurrent with the perturbation. As the trunk regains a static post-perturbation state before the induction of the paraspinal reflex, the effective trunk stiffness should be calculated based on the force and horizontal displacement at this state [10]. At the outset of the displacement perturbation, the trunk is in a state where all muscles are activated, at least at tone level, and in a stretched state. After perturbation, the muscle’s active force and stiffness remain constant or drop due to excessive muscular shortening or shortening velocity [23]. Solving the equilibrium problem while assuming that all the muscles maintain their pre-perturbation activation leads to an upper estimate K_max of the trunk stiffness. On the other hand, a lower estimate K_min is obtained by considering the shortened muscle to be inactive with no stiffness, and the others to preserve their pre-perturbation activation. The equilibrium problem is solved again, with the inactive muscles deleted from the model, followed by straightforward calculation of the trunk stiffness. These upper and lower stiffness estimates will serve in solving the tone identification problem described in the following section.

2.4. Tone Identification Problem Formulation

The tone identification relies on published experimental trunk stiffness data, collected from multidirectional trunk tests [18]. The muscle tones are treated as input parameters in the numerical trunk model. Since the position perturbations ideally start from the same sitting posture the muscle tones to be identified are accordingly presumed to be independent of the direction of perturbation. Thus, the tones are sought such that the model-generated trunk stiffness matches the experimental stiffness K_exp,j for each (the j-th) of the five directions of displacement perturbation. Given the uncertainty associated with the loss of muscle activation during the displacement perturbation, the model stiffness can only be approximated. Therefore, instead of attempting to match the experimental with the model stiffness, we will only require it to lie within the interval [K_min,j, K_max,j] of model trunk stiffness. As a result, the identification amounts to determining a domain D of “feasible” RMTT values. The boundary of D is constructed pointwise by finding it intersections with selected lines crossing it. Along a line in the space T, the RMTTs are represented by a scalar parameter α. A line crossing the domain D intersects it with a segment defined by an interval [α_min, α_max]. The bounds α_min and α_max are solutions to the problems max{α_j/K_min,j(α _j) = K_exp,j, j = 1, …, 5} and min{α_j/K_max,j(α_j) = K_exp,j, j = 1, …, 5}, respectively.

The search process involves a repetitive calculation of the trunk stiffness bounds K_min,j and K_max,j. This requires solving a series of linear equilibrium problems (Equation (1)) for the trunk subjected to small horizontal perturbations in the specified directions. For a given set β of RMTTs and a selected set of muscles assumed to remain active during the perturbation, the equilibrium problem for the j-th direction is stated as:

[K({k})]{U} = {ΔF_j}

(1)

k_m = q β_m F_m0,max/L_m m = 1, …, M_a

where

[K] denotes the stiffness matrix of the trunk system;

k_m: the stiffness of the mth retained muscle;

{k}: the M-vector of retained muscle stiffnesses;

{U}: the vector of nodal displacements (output parameter);

{ΔF_j}: the perturbation force vector representing a single horizontal force applied at T10 in the m-th direction;

q: a dimensionless coefficient taken equal to 10 [10,24];

β_m: the RMTT of the m-th muscle;

F_m0,max: the maximal force developed by the m-th muscle;

L_m: the length of the m-th muscle.

The product β_m F_0m expresses the active muscle force, that is, the muscle tone T_m in the present case of free sitting posture. The trunk flexural stiffness K_j, in the j-th direction, is deduced from

K_j= ΔFj/δ_j

(2)

For the selected set of retained active muscles, where ΔF_j is the horizontal force applied at T10 in the direction j, and δ_j the associated displacement component in the solution {U}. The calculated stiffness K_j corresponds to the upper bound K_max,j when all the muscles are assumed to remain active, and to the lower bound K_min,j when all the shortened muscles are suspended from the model. The shortened muscles are identified from the solution of the equilibrium problem, and run with all the muscles assumed to be active.

The tone identification is conducted according to two distinct formulations. In the first, it is carried out with the crude assumption that the RMTT is identical for all muscles. In an alternative formulation, in order to better capture the diversity of RMTTs, three classes of trunk muscles are defined, and a distinct tone ratio is assigned to each class.

3. Results and Analysis

The identified domain D of possible tone ratios is presented, in this section, under either of the hypotheses of (i) a unique tone ratio for all muscles, and (ii) three tone ratios associated with three classes of muscles.

• Hypothesis of a Unique RMTT

Let the unknown tone ratio be denoted β. The search space T is one dimensional and the domain D to be determined is an interval [β_min, β_max]. For each direction of perturbation j, a solution β_j is obtained for the tone ratio by matching K_max,j with K_exp,j, as indicated by the circled intercept points in Figure 3. The lower bound is β_min = max(β_j, j = 1, …, 5). Likewise, the upper bound β_max is min(β_j, j = 1, …, 5), where β_j is the tone ratio for which K_min,j = K_exp,j (squared intercept in Figure 3).

• Hypothesis of Three Distinct RMTTs

Three classes of trunk muscles are defined, and a distinct tone ratio is assigned to each class. The classification is made on the basis of muscle size and location. The first class, C1, comprises the external muscles (longissimus thoracic, iliocostalis lumborum pars thoracis, longissimus thoracis pars lumborum, iliocostalis lumborium pars lumborium, rectus abdominis and external oblique), class C2 includes the small internal muscles (thoracic multifidus and multifidus) and class C3 contains the large internal muscles (psoas major and quadratus lumborum). The RMTTs are described by three independent variables β^k, k = 1, 2, 3 (ranging between 0.05 and 0.3). The domain D is a volume in the three dimensional space T.

3.1. Shortened Muscles

The results include, for each direction j of perturbation, the list of muscles that are shortened and considered inactive in the calculation of the minimal trunk stiffness K_min,j. A muscle is considered to lose activation if it undergoes shortening during the position perturbation [10].

Table 2. Muscles of the three classes.

Classes	Muscles
C1	RA—(2)/EO—(6)/LT—(24)/Ilpt—(16)/Ltpl—(10)/Ilpl—(8)
C2	TM—(10)/M—(40)
C3	PM—(22)/QL—(36)

Muscles abbreviations: RA: rectus abdominus, PM: psoas major, QL: quadratus lumborum, EO: external oblique, LT: longissimus thoracic, Ilpt: iliocostalis lumborum pars thoracis, Ltpl: longissimus thoracis pars lumborum, Ilpl: iliocostalis lumborium pars lumborium, TM: thoracic multifidus and M: multifidus. Numbers of shortened muscles are put between parentheses when muscles for left (L) and right (R) sides are identical.

lists the number of shortened muscles in each of the 10 groups of muscles, for the five displacement directions. Since the displacements may be non-symmetrical, the numbers of shortened muscles on the right and left sides (of

Table A1. Number of shortened muscles for the five directions.

Muscle (Group)				Direction of Position Perturbation
#	Class	Label	Number of Muscles	1		2		3		4		5
				R	L	R	L	R	L	R	L	R	L
1	C1	RA	2	(1)	(1)	(1)	(1)	-	1	-	-	-	-
2	C3	PM	22	(11)	(11)	-	11	-	11	-	11	-	-
3		QL	36	(2)	(2)	-	15	-	18	(6)	(6) + 11	-	-
4	C1	EO	6	(3)	(3)	-	3	-	3	-	2	-	-
5		LT	24	-	-	-	-	-	9	(9)	(9) + 1	(10)	(10)
6		Ilpt	16	-	-	-	2	-	8	(8)	(8)	(8)	(8)
7		Ltpl	10	-	-	-	-	-	5	(5)	(5)	(5)	(5)
8		Ilpl	8	-	-	-	1	-	4	-	-	(4)	(4)
9	C2	TM	10	-	-	-	-	-	3	-	-	(5)	(5)
10		M	40	-	-	-	-	-	20	(20)	(20)	(20)	(20)
	Total		174	17	17	1	33	0	82	48	73	52	52

) are listed in separate columns.

In addition, the numbers of anatomically symmetrical muscles that are shortened during the same position perturbation are written in parentheses. The shortenings of symmetric muscles are identical in the (sagittal) perturbation directions 1 and 5. For non-sagittal directions, the muscle elongations are asymmetric, and some pairs of symmetric muscles may experience shortening on one side and stretching on the other. The numbers of muscles that undergo shortening during position perturbations in directions 3, 4, and 5 are significantly larger than those for directions 1 and 2. The stiffnesses in these two directions are smaller than those in the other three directions. The number of shortened muscles, their sizes, and the importance of their functions in the considered directions are parameters that affect trunk stiffness.

3.2. Identified Domain of Tone Ratios

3.2.1. Case of Unique Tone Ratio ${\beta }_m$

Figure 3 plots the maximal and the minimal model trunk stiffnesses as functions of tone ratio for the five directions of perturbation.

As expected, the upper and lower estimates of trunk stiffness increase with muscle tone. For each perturbation direction, the circled (squared) tone value, for which the maximum (minimum) model stiffness matches the experimental one, sets a lower (upper) bound on the RMTT (

Table 3. Tone lower and upper limits for the five directions.

Direction #	RMTTs (%)
	βmin	βmax
1	6.7	11.9
2	5.4	10.8
3	4.8	10.2
4	8.6	31
5	9.6	48

). The maximum (minimum) of these values, among the five directions, is a net lower (upper) bound for the RMTT.

The net range D of tone ratio found, common to all directions, is the interval [9.6%,10.2%]. The lowest and largest ratios are those relative to directions 5 and 3, respectively.

3.2.2. Case of Three-Tone Ratios

Figure 4 visualizes the trunk stiffness, for the five directions of position perturbation, as functions of the tone ratios for classes C1 and C2. For each direction j, the shown three surfaces describe the variation of the model stiffnesses K_min,j and K_max,j and the experimental stiffness K_exp,j for β₁ and β₂ in the data range [0.05, 0.3] and a fixed value β₃ = 0.15 (tone ratio for class C3).

The plane representing the constant experimental trunk stiffness intersects the surfaces associated with the maximum and minimum model stiffnesses at two lines. These intersections are determined numerically for the five directions (Figure 4). The projections of the intersection lines onto the (β¹, β²) plane, as depicted in Figure 5, delimit a plane section of the 3D solution domain D, feasible for the five directions. For instance, the section marked with yellow fill corresponds to a cut by the planes β³ = 0.15.

Figure 6 visualizes the plane sections (Figure 6a) of the solution domain D for eight values of class C3 tone ratio, ranging between 2% and 33.5%, and the constructed domain D (Figure 6b). Each of the outer sections is determined through extrapolation of the vertex nodes from the second and third nearest end sections.

The tone ratios of the external (C1) and the small internal (C2) muscles are found to lie, globally, in the intervals [5%, 12%] and [6%, 33%], respectively, as those of the large internal (C3) muscles are varied between 2% and 33.5%. More accurate, narrower ranges can be obtained over a section of D, relative to a selected value of class C3 muscle tone ratio.

There is a clear dependency between ${\beta }_{mi}$. The more the ${\beta }_m$ of the small internal muscles increases, the more the area of the solution domain is reduced.

4. Discussion

The proposed method aims at identifying a set of possible tone ratios for the trunk muscles based on trunk test data. By matching model-generated trunk stiffness with experimental multidirectional stiffness results, the muscular tone values have been approximated for both the superficial and the deep trunk muscles that are hardly accessible to standard, non-invasive, instrumentation.

Given the complexity of the trunk behavior, the identification is simplified by adopting a number of assumptions and approximations, in addition to those pertaining to the numerical model that was reported previously in [10]. A major assumption is that large families of muscles are assumed to share the same tone ratio. This obviously forces muscles that belong to the same class to artificially take on identical RMTTs. A gross approximation was made in the first, single-class scenario, where the tone ratios were identical for all muscles. In the second, more realistic scenario, three muscle classes were distinguished, each sharing one of three tone ratios. The significance of the identification should improve as the classification is further refined. Nevertheless, as the number of unknowns increases, several issues should arise. The identification procedure should be automated since it cannot be carried out manually. Recalling the tones to be the active forces that achieve a stable equilibrium at the pre-perturbation posture, the redundancy in the muscle recruitment mechanism becomes significant, and requires additional constraints to guard against overestimation of the multidimensional solution domain D. Stability [24] and postural constraints [2], which were not considered in the present study, are relevant additions to the formulation of the tone identification problem. Another limitation lies in the differences in the characteristics of the experiment and the numerical model. For instance, the experimental and the model trunk stiffnesses are assumed to carry the same definition [10,18]. The model trunk is considered to be clamped at the hip where the sitting subject is only partially restrained [18]. The experimental trunk stiffnesses are averaged over a sample of 15 subjects while the model stiffness relates to a specific skeletal morphology.

The primary advantage of the proposed identification method is its ability to safely provide quantitative evaluation of the tones for the deep muscles, without resorting to the existing invasive techniques [12]. The method is safe in that it causes no pain, and raises no risk of infection.

In contrast, for the superficial muscles, the tone can be measured easily by existing techniques [7,12,13]. While the proposed methodology provides an alternative means for estimating the tone of these muscles, it can also take advantage of the tone data provided by other techniques [12]. The superficial muscle tone data can be used as a reference in the validation of the proposed methodology by comparing the identification results with the direct measurement data. On the other hand, these data can be used as input for the identification of unmeasured tones for a specific subject. As a result, the identification problem is reduced to finding the tones of only a subset of muscles, such as the deep ones, for the subject. The proposed method can then be applied with further gains in accuracy and efficiency in estimating the unmeasured muscle tones [25].

Table A2. Comparison of results with published reference data.

Muscle Class	Identified RMTT Range (%)	Muscle Groups	Reference RMTT Range/(Overlap)
			McGrill 1991	Park 2014	Kiesel 2007	Lee 2011 *
C1	5–12	RA	2–5 (5)	-	-	-
		EO	4–6 (4–6)	9–12 (9–12)	-	-
		LT-Ilpt	2–21 (5–12)	8–14 (8–12)	-	18–20 * (No overlap)
		Ltpl/-Ilpl	4–7 (5–7)	-	-	-
C2	6–33	TM	-	-	-	19–32 * (19–32)
		M	-	-	9 (9)	-
C3	2–33.5	PM	-	4–8 (4–8)	-	-
		QL	-	2–5 (2–5)	-	-

* EMG measured with pre-activation.

shows published reference data [12,26,27,28] in the form of ranges of tone ratios for particular muscles/groups in addition to the gross ranges identified on the basis of three muscle classes. For class C1 and C2 muscles, the minimum ratios of 5% and 6% are significantly large compared to the minimum for class C3. In fact, the latter can numerically be set to zero while maintaining equilibrium because the C3 muscle’s role is not to maintain and stabilize vertical posture. The psoas major plays a role in flexion, and the quadratus lumborum in flexion and extension. For the class C1 muscles, the 2% bound was taken from the literature [26] and used for the extrapolation down from the 5% level data points. An inspection of the data in Table A2 reveals a large agreement between the results of the identification and the published experimental tone data, except for the data from [12], which were taken for subjects exerting pre-perturbation contractions that resulted in an increase in the muscle active forces. Nevertheless, the results from [12] have been considered here because they are the only ones found that included data on the deep, TM muscles. In [12], fine-wire electrodes were used to capture the electromyographic (EMG) activity from the right deep (multifidus/rotatores) and superficial (longissimus) muscles during sudden loading and unloading of the trunk. The EMG amplitude was calculated during 10 ms intervals, 50 ms before the onset of the perturbation, and 150 ms thereafter. The tone level was estimated proportionally to the EMG peak, based on the beginning of these records.

Considering exclusively the data relative to subjects without pre-contraction [26,27,28], all muscle classes appear to have an overlap between the identified and RMTT ranges and the reference experimental ones. The overlap reduces to a single value for some muscles, such as the RA [26] and the M [28], it covers the complete identified range for LT-Ilpt [26], and includes the entire experimental range for the PM and QL [27]. It should be noted that if experimental or modeling errors were taken into account the overlaps would be even larger.

5. Conclusions

An identification methodology is proposed to estimate the tones of the trunk muscles, including the deep ones that are difficult to reach with nonintrusive means. It combines trunk test experimental data with numerical modeling and solves a series of equilibrium problems where the muscle tone ratios are unknowns. By dividing the muscles into a few classes, assumed to each share the same tone ratio, the solution domain is found to largely overlap with the ranges of tone ratio measurements taken from the literature. The accuracy of the identification can be improved by considering more muscle classes while incorporating postural and stability constraints into the search methodology.