Muscle structure and the mechanism of force production

Skeletal muscles are made up of fiber bundles, nerves, blood vessels, cells, and other tissues. The muscle fibers within the skeletal muscles can be divided into intrafusal muscle fibers (providing proprioception) and extrafusal muscle fibers (providing contraction) [124]. The intrafusal muscle fibers, also known as muscle spindle, lie within and alongside the extrafusal muscle fibers [124] (Figure 2.3).

Figure 2.3: Structure of a typical skeletal muscle.

The extrafusal muscle fibers are organized in bundles, called fascicles (Figure 2.3). These muscle fibers vary in length from a few millimeters up to ~30 cm and in diameter from 10 to 500 µm [80]. They also have different arrangements in muscles, such as parallel to the muscle line of action, pennated, convergent, or fusiform, influencing some of the mechanical properties of muscles. For instance, muscles designed for speed tend to have parallel fibers, whereas the muscles intended for strength are typically pennated [125]. Each fiber within the fascicles is made up of myofibrils. Myofibrils, wrapped in the sarcoplasmic reticulum, are composed of sarcomeres arranged mechanically in series along the myofibril length (Figure 2.3).

The main function of sarcoplasmic reticulum is to store the calcium ions involved in muscle contraction. The sarcomere, as the basic contractile unit of the muscle, is composed of two contractile proteins: myosin and actin (Figure 2.3).

It is thought that the muscle contraction happens when an actin filament slides along the myosin filament, known as the cross-bridge theory [126]. When a motor unit activates the muscle, action potential travels through the muscle fiber causing sarcoplasmic reticulum to release calcium ions. This increase in the calcium ion concentration initiates the cross-bridge formation between the actin and myosin fil-aments leading to muscle contraction. The relaxation of the muscle, on the other hand, happens due to the re-uptake of calcium ions into the sarcoplasmic reticulum.

Both the release and re-uptake of the calcium ions are not instantaneous, causing a delay (~5-50 milliseconds) between muscle excitation (i.e., firing of the motor unit) and muscle activation (i.e., concentration of calcium ions). The relationship between the excitation of a muscle and its activation is known as activation dynamics.

In addition to the activation dynamics, the force generation capacity of the skele-tal muscles is also limited to the current length and velocity of the muscle, termed as muscle contraction dynamics [127]. The total force generated by the muscle con-sists of an active force (due to excitation of the muscle) and a passive force (due to resistance of the muscle against elongation) [128]. Both the active and passive com-ponents vary nonlinearly with the muscle length, represented by the force-length curve (Figure 2.4-A). Moreover, the force generated by the muscle varies nonlinearly with the rate of changes of the muscle’s length, represented by the force-velocity curve (Figure 2.4-B). Tendons, as the connector elements between the muscles and

Figure 2.4: Schematic length-velocity relationship of a muscle and length relationship of a tendon. (A): Active (in blue) and passive (in orange) force-length relationship, (B): force-velocity relationship and (C): tendon force-force-length re-lationship of a skeletal muscle (i.e., contraction dynamics). Forces are normalized to the maximum isometric force of the muscle [127–133].

bones, also display a nonlinear force-length relationship (Figure 2.4-C). The elonga-tion in the tendon alters the whole muscle-tendon’s length that can potentially affect the muscle force-generation capacity, i.e., contraction dynamics [127].

The CNS controls the level of muscle activation, and consequently, the active muscle forces by altering the frequency and magnitude of stimulation passed to motor neurons [80]. Three groups of motor neurons (α,β, andγ) are distinguished by the targeted fiber type, diameter, and propagation velocity. The extrafusal muscle fibers are innervated at neuromuscular junctions along their length by branches ofα (the longer and faster) orβmotor neurons. The muscle spindle, on the other hand, is innervated byβorγmotor neurons for sensing muscle length and force [80,124,134].

Although the mechanism underpinning muscle structure and force production is relatively well-established, the muscle recruitment strategy controlled by CNS is not fully understood [135]. Hence, electromyography (EMG) is conventionally used to measure the train of action potentials along the muscles as a way of exploring the muscle activation patterns and CNS muscle requirement strategy [136].

3 KNEE JOINT MODELING AND SIMULATION

