Gene Network Study Reveals How Treadmill Training May Aid Spinal Cord Recovery

Treadmill training after spinal cord injury appears to act through a specific set of molecular switches, according to a new study that combined computational gene network analysis with laboratory validation in rats. Researchers in Shenyang, China, report that five molecules—AKT1, JAK2, TNF, STAT1, and the microRNA miR-135a—change their expression in the soleus muscle of injured rats that undergo treadmill training, and that miR-135a can directly bind the messenger RNA of both JAK2 and STAT1. The findings, published in Frontiers in Neurology, do not demonstrate a cure, but they offer a detailed molecular map of how a widely used rehabilitation therapy leaves a measurable signature in muscle tissue after spinal cord injury.

Spinal cord injury, which can arise from traumatic events such as falls, vehicle collisions, and sports accidents, or from nontraumatic causes including infection, tumors, and reduced blood supply, disrupts muscle function, sensation, and autonomic control. Its global incidence has been estimated to range from 236 to 4,187 cases per million people, and survivors face complications ranging from pressure sores and muscle atrophy to respiratory problems and infections. Exercise has long been known to improve strength and physical capacity in people with chronic spinal cord injury, and treadmill training after incomplete injury can enhance lower-extremity motor function and gait recovery. In animal models, appropriately dosed treadmill training promotes functional recovery, neuromuscular regeneration, and skeletal muscle plasticity. Earlier work has linked exercise after injury to changes in microRNAs that regulate the PTEN and mTOR signaling pathways, activation of the ERK1/2 pathway that supports corticospinal tract axon outgrowth, and increased expression of brain-derived neurotrophic factor, which helps suppress spasticity and allodynia by boosting KCC2 expression. What has remained unclear is the fuller gene-level architecture underlying these benefits.

To close that gap, the team mined a publicly available microarray dataset, GSE45550, originally generated by Baligand and colleagues, in which moderate contusion injuries were produced in rats and quadrupedal locomotor training began on postoperative day 7. The dataset covered soleus muscle samples collected at 8 and 14 days after injury, both with and without treadmill training, six animals per condition, profiled on the Affymetrix Rat Genome 230 2.0 Array. Training consisted of 20-minute treadmill stepping sessions—performed three times by day 8 and twice daily for five days by day 14—with manual body-weight support and hindpaw placement assistance withdrawn as locomotor capability improved. After normalizing the raw data with the robust multi-array average method, the researchers used the limma package to identify differentially expressed genes between trained and untrained injured animals at each time point. This yielded 640 differentially expressed genes at day 8, of which 428 were upregulated and 212 downregulated, and 877 genes at day 14, split almost evenly between 439 upregulated and 438 downregulated. Because multiple-testing corrections left too few genes for network analysis, the team screened with an unadjusted p-value below 0.05 and no fold-change threshold, a compromise they later flagged as a limitation.

The merged gene list was then fed into weighted gene co-expression network analysis, a systems biology method that groups genes into modules based on correlated expression patterns across samples. After selecting a soft-thresholding power of 16—the first value at which the squared correlation coefficient reached 0.9, indicating a scale-free network topology—the algorithm partitioned the differentially expressed genes into 30 modules. Seven modules correlated strongly with the treadmill-training phenotype at correlation coefficients above 0.8, and the top four, labeled MElightcyan1, MEbisque4, MEdarkseagreen4, and MEgreen, were carried forward. Functional annotation revealed distinct biological identities: MElightcyan1 was tied to the regulation of synaptic transmission and the glutamatergic pathway, MEbisque4 to the negative regulation of actin filament polymerization, MEdarkseagreen4 to calcium ion transmembrane transport, and MEgreen to the negative regulation of DNA-templated transcription.

When the nodes of all four modules were merged, they comprised 492 genes—308 upregulated and 190 downregulated. Gene ontology analysis showed that the upregulated genes were concentrated in ion transmembrane transport and the glutamatergic synapse pathway, while the downregulated genes clustered around neural tube closure and ribosome biogenesis in eukaryotes. Pre-ranked gene set enrichment analysis, run independently at both time points, added a temporal nuance: the glutamatergic synapse gene set was significantly enriched at day 8, with a normalized enrichment score of 1.67 and an adjusted p-value of 0.0108, but neither glutamatergic synapse at day 14 nor ribosome biogenesis at either time point reached significance. This suggests the glutamatergic signal is an early response to training, fading by two weeks after injury, while the ribosome biogenesis association rests on the enrichment analysis of the merged gene set rather than on significant enrichment at either individual time point.

To find the molecules most likely to matter, the researchers built a protein–protein interaction network from the 492 genes using the STRING database at medium confidence, producing a graph of 249 nodes and 476 edges. AKT1 emerged as the most connected hub with a degree of 35, followed by JAK2 with 18, and TNF and STAT1, each with 14. Two dense clusters, module a with a score of 7.714 and module b with a score of 5.419, each contained eight nodes, and JAK2 and STAT1 were found to interact directly within module b. Cross-referencing these hubs with spinal cord injury–related genes curated in the Comparative Toxicogenomics Database and GeneCards confirmed all four as overlapping candidates. Expression profiling indicated that in the 14-day trained group, AKT1 and STAT1 were significantly increased while JAK2 and TNF were significantly decreased compared with untrained injured animals.

