Sensorimotor and pain-related alterations of the gray matter and white matter in Type 2 diabetic patients with peripheral neuropathy.

Although diabetic peripheral neuropathy (DPN) has long been considered a disease of the peripheral nervous system, recent neuroimaging studies have shown that alterations in the central nervous system may play a crucial role in its pathogenesis. Here, we used surface-based morphometry (SBM) and tract-based spatial statistics (TBSS) to investigate gray matter (GM) and white matter (WM) differences between patients with DPN (n = 67, 44 painless and 23 painful) and healthy controls (HCs; n = 88). Compared with HCs, patients with DPN exhibited GM abnormalities in the pre- and postcentral gyrus and in several deep GM nuclei (caudate, putamen, medial pallidum, thalamus, and ventral nuclear). They also exhibited altered WM tracts (corticospinal tract, spinothalamic tract, and thalamocortical projecting fibers). These findings suggest impaired motor and somatosensory pathways in DPN. Further, patients with DPN (particularly painful DPN) exhibited morphological differences in the cingulate, insula, prefrontal cortex, and thalamus, as well as impaired WM integrity in periaqueductal WM and internal and external capsules. This suggests pain-perception/modulation pathways are altered in painful DPN. Intermodal correlation analyses found that the morphological indices of the brain regions identified by the SBM analysis were significantly correlated with the fractional anisotropy of brain regions identified by the TBSS analysis, suggesting that the GM and WM alterations were tightly coupled. Overall, our study showed sensorimotor and pain-related GM and WM alterations in patients with DPN, which might be involved in the development of DPN.

INTRODUCTION

Diabetic peripheral neuropathy (DPN) is one of the most common and distressing complications in patients with diabetes mellitus (DM), affecting ~30–50% of the diabetic population worldwide. Typically, DPN manifests with positive sensory symptoms in the feet such as tingling, prickling, and pain, as well as negative symptoms such as numbness, which results in considerable disability and suffering (Feldman, Nave, Jensen, & Bennett, 2017; Selvarajah et al., 2014). For years, DPN has been considered a disease solely of the peripheral nervous system; However, evidence now suggests that alterations in the central nervous system (CNS) also play a crucial role in its pathogenesis (Selvarajah, Wilkinson, Davies, et al., 2011; Selvarajah, Wilkinson, Gandhi, et al., 2011; Selvarajah et al., 2014; West, Bannister, Dickenson, & Bennett, 2015). Using genetic engineering and tracing technology, neurophysiologists have reported that DPN is associated with altered Rac‐1‐mediated dendritic spine morphology in spinal dorsal horn neurons and P2X4 receptor‐mediated release of brain‐derived neurotrophic factor in dorsal horn microglia. These phenomena were associated with spontaneous activity and hyperexcitable responses to peripheral stimuli (Beggs, Trang, & Salter, 2012; Feldman et al., 2017). However, changes in the brain remain largely unknown. Investigating the association between DPN and any changes to specific brain regions might help us better understand the pathologic mechanisms of this disorder.

Using advanced neuroimaging techniques, studies have provided valuable insights into abnormalities of the brain parenchyma in patients with DPN. For example, a resting state functional MRI (fMRI) study found decreased functional connectivity between the thalamus and the cortex in patients suffering from peripheral neuropathic pain (Cauda et al., 2009). This finding suggested that diabetic pain could alter thalamocortical connections, causing a disruption of thalamic feedback (Cauda et al., 2009). A perfusion‐weighted imaging study demonstrated that patients with painful DPN had greater thalamic relative cerebral blood flow and longer first moment transit time than patients with painless DPN, suggesting that painful DPN is related to increased thalamic vascularity (Selvarajah, Wilkinson, Davies, et al., 2011; Selvarajah, Wilkinson, Gandhi, et al., 2011). Similarly, using magnetic resonance spectroscopy, Selvarajah et al. (2008) observed a significantly lower N‐acetyl aspartate to creatine ratio in the thalamus of patients with DPN, suggesting thalamic neuronal dysfunction. One task‐based fMRI study (Li et al. 2018) found that after several sessions applying thermal stimuli to the lower extremities, patients with DPN showed stronger brain activation in somatosensory pathways (e.g., in the right insula and left caudate nucleus) than healthy volunteers or diabetic patients without DPN. A multimodal study combining blood oxygenation level‐dependent and arterial spin labeling data (Segerdahl, Themistocleous, Fido, Bennett, & Tracey, 2018) reported altered ventrolateral periaqueductal gray (PAG) functional connectivity in patients suffering from painful DPN, and suggested this could be an index of brain‐based pain facilitation for painful DPN. While these studies were conducted primarily in the functional domain, far less attention has been paid to anatomical alterations of the brain in DPN. Of note, a voxel‐based morphometry study reported decreased gray matter (GM) volume in the supramarginal gyrus, cingulate cortex, and primary somatosensory cortex in Type‐1 patients with diabetes and DPN (Selvarajah et al., 2014), indicating that brain structural abnormalities play an important role in the pathogenesis of DPN in Type‐1 DM. Given that Type‐1 and Type‐2 DM are considered to have a disparate pathogeneses (Callaghan, Hur, et al., 2012; Callaghan, Little, et al., 2012; Feldman et al., 2017), whether GM alterations are similarly associated with Type‐2 DPN has yet to be determined. Additionally, considering the frequently reported pathological changes of axonal loss, focal areas of demyelination, and a reduction in myelinated fiber density in nerve tracts of patients with DPN (Said, Baudoin, & Toyooka, 2008; Selvarajah, Wilkinson, Davies, et al., 2011; Selvarajah, Wilkinson, Gandhi, et al., 2011), a comprehensive assessment of this disease should simultaneously consider alterations in GM and WM.

