Rate of lumbar paravertebral muscle fat infiltration versus spinal degeneration in asymptomatic populations: an age-aggregated cross-sectional simulation study

Background The spinal column including its vertebrae and disks has been well examined and extensively reported in relation to age-aggregated degeneration. In contrast, paravertebral muscles are poorly represented in describing normative degeneration. Increasing evidence points to the importance of paravertebral muscle quality in low back health, and their potential as a modifiable factor in low back pain (LBP). Studies examining normative decline of paravertebral muscles are needed to advance the field’s etiological understanding. With a novel approach and based on published data, we establish and compare decline rates of imaging features for degeneration of lumbar vertebrae and disks, versus fatty infiltration in paravertebral muscles in asymptomatic adults. Methods Our cross-sectional simulation study examined age-aggregated data from three published studies who reported on asymptomatic adults spanning 18–60 years. Prevalence rates of imaging degenerative features of the spinal column were examined via logistic regression and compared with percentage fatty infiltration in erector spinae, multifidus and psoas using synthetic data and Monte Carlo simulation with 10,000 endpoint-specific regression iterations. General linear regression models were employed to estimate marginal effects of age reported as a one-year change rate (with 95 % confidence intervals) for comparisons between all reported spinal features. Results Declines in multifidus (0.24 & 0.11 %/year), erector spinae (0.13 & 0.07 %/year), and psoas (0.04 %/year) occur at similarly slow rates to disk protrusion (0.25 %/year), annular fissure (0.15 %/year), and spondylolisthesis (0.29 %/year). Multifidus showed a trend for faster decline than erector spinae, particularly in men. Of the features examined, disk signal loss declined fastest, and psoas muscle the slowest. Conclusions Degeneration of lumbar paravertebral muscles occurs slowly in asymptomatic adults, with a tendency to be most pronounced in multifidus. Rate of decline of spinal structures represents a novel variable that warrants inclusion as a known feature of the expected degenerative cascade, and to provide a basis for comparison to diseases of the spine in research and clinical practice. Concurrent examination of spinal features using advanced imaging to improve muscle analysis would be a strong addition to the field.


