Regional Homogeneity and Multivariate Pattern Analysis of Cervical Spondylosis Neck Pain and the Modulation Effect of Treatment

Objects: We investigated brain functional alteration in patients with chronic cervical spondylosis neck pain (CSNP) compared to healthy controls (HCs) and the effect of intervention. Methods: 104 CSNP patients and 96 matched HCs were recruited. Patients received 4 weeks of treatment. Resting-state fMRI and Northwick Park Neck Pain Questionnaire (NPQ) were collected before and after treatment. Resting state regional homogeneity (rs-ReHo) and multivariate pattern analysis (MVPA) were applied to (1) investigate rs-ReHo differences between CSNP patients and controls and the effect of longitudinal treatment and (2) classify CSNP patients from HCs and predict clinical outcomes before treatment using MVPA. Results: We found that (1) CSNP patients showed decreased rs-ReHo in the left sensorimotor cortex and right temporo-parietal junction (rTPJ), and rs-ReHo at the rTPJ significantly increased after treatment; (2) rs-ReHo at rTPJ was associated with NPQ at baseline, and pre- and post-treatment rs-ReHo changes at rTPJ were associated with NPQ changes in CSNP patients; and (3) MVPA could discriminate CSNP patients from HCs with 72% accuracy and predict clinical outcomes with a mean absolute error of 19.6%. Conclusion: CSNP patients are associated with dysfunction of the rTPJ and sensorimotor area. Significance: rTPJ plays on important role in the pathophysiology and development of CSNP.


INTRODUCTION
Cervical spondylosis (CS) is a common age-related chronic disease characterized by stiffness and neck and upper back pain. Chronic CS neck pain (CSNP) has become an important public health and social issue because of its high prevalence, unsatisfactory treatment options, large medical burden, and reduction in quality of life (Binder, 2007;Hoy et al., 2010;Baime, 2016). Recent studies have demonstrated that patients with cervical spondylotic myelopathy (CSM), the later stage of CS, have cerebral cortex reorganization as a result of spinal cord compression (Holly et al., 2007;Dong et al., 2008;Holly, 2009;Duggal et al., 2010). However, the characteristics of cortical activity in CSNP patients without spinal cord compression have not been reported.
Recent studies also found improvements in brain functional activity and connectivity in patients with spinal cord decompression following intervention (Holly et al., 2007;Dong et al., 2008;Holly, 2009;Duggal et al., 2010;Tam et al., 2010;Hrabalek et al., 2014;Green et al., 2015;Hrabalek et al., 2015;Tan et al., 2015;Bhagavatula et al., 2016). Dong and colleagues found that abnormalities in the cortical sensorimotor recruitment pattern in CSM patients gradually disappeared and were accompanied by behavioral improvements after surgery (Dong et al., 2008). Effective longitudinal treatment can also modulate the aberrant functional neural activity and brain structure of chronic low back pain patients (Lorimer Moseley, 2005;Baliki et al., 2008;Seminowicz et al., 2011;Younger et al., 2011;Li et al., 2014;Ceko et al., 2015;Louw et al., 2015;Shi et al., 2015;Vrana et al., 2015;Braden et al., 2016). Nonetheless, there is a lack of research evaluating the modulation effect of effective treatment on the brain activity of CSNP patients.
In recent years, machine learning has gained popularity in the translational neuroimaging community for it allows researchers to identify brain signatures of patients with psychiatric and neurologic diseases and make accurate classifications (Arbabshirani et al., 2017). A widely applied machine learning approach, multivariate pattern analysis (MVPA) uses pattern classifiers to capture the information encoded by spatially correlated voxels in fMRI and has proven to be more sensitive to the functional organization of the cortex (e.g., general linear model; Haxby, 2012). This technique has been used to explore the neurophysiology of chronic pain disorders, including migraine (Chong et al., 2017), chronic low back pain (Ung et al., 2014) and fibromyalgia (Lopez-Sola et al., 2017).
In this study, we explored (1) the differences in resting state brain activity between CSNP patients and matched healthy controls (HCs), (2) how longitudinal treatment can modulate brain abnormalities in CSNP patients using regional homogeneity (ReHo), and (3) if MVPA approaches can classify CSNP patients and HCs and predict clinical outcomes before treatment. We hypothesized that (1) CSNP patients would be associated with altered activity in specific brain regions (e.g., sensorimotor cortex or other regions) compared with HCs, (2) symptom reduction after an effective treatment could modulate the altered brain activities of CSNP patients, and (3) CSNP patients and HCs could be discriminated using MVPA.