In this study, we performed a range of structural analyses using T1‐weighted MRI and diffusion tensor imaging (DTI) data to characterize GM and WM alterations in a cohort of Type 2 diabetic patients with peripheral neuropathy (further divided into subgroups of painful and painless DPN due to their disparate disease stages, treatment options and prognoses [Sorensen, Molyneaux, & Yue, 2006; Tesfaye et al., 2010]), compared with healthy controls (HCs). Considering the aforementioned neuroimaging studies and the fact that DPN is a predominantly sensory and autonomic neuropathy with sensory defects, neuropathic pain, and occasionally motor deficits (Said et al., 2008; Selvarajah, Wilkinson, Davies, et al., 2011; Selvarajah, Wilkinson, Gandhi, et al., 2011; Selvarajah et al., 2019), we hypothesized that patients with DPN would exhibit GM and WM structural alterations in (a) areas associated with somatic sensation (e.g., primary somatosensory cortex, thalamus, spinothalamic tract, and thalamocortical projecting fibers), (b) areas associated with pain perception and modulation (e.g., anterior and mid‐cingulate cortex [ACC; MCC], thalamus, insula and prefrontal cortex [PFC], the internal and external capsule, corpus callosum, PAG and the adjacent WM), and (c) areas associated with motor regions (e.g., primary motor cortex and the basal ganglia).

MATERIALS AND METHODS

Participants

Sixty‐seven consecutive patients with Type‐2 diabetes and a diagnosis of peripheral neuropathy were enrolled from the Department of endocrinology at Xiangya Hospital, Central South University in China. The diagnosis of DPN was made according to the American or Toronto consensus criteria (Boulton, Valensi, & Tesfaye, 2011; Pop‐Busui et al., 2017). Specifically, the presence of a nerve‐conduction abnormality and at least one symptom (e.g., stabbing, prickling, burning, or aching pain in the legs, feet, or toes) or sign (e.g., decreased or absent ankle reflexes or symmetrically decreased distal sensation) of neuropathy confirmed DPN. If nerve conduction was normal, a validated measure of small fiber neuropathy was also sufficient for a diagnosis of DPN (Tesfaye et al., 2010). Inclusion criteria were as follows: right‐handed, aged between 31 and 68 years, and stable glycemic control (HbA1c < 11%) (Selvarajah et al., 2014; Selvarajah et al., 2019). Patients were excluded from the study if they exhibited: (a) nondiabetic neuropathies or diabetic neuropathies other than DPN; (b) cerebrovascular, neurological, or psychiatric diseases; (c) an unknown hypoglycemic status; (d) prior substantial head trauma or surgery; (e) brain tumors; (f) recurrent severe hypoglycemia; (g) alcoholism or drug abuse; or (h) serious brain atrophy or any other major intracranial disease, such as epilepsy. Eighty‐eight age‐ and sex‐matched, nondiabetic, right‐handed HCs were also included in this study. The specific inclusion criteria were as follows: (a) no chronic pain conditions, analgesic medications, or alternative therapies for treatment of pain; (b) no history of brain tumors, head trauma or surgery, alcoholism, or drug abuse; (c) absence of neurological and psychiatric diseases; (d) no use of sedatives or hypnotics such as benzodiazepines for at least 1 month before the study; (e) no history of hypertension, cardiopathy, nephritis, chromatodysopia, or other major medical illness. To minimize the confounding effects of white matter hyperintensities (WMHs) on the brain structure and function (Chao et al., 2013; Papma et al., 2014; Prins & Scheltens, 2015; Ylikoski et al., 1995), the subjects those who were rated as Fazekas 2 and 3 (Fazekas, Chawluk, Alavi, Hurtig, & Zimmerman, 1987) on T2‐weighted MR images were excluded from this study (including 17 DPN patients and 6 HCs). This prospective study was approved by the Medical Research Ethics Committee of Xiangya Hospital, Central South University, and written informed consent was obtained from all participants before participation.

Assessment of clinical characteristics

All patients completed the standardized screening forms to assess their basic clinical parameters (blood pressure, height, and weight), medical history, pain symptoms, and medication use. All participants completed a structured neurological examination including standard Quantitative Sensory Testing (QST), the Neuropathy Symptom Score (NSS) questionnaire, Neuropathy Disability Score (NDS), and the Douleur Neuropathique 4 (DN4) questionnaire. Nerve conduction tests were performed at a stable skin temperature of 31°C and a room temperature of 24°C. The nerve conduction attributes (including the nerve conduction velocity and nerve amplitude) in both the motor and/or sensory nerves in the upper and lower extremities were documented. Specifically, motor nerve conduction testing was performed in the bilateral median nerve, ulnar nerve, peroneal nerve and tibial nerve, and sensory conduction testing was performed only in bilateral median nerve, ulnar nerve, and sural nerve. Medication history as well as laboratory and clinical investigation data from the clinical records were also collected. The presence of DPN was confirmed using a comprehensive upper and lower limb neurological examination, nerve conduction testing, electromyography, the battery of QST measures, and the assessment of clinical symptoms and signs. The presence of chronic neuropathic pain caused by DPN was determined at the time of the clinical assessment and was in line with the International Association for the Study of Pain (IASP) definition of neuropathic pain. Based on these clinical and neurophysiological assessments, DPN participants were further divided into two groups: (a) a painless DPN group (n = 44), comprising patients with abnormal clinical findings, painless neuropathy, and at least two abnormalities observed on neurophysiological examinations; and (b) a painful DPN group (n = 23), consisting of patients with symptoms of chronic pain lasting at least 6 months and both clinical and neurophysiological abnormalities (Boulton et al., 2011; Selvarajah et al., 2014).

MAGNETIC RESONANCE IMAGING

Data acquisition

Brain MRI scans were obtained on a Siemens Magnetom Prisma 3.0T MR scanner using a 64‐channel head coil. Routine imaging sequences, including axial T1‐weighted images, T2‐weighted images, and T2‐weighted FLAIR images were obtained for every participant to rule out gross brain pathology. For each participant, the three‐dimensional T1‐weighted MRI scans were obtained using a magnetization prepared rapid acquisition gradient‐echo sequence with the following parameters: thickness/gap = 1.0/0 mm, repetition time (TR) = 2,300 ms, echo time (TE) = 2.98 ms, inversion time = 900 ms, 176 sagittal slices, field of view (FOV) = 256 mm × 256 mm, matrix size = 256 × 256, flip angle = 9°, voxel size = 1.0 × 1.0 × 1.0 mm3, sequence scan time = 5.2 min. An echo planar imaging sequence with a simultaneous multi‐slice technique was applied to the DTI scans using the following scan parameters: TR/TE = 4,200/68 ms; slice thickness = 2.0 mm; matrix = 128 × 128; FOV = 264 × 264 mm; number of diffusion‐encoding directions = 64; b‐value = 1,000 s/mm2; number of excitations for b0 = 3, and for b1000 = 1; layer spacing = 0; and total acquisition time = 5.05 min.