The team then mapped the regulatory layer above these genes. Overrepresentation enrichment analysis predicted 14 transcription factors and 7 microRNAs acting on the network genes, yielding 310 transcription factor–microRNA–target pairs. Within this regulatory network, both STAT1 and JAK2 were predicted targets of miR-135a. The researchers then validated the five key RNAs experimentally in a new cohort of 20 adult Sprague Dawley rats, in which spinal cord contusion was produced at the ninth thoracic vertebra by clamping the exposed cord for one minute. After one modeling failure, the injured animals were divided into four groups of three—8 and 14 days after injury, each with and without treadmill training at 10 meters per minute for 20 minutes per session—alongside uninjured controls with and without training. Reverse transcription quantitative polymerase chain reaction on soleus muscle confirmed the computational predictions: AKT1 and STAT1 were significantly upregulated in both trained groups, while JAK2, TNF, and miR-135a were significantly downregulated relative to their untrained injured counterparts.

To test whether the predicted microRNA–gene interactions were real rather than statistical artifacts, the team turned to dual-luciferase reporter assays in human embryonic kidney 293T cells. They cloned the wild-type 3-prime untranslated regions of human JAK2 and STAT1, each containing the predicted miR-135a-5p binding sequence, into reporter vectors, alongside mutant versions with substitutions in the seed-binding sites. Co-transfection with a miR-135a-5p mimic significantly reduced luciferase activity from the wild-type JAK2 and STAT1 reporters compared with negative controls, while the mutant reporters showed no significant change. This pattern supports direct binding of miR-135a-5p to the predicted sites within both 3-prime untranslated regions, providing a plausible mechanistic link between the training-induced drop in miR-135a and the divergent expression changes of its two targets.

The biological logic of these targets fits established themes in spinal cord injury research, though the authors frame their own contribution as connecting them to treadmill training specifically. AKT signaling is known to support neuronal survival and muscle protein synthesis, and previous studies have reported that treadmill exercise upregulates neurotrophic factors through the PI3K/Akt pathway and that moderate-intensity treadmill exercise activates cortical Akt/mTOR signaling in ways associated with neuroplasticity. TNF-alpha, by contrast, induces apoptosis in oligodendrocytes and neurons and drives inflammatory pain and synaptic plasticity, and the TNF antagonist etanercept has reduced tissue damage and promoted myelin regeneration in preclinical rat models. STAT1 inhibition has been reported to reduce histological damage and motor dysfunction by modulating apoptosis and inflammation. The ribosome biogenesis connection is also grounded in muscle physiology: because ribosomal capacity governs protein synthesis, impaired ribosomal activity is thought to contribute to the skeletal muscle atrophy that follows spinal cord injury. The glutamatergic synapse finding echoes earlier work showing that surgically rerouted descending spinal fibers can innervate abdominal muscle and drive newly formed neuromuscular junctions to switch from cholinergic to glutamatergic signaling.

The authors are explicit that these results are preliminary in important ways. The unadjusted p-value threshold raises the risk of false positives, and because both the discovery dataset and the validation samples came from soleus muscle, the observed signals may reflect peripheral muscle adaptation rather than molecular repair within the spinal cord itself. Validation was limited to RNA expression and reporter assays, leaving protein abundance, phosphorylation, cell-specific localization, and the functional consequences of manipulating AKT1, JAK2, TNF, STAT1, and miR-135a-5p untested. Only days 8 and 14 and male rats were examined, individual behavioral scores were not correlated with the molecular measurements, and no proteomic, phosphoproteomic, metabolic, or single-cell data were integrated. Within those boundaries, the study identifies AKT1, JAK2, TNF, STAT1, and miR-135a as candidate molecular factors associated with the response to treadmill training after spinal cord injury, and points to the glutamatergic synapse and ribosome biogenesis pathways as processes the training may engage. Whether these candidates can be exploited therapeutically—through pharmacological or genetic manipulation—remains a question for larger, sex-balanced, longitudinal, multi-omic studies with functional experiments.

Subject of Research: Molecular targets of treadmill training in rat spinal cord injury identified by weighted gene co-expression network analysis

Article Title: Application of weighted gene co-expression network analysis to explore the targets of treadmill training in rats with spinal cord injury

Article References: Zhang, H., Guo, M., Meng, L., Piao, M., & Yu, H. (2026). Application of weighted gene co-expression network analysis to explore the targets of treadmill training in rats with spinal cord injury. Frontiers in Neurology, 17, Article 1846393. https://doi.org/10.3389/fneur.2026.1846393

Image Credits: AI Generated

DOI: 10.3389/fneur.2026.1846393

Keywords: spinal cord injury, treadmill training, weighted gene co-expression network analysis, AKT1, JAK2, STAT1, TNF, miR-135a, soleus muscle, glutamatergic synapse, ribosome biogenesis, neurorehabilitation