Participants
This study was approved by the Ethics Committee at the Second Affiliated Hospital of Guangzhou University of Chinese Medicine. The experiment was performed in accordance with approved guidelines. Participants were recruited from the outpatient clinics of the Second Affiliated Hospital, Guangzhou University of Traditional Chinese Medicine. This study was registered on http://www.chictr.org.cn/index.aspx (ChiCTR-IPR-14005282). Written informed consent was obtained from each participant. Subjects were recruited from July 2010 to December 2012.
The inclusion criteria were as follows: (1) a confirmed diagnosis of CSNP in accordance with the diagnostic criteria published by the Chinese Association of Rehabilitation (2010) and in accordance with the International Classification of Diseases, 10th edition (ICD-10) codes 1 ; (2) male or female 18-35 years of age; (3) no signs of musculoskeletal system abnormalities upon physical examination; (4) not receiving any treatment within the last 7 days; and (5) informed consent obtained. Patients were excluded if they met any of the following criteria: (1) history of neck trauma; (2) vertebral body or spinal canal cancer, tuberculosis, or severe osteoporosis; (3) history of neck surgery or presence of congenital malformation of the cervical vertebrae; (4) pregnant or breastfeeding (females); and (5) presence of a severe systemic disease such as tumors, diabetes mellitus, kidney disease, or digestive system disease.
Additionally, pain-free age and gender-matched individuals were recruited for this study as HCs. Each subject underwent a medical history evaluation and physical examination.

Interventions
All CSNP patients received acupuncture treatments for one month at the Second Affiliated Hospital of Guangzhou University of Chinese Medicine. The CSNP patients were assigned to one of five acupuncture treatment groups based on the acupuncture points applied. The acupoints selected were bilateral Bailao (EX-HN15) in acupuncture group 1 (AG1), bilateral EX-HN15 and Hegu (LI4) in AG2, bilateral EX-HN15 and Zhongzhu (TE3) in AG3, bilateral LI4 in AG4, and bilateral Zusanli (ST36) in AG5.
Each patient received about 30 min of acupuncture treatment approximately once every three days for a total of 10 sessions over a 4-week period. The needles were inserted to a depth of 20-30 mm, perpendicular to the surface of the skin. Electroacupuncture was performed to stimulate the relevant acupoints for 30 min using a continuous-wave, 1-Hz frequency, and a comfortable strength until deqi sensation (soreness, numbness, distention, and heaviness) was obtained.
The Northwick Park Neck Pain Questionnaire (NPQ) was administered to CSNP patients to assess neck pain and disability (Leak et al., 1994). The NPQ scores were assessed before the first treatment as baseline data and again after the final treatment by the patients themselves following thorough instruction from the researchers.

Image Acquisition
All participants underwent a structural and functional MRI scan using a 3.0T MR system (Magnetom Verio, Siemens, Germany) at the Department of Radiology at the Second Affiliated Hospital of Guangzhou University of Chinese Medicine. Cushions were used to support participants' heads within the coil to prevent excessive head motion. Participants were asked to stay still with their eyes closed. All images were acquired parallel to the anterior-commissure-posterior-commissure line with an autoalign technique. The functional data were collected using an echo-planar imaging sequence: 31 axial slices; repetition time (TR) = 2000 ms; echo time (TE) = 30 ms; slice thickness = 3.5 mm; gap = 0.35 mm; flip angle (FA) = 90 • ; matrix = 64 × 64; field of view (FOV) = 224 mm × 224 mm, 240 time points. 3D structural images were acquired using a T1-weighted MP-RAGE sequence: 176 sagittal slices; TR = 1900 ms; TE = 2.27 ms; inversion time = 900 ms; slice thickness = 1.0 mm; no gap; FA = 9 • ; matrix = 256 × 256; FOV = 256 mm × 256 mm.

