Electrical stimulation mitigates muscle degradation shift in gene expressions during 12-h mechanical ventilation – Scientific Reports

Overall characterization of changes in ABG parameters and gene expression during V12, S-V12, and E-V12

The ABG parameter values at the start and end of 12-h MV during V12, S-V12, and E-V12 are shown in Fig. 1. Although significant differences were observed in the Kruskal‒Wallis test for pH at the start, no significant differences were observed for any ABG parameters at either the start or end of 12-h MV in pairwise comparisons among V12, S-V12, and E-V12 using the Mann‒Whitney U test. Comparisons of ABG parameter values between the start and end of 12-h MV also showed no significant differences in the V12, S-V12, and E-V12 groups.

As shown in Fig. 2, the CON, V12, S-V12, and E-V12 groups were significantly clustered in the 3-dimensional PCA plot of the microarray data (Pillai’s trace = 2.4838, approx F = 19.247, df = 3, 36, p = 3.001e-11). There were 1921 differentially expressed genes (845 upregulated, 1076 downregulated) in the V12 group (V12/CON), 1556 differentially expressed genes (925 upregulated, 631 downregulated) in the S-V12 group (S-V12/V12), and 491 differentially expressed genes (236 upregulated, 255 downregulated) in the ES group (E-V12/S-V12), as shown in the Venn diagrams (Fig. 3a and b).

PCA plot of microarray expression data. The figure represents the three principal components (PC1, PC2 and PC3) of the microarray analysis data. Each circle indicates the normalized gene expression of one microarray sample, and the color of the circle indicates the treatment condition. CON, V12, S-V12, and E-V12 samples (n = 4) were significantly clustered based on their respective groups by MANOVA.

Venn diagram of genes altered by ES among the differentially expressed genes in the diaphragm in the VIDD rat model. (a) Out of the 845 upregulated genes in the V12/CON comparison, 41 genes that overlapped with the 255 downregulated genes in E-V12/S-V12 were downregulated by ES. Excluding the overlap with 1556 differentially expressed genes affected by the laparotomy procedure of ES resulted in 18 genes. (b) Out of the 1076 downregulated genes in V12/CON, 8 genes that overlapped with the 236 upregulated genes in E-V12/S-V12 were upregulated by ES. Excluding the overlap with 1556 differentially expressed genes affected by the laparotomy procedure of ES resulted in 6 genes.

Functional annotation and enrichment analysis

We utilized Metascape to conduct annotation and enrichment analyses on multiple gene sets derived from the diaphragm in the VIDD rat model. As shown in the heatmap representations (Fig. 4a), we observed significant functional overlap, with 9 out of the top 20 enriched terms corresponding to the same biological processes present in the 1921 genes (845 upregulated and 1076 downregulated genes) in the V12/CON comparison. Interestingly, for the downregulated genes (1076 genes), we identified five selectively enriched pathways, including “muscle system process” and “striated muscle cell differentiation”, both of which are associated with muscle function. In contrast, the upregulated genes (845 genes) did not exhibit significant enrichment in muscle-related GO terms. Instead, we identified seven selectively enriched clusters and pathways, including the “forkhead box protein O1 (FoxO1) signaling pathway”, “positive regulation of cytokine production”, and the “apoptotic signaling pathway”. Table 1a shows a list of genes that were duplicated among the significantly enriched GO terms in the V12/CON comparison and the fold changes corresponding to the list of overlapping genes in the microarray.

