Pain trajectories and possible predictors of a favourable course of low back pain in patients consulting musculoskeletal physicians in The Netherlands

Abstract

Background

In The Netherlands, low back pain patients can consult physicians specialized in musculoskeletal (MSK) medicine. Previous studies have reported on the characteristics of patients consulting MSK physicians, and the treatment options used. There are no studies yet reporting on the course of Low Back Pain (LBP) after treatment by musculoskeletal (MSK) physicians in The Netherlands.

Methods

In an observational cohort study MSK physicians recorded data about all low back pain patients presenting for a first consultation. At baseline they recorded age, gender, type and duration of the main complaint, and concomitant complaints. At the end of treatment they recorded the type of treatment and the number of treatment sessions. Patients were recruited to answer questionnaires at baseline, and at 6-weekly intervals during a follow-up period of six months. Patient questionnaires included information about previous medical consumption, together with PROMs measuring the level of pain and functional status. Latent Class Growth Analysis (LCGA) was used to classify patients into different groups according to their pain trajectories. Baseline variables were evaluated as predictors of a favourable trajectory using logistic regression analyses, and treatment variables were evaluated as possible confounders.

Results

A total of 1377 patients were recruited, of whom 1117 patients (81%) answered at least one follow-up measurement. LCGA identified three groups of patients with distinct pain trajectories. A first group (N = 226) with high pain levels showed no improvement, a second group (N = 578) with high pain levels showed strong improvement, and a third group (N = 313) with mild pain levels showed moderate improvement. The two groups of patients presenting with high baseline pain scores were compared, and a multivariable model was constructed with possible predictors of a favourable course. Male gender, previous specialist visit, previous pain clinic visit, having work, a shorter duration of the current episode, and a longer time since the complaints first started were predictors of a favourable course. The multivariable model showed a moderate area under the curve (0.68) and a low explained variance (0.09).

Conclusions

In low back pain patients treated by musculoskeletal physicians in The Netherlands three different pain trajectories were identified. Baseline variables were of limited value in predicting a favourable course.

Background

Low back pain (LBP) is a major health problem, with a point prevalence in Western countries of 12–30% [1]. In 2010, in The Netherlands, total costs were estimated to be €3.5 billion, including both direct costs, related to the consumption of medical care, and indirect costs, related to loss of productivity and disability pensions [1]. Because in most cases the mechanism of LBP is not known, there is no intervention that can be directed at the cause of the pain, and while many interventions are available, none has shown to be superior [2]. Although the course of low back pain has long been considered favourable, recurrences are common [3, 4], and many patients (65%) still reported pain 1 year after onset [5]. Considering the recurrent character of LBP, recent research has increasingly focused on identifying LBP trajectories [6, 7]. Distinct clusters of pain trajectories were identified [8], and over the course of their LBP, patients showed consistent cluster membership [6]. Rather than studying prognostic factors based upon outcome measurements at one follow-up moment, it may therefore be more informative to follow patients for longer periods of time, and to identify prognostic factors that predict the trajectory of LBP. This knowledge can potentially be used in outcome research, studying whether interventions can influence patients pain trajectories, rather than offering momentary improvement [9]. Measuring pain trajectories has become easier with the development of automated systems distributing patient reported questionnaires over the internet, or by using text messages. A recent study by Ailliet et al., for example, used text messages to study the pain trajectories of patients with low back and neck pain in patients consulting chiropractors in The Netherlands and Belgium [10]. A review of LBP trajectory studies by Kongsted et al. stated that it would be highly relevant to investigate whether LBP trajectories identify phenotypes of LBP that benefit from different care pathways [11].

In The Netherlands, among the various diagnostic and treatment possibilities, patients can consult physicians specialized in Musculoskeletal Medicine (MSK). MSK physicians are generally consulted by patients with musculoskeletal pain, and about half of the patients consult MSK physicians because of LBP [12]. MSK physicians use an array of diagnostic and treatment options, almost invariably including a form of Spinal Manipulative Therapy (SMT). The aim of our study was to assess whether different pain trajectories could be identified in LBP patients after consulting MSK physicians, and to identify possible predictors of a favourable course.