Data Preprocessing
Data preprocessing was carried out using Data Processing Assistant for Resting-State fMRI (DPARSF 2 ), which is based on Statistical Parametric Mapping (SPM 3 ) and the toolbox for Data Processing and Analysis of Brain Imaging (DPABI 4 ). The preprocessing procedures included: time alignment across slices, motion correction, within-subject registration between T1 and EPI images, T1 segmentation, and normalization of the BOLD fMRI datasets to register them to Montreal Neurologic Institute (MNI) space with voxels resampled at 3 mm × 3 mm × 3 mm.
Participants were excluded if head motion exceeded 1.5 mm and 1.5 • . Finally, the waveform of each voxel was passed through a band-pass filter (0.01-0.08 Hz) to reduce the effects of lowfrequency drift and high-frequency physiological noise. Finally, several nuisance signals were regressed out from each voxel's time course, including 24-parameter head-motion profiles, mean white matter (WM) and cerebrospinal fluid (CSF) time series.

Regional Homogeneity (ReHo)
ReHo was calculated by the Kendall's coefficient of concordance (KCC) (Kendall, 1990) using the REST toolkit 5 . Individual resting state ReHo (rs-ReHo) maps were generated by assigning each voxel a value corresponding to the KCC of its time series with its 26 nearest neighboring voxels (Zang et al., 2004). The individual rs-ReHo maps were standardized by their own mean KCC. Then, a Gaussian kernel with a full-width at half-maximum of 4 mm was used to smooth the images in order to reduce noise and residual differences.
Group analysis was calculated with a random effects model using the Data Processing and Analysis of Brain Imaging toolbox (DPABI 6 ). We first compared the rs-ReHo differences between CSNP patients and HCs using two-sample t-tests. Then, a partial correlation analysis was used to investigate the relationship between the rs-ReHo of these aberrant brain regions and NPQ scores in the CSNP patients, controlling for age and gender as covariates using SPSS.
In order to explore whether clinical improvement after an effective treatment could affect the brain's intrinsic neural activity in those aforementioned aberrant brain regions, we defined all brain regions that showed significant rs-ReHo changes between CSNP patients and HCs as regions of interest (ROI) and compared the pre-and post-treatment rs-ReHo changes in CSNP patients using a paired t-test within the ROIs as well as the whole brain. Then, partial correlation analysis was used to investigate the relationship between pre-and post-treatment rs-ReHo changes within the ROIs and corresponding NPQ score changes in the CSNP patients, controlling for age and gender as covariates using SPSS. A threshold of voxel-wise p < 0.005 and p < 0.05 family-wise error corrected at cluster level was applied for whole brain analyses, and a small volume correction method (voxel-wise p < 0.005 and p < 0.05 FWE corrected) was used for all ROIs (derived from a comparison between the CSNP patients and healthy controls before and after acupuncture treatments).

Multivariate Pattern Analysis (MVPA)
The ROIs defined in the ReHo analyses were used for the MVPA analyses with two objectives: (1) classify CSNP patients and HCs and (2) predict treatment responses using baseline (pretreatment) ReHo values. To avoid the risk of overfitting, all analyses were based on fivefold cross-validation (CV) (Pereira et al., 2009).
In the first step, machine learning models were trained to classify CSNP patients and HCs using ReHo values within the ROIs. A support vector machine (SVM) classifier was used and the implementation of SVM was based on LIBSVM (Chang and Lin, 2011). To quantify the performance of the SVM classifier, classification accuracy, sensitivity and specificity were calculated. Here, sensitivity and specificity represent the proportion of patients correctly classified and the proportion of HCs correctly classified, respectively. To further assess the performance of the classifier and evaluate the significance of classification accuracy, we ran permutation testing. In each test, we randomly permuted the class labels of the data prior to training. Fivefold CV was then performed on the permuted dataset and the procedure was repeated 10,000 times. If the classifier trained on real class labels had an accuracy exceeding the 95% confidence interval generated from the accuracies of the classifiers trained on randomly relabeled class labels, this classifier was considered to be well-performing.
In the second step, we aimed to predict treatment responses (percentage change of NPQ scores) using baseline ReHo values. Here, we used a novel fMRI prediction approach, namely sliced inverse regression (SIR), to extract the most predictive information from ROIs and subsequently predict responses. SIR is powerful because it can work efficiently regardless of the linearity of the relationship between fMRI data and their labels. The details of this implementation of SIR-based fMRI prediction can be found in Tu et al. (2017). Briefly, in the first step, we performed SIR on the features ∈ R n×p and corresponding behavioral data ψ ∈ R n (i.e., percentage changes of NPQ), resulting in the effective dimension reduction (e.d.r.) direction β ∈ R p×K and SIR-derived feature set S = β ∈ R n×K . In the second step, we used a support vector regression (SVR) ψ = C(S t α; t = 1, · · · , n) with α ∈ R K to estimate the predicted value and the contribution of each variable in the SIR reduced space to predict ψ. In the third step, we transformed α back to the original space by multiplying β, resulting in the weights for all features θ = βα. To quantify the performance of prediction, we used the prediction-outcome correlation which is defined as the correlation between the actual and predicted differences of NPQ scores, as well as mean absolute error (MAE) (Tu et al., 2016), which is defined as where ψ n andψ n denote, respectively, the actual and predicted differences of NPQ scores (1-Post NPQ/Pre NPQ) for patient n, and N is the number of CSNP patients. The number of slices for SIR was set to 5 and the number of directions was set to 4. Performance measures were assessed using permutation testing (as mentioned in the classification session).