Background
Low back pain (LBP) causes more global disability and lost healthy years than any other condition [1]. It is a common health problem forecast to have a wider societal impact alongside an increasingly ageing population [2]. Identifying modifiable risk factors associated with LBP is important in influencing the disease, and necessary in understanding its etiology to develop targeted and effective interventions. While LBP is a multifactorial condition, lumbar paravertebral muscles are receiving increased attention for their therapeutic, diagnostic and prognostic (referred to as 'theranostic') potential.
Imaging features of spinal degeneration that are associated with LBP are also prevalent in asymptomatic individuals [3]. Determining the natural history of ageing is fundamental in understanding the spine in both health and disease. Despite an established association between paravertebral muscle degeneration and LBP [4][5][6][7][8][9][10], only a modest literature describes the role of these muscles in normal ageing [11,12]. The importance of paravertebral muscles in optimizing back health is increasingly acknowledged and investigations determining normative degenerative change and muscles' capacity to influence the course of recovery of LBP are needed.
In their meta-analysis of imaging studies, Brinjikji and colleagues [3] describe the prevalence of spinal column degenerative features in asymptomatic adult spines according to decade, but without reference to muscle. This omission likely reflects the paucity of available literature in describing normal (or asymptomatic) ageing of paravertebral muscles [11][12][13], and thus a relative lack of recognition towards non-invasive and reproducible quantitative measurement of soft-tissues exists in clinical practice. Traditional studies examining spinal muscle quality in LBP have employed qualitative grading based on the method of Kjaer et al. [7], or various semiquantitative and quantitative methods that determine fat proportion within a cross-sectional area and defined region of interest based on counting pixel number and signal intensity [4,14,15]. As an advance from this practice, multi-echo imaging techniques like Dixon [16,17] and proton-density methods [18] are, or should be, considered for contemporary studies as quantification of skeletal muscle fat content based on these imaging sequences is shown to have superior accuracy, visualization, speed of acquisition and utility [5,[19][20][21][22][23]. Comparing the wealth of literature reporting the natural history of degeneration of the vertebra and disks, with the recently promoted quantification studies describing muscle change is therefore somewhat difficult.
Fatty infiltration (FI) in skeletal muscle is an accepted feature of declining muscle structure and quality [22,24]. While studies have identified increasing paravertebral FI with age in volunteers with [6,[25][26][27] and without [11][12][13] LBP, its etiology in normative decline, and relationship to other features of spinal degeneration, is poorly understood. Divergent from an earlier understanding, Hodges et al. [28] describe muscle FI as a feature of structural remodeling rather than muscle atrophy in the ovine multifidus (MF) after experimental and controlled injury involving the intervertebral disc. Human lumbar paravertebral muscle fat content has been associated with spinal features of degeneration (e.g. facet joint osteoarthritis [6], spondylolisthesis [6], disc narrowing [6,9], and type 2 Modic change [9]). However, further information regarding the relationship between normative declining muscle quality and other features of spinal degeneration is necessary to better understand their combined impact, causation and potential for change with clinicallyrelevant and more informed interventions.
In this cross-sectional study, we use a novel simulation approach to compare age-related degenerative imaging features of the spinal column and paravertebral muscle FI of asymptomatic volunteers as documented in three published studies: Brinjikji et al. [3] describe age-specific prevalence rates of spinal column degeneration based on a meta-analysis of 33 imaging studies reporting across seven decades of life   [3]; Valentin et al. [12] report percentage of FI for psoas major (psoas), MF and erector spinae (ES) comparing two age-groups (18-25 with 45-60 years); and Crawford et al. [11] report FI in MF and ES across four decades (20-29, 30-39, 40-49, 50-60 years). Both muscle studies employed a standard multi-echo MRI sequencing that has been validated against biopsy and across magnetic fields [29]. Our primary aim was to determine a yearly rate of decline of muscle tissue (defined as increased percent FI) as compared with reported degenerative changes of the spinal column, given the latter typically benchmarks the natural history of age-related spinal change.

Methods
To enhance comparability between the studies, we included published data for each study that represented subjects between 18 and 60 years. Based on published age-aggregated information (Table 1), we determined the marginal effect of age with corresponding 95 % confidence intervals (CI) for each degenerative feature as follows: For the Brinjikji et al. [3] study that reported prevalence rates for eight degenerative features, a sample of observations for each of the degenerative features corresponding to the summed total of subjects in the first four age groups from their study was created. All statistical analyses, including simulations, were carried out in Stata 14 (StataCorp, College Station, TX).
More specifically, we created samples of n = 1603 (disk degeneration), n = 613 (disk signal loss), n = 572 (disk height loss), n = 430 (disk bulge), n = 1408 (disk protrusion), n = 996 (annular fissure), n = 849 (facet degeneration), and n = 284 (spondylolisthesis). Subjects within each sample were distributed across the four age groups of 20-50s according to the sample of subjects reported by Brinjikji et al. [3] in the respective age groups; reported prevalence of degeneration per age group (%) was then used to assign corresponding proportions of 1 (degenerative) and 0 (not degenerative) to subjects within each age group. For example, synthetic data for disk degeneration comprised a total sample of 1603 subjects (273 + 604 + 415 + 311). In age group 20, 101 of 273 subjects were assigned values 1 (37 % degenerative) and 172 subjects were assigned values 0 (63 % not degenerative). Similarly, 314 of 604 (52 %), 282 of 415 (68 %) and 294 of 311 (80 %) subjects in age groups 30, 40 and 50 respectively were classified as having a disk degeneration (Table 1). Logistic regression analyses were then used to estimate prevalence rates from the synthetic data and to determine the marginal effects of age together with CIs for each degenerative feature. For two degenerative features (facet degeneration and spondylolisthesis), no information on the number of subjects in the age-groups 20-29 and 30-39 was available; we assumed a sample size of 100 in each of these two groups. By way of justification, this sample size is roughly twice the size of any reported age-specific sample in the spondylolisthesis group, and of four out of five reported age-specific samples in the facet degeneration group. This yielded CIs that were narrower than they would have been had we assumed sample sizes similar to the reported ones. Consequently, choosing higher sample sizes made it more challenging to demonstrate that marginal effects of age were not different between end-points, i.e. showed a similar rate of decline.
Similarly, data from Valentin et al. [12] and Crawford et al. [11] were simulated from information on their published age-aggregated data ( Table 1). The original age-and sex-specific sample sizes, means and standard deviations for each muscle's FI% were used in a Monte Carlo simulation where 10,000 samples of age-specific normal random variates with sample size, mean and standard deviation were drawn equivalent to the original studies. General linear regression models (Gaussian) were then used to estimate marginal effects of age with corresponding CIs. Figure 1 depicts the general flow of the simulation for each endpoint of the Crawford/Valentin muscle studies.
For the Crawford study [11], we created samples of n = 80 subjects for each endpoint (i.e. MF + ES). Subjects within each sample were distributed across the four age groups 20-29 to 50-60 according to the Table 1 Reference data used for age-specific prevalence estimates (%) for degenerative features of the spinal column (Brinjikji et al. [3]) and paravertebral muscle fat (%; ± standard deviation; Crawford et al. [11] and Valentin et al. [12] MF + ES multifidus and erector spinae combined, n total sample size a Sample size of 100 assumed in simulation sample of subjects reported by Crawford et al. [11] in each age group (n = 20) and reported means and standard deviations within each age group (Table 1) were used to create corresponding normal random variates using Stata's rnormal command. Similarly, we created samples of n = 24 subjects for each endpoint in the Valentin study [12] (psoas, MF + ES) and subjects within each sample were distributed across the two age groups 18-25 and 46-60 according to the sample of subjects reported by Valentin et al. [12] in each age group (n-12). Again, reported means and standard deviations within each age group (Table 1) were used to create corresponding normal random variates. For each endpoint reported in the Crawford [11] and Valentin [12] studies, the process of data generation and a corresponding regression analysis was repeated 10,000 times in order reduce the effect of noise in the simulated data. For example, the total Crawford MF + ES sample of 80 subjects (20 + 20 + 20 + 20) in Table 1 was simulated by generating 20 normal random MF + ES variates with mean 12.4 and standard deviation 3.7 for age group 20-29 using the Stata rnormal (12.4, 3.7) command. Similarly, random MF + ES variates with mean (standard deviation) 15.9 (5.0), 15.7 (4.5) and 18.1 (7.8) were generated for the remaining age groups 30-39, 40-49, 50-60. The generated MF + ES variates were then regressed on age. The MF + ES simulation process, meaning sample generation and regression, was iterated 10,000 times. Where appropriate, a similar simulation procedure was used to estimate gender-specific rates of decline.
Likelihood ratio test was used to make comparison between saturated models using original age-group factor variables, with models employing a linear age predictor derived from the mid-points of the original age-groups. All tests showed that the linearity assumption was not violated. Hence, we report marginal effects of a 1-year change in age (and CIs) for all imaging-identified features of spinal degeneration, including increasing FI in paravertebral muscles, marginal effects per year were used to depict a relative rate of decline per variable.   [11] and Valentin et al. [12] data. Footnote: epj corresponds to a random sample of endpoint variates in the j-th age group. muj corresponds to the mean value of the endpoint in the j-th age group. sdj corresponds to the standard deviation of the endpoint in the j-th age group. sj corresponds to the sample size of the endpoint in the j-th age group  (Table 3). Tables 2 and 3 give additional and more detailed information on the estimated relative decline rates based on synthetic data and Monte Carlo simulation, including number of iterations, sample size, range and standard errors.

Discussion
Our study reports a novel method for comparing yearly rate of decline of accepted features of spinal degeneration as determined through imaging studies involving asymptomatic participants. Providing reference values for normative rate of degeneration that includes all structures of the spinal column including muscle represents a new direction for examining spinal tissues, and may have theranostic value in managing the effects of both ageing and LBP.
While muscle tissues are alluded to in the degenerative cascade [13,30], literature to date has focused on the natural history of degeneration of the vertebra and disks for describing age-aggregated decline of the spinal column. The yearly rate of increasing paravertebral muscle fat content determined in our study was low and in agreement in both muscle studies used, which suggests a relatively slow decline to lumbar muscle quality into healthy adulthood. Whether this trend continues into older age (>60 years) for both genders, and other regional muscles, requires further investigation.  Two population-based studies report longitudinal change to paravertebral muscle quality determined from MRI in participants aged 40 plus [15,31] and offer a comparison for actual decline to our simulated data. Hebert et al. [31] reported 28.8 % MF fat content at age 40, 28.7 % at 45, and 31.6 % at 49 years of age for their group. While these values indicate a non-linear increase in fat content, the yearly rate of decline calculated between 40 and 49 age time-points is calculated to be 0.31 %. This rate is slightly higher than shown for MF in our two studies (Crawford = 0.24 % and Valentin = 0.11 %), and may be explained by their sample including a proportion of individuals with LBP and/or older age. Fortin et al. [15] reported reduction of functional cross-sectional area (FCSA; fatfree muscle tissue) for MF/ES of 29/42 % at L3/4, and 37/ 78 % at L5/S1 over their 15 year follow-up of male twins. Assuming linear decline, this suggests yearly rates of 1.9/ 2.8 % and 2.5/5.2 %, respectively, which are higher again and may reflect their older average age, sex differences [11], a partly symptomatic sample, and/or the differing variable. Interestingly, Fortin et al. [15] reported FCSA/ CSA percent reduction for MF/ES to be 0.21/0.30 % and 0.21/0.40 % at L3/4 and L5/S1, respectively. Whether this variable relates best to change to FI requires further substantiation in comparative studies. However, that two population cohorts showed faster paravertebral muscle decline than for the asymptomatic cohorts used in our modelling, suggests potential for differences in rate secondary to painful (or possibly non-painful) symptoms that at least warrants further investigation. Furthermore, the trends for rate differences noted between paravertebral muscles in ours and the Fortin et al. [15] study indicate a non-uniform decline between muscles that may be important in understanding the etiology of FI and any implications for targeted (and likely muscle-differential) therapies. The relative rate of increasing fat content in paravertebral muscles compared best with that of disk protrusion, spondylolisthesis and annular fissure. While caution should be exercised in making comments regarding common causation between these features, we postulate that an etiological relationship exists. Similar decline rates for spondylolisthesis and muscle quality can be rationalized on a neuro-mechanical basis. Paravertebral muscle denervation is present in asymptomatic individuals [32], with MF purported to be morphologically susceptible to effects of neural stretching after disk height loss and subsequent listhesis resulting in asymptomatic denervation [33]. Kalichman et al. [6] and Teichtahl et al. [9] reported a relationship between disk height loss and paravertebral muscle fat in patients, with the former also showing correlation to spondylolisthesis. Disuse-related muscle atrophy relates to deconditioning, local tensile unload and altered muscle recruitment [4,34,35], which Hodges et al. [28] purported was based on changes secondary to structural remodeling rather than atrophy. It seems probable that the faster declining loss of disk height shown in our study precedes declining muscle quality and spondylolisthesis within the normative ageing cascade. Furthermore, we speculate that our results showing descending rates of decline starting from disk signal loss, to disk degeneration, then disk height loss, facet degeneration and disk bulge, preceding disk protrusion, annular fissure, spondylolisthesis and paravertebral muscle tissue decline, are a representation of Kirkaldy-Willis et al.'s [30] degenerative cascade. Longitudinal investigations are required to substantiate this.
The relationship between muscle decline and disk protrusion or annular fissure is difficult to explain without details describing criteria determining each feature from the primary sources referenced by Brinjikji et al. [3], particularly in light of bulging disks being included as a separate category. Both disk protrusion and annular fissure predominate in the posterolateral and posterior disk, respectively and represent disruption to the nociceptive annular layers [36]. Whether the proximity of the posterior disk to the paravertebral muscle tissues and likelihood for exposure to common inflammatory mediators has a bearing on any etiological association to degenerative decline is speculation. However, several animal studies from Hodges and colleagues [28,35,37] describe rapid atrophy or re-modelling after experimental disk injury, so an association probably exists.
Age-related change to skeletal muscle quality differs between sexes, wherein, men lose more muscle with ageing, yet women suffer greater functional consequences [24]. Our trend revealing lumbar MF declining faster and ES slower in men than women may reflect sexdependent degeneration of paravertebral muscle quality whose influence is likely multifactorial. Furthermore, the kinematics of MF and ES based on their different morphology and architecture may point to a need for differential therapeutic strategies in optimizing muscle quality and function.
Comparisons between both muscle studies we used showed consistently lower decline rates in the Valentin et al. [12] parameters, despite reporting higher fat content compared to the Crawford et al. [11] study; sample differences including demographics and activity levels offer a probable rationale for the disparity. While our study references asymptomatic cohorts, we cannot assume that samples included lifetime pain-free subjects, and as such, previous common (or ongoing and underlying) inflammatory mechanisms to decline cannot be discounted.
Our study should be interpreted in light of its limitations. The cross-sectional nature limits its generalizability wherein concordance with secular changes cannot be assumed. Longitudinal studies that concurrently examine degenerative features of the spinal column and paravertebral musculature would offer improvement but are unrealistic in determining rate of decline over the wide age-range described here. Examining normative imaging datasets derived using advanced sequencing like proton density fatfraction and Dixon fat/water MRI would be warranted in enabling concurrent visualization of both spinal column and adjacent musculature [16][17][18]. While we consider our simulation design to be statistically sound, drawing comparisons between prevalence rates for spinal degenerative features and actual paravertebral muscle FI is a new approach in determining rate of decline. The robustness of this method might be best tested with comparisons to longitudinal studies that employ methods able to reliably detect serial yearly changes of the small magnitude defined by our study. This is not limited to the lumbar spine region.

Conclusions
In conclusion, lumbar paravertebral muscle degeneration as defined by increasing FI appears to occur relatively slowly in asymptomatic adults aged 18 to 60 years in comparison to degeneration of the spinal column. Muscle decline warrants inclusion as a feature of the normative degenerative cascade. Further investigation is necessary to determine the functional significance of the rate of decline of paravertebral muscle tissues, and by comparison to people with LBP and other musculoskeletal disease.