Study outline

The main analytical procedures of the present study can be outlined as follows: First, surface‐based morphometry (SBM) analysis was conducted to compare the cortical thickness and cortical surface area between pairs of groups. Second, volumetric and shape analyses were performed to identify differences between groups in several deep GM nuclei, such as the thalamus and caudate. Third, tract‐based spatial statistics (TBSS) analysis of four diffusion parameters (fractional anisotropy [FA], mean diffusivity [MD], axial diffusivity [AD], and radial diffusivity [RD]) was conducted to determine whether patients and controls differed in WM microstructure. Four, intermodal correlation analyses were performed to investigate the relationships between the cortical morphological indices (including cortical thickness and cortical surface area) and FA.

Preprocessing for SBM

Each T1‐weighted MRI scan was processed using the Freesurfer package (version 6.0.0) to obtain a three‐dimensional model of the cortical surface for cortical thickness and surface area measurements. The technical details of these processing procedures can be found in previous publications (Wang et al., 2018; Zhang et al., 2014, 2017). Briefly, the T1‐weighted MRI scan was intensity‐corrected and skull‐stripped. The resultant volume was segmented with a connected‐components algorithm to obtain a voxel‐based GM/WM boundary for each hemisphere, which was then tessellated to produce a triangle‐based GM/WM boundary surface. This initial triangle‐based GM/WM surface was further refined by an automated topological correction procedure, resulting in a topologically correct GM/WM surface, which hereafter is referred as the white surface. Subsequently, the white surface was deformed outward using a deformable surface algorithm to identify the pial surface. For all participants, both the white surface and the pial surface were visually inspected, and if necessary, manually corrected according to the software guidelines. After the generation of these surface models, cortical thickness and cortical surface area maps were obtained. Specifically, cortical thickness was calculated as the distance between the white surface and the pial surface using the T‐average algorithm; cortical surface area for each vertex was calculated as one third of the area of each triangle to each of its vertices. Before statistical analysis, all the cortical maps were resampled onto the standard space (fsaverage template) and further smoothed with a 20‐mm heat kernel. The preprocessing for deep GM nuclei volumetric analysis can be found in Supporting Information.

Preprocessing for deep GM nuclei shape analysis

For vertex‐based shape analysis, deep GM nuclei, including the thalamus, hippocampus, amygdala, putamen, globus pallidus, caudate, and nucleus accumbens (NAc), were segmented using the FIRST procedure in FSL (Du et al., 2011; Zhang et al., 2013). For each deep GM nucleus, FIRST created a surface mesh using a deformable mesh model. The surface mesh was composed of a set of vertices and triangles. The number of vertices for each structure was fixed so that the spatial location of corresponding vertices could be compared across participants. Although the vertices were in correspondence, the surface mesh resides in the native space of each participant, and therefore may have an arbitrary orientation and/or position. To remove potential pose differences (rotation and translation) for each structure, the surface mesh of each participant was aligned to the mean surface in standard space by minimizing the sum‐of‐squares difference between a participant's surface and the mean surface. Before statistical analysis, all segmentations were visually inspected for accuracy.

Preprocessing for TBSS

DTI data were processed using FSL (Castellano et al., 2016; Guimaraes et al., 2013). The DTI scans for each participant were first corrected for head motion and eddy currents using affine registration and then brain‐extracted using the Brain Extraction Tool. The FA, MD, AD, and RD maps were obtained by fitting a diffusion tensor model to each voxel. Voxel‐wise analysis of these diffusion maps was carried out using the TBSS procedure in FSL. The FA maps for all participants were aligned to Montreal Neurological Institute (MNI) space using the nonlinear registration tool (FNIRT). The mean image of all aligned FA images was created and thinned to obtain a mean FA skeleton with an FA threshold value greater than 0.2. The aligned FA map of each participant was then projected onto the skeleton, and the resulting data were entered into voxel‐wise statistics. By using the protocol of non‐FA images in TBSS, the MD, AD, and RD images were also aligned to MNI space and projected onto the mean FA skeleton.

STATISTICAL ANALYSIS

Demographic and clinical variables

Statistical analyses were performed using SPSS 18.0 (SPSS Inc., Chicago, IL). Participant clinical characteristics are presented as the mean ± standard deviation (SD) with a normal distribution, median and interquartile range with a nonnormal distribution, or distribution of the frequency with qualitative data. Group differences in QST z‐scores were assessed using one‐way ANOVA analysis and all other data with Kruskal–Wallis nonparametric tests. Categorical data were analyzed with a Chi‐squared test. The significance level for the comparisons was set at p < .05.

SBM analysis