Methods

Study design and participants

We conducted a prospective cohort study with a follow-up period of six months.

All MSK physicians registered with the Dutch Association for Musculoskeletal Medicine were invited to participate in our study. Participating physicians were instructed to register all patients who presented for the first time in an MSK practice through a web-based registry. Patients were invited to take part in the study when making a first appointment. Inclusion criteria were LBP, age ≥ 18, and sufficient mastery of the Dutch language to answer questionnaires in Dutch. If patients gave informed consent, the treating physician entered email addresses of the recruited patients in the web-based registry. Thereafter, a specially designed computer program (Readmail) was used to automatically distribute invitations to patients by email to fill in web-based questionnaires.

Study procedure

Both the treating physicians and the individual patients provided data via web-based registries. The treating physicians recorded data at baseline and at the end of the treatment. Study procedures were explained to participating physicians during specially organized information sessions. In addition, a research assistant visited all participating practices to explain the procedures. Practices that agreed to participate at a later stage were informed by telephone. Instructions were to ask all consecutive patients presented for a first consultation to participate in the study. Recruited patients received invitations to fill in web-based questionnaires within three weeks before the first consultation, and at six weekly intervals during the ensuing six months. When patients did not respond, a maximum of three reminders were sent within a period of two weeks. Both the invitation email and the web-based questionnaires contained links to a leaflet with information about the study.

Measurement

At baseline physicians registered data about age, gender, type and duration of the main complaint, and the existence of concomitant complaints. Complaints were registered according to the International Classification of Primary Care (ICPC) [13]. At the end of treatment, data were registered about the number of treatment sessions, the type of treatment used, and further referrals.

At baseline, patients were asked to indicate whether their main complaint was low back pain, neck pain or any other complaint. This question was supported by text and manikins, explaining which area was considered to cover neck pain or low back pain. For other complaints, patients could explain these in text. Patients were asked to indicate whether their pain radiated to the legs or arms, and whether they had numbness or pins and needles in their legs or arms. Patients were also asked about the duration of the current episode, the time since the first episode, educational level, work status, previous specialist consultations, and previous treatments. The effect of previous treatments was measured on an ordinal scale, with four possible answer categories; 1.strong improvement, 2.little improvement, 3.unchanged, 4.deteriorated. Furthermore, all patients were asked to answer a set of Patient Reported Outcome Measures (PROMs), including a Numerical Rating Scale (NRS) for pain severity, the SF6D [14], and the Fear Avoidance Beliefs Questionnaire (FABQ) [15, 16]. Patients who indicated LBP as their main complaint were asked to answer the Oswestry Disability Index (ODI) [17, 18]. The SF6D is a short version of the SF36, measuring health related quality of life. Scores range from 0–1, with higher scores indicating lower quality of life. The FABQ consists of 16 items, and measures pain related fear in LBP patients. Higher scores indicate more pain related fear. The ODI consists of 10 items with scores ranging from 0–50, with higher scores indicating more disability because of LBP. At all follow-up points patients were asked to answer the same PROMs, except for the FABQ. A question about the Global Perceived Effect (GPE) of treatment was added.

Statistical analyses

Identification of pain trajectories