RESULTS
A total of 104 right-handed CSNP patients (45 males and 59 females, age 24.90 ± 1.98 years, neck pain duration 36.49 ± 24.93 months) and 96 HCs (50 males and 46 females, age 24.80 ± 1.52 years) were recruited. There were no significant differences in age and gender between CSNP patients and HCs (P > 0.05).
Compared with baseline scores, all treatment groups (AG1, AG2, AG3, AG4, and AG5) showed improvement in NPQ scores (P < 0.01). We found no significant differences among AG1, AG2, AG3, AG4, and AG5 in NPQ score improvement (P > 0.05). The average NPQ score across all patients was 24.1 ± 9.2 before treatment and 12.7 ± 9.0 at the end of treatment. There was a significant decrease in NPQ scores after 1 month of treatment (P < 0.01).
Based on these clinical findings, we merged all acupuncture treatment groups (AAG) in the rs-ReHo analyses to investigate how an effective treatment can modulate rs-ReHo.

ReHo Results
To explore the neural pathophysiology of neck pain, we first compared all CSNP patients to HCs using a two-sample t-test. Compared to HCs, CSNP patients showed decreased rs-ReHo in the right temporo-parietal junction (rTPJ) and left sensorimotor cortex (left post-central gyrus and precentral gyrus) (Figure 1 and Table 1). Regression analysis indicated that the rs-ReHo in the rTPJ was positively associated with NPQ scores at baseline in CSNP patients (Figure 1).
We then compared post-and pre-treatment rs-ReHo differences using paired t-tests. Results revealed that after longitudinal clinical treatment, CSNP patients showed increased rs-ReHo in the rTPJ (peak 60, −48, 6, 13 voxels). There were no significant differences between the pre-and post-treatment rs-ReHo in other ROIs and non-ROIs at the threshold we set. Regression analysis indicated that the pre-and post-treatment rs-ReHo changes at the rTPJ were positively associated with corresponding NPQ changes in CSNP patients.

MVPA Results
Classification accuracy for separating CSNP patients and HCs is summarized in Figure 2, upper panel. Using the ReHo values within the rTPJ, left postcentral gyrus, and precentral gyrus, the classifier achieved an accuracy of 72.0% (76.0% for CSNP and 67.7% for HCs; p < 0.0001, permutation testing).
Similar to a previous study (Doehrmann et al., 2013), two subjects were removed as outliers from the prediction analysis since their NPQ changes were more than 3 SDs higher than the average of our samples. The correlation between actual and predicted treatment response was 0.30 (p = 0.002; Figure 2, lower panel). We obtained an MAE of 19.6% for predicting treatment responses (p = 0.004, permutation testing).