As briefly explained in chapter 2, the mechanical responses of the knee joint load-bearing tissues (i.e., cartilages and menisci) are attributed to the intricate combi-nation of macroscale and microscale biomechanical and biological factors. Sev-eral studies have provided convincing evidence that altered joint kinematics and kinetics, and consequently, altered tissue-level mechanical responses can effectu-ally lead to undesirable remodeling in the structure, composition, and metabolism of joint load-bearing tissue, and as a result, initiate MS disorders such as KOA [2–5,7,9,31,137–139]. Hence, joint realignment strategies (either non-invasive or sur-gical based), as well as rehabilitation protocols, have been developed to favorably re-distribute and manipulate joint kinematics and kinetics in the interest of prevent-ing or deceleratprevent-ing KOA and increasprevent-ing the joint’s functionality [1, 10, 16, 140, 141].

Yet, reliable assessment of these interventions requires a thorough knowledge of associated joint-level and tissue-level mechanical responses.

In vivo and in situ measurements of knee joint mechanics, such as JCF, con-tact area, and concon-tact pressure [4, 17–20, 22, 142, 143], are either limited to specific subjects (e.g., those with instrumented implants) or require highly invasive proce-dures. More importantly, experimental approaches are unable to measure crucial mechanical quantities such as stress, strain, or fluid flow within the tissue. Alterna-tively, computational methods have become a tool of choice to investigate multiscale biomechanics of the joint tissue. MS models, on the basis of rigid-body dynam-ics, have been developed and used to predict joint-level quantities such as muscle forces, JCF, and contact pressure. At lower spatial scales (e.g., tissue or cell level scales), continuum mechanics along with numerical approaches (predominantly the FE method) have been extensively utilized to describe the mechanobiology of the tissue. Nonetheless, rather few multiscale models are developed to simulate both joint-level and tissue-level joint mechanics while taking into account different subject-specific aspects such as muscle synergy, joint geometry, kinematics, and ki-netics. The most conventional modeling approaches with their potential advantages and limitations are discussed in this chapter.

3.1 MUSCULOSKELETAL MODELING

MS analyses provide multibody dynamic simulations on the basis of an integrated model of anatomy and the neuromusculoskeletal system. To this end, a multibody system (i.e., the MS model) composed of rigid or partly elastic bodies along with joint drivers (e.g., muscles), occasionally with passive resistances (e.g., ligaments), are defined according to the anatomical dimensions of the subject of interest. The MS model of a participant can be created by either scaling of a generic MS model [23], constructing a subject-specific MS model (e.g., using subjects MRIs) [30, 144], or a combination of both methods [24].

In the scaling approach, mass properties (i.e., mass and inertia tensor) and the body segments and length-dependent muscle properties of the MS model are lin-early scaled [23]. The scaling is executed based on relative distances between pairs

of markers attached to the subject (e.g., acquired from the motion capture system) and the corresponding marker locations in the MS model [23]. Cautions (such as precise placement of the model’s markers consistent with those attached to the sub-ject) are required when scaling an MS model to avoid possible scaling errors which may considerably influence the estimated kinematics and kinetics [145]. Although the scaling approach accommodates variations in subject size, subject-specific MS modeling (based on segmentation of the subject’s medical images) [30, 144,146] may be required to account for specific variations in the MS geometries, such as torsional deformities in individuals with cerebral palsy [146].

Nonetheless, MS models in general can provide two different analysis ap-proaches: forward dynamics (FD) and inverse dynamics (ID). FD analyses are in-tended to estimate movements of the body and the forces applied from the body to the environment (e.g., from the foot to the ground) for a given set of muscle activations. Nevertheless, estimation of physiologically realistic muscle activation patterns in FD problems requires sophisticated predictive algorithms, and conse-quently, is computationally expensive [147–150]. Conventionally, FD simulations are utilized to explore the principles underlying the locomotion and muscle syn-ergy [147–151]. It worth mentioning that the FD analysis was not used in this thesis, and the developed workflows (i.e., MS analyses) were based on the ID approach.

Hence, the FD approach is described briefly here.

In the second approach of the MS analyses, i.e., ID, the equations of motion are employed to describe the movement of the body segments, the external loads (i.e., forces and moments) applied on the joints, and finally, to estimate the mus-cle activations/forces required to counterbalance those external loads on the joint during a specific task. By using the equations of motion (i.e., Newton’s second law of motion), the force applied on a body (~F) is related to the rate of change in the momentum (~P) of the body over time (t). Hence, we can write:

~F= d~P

dt (3.1)

For an object with constant mass (m), Equation (3.1) can be re-stated as:

~F= d(m~v)

dt = md(~v)

dt (3.2)

where~vis the translational velocity of the object. Similarly, the net external moment on a body (~τ) is related to its moment of inertia relative to its center of mass (Ic) and its angular velocity (~ω):

~τ= d(Icω)~

dt (3.3)

Accordingly, ID is specifically used to determine forces and moments at each joint (i.e.,~Fand~τin Equations 3.2 and 3.3) responsible for a given body movements (i.e.,~v andω~ in Equations 3.2 and 3.3). Using inverse kinematics (IK), the body movement required by Equations 3.2 and 3.3 are obtained according to the experi-mental data (Figure 3.1). Several measurements techniques may be utilized to mea-sure body motion, such as radiography based devices (e.g., dual fluoroscopy) [152]

or reflective markers [153]. Radiography based motion capture may provide pre-cise motion data for a specific region of interest; however, reflective-marker based

Figure 3.1: (A): Inverse kinematics with which joint angles are determined while minimizing the error between the experimental markers (blue) and the correspond-ing model markers (pink). (B): Inverse dynamics with which the external moments (represented in local coordinate systems) on the joints are calculated.

approaches are more extensively used due to their lower costs, ease of use, and cap-turing a much wider field of view [154]. In this respect, a calibrated set of video cameras is used to capture the 3-dimensional (3D) location of the reflective mark-ers attached to the body. The same number of markmark-ers are also defined on the MS model of the subject in the corresponding anatomical locations. In some cases, tracking of an enforced joint angle may be of interest in addition to following the measured marker trajectories.

Finally, an optimization algorithm is utilized to estimate the joint coordinates (i.e., IK) of the MS model at which the experimentally measured joint angles and markers best match with the corresponding joint angles and markers on the MS model, separately at each time point of the motion. Although there are several op-timization algorithms [23, 155], the weighted least square opop-timization method has been most widely used in which the following objective function (fobj) is minimized separately at each time point [23]:

fobj=

∑

wi(~r^s_i −~r^m_i )²+

∑

wj(θ^s_j−θ^m_j )² (3.4) where~r^s_i and~r^m_i are the position of thei^th marker respectively on the subject and on the model and θ^s_j andθ^m_j are the values for the j^th joint of the subject and the model, respectively. wi andwj are factors to favorably weight the markers and joint angles, separately.

Following the IK, joint velocities and accelerations required by ID (i.e., Equations 3.2 and 3.3), and consequently, the external moments applied to each joint’s DoF can be calculated (Figure 3.1). Next, these external moments should be resolved into the muscle forces/activations, as follows:

τ_i =

∑

f_ir_j,i (3.5)

whereτ_i is the external moment applied to thei^thDoF of the model (i.e., obtained from ID, f_j is the force that the j^th muscle of the model should generate, r_(j,i) is

the moment arm of the j^th muscle around the rotation axis of the i^th DoF of the model, and N is the total number of muscles acting on the i^th DoF of the model.

The fj in Equation (3.5), which will be discussed later in this section, is generally a function of the muscle activation level, the maximum isometric force of the muscle, and current length, and contraction velocity of the muscle-tendon unit (replicating muscle contraction dynamics).

Muscle forces, and consequently, muscle activations, are the unknowns in the Equation (3.5). However, Equation (3.5) is underdetermined since almost all the DoFs of the body joints are driven by more than one muscle (e.g., agonist, antago-nist, etc.) and, more importantly, the muscle recruitment strategy controlled by the CNS is far from being fully understood [135]. Accordingly, Equation (3.5) should be solved by considering two aspects. First, to resolve the external moments into the muscle forces, and second, to respect the muscle activation dynamics and con-traction dynamics in order to estimate physiologically feasible muscle forces. As a result, the calculation of the muscle forces/activations (Equation 3.5) relies on optimization-based methods, i.e., the assumed criterion of the recruitment strategy of the CNS [23, 24, 27, 135, 156] along with muscle activation/contraction dynamics.

The employed optimization algorithms may estimate muscle activations according to the states of the MS model (e.g., static-optimization-based MS models), or alter-natively they might be assisted with extra experimental inputs such as measured EMGs (i.e., EMG-informed MS models). These two approaches will be discussed in more detail in the following sections.