Sarcomere length (SL) instability and SL non-uniformity have been used to explain fundamental properties of skeletal muscles, such as creep, force depression following active muscle shortening and residual force enhancement following active stretching of muscles. Regarding residual force enhancement, it has been argued that active muscle stretching causes SL instability, thereby increasing SL non-uniformity. However, we recently showed that SL non-uniformity is not increased by active muscle stretching, but it remains unclear if SL stability is affected by active stretching. Here, we used single myofibrils of rabbit psoas muscle and measured SL non-uniformity and SL instability during isometric contractions and for isometric contractions following active stretching at average SLs corresponding to the descending limb of the force–length relationship. We defined isometric contractions as contractions during which mean SL remained constant. SL instability was quantified by the rate of change of individual SLs over the course of steady-state isometric force and SL non-uniformity was defined as deviations of SLs from the mean SL at an instant of time. We found that whereas the mean SL remained constant during isometric contraction, by definition, individual SLs did not. SLs were more stable in the force-enhanced, isometric state following active stretching compared with the isometric reference state. We also found that SL instability was not correlated with the rate of change of SL non-uniformity. Also, SL non-uniformity was not different in the isometric and the post-stretch isometric contractions. We conclude that since SL is more stable but similarly non-uniform in the force-enhanced compared with the corresponding isometric reference contraction, it appears unlikely that either SL instability or SL non-uniformity contribute to the residual force enhancement property of skeletal muscle.
Skeletal muscles are made of billions of micrometer-sized sarcomeres that are organized within a complex structural network of connective tissue. Sarcomeres are the basic building blocks of skeletal muscles and are made of actin (thin) and myosin (thick) contractile proteins, as well as a series of structural proteins, such as titin, nebulin and desmin. For decades, the functional properties of skeletal muscles have been treated like an up-scaled version of a single sarcomere or half-sarcomere (Huxley, 1957; Huxley and Simmons, 1971; Huxley and Tideswell, 1996; Lombardi and Piazzesi, 1990). One such property is the sarcomere force–length relationship. According to the sliding filament (Gordon et al., 1966) and cross-bridge (Huxley, 1957) theories, the maximal isometric force a sarcomere can produce depends on the amount of overlap between the thin and thick myofilaments, and can be described by an inverted ‘U-shaped’ force–length curve. This sarcomere force–length relationship can be theoretically derived from the lengths of the thin and thick filaments (Gordon et al., 1966; Herzog et al., 1992) and has been validated in amphibian (Gordon et al., 1966) and mammalian muscles (Edman, 2005).
The sarcomere force–length relationship was originally derived from isometric contractions of single intact frog muscle fibres, during which the forces and lengths of sarcomeres were assumed implicitly to stay at constant values (or in a stable condition). Stability of sarcomere forces and lengths has not been challenged for the so-called ascending and plateau regions of the force-length relationship, but has been said to not exist for sarcomeres working on the descending limb of the force–length relationship. Specifically, Hill (1953) introduced the concept of muscle force and segment length instability on the descending limb of the force–length relationship and stated that for this condition: ‘… the system would be unstable and a long-continuous creep would be observed’. The notion of sarcomere force and sarcomere length (SL) instability was further extended by observations of ‘creep’, which is a slow continuous increase in force observed during isometric contractions, and was related to SL non-uniformity by Gordon et al. (1966), who explained creep explicitly by stating: ‘There can be little doubt that this [the creep behaviour] is due to the progressive development of irregularities of striation spacing [SL non-uniformity] during the tetanus which is to be expected because of the instability represented by the negative slope of this part [the descending limb] of the length–tension relation’.
SL instability and SL non-uniformity have become an accepted part of muscle contraction, and have been used to explain a series of experimentally observed results, including creep (mentioned above), force depression following muscle shortening (Morgan et al., 2000) and the residual force enhancement property following stretching of an active muscle (Edman and Tsuchiya, 1996; Morgan, 1990, 1994). Of particular interest is the residual force enhancement property, because it cannot be explained within the framework of the classic cross-bridge theory (Herzog, 2018b; Walcott and Herzog, 2008). Residual force enhancement is defined as the additional isometric steady-state force obtained following elongation of an activated muscle (eccentric action) compared with the corresponding (same length, same activation) steady-state force obtained in an isometric contraction (Edman et al., 1982). Residual force enhancement is explained by the theory of SL non-uniformity (Morgan, 1990, 1994) as being caused by the development of SL non-uniformity during stretching of an active muscle on the ‘unstable’ descending limb of the force–length relationship (Morgan, 1990, 1994). It has been shown more recently (Johnston et al., 2016; Joumaa et al., 2008b, 2018) that SL non-uniformity following active stretching of a muscle was comparable to that during an isometric contraction on the descending limb of the force–length relationship. However, SLs barely changed following active stretching of a muscle, suggesting that the SLs were stable (Joumaa et al., 2008b). Owing to the lack of systematic studies into SL stability and SL non-uniformity following an active stretch and during isometric contractions, the intricate causal relationship between SL stability and SL non-uniformity remains unexplained.
Therefore, the purpose of this study was to measure SL instability and SL non-uniformity within the same myofibril during isometric contractions and during isometric contractions following active stretching at average SLs corresponding to the descending limb of the force–length relationship. We defined SL instability as the rate of change of individual SL over the course of steady-state isometric contractions and SL non-uniformity as deviations of SLs from the mean SL at an instant of time. This aim was accomplished using single myofibril preparations (Bartoo et al., 1993; Joumaa et al., 2008a; Leonard and Herzog, 2010; Rassier et al., 2003) where sarcomeres are arranged strictly in series and individual SLs can be measured accurately. We hypothesized that SLs are more stable during isometric contractions than during isometric contractions following active stretch and that SL instability can explain the occurrence of SL non-uniformity.
MATERIALS AND METHODS
This study involves additional analyses of previously published data (Johnston et al., 2016). While the experimental setup was the same as in the previous paper, the data analysis was new and targeted to unravel the relationship between SL instability and SL non-uniformity during isometric contractions, and isometric contractions following active stretching.
Myofibril preparation and experimental set-up
All aspects of animal care and experimental protocol were approved by the Life and Environmental Sciences Animal Care Committee of the University of Calgary. The rigor, relaxing and activating solutions used here were as described previously (Johnston et al., 2016) and consisted of: (1) rigor solution: 50 mmol l−1 Tris, 100 mmol l−1sodium chloride, 2 mmol l−1 potassium chloride, 2 mmol l−1 magnesium chloride and 10 mmol l−1 ethylene glycol bis(2-aminoethyl ether)-N,N,N'N'-tetraacetic acid (EGTA), pH 7.0; (2) relaxing solution: 10 mmol l−1 3-(N morpholino)propanesulfonic acid (MOPS), 64.4 mmol l−1 potassium proprionate, 9.45 mmol l−1 sodium sulfate, 5.23 mmol l−1 magnesium proprionate, 2 mmol l−1 potassium EGTA, 7 mmol l−1 adenosine triphosphate (ATP) and 10 mmol l−1 creatine phosphate, pH 7.0; (3) activating solution: 10 mmol l−1 MOPS , 45.1 mmol l−1 potassium proprionate, 5.21 mmol l−1 magnesium proprionate, 9.27 mmol l−1 sodium sulfate, 1 mmol l−1 sodium EGTA, 7 mmol l−1 ATP, 10 mmol l−1 creatine phosphate, 0.75 mmol l−1 calcium chloride at pCa=3.12 and pH 7.0.
The Life and Environmental Sciences Animal Care Committee (LESACC) of the University of Calgary approved the ethics protocol. Eight 6-month-old female New Zealand white rabbits were used in this study. The animals were euthanized by an intravenous injection of sodium pentobarbital prior to tissue harvesting. Strips of psoas muscle were harvested from freshly euthanized animals and tied to wooden sticks to preserve the in situ fibre length. Samples were stored in a rigor solution containing glycerol (50:50 v/v, pH 7.0) and protease inhibitors (Complete, Roche Diagnostics, Quebec, Canada) and were stored at −20°C for 10–20 days.
On the day of experiments, muscle strips were removed from the freezer, homogenized and placed in a relaxing solution at 4°C. Individual myofibrils (N=12) containing 5–23 sarcomeres and with a good striation pattern were selected and attached to a glass needle at one end and one of the cantilever pairs at the opposite end (Fig. 1A). The glass needle was attached to a radial piezotube motor (Boston Piezo-Optics Inc., Bellingham, MA) that allowed for automated myofibril length changes while the pair of cantilevers (stiffness, 132 nN μm−1) served to quantify force production (Johnston et al., 2016; Joumaa et al., 2007, 2008b; Leonard et al., 2010a; Powers et al., 2014; Rassier et al., 2003).
Experiments were conducted using an inverted phase-contrast microscope (Zeiss Axiovert 200M, Germany) equipped with a 100×/1.3 NA oil immersion objective and a 2.5× Optovar. A complementary metal-oxide semiconductor (CMOS) camera (Rolera Bolt, Quantitative Imaging Corp., Surrey, Canada) attached to the microscope was used to record all of the experiments at 30 Hz. As the optical resolution of the microscope system is ∼0.25 µm, the pixel size was set at 87 nm pixel−1 to satisfy the Nyquist sampling criterion.
Automated length changes of myofibrils were achieved by a custom-written LABVIEW (National Instruments Corp., Austin, TX, USA) program that controlled the radial piezotube motor, and thus the position of the attached glass needle. The length of myofibrils was controlled throughout the experiment, but the experimental protocol is described in terms of the mean SL for individual myofibrils (Fig. 1B). Briefly, the myofibrils were set to a starting length of 2.5 μm per sarcomere, which represented the ‘passive short’ (PS) state. They were then passively stretched (0.1 μm s−1 per sarcomere) to 3.2 μm per sarcomere to exhibit the ‘passive long’ (PL) state. At this length, the myofibrils were activated using a jetted fluid (Jacuzzi) technique, which delivered a calcium-based activating solution uniformly to the myofibrils in order to elicit the ‘active long’ (AL) state (Johnston et al., 2016). Owing to the force production and the compliance of the cantilever, the myofibrils shortened to an average SL of 2.8 μm. Following activation, myofibrils were rapidly shortened (0.8 μm s−1 per sarcomere) to 2.2 μm per sarcomere and held for ∼8 s in the ‘active short’ (AS) state. Myofibrils were then stretched (0.1 μm s−1 per sarcomere) back to the same motor position as in the AL state, finishing in the ‘active post-stretch’ (APS) state. All experiments were performed at room temperature (21°C). The hold time for every steady state varied between 8 and 13 s, except for the steady state of PS (5 s, see Results). The distribution of mean SLs in individual myofibrils are presented in Fig. S1.
Analysis and definitions
The centroids of each sarcomere A-band within a myofibril were identified, and the deflection of the cantilever holding the myofibril relative to the free reference cantilever was tracked using a custom-written MATLAB (MathWorks, Natick, MA, USA) code. Individual SLs (N=137) were defined as the distance between the centroids of adjacent A-bands. Force was calculated by multiplying the deflection of the cantilever from its reference position by its stiffness (132 nN μm−1), and was then normalized to myofibril cross-sectional area using the myofibril diameter at length equivalent to 2.5 μm per sarcomere and expressed in units of stress (nN μm−2). Time-history graphs for individual and mean SLs and stress were generated for each myofibril through video analyses completed at a nominal frequency of 3 Hz. From the mean SLs, the five aforementioned steady states (PS, PL, AL, AS, APS) were defined for each myofibril. Transient SL changes were observed during all myofibrillar length changes and upon fixed-end activation of the myofibrils (Fig. 2A). We only measured the isometric (in terms of mean SL), steady-state phase between 700 ms after the mean SL stabilized to 700 ms prior to the start of the next transient phase, respectively (Fig. 2B). Mean SL and stress were determined for each myofibril and for each steady-state condition.
SL instability was defined as the transient length change of sarcomeres during the steady-state conditions, and was quantified as the absolute slope of the best-fitting linear approximation to the individual SL-time curves (nm s−1). A slope of zero (0 nm s−1; that is, no net SL change during the steady-state phase) was considered a perfectly stable sarcomere; while an increase in slope (greater SL change) was associated with an increased instability.
SL dispersion (non-uniformity)
SL dispersion was quantified during each of the five steady-state conditions (PS, PL, AL, AS, APS) and was expressed in terms of the standard deviation (SD) and the coefficient of variation (CV=s.d. mean−1) from the mean SL. The rate of change of SL dispersion was quantified as the slope of the best-fitting straight line to the coefficient of variation vs time data.
Correlation between SL instability and SL non-uniformity
Previous studies found increased SL non-uniformity upon Ca2+ activation (Joumaa et al., 2008b, 2018). In order to investigate if SL non-uniformity was associated with SL instability, the sum of SL instabilities during the steady states of AL, AS and APS of activated myofibrils normalized by the total number of sarcomeres present in individual myofibrils was compared, respectively, with the corresponding coefficient of variation of SLs at the start of the steady state of isometric contractions as well as with the corresponding rate of change of the coefficient of variation of SLs to reveal any correlation between these variables.
Correlation between SL instability and sarcomere position along a myofibril
The myofibrils contained 5–23 sarcomeres in series. In order to investigate if sarcomeres located near the cantilevers have higher SL instability than sarcomeres in the middle of a myofibril, the position of each sarcomere along a myofibril was defined relative to the total myofibril length, with a value of 1.0 assigned to the end sarcomeres and 0.0 assigned to the middle (odd number of total sarcomeres) or two middle (even number of total sarcomeres) sarcomeres. SL instability was compared with sarcomere position along a myofibril for any correlation.
Statistical analyses were performed using SPSS (version 25, SPSS Inc.). The means of sarcomere instability, s.d. and CV of SLs were analysed for the five steady-state conditions using a generalized estimating equation (GEE, under Genlin procedures in SPSS) to take into account the correlated nature of the observations and the unbalanced study design. The dependence of sarcomere instability on the position of sarcomeres along the myofibril was analysed with GEE for the PS, PL, AL, AS, and APS steady states. The relationship between SL instability and SL non-uniformity as well as between SL instability and rate of change of SL non-uniformity was, respectively, analysed with non-parametric Spearman's rank-order correlation in AL, AS and APS steady states. All statistical tests were performed using a two-sided approach, with a type I error, α, set at 0.05. Multiple comparisons were accounted for through Bonferroni adjusted P-values. Since the data analysed here were interval data, normality was not needed for the GEE approach. Also, since all measurements were repeated measures, no estimate of within-group variation was required. Unless otherwise stated, results were expressed as estimated marginal means (EMM)±1 s.e. (N=8), which were estimated by using GEE.
The average cross-sectional area of the myofibrils was 1.5±0.1 µm2 at a mean SL of 2.5 µm. A summary of the SLs and stresses for the 137 sarcomeres in each of the five steady states is provided in Table 1. The mean duration of the steady state conditions for the PS, PL, AL, AS and APS states were 4.8±1.5 s, 10.6±2.9 s, 11.5±3.7 s, 8.2±0.7 s and 12.9±2.7 s, respectively. Myofibril forces were enhanced by 67% on average during isometric contractions following active stretch compared with regular isometric contractions (Table 1, Fig. 2).
A representative example of length changes as a function of time for sarcomeres from a single myofibril is shown in Fig. 2. Although the mean SL was kept constant during the steady-state conditions, only a few sarcomeres remained at a perfectly constant length. In the active steady states, most sarcomeres continuously shortened or elongated slowly, especially during the AL state.
Prior to activation (PS and PL states), SL instability showed baseline values of 11.1±1.1 nm s−1 (Fig. 3A). When activated at 3.2 µm (AL state), SL instability increased to 67.5±14.3 nm s−1 (Fig. 3A). Following shortening and stretching of the active myofibrils (AS and APS states), SLs were more stable than in the isometric state (AL) with mean instability values of 28.8±5.1 nm s−1 and 39.4±7.0 nm s−1, respectively.
SL dispersions were described in terms of standard deviations and coefficients of variation (Fig. 4). For the passive conditions (PS, PL), the CV was approximately 8%. For the active conditions, SL dispersion ranged from 26 to 30% with no significant difference between the AL, AS and APS conditions (Fig. 4B). During steady states, sarcomeres in the AL state increased in length dispersion at a rate of 0.7% s−1. Although the AS and APS states had slower average rates of change of CV at 0.2% s−1, they were not significantly different from the AL state.
The sum of SL instabilities was significantly correlated with the corresponding coefficient of variation of SLs at the start of a steady state for the AL (r=0.79) and AS (r=0.73) states (P<0.01, Fig. 5A,B), but not for the APS state (P=0.191, Fig. 5C). The rate of change of SL non-uniformity was not correlated with the sum of SL instabilities for any of the steady states (AL, AS or APS; P>0.05, Fig. 6). Finally, there was no relationship between the amount of SL instability and the position of sarcomeres along myofibrils for any of the five steady-state conditions.
This study was motivated by previous findings showing that approximately 70% of sarcomeres did not change lengths, and thus showed length stability, in isometric contractions after active stretching (Joumaa et al., 2008b). Here, we used a novel experimental protocol that allowed direct comparison of SL instability and SL non-uniformity in isometric contractions and isometric contractions following active stretching within the same myofibrils. The 8-s isometric contractions following active shortening (AS) that were inserted between the AL and APS steady states allowed the aforementioned comparison and has been shown to not influence the residual force enhancement in the APS state (Fortuna et al., 2016; Rassier and Herzog, 2004). Stresses in the isometric reference contractions (35-151 nN µm−1) were comparable to values observed previously (95–160 nN µm−1) (Herzog et al., 2010; Joumaa and Herzog, 2010; Joumaa et al., 2008b; Leonard and Herzog, 2010; Pavlov et al., 2009).
The primary results of this study were: (1) SL instability was greater for isometric contractions compared with isometric contractions following active muscle stretching (Fig. 3); (2) SL instability was correlated with the initial establishment of SL non-uniformity only for the AL and AS steady states, but not for the APS steady state (Fig. 5), nor was it correlated with the rate of change of SL non-uniformity during the active steady states (AL, AS and APS; Fig. 6). These results have important implications for our understanding of muscle function. The first result indicates that SL stability is enhanced by active muscle stretching; the second result shows that the contraction dynamics for the isometric reference contractions were different from the force-enhanced isometric contractions, and that SL instability did not lead to increased SL non-uniformity over time.
Hill (1953) was among the first who stated that muscle force and muscle segmental lengths were unstable on the descending limb of the force-length relationship. His argument was based on the negative slope of the descending limb of the force–length relationship and implying negative muscle stiffness. However, the force–length relationship was obtained through a series of independent isometric contractions over a range of muscle lengths. It represents a static property, rather than a dynamic response evoked during muscle stretching. When an active muscle is stretched, force always increases regardless of its initial length, thus exhibiting positive stiffness (except for cases when stretching is done at very high speeds and the muscle is thought to undergo so-called ‘slippage’ (Bagni et al., 2005; Fukutani et al., 2019; Griffiths et al., 1980). It is also because of this positive stiffness that steady-state forces can be obtained easily for isometric contractions on the descending limb of the force–length relationship (Herzog and Leonard, 2002). Arguably the most puzzling and least understood mechanical property of skeletal muscle contraction is the residual force enhancement (Abbott and Aubert, 1952; Edman et al., 1982). Residual force enhancement increases with increasing stretch magnitude (Bullimore et al., 2007; Edman et al., 1982), but is virtually independent of stretch speed (Lee and Herzog, 2003; Morgan et al., 2000). It is associated with a decrease in metabolic energy cost (Joumaa and Herzog, 2013), can be abolished instantaneously by deactivating muscles (Abbott and Aubert, 1952; Julian and Morgan, 1979) and is associated with a passive structural component, the so-called passive force enhancement (Herzog and Leonard, 2002). Residual force enhancement has been observed on all structural levels of muscle, from single sarcomeres (Leonard et al., 2010b) to single fibres (Edman et al., 1982; Powers et al., 2017), isolated muscles (Herzog and Leonard, 2002) and in vivo human skeletal muscles measured by an increase in force (Lee and Herzog, 2002) or a decrease in activation for the same force (de Brito Fontana et al., 2018; Oskouei and Herzog, 2005).
SL instability, and associated development of SL non-uniformity, has been used as an explanation for residual force enhancement (Edman and Tsuchiya, 1996; Morgan, 1990, 1994). In the SL instability/non-uniformity theory, residual force enhancement is achieved ‘… by extremely rapid and uncontrolled lengthening of individual sarcomeres’, and this rapid, uncontrolled lengthening was thought to be caused by ‘…the instability inherent in the descending limb of the length-tension curve’ (Morgan, 1994). In the current study, the average residual force enhancement was measured to be approximately 67% [Table 1, see also Johnston et al. (2016) for more detail]. But this substantial amount of additional force following active stretching of the myofibrils was achieved with no difference in SL non-uniformity (Fig. 4B), and a significantly increased stability of SLs (Fig. 3A) in the force-enhanced compared with the corresponding isometric reference state. Our results are consistent with previous studies which found that active stretching of muscle had a stabilizing effect on the SLs (Edman et al., 1982; Telley et al., 2006a, 2006b). Also, we found that SL instability was not correlated with the initial SL non-uniformity, or the rate of change of SL non-uniformity for the force-enhanced isometric contractions (Figs 5 and 6); thus, further strengthening the idea that SL instability does not necessarily lead to SL non-uniformity, and vice versa. Therefore, the development of residual force enhancement involves more complicated processes than those suggested by the SL non-uniformity theory.
In the current study, the observation that sarcomeres are more stable following an active stretch than during an isometric reference contraction at the same SL is novel, and deserves particular attention. First, this result is in direct contradiction to the underpinnings of the SL non-uniformity theory. Second, it suggests that resistance to SL changes is increased following active stretching compared with the corresponding isometric reference contraction. This result might be explained by the engagement of structural braces, such as titin proteins (Herzog, 2018a; Herzog et al., 2015), that provide stability within and between sarcomeres in a myofibril. It has been suggested that such bracing could be accomplished either by an increase in titin stiffness upon muscle activation and stretch, or by a shortening of the effective spring length of titin by binding to other rigid structures of the sarcomere. Both these mechanisms have been proposed (Forcinito et al., 1998; Noble, 1992; Rode et al., 2009), but the molecular details still need verification (Herzog, 2016, 2017, 2018a).
Attention should be paid to the definition of stability or instability that we introduced here. Nominally, we defined instability by the sarcomere length change over several seconds of steady-state force production and constant average SL. Technically, this is not a mechanical definition of (in)stability but we chose it in order to compare our results with the literature. The results presented here were obtained from single myofibril preparations, and may not hold at higher structural levels, such as in whole muscle where transverse connections between sarcomeres, fibres, fascicles and muscles may influence sarcomere dynamics (Lichtwark et al., 2018; Llewellyn et al., 2008; Moo and Herzog, 2018; Moo et al., 2016, 2017a,b).
Based on the results of this study, we conclude that: (1) SLs are more stable in the force-enhanced compared with the isometric reference state, and (2) SL instability does not necessarily correlate with SL non-uniformity and they should therefore not be treated as the same phenomena. These results have profound implications in how we interpret basic mechanical properties of sarcomeres and muscles and place a need for a causal explanation of residual force enhancement, and the reconciliation of this property to the overall function of skeletal muscles.
We would like to thank Dr Tim Leonard for help with data collection and for providing the nanofabricated cantilevers built at the NanoScale Facility of Cornell University (Ithaca, NY), which is supported by the National Science Foundation. We would also like to thank Dr Tak-Shing Fung for help in statistical analysis and Dr Erica Buckeridge for help with data analysis.
Conceptualization: K.J., E.K.M., A.J., W.H.; Methodology: K.J., E.K.M., A.J., W.H.; Software: K.J., E.K.M., A.J., W.H.; Validation: K.J., E.K.M., A.J., W.H.; Formal analysis: K.J., E.K.M., A.J., W.H.; Investigation: K.J., E.K.M., A.J., W.H.; Resources: W.H.; Data curation: K.J., E.K.M., A.J., W.H.; Writing - original draft: K.J., E.K.M., A.J., W.H.; Writing - review & editing: K.J., E.K.M., A.J.; Visualization: K.J., E.K.M., A.J., W.H.; Supervision: W.H.; Project administration: W.H.; Funding acquisition: W.H.
This research was supported by an Alberta Innovates – Health Solutions (AI-HS) postdoctoral fellowship (grant number: 10013510), a Canadian Institutes of Health Research (CIHR) postdoctoral fellowship (grant number: 10010953), Alberta Innovates – Technology Futures (AITF), the Natural Sciences and Engineering Research Council of Canada (NSERC) Collaborative Research Training Experience (CREATE) programme, the Canada Research Chair Program for Molecular and Cellular Biomechanics and The Killam Foundation.
All data and the MATLAB source code used for obtaining SLs and myofibril forces have been uploaded onto ResearchGate repository and are accessible at: https://www.researchgate.net/publication/335777136_experimental_data_and_Matlab_code_used_in_Johnston_Moo_et_al_JEB_2019.
The authors declare no competing or financial interests.