DISCUSSION
In this study, we investigated rs-ReHo differences between CSNP patients and HCs as well as how symptom relief after treatment can modulate rs-ReHo. We found that CSNP patients were associated with reduced rs-ReHo in brain regions including the rTPJ and left sensorimotor cortex as compared with HCs. In addition, we found that effective longitudinal treatment (acupuncture treatment) can normalize (increase) decreased rs-ReHo in CSNP patients. Moreover, rs-ReHo in the rTPJ was associated with neck pain intensity, as measured by NPQ scores  at baseline. Following treatment, rs-ReHo changes at the rTPJ were associated with alleviation of neck pain intensity (i.e., decrease in NPQ scores) in CSNP patients. We applied MVPA for discriminating CSNP patients from HCs with an accuracy of 72% as well as predicting treatment responses using baseline data with an error of 19.6%. ReHo has been successfully used to detect local abnormalities in subjects with different disorders such as chronic pain Ke et al., 2015;Huang et al., 2016;Zhang et al., 2016) and psychiatric disorders (Liang et al., 2013;Gao et al., 2014;Wang et al., 2015;Yang et al., 2015Yang et al., , 2016. Compared to other resting state functional connectivity (rs-FC) methods such as seed-based rs-FC analysis or independent component analysis, ReHo can target key regions involved in the pathophysiology of neck pain rather than elucidating the connectivity between and among different regions (Yu et al., 2014b).
Neck pain accompanied by neck stiffness and discomfort is the earliest symptom of CS. Patients who do not receive appropriate treatment in the initial stage of CS later develop a complex series of clinical symptoms including spinal cord compression (i.e., CSM). Recently, studies have investigated brain functional changes in CSM (Kowalczyk et al., 2012;Zhou et al., 2014;Zhou F. Q. et al., 2015) and found patients with CSM are associated with cortical reorganization in the primary motor cortex (Holly et al., 2007;Duggal et al., 2010;Tam et al., 2010;Kowalczyk et al., 2012;Bhagavatula et al., 2016), primary somatosensory cortex Bhagavatula et al., 2016), and sensorimotor cortex (Dong et al., 2008;Zhou et al., 2014;Tan et al., 2015;Zhou F. Q. et al., 2015).
Consistent with previous findings, we found CSNP patients were associated with reduced rs-ReHo in the sensorimotor cortex compared to HCs. This result is also consistent with previous studies that report alterations in sensorimotor cortices in conjunction with clinical deterioration in CSM patients (Dong et al., 2008;Tan et al., 2015). Some studies have found that sensorimotor cortices normalized with the recovery of function following surgical spinal cord decompression in CSM patients (Holly et al., 2007;Dong et al., 2008;Duggal et al., 2010;Tam et al., 2010;Green et al., 2015;Tan et al., 2015;Bhagavatula et al., 2016). In our study, we did not find any significant rs-ReHo changes in the sensorimotor cortex after treatment. We speculate these differences might be due to the different stages of CS investigated in the two studies. Patients with CSM often show motor and sensory deficits, while CSNP patients often show discomfort of the neck including chronic neck pain and stiffness.
In the present study, decreased ReHo in the sensorimotor cortex was mainly localized to S1. Previous studies have indicated that S1 plays a crucial role in nociceptive processing pathways in the human brain (Bushnell et al., 1999;Ploner et al., 1999). In our previous study we also found that chronic low back pain patients have increased cortical thickness and increased volume in S1 and decreased voxel-by-voxel functional connectivity during low pain conditions at S1 compared to HCs . Our results provide further evidence of rs-ReHo changes in the sensorimotor cortex in CSNP patients, and we speculate this may reflect aberrant central processing of neck pain afferentiation.
We also found reduced rs-ReHo in the rTPJ, including the posterior end of the right superior temporal gyrus and right inferior parietal lobule, in CSNP patients compared to HCs. Our findings were in line with previous studies that reported abnormal structure and function in the TPJ for various chronic pain conditions. It is worth noting that the literature sometimes uses the superior temporal gyrus and inferior parietal lobule instead of the TPJ. For instance, Friebel et al. (2011) used coordinate-based meta-analysis to find that patients with chronic persistent neuropathic pain were associated with activation of the right inferior parietal lobe, a component of the TPJ. The TPJ is an approximate term that generally refers to the regions including the posterior end of the superior temporal gyrus, the inferior parietal lobule, and the lateral occipital gyrus in the human brain (Mars et al., 2012;Krall et al., 2015). Studies have suggested that the TPJ is an area in the mirror neuron system, which maintains a coherent representation of our own, as well as others' , body (Iacoboni and Dapretto, 2006;Tsakiris et al., 2008;Rizzolatti and Sinigaglia, 2010) and is also involved in imitative motor learning and control (Buccino and Riggio, 2006;Fogassi, 2011). The mirror neuron system also participates in pain perception (Hojat and Cohen, 2012) and includes the brain regions associated with pain experience, such as the cingulate and insular cortices (Ramachandran and Oberman, 2006). Studies have suggested that the TPJ receives information from the thalamus and limbic system and projects to many other brain areas associated with pain perception and modulation, such as the lateral anterior prefrontal cortex, ventral prefrontal cortex, anterior insula, posterior cingulate, and anterior medial prefrontal cortex (Mars et al., 2012). In our study, we found that effective longitudinal acupuncture treatments can increase reduced rs-ReHo in the rTPJ, and this normalization is associated with a reduction in neck pain intensity in CSNP patients. More studies are needed to explore the role of the TPJ in the development of CSNP.
Machine learning techniques have been widely applied in translational neuroimaging studies to provide a basis for identifying neuropathological features of different diseases and show potential clinical utility beyond current clinical diagnostic categories (Pereira et al., 2009;Arribas et al., 2010;Orru et al., 2012;Zeng et al., 2012;Abi-Dargham and Horga, 2016;Woo et al., 2017). The earliest work identified brain signatures for neurological diseases, particularly Alzheimer's Disease and related dementia (Kippenhan et al., 1992). The studies have been extended to other neurological and psychiatric conditions, such as Parkinson's disease (Tang et al., 2010), major depression (Zeng et al., 2012), schizophrenia (Shen et al., 2010), and ADHD (Hart et al., 2014). Most studies have focused on identifying brain signatures for discriminating patients from HCs and consequently established a meaningful neurophysiological basis for the disorder of interest. Another machine learning model uses a patient's brain measurements to predict treatment outcomes. Most focus on depression and anxiety disorders (Doehrmann et al., 2013;Hahn et al., 2015;Whitfield-Gabrieli et al., 2016) as well as schizophrenia (Sarpal et al., 2016) and Parkinson's disease (Ye et al., 2016).
To the best of our knowledge, the present study is the first attempt to classify CSNP patients from healthy controls using resting-state fMRI and a state-of-the-art machine learning approach. Therefore, it provides a possible method of using objective measures to diagnose CSNP patients. More importantly, we used baseline resting-state fMRI to predict clinical symptom changes for a longitudinal treatment with an error of 19.6% and prediction-outcome correlation of 0.30. The present approach offers a useful tool for clinicians when determining the effectiveness of acupuncture treatment before executing treatment procedures.