Functional enrichment analysis for genes altered by ES among differentially expressed genes in the diaphragm in the VIDD rat model. Pathway and process enrichment analyses were performed using Metascape (http://metascape.org, accessed on 9 June 2023) for (a) the 1921 genes (845 upregulated and 1076 downregulated genes) in V12/CON, (b) the 491 genes (236 upregulated and 255 downregulated genes) in E-V12/S-V12, (c) the 1556 genes (925 upregulated and 631 downregulated genes) in S-V12/V12, and (d) the 18 and 6 genes as illustrated in the Venn diagrams in Fig. 3a and b. The heatmaps display the results of enrichment analysis of all statistically enriched ontology terms (GO biological processes, KEGG pathways, Reactome Gene Sets, WikiPathways), with the color representing the − log10(p value). The thresholds used were a 0.3 kappa score and a similarity score greater than 0.3 (a–c). Darker shades of orange indicate greater significance for the term, while a gray color indicates a lack of significance. The dendrograms show the closeness of the enriched term clusters, with one per row. The grouped colors represent the clustering of enriched GO terms, and the heatmap visualizes the GO results of Metascape using heatmapper (http://www.heatmapper.ca/expression/) (d).

In the E-V12/S-V12 comparison, selectively enriched clusters were observed for both the 236 upregulated and 255 downregulated genes (Fig. 4b). For the 236 upregulated genes, biological processes related to oxygen and mitochondria were enriched, including “protein import into mitochondrial matrix”, “oxygen transport”, and “mitochondrial electron transport (ubiquinol to cytochrome c)”. For the 255 downregulated genes, selective enrichment was found in processes related to muscles and hormones, such as “regulation of muscle system process”, “cellular response to hormone stimulus”, “regulation of intracellular steroid hormone receptor signaling pathway”, and “insulin signaling”. Table 1b shows a list of genes that were duplicated among the significantly enriched GO terms in the E-V12/S-V12 comparison and the fold changes corresponding to the list of overlapping genes in the microarray.

In the S-V12/V12 comparison, a total of 1556 genes (925 upregulated and 631 downregulated genes) were analyzed for enrichment (Fig. 4c). For the upregulated genes and downregulated genes, there were 9 and 8 selectively enriched clusters, respectively. The 925 upregulated genes showed enrichment in clusters related to stress, tissue damage, and apoptosis, such as “cellular responses to stress”, “regulation of epithelial cell apoptotic process”, “response to wounding”, and “cellular response to external stimulus”. In contrast, the 631 downregulated genes were selectively enriched in clusters related to mitochondrial function and oxidation, including “mitochondrial cytochrome c oxidase assembly” and “mitochondrial translation elongation.” Table 1c shows a list of genes that were duplicated among significantly enriched GO terms in the S-V12/V12 comparison and the fold changes corresponding to the list of overlapping genes in the microarray.

Gene changes related to inflammatory cytokines, stress, and skeletal muscle occurred after 12-h MV, sham operation, and ES

We investigated genes associated with inflammatory cytokines, stress, and skeletal muscle according to the results of enrichment analyses. Table 2 summarizes the genes that were selected from selectively enriched clusters in the 12-h MV (V12/CON), sham operation (S-V12/V12), and ES (E-V12/S-V12) conditions, as shown in Table 1a to c. After 12 h of MV management, the expression of FoxO1 was increased by 37.8-fold, and the expression of Ppargc1a (the gene symbol of PGC1-α) was decreased by 37.1-fold. Genes associated with muscle atrophy, such as myostatin (Mstn), tripartite motif-containing 63 (Trim63), and f-box protein 32 (Fbxo32), and genes associated with inflammatory cytokines and stress, such as mitogen-activated protein kinase 14 (Mapk14), sirtuin 1 (Sirt1), and CCAAT/enhancer binding protein beta (Cebpb), were upregulated, whereas genes related to muscle synthesis, regeneration, and contraction, such as tripartite motif-containing 72 (Trim72), calsequestrin 2 (Casq2), and myogenic differentiation 1 (Myod1), were downregulated. The expression of FoxO1 and PGC1-α, which was particularly changed after 12 h of MV, was altered by no more than twofold by additional sham operation or ES (FoxO1 expression increased by 1.9-fold in the sham operation group and decreased by 1.9-fold in the ES group; PGC1-α expression decreased by 1.1-fold in the sham operation group and by 1.6-fold in the ES group). Sham operation resulted in a 7.8-fold increase in the expression of an energy metabolism-related gene, solute carrier family 2 member 1 (Slc2a1), and twofold increases in the expression of stress response genes, such as forkhead box O3 (FoxO3), Mapk14, and nuclear factor erythroid 2-related factor 2 (Nfe212). Of the 4 genes that showed twofold or greater decreases in expression under ES, the gene that showed opposite changes in expression under 12-h MV was peroxisome proliferator-activated receptor delta (Ppard). Ppard expression was increased by 1.7-fold by 12-h MV, increased by 4.7-fold by the sham operation, and decreased by 2.2-fold by ES.

Genes for which ES counteracted expression changes during 12-h MV

For the ES group, laparotomy was required to insert electrodes directly into the diaphragm to ensure contraction of the diaphragm muscle. Therefore, the influence of laparotomy was accounted for to identify genes for which ES counteracted expression changes during 12-h MV.

Of the 845 genes whose expression increased by more than twofold in the 12-h MV group and the 255 genes whose expression decreased by more than twofold in the ES group, 41 genes were duplicated. The 23 genes that were upregulated or downregulated by more than twofold in the sham operation group were removed, and the remaining 18 genes were extracted as those whose expression was increased by 12-h MV but partially restored by ES (Fig. 3a). Of the 1076 genes that showed a twofold or greater decrease in expression after 12-h MV and the 236 genes that showed a twofold or greater increase in expression after ES, 8 genes were duplicates. The genes whose expression decreased by twofold or greater after the sham operation were removed, after which the number of genes whose downregulation due to 12-h MV was at least partially counteracted by ES was 6 (Fig. 3b). To summarize the Venn diagram analysis, a total of 24 of the genes that were altered by 12-h MV (18 genes upregulated by 12-h MV but downregulated by ES; 6 genes downregulated by 12-h MV but upregulated by ES) were not affected by the sham operation but were affected by ES.

For further analysis, we focused on 1100 genes (including 42 duplicates), consisting of 845 upregulated genes in V12/CON and 255 downregulated genes in E-V12/S-V12 (Fig. 3a), and 1312 genes (including 8 duplicates), consisting of 1076 downregulated genes in V12/CON and 236 upregulated genes in E-V12/S-V12 (Fig. 3b). We analyzed the commonly enriched pathways for the 1100 genes and the 1312 genes. Consequently, 18 and 6 overlapping genes whose expression changes were counteracted by ES in the VIDD rat model were found in the five identified clustering GOs (Fig. 4d). Details of the 24 genes (overlap of 18 and 6 genes) are summarized in Table 3a. For the clustered GOs, the classifications are summarized in Table 3b, with colors representing the clustering of enriched GO terms.

Among the genes included in the selectively enriched GO terms, the gene associated with stress was periodic circadian clock 1 (Per1), those associated with cytokines were ciliary neurotrophic factor receptor (Cntfr) and Per1, and those associated with skeletal muscle were Cntfr and patatin-like phospholipase domain-containing 7 (Pnpla7). Other changes in gene expression were observed for genes associated with the apoptotic signaling pathway, catabolic process, hormones, and muscle. The six genes that were downregulated by 12-h MV and upregulated by ES were related to internal respiration (mitochondrial respiration and oxygen transport). The gene related to mitochondrial respiration was NADH:ubiquinone oxidoreductase subunit A2 (Ndufa2), and the genes related to oxygen transport were hemoglobin genes.