Our study population consisted of all LBP patients who completed the baseline questionnaire. For the analyses of pain trajectories, patients were selected who completed the baseline questionnaire and at least one of the follow-up questionnaires. Latent Class Growth Analyses (LCGA) were used to explore whether subgroups of patients following distinct pain trajectories could be identified, using the NRS for pain scores [19]. Several LCGA models were evaluated with different numbers of trajectories, allowing linear or quadratic pain trajectories, and allowing more or less heterogeneity in pain trajectories within subgroups. For competing models the posterior probabilities were compared (average, range and standard deviation of the latent class probability per class). A final model was selected based upon model fit and considerations of interpretability and clinical practicality, proportion of patients in each class, and the final outcome after a follow-up of six months. Model fit was evaluated with the Vuong-Lo-Mendell-Rubin likelihood ratio test (LMR-LRT) and the Bayesian Information Criterion (BIC) [20]. The BIC considers both the likelihood of the model as well as the number of parameters in the model, with lower values showing better model fit. The LMR-LRT provides a p-value. A significant p-value indicates that a model with k classes fits better than a model with k-1 classes. LCGA was carried out using Mplus (Version 7) [20, 21].

Identifying possible predictors

Descriptive analyses of baseline variables were carried out for the complete population included in the analyses, and for each group of patients with a distinct pain trajectory separately. For the patients that presented with high baseline pain scores two distinct pain trajectories were identified (see results). One trajectory identified a group of patients who did not improve (class 1), and one trajectory identified a group of patients who improved (class 2). Logistic regression analyses were conducted to study the univariate relationship between baseline variables and the dependent variable (i.e. high baseline pain and not improved (class 1) vs high baseline pain and improved (class 2)).

A backward selection procedure was carried out on cases with complete data on all variables to construct a multivariable model, based upon a p-value of 0.157 (Akaike criterion). Treatment variables were considered as possible confounders instead of predictors. Although not known at baseline, they could possibly influence the outcome. To evaluate the influence of treatment variables the backwards selection was therefore carried out twice, with and without including treatment variables (type of treatment, number of treatment sessions, and adjuvant McKenzie treatment). The relationship between continuous predictors and the outcome was tested for linearity, and non-linear variables were entered as splines. The fit of the final model was evaluated with the loglikelihood and the Hosmer Lemeshow test [22]. Discriminative properties of the model were evaluated by calculating the Area Under the Curve (AUC), and Nagelkerke R² [23] was used as an overall measure to quantify the amount of information in the model. Bootstrapping was used for internal validation [24]. Descriptive analyses and univariate analyses were carried out using SPSS 22, except for the univariate analyses of non-linear variables. Linearity, univariate analyses of non-linear variables and internal bootstrap validation was carried out with the R-package rms (version 5.1–2). In the multivariable analyses, backwards selection was carried out with the R-package pfsmi [25].

Missing data and evaluation of loss to follow-up