Limitations
It is worth noting that the aim of this study was to investigate the brain's neural plasticity in CSNP patients, how an effective treatment can modulate altered neural plasticity, and the association between ReHo changes and clinical symptom changes. This paradigm has been used in previous studies from our group  and from other groups (Zhou and Jin, 2008;Rodriguez-Raecke et al., 2009;Sato et al., 2009;Fang et al., 2012). Although our results may provide some novel insights into the mechanisms underlying effective treatment, the lack of sham treatment as a control has significantly limited our ability to infer an underlying longitudinal acupuncture treatment mechanism. Additional controls should be applied in future studies to elucidate the specific effect of the corresponding treatment.

CONCLUSION
We found that CSNP patients are associated with increased regional coherences at the right temporo-parietal junction and left sensorimotor cortex. Additionally, the alleviation of symptoms after treatment can normalize the reduced rs-ReHo in CSNP patients, and this normalized rs-ReHo value is also associated with a reduction in neck pain intensity. We were also able to use MVPA to discriminate CSNP patients from HCs with an accuracy of 72% as well as predict treatment responses using baseline data with an error of 19.6%. Our results demonstrate the crucial role of the rTPJ in the pathophysiology and development of non-specific neck pain.

AUTHOR CONTRIBUTIONS
BL and JL proposed the study and were the guarantors. JC, ZW, YT, KJ, JP, CL, and JK analyzed the data and/or wrote the original and revised drafts of the manuscript. All authors contributed to the design and/or interpretation of the study.