Vertex‐wise contrasts of cortical thickness and cortical surface area were performed for each pair of groups using the SurfStat package (http://www.math.mcgill.ca/keith/surfstat/). Specifically, each contrast was entered into a vertex‐wise general linear model (GLM) with group, sex, age, and intracranial volume (ICV) as covariates. Results were first thresholded at vertex‐wise p < .005 and then corrected for multiple comparisons at the cluster level using random field theory (RFT). The significance level was set at corrected p < .05.

Deep GM nuclei shape analysis

Contrasts of the spatial locations for each vertex in each deep GM nucleus were performed in MNI space for each pair of groups using GLM in FSL. Specifically, each contrast was entered into a vertex‐wise GLM with group, sex, age, and ICV as covariates. A nonparametric permutation test with 5,000 repetitions was then used to perform the group‐wise statistics. The results were corrected for multiple comparisons using threshold‐free cluster enhancement (TFCE). The significance level was set at corrected p < .05. The statistical analysis for the deep GM nuclei volumetric analysis can be found in Supporting Information.

TBSS analysis

Voxel‐wise contrasts for FA, MD, AD, and RD maps were performed for each pair of groups using GLM in FSL. Specifically, each contrast was entered into a voxel‐wise GLM with group, sex, and age as covariates. A nonparametric permutation test with 5,000 repetitions was then used to perform the group‐wise statistics. The results were corrected for multiple comparisons using TFCE. The significance level was set at corrected p < .05.

Correlations between imaging parameters and neurological test scores

In the DPN group, vertex‐wise correlation analyses were performed to investigate the relationships between some neurological test scores (NSS, NDS, and DN4) and the cortical thickness/cortical surface area of the significant regions intersection that was revealed by the painless DPN versus HCs and painful DPN versus HCs contrasts in the SBM analysis. Specifically, the vertex‐wise GLM was fitted with the individual neurological test score as a covariate. Results were first thresholded at vertex‐wise p < .005 and then corrected for multiple comparisons at the cluster level using RFT. The significance level was set at corrected p < .05. Using GLM in FSL, vertex‐wise correlation analyses were performed to examine the relationship between deep GM nuclei shapes and the neurological parameters. Specifically, each contrast was entered into a vertex‐wise GLM with each neurological test score as a covariate. A nonparametric permutation test with 5,000 repetitions was then used to obtain the vertex‐wise statistics. The results were corrected for multiple comparisons using TFCE. The significance level was set at corrected p < .05.

In the DPN group, voxel‐wise correlation analyses were performed to investigate the relationships between the neurological test scores and the FA/MD/RD of the intersection of the significant regions that was revealed by the painless DPN versus HCs and painful DPN versus HCs contrasts in the TBSS analysis. Specifically, the voxel‐wise GLM was fitted with the individual neurological test score as a covariate. A nonparametric permutation test with 5,000 repetitions was then used to perform the group‐wise statistics. The results were corrected for multiple comparisons using TFCE. The significance level was set at corrected p < .05.

Exploratory intermodal correlation analysis

Using the pooled data from the three groups, correlation analyses were performed to investigate the relationships between the cortical morphological indices and FA. Specifically, for each cortical morphological index, the union set of all the brain regions showing a significant between‐group difference was saved as a binary mask. For the corresponding morphological index, the mean value for all the vertices within the mask was extracted. Likewise, the mean FA of the union set for all the brain regions showing a significant between‐group difference in FA was extracted. Then, Pearson's correlation coefficient was calculated between the mean FA value and the mean thickness and surface area values. The significance level was set at p < .05.

RESULTS

Demographic characteristics and neurophysiological assessments

Demographic characteristics and results of the neurophysiological assessments for patients with painless DPN, patients with painful DPN, and HCs are presented in Table 1. Participants were matched for age (p = .14), sex (p = .36), HbA1c level (p = .69), body mass index (BMI; p = .67), duration of diabetes (p = .08), diabetic retinopathy (p = .94), nerve conduction velocity, and QST parameters. NSS, NDS, and DN4 scores were significantly higher in the painful DPN group than in the painless DPN group (NSS: p = .001; NDS: p < .001; DN4: p < .001). The duration of painful DPN was 27.41 ± 20.53 months (range: 6–72 months), with data missing for seven patients. In the painful DPN group, three patients were treated with the painkiller gabapentin, and no analgesic treatment was given in 12 patients. The detailed data for neuropathic pain treatment in the remaining eight patients was missing.

Table 1

 Clinical characteristics and neurophysiological assessments of study subjects

Characteristics	HCs	Painless DPN	Painful DPN	p value
	(n = 88)	(n = 44)	(n = 23)
Age (years), mean ± SD	55.58 ± 0.83	54.07 ± 1.15	58.74 ± 1.91	.14
Gender, n
Male	56 (63.6)	28 (63.6)	11 (47.8)	.36
Female	32 (36.4)	16 (36.4)	12 (52.2)
HbA1c (%) [mmol/mol]		9.44 ± 1.93	9.20 ± 2.52	.69
BMI		23.71 ± 3.35	23.27 ± 3.03	.67
Duration of diabetes		7.43 ± 5.22	10.19 ± 5.42	.08
Diabetic retinopathy		13/44 (29.5)	7/23 (30.4)	.94
Motor nerve fiber‐nerve conduction velocity
Right PN_MCV (m/s)		42.95 ± 5.67	42.65 ± 4.93	.85
Left PN_MCV (m/s)		44.13 ± 3.73	42.07 ± 6.91	.22
Right TN_MCV (m/s)		41.91 ± 4.61	40.17 ± 7.08	.4
Left TN_MCV (m/s)		41.90 ± 5.11	39.82 ± 6.23	.29
Right MN_MCV (m/s)		54.74 ± 5.35	54.73 ± 5.59	.17
Left MN_MCV (m/s)		54.93 ± 4.33	52.97 ± 4.94	.99
Right UN_MCV (m/s)		52.98 ± 6.34	51.22 ± 9.37	.47
Left UN_MCV (m/s)		53.35 ± 5.90	51.68 ± 5.93	.38
Motor nerve fiber‐nerve amplitude
Right PN (mv)		5.28 ± 3.20	4.45 ± 2.54	.36
Left PN (mv)		5.36 ± 3.22	3.95 ± 2.88	.14
Right TN (mv)		15.27 ± 6.30	10.94 ± 5.14	.02a
Left TN (mv)		14.85 ± 5.44	13.36 ± 7.61	.45
Right MN (mv)		12.71 ± 4.50	11.87 ± 4.12	.54
Left MN (mv)		11.95 ± 4.09	12.73 ± 2.96	.50
Right UN (mv)		13.56 ± 3.64	13.05 ± 4.12	.67
Left UN (mv)		13.34 ± 3.11	11.51 ± 2.45	.046a
Sensory nerve fiber‐nerve conduction velocity
Right SN_SCV (m/s)		44.31 ± 3.53	41.16 ± 7.27	.14
Left SN_SCV (m/s)		41.3 ± 9.77	42.00 ± 5.02	.80
Right MN_SCV (m/s)		50.97 ± 6.13	51.38 ± 6.41	.84
Left MN_SCV (m/s)		52.85 ± 7.39	52.48 ± 7.68	.88
Right UN_SCV (m/s)		51.92 ± 5.26	51.69 ± 7.62	.92
Left UN_SCV (m/s)		51.29 ± 5.64	51.17 ± 7.77	.96
Sensory nerve fiber‐nerve amplitude
Right sural nerve (uv)		11.68 ± 7.96	8.12 ± 4.29	.16
Left sural nerve (uv)		10.42 ± 8.75	8.55 ± 5.20	.45
Right MN (uv)		16.46 ± 12.48	12.04 ± 11.16	.26
Left MN (uv)		14.83 ± 11.59	11.25 ± 7.02	.27
Right UN (uv)		12.15 ± 15.05	9.83 ± 13.93	.64
Left UN (uv)		10.13 ± 7.05	8.67 ± 5.54	.53
Quantitative sensory testing
CDT		26.43 ± 3.24	25.30 ± 6.11	.51
WDT		43.22 ± 3.41	44.71 ± 3.26	.24
Z‐score for CDT		−0.02 ± 1.34	−1.07 ± 3.59	.31
Z‐score for WDT		1.36 ± 1.32	1.88 ± 1.57	.29
NSS		3.28 ± 2.75	6.17 ± 2.57	.001a
NDS		1.22 ± 1.50	2.94 ± 1.51	<.001a
DN4		1.31 ± 1.91	4.29 ± 2.26	<.001a

Note: Data in parentheses are percentages. The lower limiting value for PNCV is 40 m/s, for TNCV is 40 m/s, for SNCV is 45 m/s. The lower limiting value for PN amplitude is 3 mv in subjects under the age of 60, for TN amplitude is 7 mv, for SN amplitude is 8 uv in male and 6 uv in femle. The normal range for Z score is from −1.96 to 1.96.

Abbreviations: BMI, body mass index; CDT, cold detection threshold; DN4, douleur neuropathique en 4 questionnaire; DPN, diabetic peripheral neuropathy; HCs, healthy controls; MN_MCV, median nerve motor conduction velocity; MN_SCV, median nerve sensory conduction velocity; NDS, neuropathy disability score; NSS, neuropathy symptom score; PN_MCV, peroneal nerve motor conduction velocity; SD, standard deviation; SN_SCV, sural nerve sensory conduction velocity; TN_MCV, tibial nerve motor conduction velocity; UN_MCV, ulnar nerve motor conduction velocity; UN_SCV, ulnar nerve sensory conduction velocity; WDT, warm detection threshold.

Denotes differences were statistically significant.

Cortical thickness analysis

Compared with HCs, patients with painless DPN exhibited thinner cortex across the brain, including the bilateral insula, cingulate cortex, right precentral gyrus and paracentral lobule (for details, see Figure 1a and Table 2; p < .05, RFT corrected). Compared with HCs, patients with painful DPN also exhibited thinner cortex across the brain, including the right insula, ACC, MCC, and right PFC (for details, see Figure 1a and Table 2) (p < .05, RFT corrected). No significant differences in cortical thickness were observed between the painless and painful DPN groups.

Figure 1

Pair‐wise intergroup comparisons of cortical thickness and cortical surface area. Brain regions with altered cortical thickness (a) or cortical surface area (b) between the indicated groups are shown. Significance is represented on a log (p value) scale, in which positive values (warm colors) are assigned to the painless DPN > HCs, painful DPN > HCs, and painful DPN > painless DPN clusters and negative values (cold colors) are assigned to the painless DPN < HCs, painful DPN < HCs, and painful DPN < painless DPN clusters. LH, left hemisphere; RH, right hemisphere

Table 2

 Clusters showing cortical thickness differences between painless DPN, painful DPN and HCs

Between‐group comparisons	Clusters	Brain regions	Talaraich coordinates			Side	Cluster size (vertices)	Peak T values	p value
			X	Y	Z
Painless DPN < HCs	C′1	ATL and insula	−31	24	5	L	17,338	5.58	<0.001
	C′2	CC	−11	−14	45	L	10,657	5.22	<0.001
	C′3	IPL and LTL	68	−42	−2	R	9,319	4.62	<0.001
	C′4	MTL (Para_Hippo_G and FFG)	34	−41	−15	R	5,701	4.64	<0.001
	C′5	Para_CL and CC	9	−24	47	R	5,827	3.80	<0.001
	C′6	Insula and IFG	37	0	14	R	4,571	4.77	<0.001
	C′7	IPL	−53	−48	48	L	3,147	3.45	0.003
	C′8	SFG	−22	11	51	L	3,682	3.46	0.005
	C′9	Pre_CUN	9	−63	22	R	2,524	4.49	0.018
Painful DPN < HCs	C′10	CC	−9	0	42	L	5,989	4.49	<0.001
	C′11	Pre_FL and insula	35	8	10	R	4,754	5.30	<0.001
	C′12	ACC and Mid_CC	3	−16	29	R	5,074	3.77	<0.001
	C′13	MTL (Para_Hippo_G and FFG)	−35	−39	−14	L	2,459	3.99	0.004
	C′14	MTL (Para_Hippo_G and FFG)	21	−12	−24	R	2,674	6.16	0.005

Abbreviations: ACC, anterior cingulated cortex; ATL, anterior temporal lobe; CC, cingulated cortex; DPN, diabetic peripheral neuropathy; HCs, healthy controls; IFG, inferior frontal gyrus; IPL, inferior parietal lobule; L left; LTL, lateral temporal lobe; Mid_CC, middle cingulated cortex; MTL, medial temporal lobe; Para_CL, paracentral lobule; Para_Hippo_G and FFG, parahippocampal gyrus and fusiform gyrus; SFG, superior frontal gyrus; Pre_CUN, precuneus; Pre_FL, prefrontal lobe; R right.

Cortical surface area analysis

Cortical area in the right insula, left posterior cingulate cortex (PCC), and precuneus was greater in painless DPN groups than in the HCs (p < .05, RFT corrected). Patients with painful DPN showed the same trend, but in a wider range of brain regions, mainly in the right insula, bilateral precentral and postcentral gyrus, bilateral paracentral lobule, and PCC (for details, see Figure 1b and Table 3; p < .05, RFT corrected). No significant differences were observed between the painless and painful DPN groups.

Table 3

 Clusters showing cortical surface area differences between painless DPN, painful DPN and HCs

Between‐group comparisons	Clusters	Brain regions	Talaraich coordinates			Side	Cluster size (vertices)	Peak T values	p value
			X	Y	Z
Painless DPN > HCs	C″1	Insula	32	26	1	R	2,414	3.79	0.002
	C″2	PCC	−4	−42	24	L	1902	3.45	0.006
	C″3	Para_CL, Pre_CG and PoCG	−4	−31	57	L	4,918	4.75	<0.001
	C″4	TP	−48	17	−21	L	4,342	4.18	0.001
Painful DPN > HCs	C″5	Insula	37	13	2	R	2085	4.64	0.005
	C″6	Para_CL, Pre_CG and PoCG	1	−29	71	R	2,293	3.66	0.024
	C″7	PoCG	−35	−29	68	L	1,788	3.60	0.030
	C″8	PCC	−2	−32	36	L	1,301	3.32	0.046

Abbreviations: DPN, diabetic peripheral neuropathy; HCs, healthy controls; L, left; Para_CL, paracentral lobule; PCC, posterior cingulated cortex; Po_CG, postcentral gyrus; Pre_CG, precentral gyrus; TP, temporal lobe; R, right.

Compared with HCs, patients with painless DPN exhibited substantial contractions in the medial posterior segment of the bilateral thalamus, left caudate, and anterior segment of the left putamen (p < .05, TFCE corrected) (Figure 2a). Similarly, the medial posterior segment of the bilateral caudate, bilateral thalamus, left NAc, and the left putamen were substantially more contracted in patients with painful DPN than in HCs. Further, these patients also showed a substantial expansion in the lateral anterior segment of the bilateral pallidum, left thalamus, and the lateral segment of the left caudate (p < .05, TFCE corrected). Overall, patients with painful DPN showed a wider range of regional contractions in deep GM nuclei than did patients with painless DPN, although this was only a trend (for details, see Figure 2b).

Figure 2

Pair‐wise intergroup comparison of deep gray matter nuclei shape. (a) Differences in the shape of deep gray matter nuclei for painless DPN versus HCs. (b) Differences in the shape of deep gray matter nuclei or painful DPN versus HCs. Colored areas denote regions where a significant difference in deep‐gray‐matter nucleus shape was observed between the indicated groups. Red hues are assigned to the painless DPN > HCs and painful DPN > HCs clusters and blue hues are assigned to the painless DPN < HCs and painful DPN < HCs clusters. S, superior; I, inferior; P, posterior; A, anterior; L, left; R, right; Thal, thalamus; Puta, putamen; Caud, caudate; Pall, globus pallidus; Accu, nucleus accumbens; Puta, putamen; Hipp, hippocampus

Compared with HCs, patients with DPN (both groups) showed lower FA and higher MD and RD throughout WM tracts, primarily including the internal and external capsules, corpus callosum, corticospinal tract, spinothalamic tract, bilateral thalamic projecting fibers, bilateral fronto‐occipital fasciculus, bilateral inferior longitudinal fasciculus, and bilateral uncinate fasciculus (p < .05, TFCE corrected). Additionally, lower FA and higher RD in periaqueductal WM were only detected in patients with painful DPN (p < .05, TFCE corrected; Figure 3).

Figure 3

Between‐group differences in DTI parameters. Compared with HCs, lower FA (a) and higher MD (b) and RD (c) were observed throughout the WM tract in painless and painful DPN. Colored areas denote regions where a significant difference in DTI parameters was observed between the indicated groups. Red hues are assigned to the painless DPN > HCs and painful DPN > HCs clusters and blue hues are assigned to the painless DPN < HCs and painful DPN < HCs clusters. L, left; R, right

Correlations between neurological test scores and imaging parameters

We found that neurological parameter scores correlated with the shape of several deep nuclei. NDS scores correlated significantly with the shape of the left caudate, both negatively (lateral segment) and positively (medial posterior segment; both p < .05, TFCE corrected). They were also positively correlated with the shape of the medial posterior segment of the right putamen (p < .05, TFCE corrected; Figure 4a). NSS scores were negatively correlated with the shapes of the lateral segment of the left caudate and the anterior lateral segment of the left putamen (both p < .05, TFCE corrected). They also significantly correlated with the shape of the right putamen, both negatively (anterior lateral segment) and positively (medial posterior segment; p < .05, TFCE corrected) (Figure 4b). No correlations were observed between DN4 scores and the shapes of any deep nuclei.

Figure 4

Correlations of deep gray matter nuclei shape with NDS (a) and NSS (b) scores. Colored areas denote regions where a significant correlation was observed. Red hues denote positive correlation and blue hues denote negative correlation. L, left; R, right; Caud, caudate; Puta, putamen

Significant correlations were also found between neurological parameter scores and other imaging parameters. NDS (p = .002, TFCE corrected) and NSS (p = .018, TFCE corrected) scores negatively correlated with FA in the bilateral thalamus (Figure 5a,b), and NDS scores positively correlated with MD (p = .003, TFCE corrected) and RD (p = .002, TFCE corrected) in the bilateral thalamus (Figure 5c,d).

Figure 5

Correlations between DTI parameters and neurological test scores. Significant negative correlations were observed between FA in bilateral thalamus and NDS (a) and NSS (b) scores. Significant positive correlations were observed between NDS score and both MD (c) and RD (d) in bilateral thalamus. L, left; R, right

No significant correlations were observed between cortical thickness or cortical surface area, and NSS, NDS, or DN4 scores.

Exploratory intermodal correlation analyses

Using the pooled data from the three groups, we found a positive correlation between the mean FA of significant clusters identified by the TBSS analysis and the mean cortical thickness of significant clusters identified by the SBM analysis (r = .453, p < .001; Figure 6a). Meanwhile, the mean FA of significant clusters identified by the TBSS analysis was found to be negatively correlated with the mean cortical surface area of significant clusters identified by the SBM analysis (r = −.315, p < .001; Figure 6b).

Figure 6

Correlations between cortical morphological indices and FA. (a) A significant positive correlation was observed between the mean FA of significant clusters identified by the TBSS analysis and the mean cortical thickness of significant clusters identified by the SBM analysis. (b) A significant negative correlation was observed between the mean FA of significant clusters identified by the TBSS analysis and the mean cortical surface area of significant clusters identified by the SBM analysis

DISCUSSION

To the best of our knowledge, this is the first study to comprehensively assess both GM and WM alterations in Type 2 diabetic patients with peripheral neuropathy. Compared with HCs, patients with DPN exhibited significant differences in the morphology of the cerebral cortex and deep GM nuclei (e.g., pre‐ and postcentral gyrus, paracentral lobule, insula, ACC, and thalamus), as well as in the microstructure of the subcortical WM (e.g., corticospinal tract, spinothalamic tract, and thalamocortical projecting fibers). These findings suggest impaired somatosensory‐, motor‐, and pain‐related pathways in painless and painful DPN. Additionally, exploratory intermodal correlation analyses of pooled data (patients and HCs) showed a significant positive correlation between cortical thickness and FA, as well as a significant negative correlation between cortical surface area and FA, indicating that the changes in GM and WM were tightly coupled. These findings indicate that alterations in GM and WM might contribute to the pathogenesis of DPN.

Sensorimotor and pain related GM and WM alterations

In the present study, we observed that DPN was associated with GM alterations in the precentral gyrus, postcentral gyrus, paracentral lobule, and deep GM nuclei (i.e., caudate, putamen, medial pallidum, and thalamus). We also observed altered WM in the corticospinal tract, internal capsule, spinothalamic tract, and thalamocortical projecting fibers. Physiologically, the somatosensory cortex, together with the thalamus, spinothalamic tract, and thalamocortical projecting fibers, make up an important component of the somatosensory pathway that governs the central integration and modulation of various peripheral sensations (Li et al., 2018). The precentral gyrus, paracentral lobule, caudate/putamen (known as the neostriatum), medial pallidum, thalamus, and corticospinal tract/internal capsule constitute a neural pathway that functions in coordinating voluntary movement. As such, the structural differences that we observed in these brain regions suggest that motor and somatosensory systems are impaired in patients with DPN. The structural alterations were supported by the results of our supplementary independent component analysis using resting state fMRI data, which showed lower intra‐network connectivity in the bilateral pre−/postcentral gyrus within the sensorimotor network ( Figure S2). Although the alterations in somatosensory and motor systems could represent neuroplastic changes in patients with DPN, the neural mechanisms underlying the alterations might differ depending on the system. Specifically, in DPN changes in higher regions of the somatosensory system might result from noxious afferent stimuli coming from the injured peripheral sensory nerves (Feldman et al., 2017). Physiologically, the peripheral sensory nerves are susceptible to injury owing to insufficient protection of the blood‐nerve barrier in the sensory neurons (such as dorsal root ganglion neurons) and the lack of protection and nutrient supplementation by Schwann cells in unmyelinated axons (known as C fibers; Feldman et al., 2017). C fibers carry information about temperature and noxious stimuli, such as mechanical forces and injurious chemicals, to the CNS (Feldman, Hughes, & Willison, 2015). In patients with DPN, they have been shown to be impaired due to exposure to systemic metabolic and hypoxic stressors (Feldman et al., 2017). Thus, peripheral noxious stimuli might induce axonal degeneration and aberrant activity (ectopic discharges), which in turn could result in the structural alterations we observed in brain regions connected to the somatosensory pathway. In contrast, the structural alterations in the motor system may reflect a protective strategy in patients with chronic pain to prevent or decrease pain sensation (Lund, Donga, Widmer, & Stohler, 1991). As a type of chronic pain, diabetic neuropathic pain may give rise to abnormalities in muscular activity. This could lead to restricted movement of an affected muscle by redistributing load and function, which in turn would protect the peripheral nervous system from further injury and support rehabilitation.

Our finding that structural alterations in the cingulate cortex, insular cortex, prefrontal lobe, and thalamus, as well as the impaired integrity of periaqueductal WM and the external capsule primarily occurred in patients with painful DPN suggest focal abnormalities of pain perception and modulation of pain‐related pathways (Davis & Moayedi, 2013; Ossipov, Morimura, & Porreca, 2014). This finding was partially supported by the results of our supplementary independent component analysis, which showed lower intra‐network connectivity in a number of brain regions within the salience network (a network frequently involved in processing the salient attention‐demanding stimuli or experiences such as pain (van Ettinger‐Veenstra et al., 2019; Hemington, Wu, Kucyi, Inman, & Davis, 2016)), including bilateral insula and the right ACC (Figure S3). Specifically, the ACC has been reported to be a pain‐intensity encoding region in both healthy volunteers and in animal models (Apkarian, Baliki, & Geha, 2009; Duerden & Albanese, 2013; Farrell, Laird, & Egan, 2005; Senapati et al., 2005). In addition, recent metabolic studies indicated that chronic knee osteoarthritis pain and fibromyalgia were associated with reduced levels of the inhibitory neurotransmitter GABA in the ACC, which was thought to reflect a disinhibition mechanism of this region in pain perception (Foerster et al., 2015; Reckziegel, Raschke, Cottam, & Auer, 2016). Further, the PFC has been shown to exert inhibitory control on the medial thalamus‐ACC loop (Cauda et al., 2009). Therefore, we speculate that the observed alterations in ACC and PFC may imply a dysregulation of the PFC‐ACC‐medial thalamus loop in DPN patients, leading to an increased level of perceived unpleasantness associated with allodynia. In contrast, the morphological alterations of the insula, together with the impaired WM integrity of the adjacent external capsule, have been implicated in both the sensory–discriminative and the affective–motivational aspects of pain processing, with the posterior region involved in augmented pain perception and the anterior region involved in pain‐related affective processing (Harris et al., 2009; Watson, 2016). Additionally, the observed WM microstructural changes around the PAG, a vital region involved in the inhibition of nociceptive inputs through its interactions with ascending and descending projections from numerous sites (Lewis & Gebhart, 1977; Ossipov et al., 2014), suggest that DPN might be associated with impaired PAG‐mediated analgesia, which would lead to exaggerated pain sensation. In fact, the insula, cingulate cortex, PFC, and PAG are linked and activated when processing painful stimuli, and are thus regarded as key regions of a network termed the pain matrix (Davis & Moayedi, 2013; Ohara, Crone, Weiss, Kim, & Lenz, 2008; Wilkinson et al., 2013). It is conceivable that the observed GM or WM alterations in these regions underlie the impairments in multidimensional aspects of pain processing that are found in DPN, including abnormal encoding of pain intensity and anticipation, abnormal pain modulation, and the experience of negative emotions, such as anxiety.

Potential neural mechanisms underlying the GM and WM alterations in DPN

In the present study, two genetically and phenotypically independent morphological indices, namely cortical thickness and cortical surface area, were adopted to assess the alterations in cortical structure (Hogstrom, Westlye, Walhovd, & Fjell, 2013; Panizzon et al., 2009; Rakic, 1988; Rimol et al., 2012; Zhang et al., 2011; Zhang et al., 2018). Because each index measures a unique aspect of cortical morphology, analyzing them separately might help uncover the neural mechanisms involved in the development of DPN. According to the radial unit hypothesis, cortical thickness is mediated by the horizontal layers in cortical columns, whereas cortical surface area is primarily determined by the number of cortical columns perpendicular to the cortical surface (Ecker et al., 2013; Jiang et al., 2016; Rakic, 1988). Thus, the alterations in cortical thickness and cortical surface area that we observed might respectively be driven by cellular events perpendicular (such as demyelination and accelerated cell death) and tangential (including gliogenesis, synaptogenesis, and intracortical myelination) to the cortical surface (Hill et al., 2010; Hogstrom et al., 2013; Jubault et al., 2011; Toth et al., 2017; Van Essen, 1997).

Our finding of concurrent contraction and expansion within or in separate deep GM nuclei of patients with DPN, together with the positive and negative correlations observed between their shapes and the neurological test scores, indicated that different deep GM nuclei or sub‐regions within a specific nucleus exhibited divergent patterns of neuroplastic changes in response to long‐term peripheral noxious stimuli. This might occur as a result of functional differences among the different deep GM nuclei or within sub‐regions of a specific nucleus. However, the exact mechanism underlying this phenomenon remains unknown, and future parcellation‐based studies are needed to further pursue this issue.

Using TBSS, we identified a pattern of lower FA and higher MD and RD throughout the altered WM tracts. The lower FA and higher MD might indicate abnormal myelin structure, although there are indications that these measures alone are insufficient (Hasan, 2006; Wu et al., 2017). In addition, experimental evidence indicates that an increase in RD might reflect myelin damage (Wheeler‐Kingshott & Cercignani, 2009). Although not significant in this study, AD measures diffusivity along the main fiber direction and is thought to be a marker of axonal injury (Metwalli et al., 2010). Pathologically, several previous studies have demonstrated that DPN is characterized by distal axonal loss, a reduction in myelinated fiber density, and focal areas of demyelination (Selvarajah, Wilkinson, Davies, et al., 2011; Selvarajah, Wilkinson, Gandhi, et al., 2011). Finding this pattern of changes in the DTI parameters, together with correlation in the thalamus between DTI parameters and neurological parameters, suggests the destruction of neuronal fiber/myelin structural integrity, such as demyelination, and continuous WM integrity impairment as the disease progresses.

Relationship between cortical morphological indices and FA

We observed significant positive correlations between cortical thickness and FA, suggesting that the GM and WM alterations were tightly coupled. Considering that WM tracts play a key role in metabolite provision, signal transduction, and neurotrophic support to GM (Han et al., 2017; Waller et al., 2016), the alterations in WM microstructure—evidenced by lower FA in the patients—likely resulted in substantial neuronal dystrophy and apoptosis, which manifested macroscopically as cortical thinning. In addition, we observed a significant negative correlation between cortical surface area and FA. This is consistent with a pervious multimodal MRI study on patients with medication‐naïve chronic schizophrenia (Liu et al., 2014), and is supported by a tension‐based morphogenesis hypothesis proposing that the tension of the underlying WM tracts plays a crucial role in shaping the morphology of the cortical surface (Van Essen, 1997). In fact, in patients with DPN, abnormalities in the cerebral cortex have been thought to represent maladaptive neuroplastic changes secondary to the spino‐thalamic tract and thalamo‐cortical projection alterations caused by the peripheral neuropathy (Feldman et al., 2017; Iadarola et al., 1995; Tan et al., 2012). As such, the cortical morphological alterations in patients with DPN seem to be closely related to the WM microstructural abnormalities that occur as a result of aberrant sensory input from the injured peripheral nerve of the extremities in this disease. However, the exact causal mechanisms underlying the relationships between the GM and WM alterations remain unclear and need further investigation.

Limitations

Our study is novel in that we were able to assess several pathways that might be affected in patients with DPN by combining GM and WM analysis. Our results point to a potential multifactorial mechanism for DPN, which may involve both GM and WM of the brain. However, the study does have limitations. First, our sample size was modest, especially for the subgroups of patients with DPN, which may obscure some subtle neuroimaging findings in the subgroups. Nevertheless, our study may allow hypothesis generation for future large‐scale multicenter studies to identify the specific changes in the brains of patients with DPN. Second, this was a single center study, which may limit its generalization to diabetic patients outside our institution. Third, we were limited by the cross‐sectional design of this study, which could not determine the causes behind the findings. We were also limited by lack of longitudinal follow‐up data for these patients. Fourth, spinal MRI scans were not obtained, which would have allowed us to assess the alterations in the spinal cords of the patients. However, this was beyond the scope of the study. Future research combining spinal and brain MRI data would be necessary to interrogate the entire neural axis and fully understand the pathogenesis of DPN. In addition, the absence of detailed quantitative assessment of pain intensity weakens the interpretability of our findings.

CONCLUSIONS

In the present study, we revealed significant GM and WM alterations in some key brain regions of the ascending spinal‐cortical somatosensory pathway, the descending motor pathway, and the pain perception and modulation pathway in patients with DPN. Our study provided further evidence that CNS abnormalities contribute to the pathogenesis of DPN.

CONFLICT OF INTEREST

The authors report no competing interest.

Data S1: Figure S1. Pair‐wise intergroup comparison of deep‐gray‐matter nuclei volume, and correlations between deep‐gray‐matter nuclei volume and the neurological parameters. Between‐group differences in deep‐gray‐matter nuclei volume (a). Colored bars denote regions where a significant difference in deep‐gray‐matter nucleus volume was observed between the indicated groups. Correlations between deep‐gray‐matter nuclei volume and NDS, DN4 and NSS, respectively (b). L left; R right; Accu nucleus accumbens; Amyg amygdala; Hipp hippocampus; Caud caudate; Puta putamen; Pall globus pallidus; Thal thalamus.

Figure S2. Pair‐wise intergroup comparison in component of sensorimotor network. (a) Differences in brain regions within sensorimotor network for painless DPN patients versus HCs. (b) Differences in brain regions within sensorimotor network for painful DPN patients versus HCs. Colored areas denote regions where a significant difference in brain regions within sensorimotor network was observed between the indicated groups. Blueness is assigned to the painless DPN < HCs or painful DPN < HCs clusters.

Figure S3. Pair‐wise intergroup comparison in component of salience network. (a) Differences in brain regions within salience network for painless DPN patients versus HCs. (b) Differences in brain regions within salience network for painful DPN patients versus HCs. Colored areas denote regions where a significant difference in brain regions within salience network was observed between the indicated groups. Blueness is assigned to the painless DPN < HCs or painful DPN < HCs clusters. P, posterior; A, anterior; L, left; R, right.

Click here for additional data file.

Appendix S1: Supporting Information

Click here for additional data file.

Appendix S2: Supporting Information

Click here for additional data file.