The relationship between complete predictors and predictors with missing values (> 20%) was studied with univariate logistic regression analyses. In this analysis, significant relationships between predictors and the variable being either missing or not missing support the assumption that missing values are probably missing at random (MAR). When including all potential predictors in the model, the percentage of missing cases was around 40%, which required 40 multiple imputed datasets. Multivariable analyses were conducted in each dataset, and results were pooled according to Rubin`s rules [26]. To evaluate the influence of loss to follow-up, the group of patients who only answered the baseline questionnaires was compared with the group of patients answering at least one follow-up measurement. Differences in baseline characteristics between these two groups were studied with logistic regression. Multiple imputation and evaluation of missing data were carried out using SPSS 22.

Results

Study population

Data was collected from February 2014 until February 2016. A discrepancy was found between the number of patients classified with low back pain by the physician, and the self-classification by the patient through the web-based questionnaire. Frequently, patients classified themselves as other, but indicated complaints in text that would classify as low back pain. It was therefore decided to use the classification of the physician to select LBP patients. In the web-based registry MSK physicians recorded 2026 patients with LBP. Of these patients 1664 were recruited for our study. Our study population consisted of 1377 patients who answered the baseline questionnaire. A total of 1117 patients (81%) answered at least one of the follow up measures next to the baseline questionnaire and were included in the LCGA and prediction analyses. Of these 1117 patients, 93% answered the first follow-up questionnaire, 74% the second, 58% the third, and 43% the fourth. Although percentages missing increased over time, patients frequently missed intermediate questionnaires. Baseline characteristics of the whole sample, the patients who answered at least one follow-up questionnaire and were included in the analyses, and the patients who were lost to follow-up are presented in Table1. Although 19 practices participated in the study, the LBP patients were recruited by 16 practices, and the number of LBP patients recruited by the various practices varied from 1–285.

Table 1 Sample characteristics and evaluation of loss-to follow-up

Missing data loss to follow-up

Additional file 1: Table S1 shows the handling of predictor variables, including the percentages of missing values. Only one variable, Baseline ODI, showed a high percentage of missing values (25.6%), which could be explained by the tailored distribution of the ODI to patients who had self-classified as LBP patient. Because not all LBP patients classified themselves as such, not all LBP patients received the ODI. Because the percentage of missing Baseline ODI values was > 20, it was decided to impute the baseline ODI.

Evaluation of loss to follow-up is presented in Table 1. It shows that baseline scores on the ODI, SF6D, and NRS did not differ significantly between patients who only answered the baseline questionnaire and the patients included in the analyses. Female patients (p = 0.042), older patients (p = 0.008), and patients treated effectively by a chiropractor (p = 0.034) were significantly more inclined to remain in our study. Patients treated effectively at a pain clinic were significantly less inclined to remain in our study (p = 0.046).

Defining subgroups with distinct pain trajectories

Model fit characteristics are presented in Additional file 1: Table S2. Although a four-class quadratic model without fixed variance showed slightly better fit (BIC 17,836 versus 17,842), it was decided to choose the three-class quadratic model without fixed variance as the preferred model, because of its better interpretability and practicality. Posterior probabilities of the three class model were slightly better than the posterior probabilities of the four class model (mean posterior probabilities ranged from 0.769–0.854 for the three class model versus 0.717–0.848 for the four class model). Furthermore, the four-class model included a small group of patients (7.0%) who showed a strong improvement in the first three months, and a return to previous pain levels in the subsequent months. Compared to the three-class model, this for a large part only changed the proportion of patients in the group that started with high pain levels and showed no improvement, suggesting that this group consisted of patients who remained unchanged during the study combined with patients who improved strongly, only to deteriorate again. Figure 1A-D shows the mean estimated trajectories (Fig. 1A), and the three separate trajectories together with the individual trajectories of all patients for class1 (Fig. 1B), class2 (Fig. 1C), and class 3 (Fig. 1D). The course of the average NRS scores of the three groups is presented in Table 2, together with the mean changes in the ODI scores. In the three-class model, a group of 226 patients started with high NRS scores at baseline and showed hardly any change during the follow-up period (mean NRS-pain changed from 6.9 to 6.7). In this group the mean ODI score changed from 24.8 to 19.4. A group of 578 patients started with high baseline scores and showed considerable improvement (mean NRS changed from 7.0 to 1.8). In this group, the mean ODI score changed from 26.4 to 6.1. A group of 313 patients started with lower baseline scores and showed moderate, but clinically relevant improvement (mean NRS changed from 3.5 to 2.2). In this group the mean ODI score changed from 15.5 to 8.0.

Figure 1

A–D Mean estimated trajectories (A) and separate trajectories per class, with individual trajectories of all patients (B–D). The X-axis represents the measurement moments (1 = baseline, 2 = 6 weeks, 3 = 12 weeks, 4 = 18 weeks, and 5 = 26 weeks follow-up). The Y-axis represents the NRS mean group scores (red line: not improved patients (26.5%), blue line: improved patients (46.0%), green line: patients with low baseline pain levels (27.5%)). In (A), squared lines represent sample trajectories, triangular lines represent estimated trajectories. B–D show the estimated trajectory with the individual trajectories of all individual patients for not improved patients (B), improved patients (C), and patients with low baseline NRS (D) separately. Grey lines represent the individual trajectories of all patients

Table 2 Mean NRS and ODI scores of the three classes

Predictors

In Table 3 the baseline characteristics for the three groups of patients with distinct pain trajectories identified with LCGA are presented. Because the group of patients with low baseline pain scores could be identified by the baseline NRS scores, our main interest was to identify predictors that distinguish patients with high NRS scores who showed a favourable course from patients with high NRS scores who did not show a favourable outcome. Further analyses therefore focused on the two subgroups that started with high pain scores, i.e. the group of patients that was considered to be improved and the group of patients that was considered to be not improved. Baseline variables were evaluated as possible predictors of a favourable course. For all baseline variables, univariate odds ratio’s for improvement are presented. The relationship between the continuous predictors and group membership was shown to be linear for all continuous variables except for the duration of the current episode, which was further analysed as a spline variable. In this spline variable, two stages could be identified. If the duration of the current episode was less than four years, the odds of a favourable outcome decreased with a longer duration. If the duration of the current episode was longer than four years, the odds of a favourable outcome hardly changed. Male gender, previous specialist visit, previous surgical treatment, and having work were associated with a favourable course. Previous consultation with a neurologist or an orthopaedic surgeon, no effect of previous treatments and concomitant complaints were associated with a non-favourable outcome. Other baseline variables did not show a significant association with the outcome.

Table 3 Characteristics of three groups of patients with different pain trajectories and univariate probability for improvement

Multivariable model

A multivariable model was constructed, and the model without treatment variables is presented in Table 4. In this model male gender, previous specialist visit, previous pain clinic treatment, having work, a shorter duration of the current episode, and a longer time since the complaints first started were predictive of a favourable course. No effect of previous chiropractic treatment was predictive of a non-favourable course, both compared with patients reporting a positive effect of previous chiropractic treatment, or patients not previously treated by a chiropractor. The fitted model showed an AUC of 0.671, with a non-significant Hosmer and Lemeshow test (0.732), supporting model fit, and the amount of information in the model (R²) was 0.10. Bootstrap validation resulted in a corrected R² of 0.06. Adding treatment variables in the backwards selection analysis resulted in a model which included the number of treatments (OR 0.91), but no other treatment variables such as the type of SMT administered or adjuvant treatment with McKenzie. The fit of this model was practically the same as the fit of the model without treatment variables (AUC 0.677, Hosmer and Lemeshow 0.734, R² 0.10, and bootstrap R² 0.06).

Table 4 Multivariable model

Discussion

Studying the clinical course of low back pain in patients consulting MSK physicians in The Netherlands with Latent Class Growth Analyses, three distinct pain trajectories were identified. More than half of all the patients (52%) presented with high baseline pain scores and showed considerable improvement during six months follow-up. A second group of patients with high baseline pain scores (20% of all patients) showed no improvement. A third group, with moderate baseline pain scores (28%) showed slight, but clinically relevant improvement. The multivariable model presented showed a moderate discrimination, with an AUC of 6.77, and poor calibration, with an R² of 0.10. Therefore one can question its usefulness in clinical practice. Apparently, with the baseline data collected in our study, it is hard to predict which patient might improve after consulting an MSK physician.

Previous studies almost invariably reported similar clusters of pain trajectories, generally including clusters with persistent high pain, clusters with more or less persistent moderate or low pain, clusters showing improvement, and clusters with a fluctuating pattern. The proportion of patients in each cluster, however, differed, possibly because of variations in the study designs. Patients were recruited in General Practice [6,7,8, 27], at chiropractic clinics [10, 28], combined in General Practice and Chiropractic clinics [29], combined in General Practice and physiotherapy practices [30], or in a population based survey [31]. Also, studies varied in recruiting patients: i.e. only patients with chronic [30], only acute [27], or a mix of both chronic and acute LBP [6,7,8, 10, 28, 29, 31]. Moreover, follow-up periods varied from 12 weeks [27] to one year [10, 29,30,31], and follow-up measurements varied form weekly text messages [10, 28, 29] to monthly questionnaires [6,7,8]. The population based study was the only study in which a cluster showing improvement was not reported [31]. And the only study recruiting acute LBP patients showed high percentages of recovery [27]. The clusters presented in our study are well comparable to those reported in other studies, with the exception of a cluster representing a fluctuating pattern. In our four class model a small cluster was added that would in other studies have qualified as fluctuating. We considered this cluster merely a subgroup of the consistent high pain cluster, eventually showing no improvement after six months follow-up, and therefor chose to use the three class model.

Most trajectory studies reported variables that were associated with group membership. Although varying variables were reported, only higher pain intensity [6,7,8, 27], longer duration [6, 7, 27, 29, 32], and more physical disability [8] were more or less consistently associated with a more severe trajectory. The same variables were reported in other studies to be associated with a worse prognosis in LBP patients, together with unemployment [33, 34]. Similarly, in our study the duration of the current episode and unemployment were both associated with a lower probability of improvement, but baseline disability was not associated with the outcome. In our univariate analysis, ineffective previous treatments were consistently predictive of an unfavourable course. Of these previous treatments reported in our study, only chiropractic treatment ended up in our prediction model.

Our study used pain trajectories to identify different courses of LBP after MSK treatment, rather than use a singular outcome measurement. Groups of patients with different trajectories were compared to assess whether baseline variables could predict a possible favourable course. We would consider that the different trajectories thus identified represent a more relevant estimation of whether patients have improved or not. A challenging question remains to what extent the clinical course represented the natural recovery or the consequence of the treatment administered. A multivariable model including treatment variables only included the number of treatments sessions, but not the type of SMT, or the adjuvant use of McKenzie treatment. No conclusions can be drawn from these findings because patients were not randomized, and physicians were at liberty to choose the type of treatment. The higher number of treatment sessions in the group of unimproved patients, for example, may well have been influenced by elongating the treatment in patients who reported no improvement to their treating physician.

Strengths and limitations

A strength of our study was the web-based data-collection, which enabled us to follow a large number of patients at a relatively low cost. In this way data could be collected from patients who consulted an MSK physician, with questionnaires that were tailored to their main complaint. A weakness was the difficulty to identify low back pain patients before consulting the physician. Our solution using web-based self-classification, aided with manikins, appeared to lead to a high proportion of miss-classification. Because we tailored the distribution of PROMs to the main complaint as reported by the patient, this miss-classification led to a high percentage of missing baseline ODI values. We therefore chose to use the physician’s diagnosis to identify patients with LBP, and we imputed the baseline ODI. Another weakness of our study set-up was the high proportion of patients that discontinued their participation. The response rate gradually diminished during the follow-up period. Out of the 1117 patients included in the baseline population 93% responded after 6 weeks, 74% after 12 weeks, 58% after 18 weeks, and 43% at six months. We found that some baseline variables were related to loss to follow-up which made the MAR assumption more plausible, supporting multiple imputation of missing values. However, loss to follow-up may have had an impact on the long-term course of the pain trajectories identified. Another limitation is the inclusion of all patients completing at least one follow-up measurement in our analyses. The fact that patients frequently missed intermediate questionnaires, however, supports the assumption that these measurements are missing at random.

Conclusion

In patients with low back pain, three different clinical courses were identified in the six months after consulting an MSK physician in the Netherlands. A large group of patients presented with high baseline pain scores, and showed improvement. In patients with a high pain score at baseline, a multivariable prediction model showed a number of predictors of a favourable course. In this model, male gender, longer time since the complaints first started, shorter duration of the current episode of pain, previous specialist visit, previous pain clinic visit, effective treatment by a chiropractor, or no previous chiropractic treatment, and having work were predictors of a favourable course. The prediction model, however, showed a low AUC and a low amount of information. It is a continuing challenge to identify predictors of a favourable outcome in LBP patients.