Skeletal Muscle Immunometabolism in Women With Polycystic Ovary Syndrome: A Meta-Analysis

Polycystic ovary syndrome (PCOS) is an endocrine and metabolic disorder affecting up to 15% of women at reproductive age. The main features of PCOS are hyperandrogenism and irregular menstrual cycles together with metabolic dysfunctions including hyperinsulinemia and insulin resistance and a 4-fold increased risk of developing type 2 diabetes. Despite the high prevalence the pathophysiology of the syndrome is unclear. Insulin resistance in women with PCOS likely affect the skeletal muscle and recently it was demonstrated that changes in DNA methylation affects the gene expression in skeletal muscle that in part can explain their metabolic abnormalities. The objective of this work was to combine gene expression array data from different datasets to improve statistical power and thereby identify novel biomarkers that can be further explored. In this narrative review, we performed a meta-analysis of skeletal muscle arrays available from Gene Expression Omnibus and from publications. The eligibility criteria were published articles in English, and baseline (no treatment) skeletal muscle samples from women with PCOS and controls. The R package Metafor was used for integration of the datasets. One hundred and fourteen unique transcripts were differentially expressed in skeletal muscle from women with PCOS vs. controls (q < 0.05), 87% of these transcripts have not been previously identified as altered in PCOS muscle. ING2, CDKAL1, and AKTIP had the largest differential increase in expression, and TSHZ2, FKBP2, and OCEL1 had the largest decrease in expression. Two genes, IRX3 and CDKAL1 were consistently upregulated (q < 0.05) in the individual analyses and meta-analysis. Based on the meta-analysis, we identified several dysregulated immunometabolic pathways as a part of the molecular mechanisms of insulin resistance in the skeletal muscle of women with PCOS. The transcriptomic data need to be verified by functional analyses as well as proteomics to advance our understanding of PCOS specific insulin resistance in skeletal muscle.

Polycystic ovary syndrome (PCOS) is an endocrine and metabolic disorder affecting up to 15% of women at reproductive age. The main features of PCOS are hyperandrogenism and irregular menstrual cycles together with metabolic dysfunctions including hyperinsulinemia and insulin resistance and a 4-fold increased risk of developing type 2 diabetes. Despite the high prevalence the pathophysiology of the syndrome is unclear. Insulin resistance in women with PCOS likely affect the skeletal muscle and recently it was demonstrated that changes in DNA methylation affects the gene expression in skeletal muscle that in part can explain their metabolic abnormalities. The objective of this work was to combine gene expression array data from different datasets to improve statistical power and thereby identify novel biomarkers that can be further explored. In this narrative review, we performed a meta-analysis of skeletal muscle arrays available from Gene Expression Omnibus and from publications. The eligibility criteria were published articles in English, and baseline (no treatment) skeletal muscle samples from women with PCOS and controls. The R package Metafor was used for integration of the datasets. One hundred and fourteen unique transcripts were differentially expressed in skeletal muscle from women with PCOS vs. controls (q < 0.05), 87% of these transcripts have not been previously identified as altered in PCOS muscle. ING2, CDKAL1, and AKTIP had the largest differential increase in expression, and TSHZ2, FKBP2, and OCEL1 had the largest decrease in expression. Two genes, IRX3 and CDKAL1 were consistently upregulated (q < 0.05) in the individual analyses and meta-analysis. Based on the meta-analysis, we identified several dysregulated immunometabolic pathways as a part of the molecular mechanisms of insulin resistance in the skeletal muscle of women with PCOS. The transcriptomic data need to be verified by functional analyses as well as proteomics to advance our understanding of PCOS specific insulin resistance in skeletal muscle.
Keywords: PCOS, transcriptomics, gene expression, immunometabolism, meta-analysis, skeletal muscle INTRODUCTION Polycystic ovary syndrome (PCOS) is a common endocrine disorder among women of reproductive age, which is characterized by biochemical or clinical hyperandrogenism, irregular cyclicity and polycystic ovarian morphology (Norman et al., 2007;Teede et al., 2018). Metabolic dysfunction is evident in women with PCOS, and manifests as impaired glucose homeostasis, dyslipidemia and abdominal obesity (Azziz et al., 2016). Women with PCOS have epigenetic and transcriptional changes in skeletal muscle that, in part, can explain the metabolic abnormalities seen in these women (Skov et al., 2007;Nilsson et al., 2018). In addition, defects in early insulin signaling in combination with fibrosis in the skeletal muscle, likely contribute to insulin resistance in women with PCOS (Stepto et al., 2020).
There is an increasing body of publications supporting a role for both innate and adaptive immunity in response to changes in metabolic status. This new field of immunometabolism builds on evidence for activation of immune-derived signals in metabolically relevant tissues and explores how immune cells support tissues to adapt to environmental challenges (Man et al., 2017). Adipose tissue is one of the most explored tissues in the field, and the expression and activation of various immune cell types and anti-inflammatory cytokines have been studied in both lean and obese states (Ferrante, 2013;Man et al., 2017). Skeletal muscle is another key metabolically active organ considering that it is the most effective organ for insulin-stimulated glucose uptake in the body (∼80%) (Thiebaud et al., 1982). It also plays an important role in the development and progression of type 2 diabetes (Petersen and Shulman, 2018). However, there is scarce knowledge on how immune cells support the metabolic function in the skeletal muscle and on the role of inflammation in modulating skeletal muscle metabolism. The possible impact of immune cells on whole body glucose uptake and insulin sensitivity is a relatively recent appreciation.
Insulin resistance is a central feature of PCOS, with 30-40% of women with the syndrome having glucose intolerance, and 10% of them develop type 2 diabetes before the age of 40 (Legro et al., 1999;Rubin et al., 2017). Clinical studies exploring the molecular pathways of PCOS-insulin resistance in the skeletal muscle are limited. We have previously performed a wide-scale transcriptomic analysis in the skeletal muscle of women with PCOS and reported that many enriched pathways were involved in immune function or immune diseases (Nilsson et al., 2018). Due to the limited sample size of our study and the need to better understand the skeletal muscle metabolism in PCOS, we aim here to perform an integrated meta-analysis of the three available gene expression arrays in skeletal muscle of overweight/obese women with impaired glucose homeostasis. We further explore whether immunometabolism pathways are dysregulated in the skeletal muscle of these women.

Selection of Studies
We collected array studies including skeletal muscle from women with and without PCOS published between January 1999 and August 21, 2020, by performing a computerized search using PubMed, Omnibus GEO and Array Express. No review protocol exists but the following key words were used: human, skeletal muscle, polycystic ovary syndrome, PCOS, array, and gene expression. The eligibility criteria were; published articles, published in English, and baseline samples (no treatment). Data was collected from the GEO database (Skov et al., 2007(Skov et al., , 2008 and from our own publication (Nilsson et al., 2018). The analysis is reported according to the PRISMA checklist (Moher et al., 2009).

Meta-Analysis
The datasets were collected from two different microarray platforms, namely Affymetrix and Illumina. Affymetrix data was downloaded from NCBI GEO (GSE8157 and GSE6798), and normalized with methods described in Skov et al. (2007). The CEL-files and additional files were downloaded using the GEOquery R package. In order to check which samples were present in both GSE8157 and GSE6798, since some samples was used in both studies (Skov et al., 2007(Skov et al., , 2008 but with different names, we used the R package DupChecker to compare the CEL-files. The Illumina data was background corrected and normalized using NormExp method in the R package limma (Shi et al., 2010). A fixed-effect meta-analysis model, with inversevariance method, was used to estimate the possible associations between gene expression and PCOS. The Q statistic from the meta-analysis was used to identify the presence of heterogeneity. To account for multiple testing and control for false positives, the p-values from the meta-analysis were adjusted, using the Benjamini-Hochberg procedure. The robust genes, showing statistically significant association with PCOS, were identified as genes with false discovery rate (FDR) q < 0.05, and are presented as effect size. A heat map for all significant genes (q < 0.05) and a forest plot for the two most robust genes were produced. All analyses and plots were performed using the R statistical programming software [Version 4.0.2]. The metaanalysis was conducted using the rma() function in the Metafor package (Viechtbauer, 2010), the heatmap was produced using the heatmap.2() function in the gplot package and the FDR was estimated using the stats package.

Meta-Analysis of mRNA Expression in Skeletal Muscle From Women With PCOS and Controls
A total of 24 studies were identified, 12 were assessed for eligibility, and three studies were included in the meta-analysis as they fulfilled the eligibility criteria; published articles, published in English, and baseline samples (no treatment) (Skov et al., 2007(Skov et al., , 2008Nilsson et al., 2018). A PRISMA flow diagram is included in Supplementary Figure 1 (Prisma checklist in Supplementary Table 1). Only 3 samples in GSE8157 were exclusive for this study, while the others were also included in GSE6798. This generated a total number of 27 controls and 36 PCOS cases in the meta-analysis. Table 1 shows the most relevant clinical characteristics of the patients from each cohort. All women were over-weight or obese (mean BMI > 30 kg/m 2 ), and premenopausal. Women with PCOS were hyperandrogenic with metabolic aberration and insulin resistance. PCOS was diagnosed according to 2 criteria; irregular periods with cycle length > 35 days, and biochemical or clinical signs of hyperandrogenism in Skov et al. (2007Skov et al. ( , 2008, and according to Rotterdam criteria in Nilsson et al. (2018) fulfilling 2 of 3 criteria; ultrasound-verified polycystic ovaries, oligomenorrhea/amenorrhea, and/or clinical signs of hyperandrogenism. In total, of the 8,351 analyzed transcripts, 1,228 were differentially expressed in skeletal muscle from women with PCOS vs. controls (P < 0.05). After correction for multiple testing by FDR (5%, q < 0.05), 114 unique transcripts remained significant ( Figure 1A, Table 2). In total, 82% of the transcripts were upregulated and 18% were downregulated in women with PCOS. The three genes with the largest differential increase in expression were ING2, CDKAL1, and AKTIP, and the three genes with the largest differential decrease in expression were TSHZ2, FKBP2, and OCEL1. Ninety-nine of the 114 identified transcripts have not been previously shown to be altered in skeletal muscle arrays from women with PCOS (Skov et al., 2007(Skov et al., , 2008Nilsson et al., 2018) (marked in bold in Table 2).
In Nilsson et al. (2018), we have previously investigated the overlap between the differentially expressed genes identified in that study (p < 0.05) and those found in Skov et al. (2007) (q < 0.1), and found four transcripts (i.e., RAPH1, INO80D, PPIE, and AKTIP) that were upregulated in both studies. Two of these genes, RAPH1 and AKTIP were significant at q < 0.05 (effect size 1.16 and 1.17, respectively) and PPIE at q < 0.1 (effect size 0.86) in the integrated meta-analysis, while INO80D was not significantly altered. The meta-analysis findings were also compared with those obtained by individual analyses in both datasets (Supplementary Table 2) to evaluate bias and reproducibility across the microarray studies. As a result, 34 genes identified in the Illumina data set were also highlighted by our meta-analysis, whereas 20 genes were shared with the Affymetrix data set ( Figure 1B). Two of these genes, IRX3 (effect size 1.09) and CDKAL1 (effect size 1.27) were consistently upregulated (q < 0.05), by the three studies (Affymetrix, Illumina and integrated meta-analysis) ( Figure 1C).

Gene Set Enrichment Analysis and Pathway Analysis
Next, we tested whether sets of biologically-related genes were altered in women with PCOS compared with controls. According to the GSEA, there were significant (q <0.05) gene expression alterations in 9 downregulated pathways in women with PCOS ( Table 3). The majority of these downregulated pathways are involved in inflammatory/immune response while myogenesis and MYC signaling, a transcription factor for growth-related genes for glycolysis and mitochondrial metabolism in T cells (Wang et al., 2011), were upregulated (p < 0.05). Finally, GO analysis including the 114 differently expressed transcripts (q < 0.05) showed 205 enriched signaling pathways (Figure 2A,  Supplementary Table 3). Genes involved in metabolic and carbohydrate biosynthetic processes, positive regulation of I-kappaB kinase/NF-kappaB signaling, and negative regulation of cell proliferation and positive regulation of cell death appeared in the enriched pathways ( Figure 2B). Muscle system processes and calcium ion transport were also enriched together with the androgen receptor signaling pathway (Figure 2B,  Supplementary Table 3). All three genes in the androgen receptor signaling pathway; RAN, RB1, and CTNNB1, were FIGURE 1 | Comparison of individual analysis in dataset 1 (Skov et al., 2007(Skov et al., , 2008 and 2 (Nilsson et al., 2018). (A) Heat map of the 114 differentially expressed genes (q < 0.05), (B) Venn diagram of the individual analysis; Illumina, Affymetrix, and integrated meta-analysis (q < 0.05). (C) Forest plot of the effect size of the two genes, IRX3, and CDKAL1, consistently upregulated in the three studies.
significantly upregulated in women with PCOS compared to controls (q < 0.05, Table 2).

DISCUSSION
Taking an integrated meta-analysis approach, this study showed an overall downregulation of transcripts involved in the immune system and several alterations in metabolic pathways in skeletal muscle of women with PCOS. Previous data suggests that insulin signaling defects in skeletal muscle partly account for the insulin resistance in obese women with PCOS. However, our data indicate that additional molecular mechanisms are involved in the PCOS-insulin resistance in skeletal muscle. Fibrosis might be one key player (Stepto et al., 2020), together with immunometabolic alterations as supported by our data.
Fibrosis in skeletal muscle can impair muscle function and affect muscle fiber regeneration after injury (Delaney et al., 2017;Mahdy, 2019). Tissue fibrosis is often initiated and maintained through TGF-beta signaling and has been suggested to be associated with insulin resistance and steatosis in PCOS (Petta et al., 2017;Stepto et al., 2019). This association is supported by a recent paper, investigating TGF-beta ligand induced fibrosis in obese women with PCOS (Stepto et al., 2020). Moreover, we have previously shown an enrichment of extracellular matrix signaling pathways and a decreased collagen gene expression in our previous study (Nilsson et al., 2018), supporting a link between dysregulated TGF-beta signaling and extracellular matrix deposition with insulin resistance in PCOS. Surprisingly, fibrosis, extracellular matrix-, or TGF-beta signaling pathways were not enriched in women with PCOS in this meta-analysis. However, the most upregulated gene in PCOS muscle, inhibitor of growth family member 2 (ING2), promotes TGF-beta-induced transcription and cellular responses (Sarker et al., 2008), which leave the issue open for further investigation.
We have previously shown that the majority of the significant gene expression pathways (GSEA) in women with PCOS were involved in immune responses or were related to immune diseases (Nilsson et al., 2018). A distinct pattern was a downregulation of human leukocyte antigen (HLA) genes in women with PCOS. This finding was not validated in this metaanalysis, although two HLA genes were downregulated (p < 0.05), but they did not survive FDR < 5%. However, here we show a dysregulated expression of many genes involved in immune pathways, including a downregulation of gene sets associated with inflammatory response, interferon alpha and gamma response and the proinflammatory cytokines interleukin-2, interleukin-6 and TNF-alpha. Among the dysregulated genes, ADK and BST2 are thought to play a role in the immune system, identified as having anti-inflammatory and anti-viral actions, respectively (Douglas et al., 2010;Boison, 2013). BST2 was one of the top downregulated genes, and is known to be an interferon-inducible gene, and to play a role in innate immunity (Douglas et al., 2010). FKBP2 was also one of the top downregulated genes, which plays a role in immunoregulation by binding the immunosuppressive compound rapamycin that inhibits mTOR signaling (Hendrickson et al., 1993). This is in line with Transcripts not previously shown to be altered in PCOS muscle arrays (Skov et al., 2007(Skov et al., , 2008Nilsson et al., 2018) are marked in bold. ES, Enrichment score.
the downregulation of enriched gene sets associated with interferon alpha and gamma responses. Looking into the most upregulated genes, ADK is a phosphotransferase that converts the purine ribonucleotide adenosine into 5 ′ -adenosinemonophosphate (Boison, 2013). This reaction controls largely the adenosine tone, and consequently affects a wide range of functions. One of the functions that ADK and adenosine exert, is the regulation of the immune system function. Adenosine is evidenced to be generated as a result of stress response in damaged tissues and to display immunosuppressive properties (Hasko and Cronstein, 2004). A downregulation of inflammatory pathways in skeletal muscle is in contrast to the chronic low-grade pro-inflammatory state in adipose tissue (Dimitriadis et al., 2016). There is a lack of conclusive evidence, but androgens likely have direct effects on the adipose tissue resident immune cells (Huang et al., 2013). Changes in adipocyte function impair the secretion of adipokines, leading to decreased adiponectin secretion (Manneras-Holm et al., 2011), which promotes susceptibility to low grade inflammation. Despite the different immune responses in skeletal muscle and adipose tissue, it is evident that immune function is impaired in women with PCOS. The enrichment score (ES) was calculated for each gene set, the primary outcome of GSEA. NES, Enrichment score normalized for differences in gene set size.
Inflammatory and metabolic pathways in skeletal muscle are closely tied to cell signaling and differentiation, which leads to tissue specific adaptations in response to different stimuli e.g., androgens, obesity and exercise. The GSEA showed an upregulated myogenesis and the GO showed an enrichment of genes involved in negative regulation of cell proliferation, positive regulation of cell death, and negative regulation of response to stimuli, supporting skeletal muscle specific adaptations in response to PCOS. Furthermore, genes involved in MYC signaling were the most upregulated enriched gene set in the GSEA. MYC is a transcription factor that drives metabolic reprogramming in activated, primary T lymphocytes, and this may be one of the key targets responsible for metabolic reprogramming in response to PCOS conditions (Wang et al., 2011). As expected, many metabolic pathways were found to be enriched in the GO, with a negative regulation of the majority of the metabolic processes. ALDOA, encoding aldolase, was upregulated in women with PCOS and seems to be a key gene involved in many carbohydrate metabolic processes. Aldolase is an insulin stimulated glycolytic enzyme that catalyzes the conversion of fructose-1,6-bisphosphate during glycolysis, but it is also involved in actin remodeling, possibly modulating both metabolic and structural tissue adaptations (Hu et al., 2016).
Lastly, two genes were consistently upregulated in all three studies, IRX3 and CDKAL1. Both genes have been associated with body mass index and type 2 diabetes in genome-wide association studies (Steinthorsdottir et al., 2007;Ragvin et al., 2010;Yengo et al., 2018). IRX3 is a transcription and neuronal progenitor factor, and acts as a regulator of energy metabolism. It is shown to be controlled by a non-coding region of the fat mass and obesity-associated (FTO) gene, which is one of the strongest obesity-associated genes found in humans (Frayling et al., 2007). Many studies have investigated the role of BMIassociated variants in women with PCOS (Jones and Goodarzi, 2016). The results are controversial, but a meta-analysis showed that FTO has an increased effect on the BMI of women with PCOS, but it is not a PCOS susceptibility locus (Wojciechowski et al., 2012). Moreover, data for CDKAL1 and FTO suggest that polymorphisms of these genes are not associated with insulin resistance or insulin secretory capacity in Asian women with PCOS (Kim et al., 2012).
The interpretation of this meta-analysis is limited by a fairly low number of women with PCOS and controls (n < 40/group), the inclusion of a single ethnicity (Caucasian), and a mix of PCOS phenotypes (Rotterdam criteria). Further, our findings cannot distinguish whether the identified dysregulated FIGURE 2 | Functional annotation clustering determined using DAVID Bioinformatics Resources with respect to the 114 differently expressed genes (q < 0.05) in the skeletal muscle from women with PCOS compared with controls. (A) The representative groups with an enrichment score of 1.3 or above are presented. The x-axis represents the number of genes, while the y-axis represents the ontology categories. (B) Selected GO terms of biological processes enriched by the differently expressed genes (P < 0.05).
pathways are responses to a stimulus or the causal effectors. Another limitation is the fact that we cannot identify the gene expression changes within specific target cells i.e., myocytes, as the muscle biopsy is a bulk of many different cell types and structures, e.g., immune cells, vessels and connective tissue.

CONCLUSION
Clinical studies and transcriptomics data clearly demonstrate that molecular dysfunction in skeletal muscle contributes to insulin resistance in women with PCOS. Intracellular insulinsignaling pathways, mitochondrial function and fat oxidation, with a possible contribution of reduced adiponectin levels, have all been in focus for the mechanism of insulin resistance in PCOS (Stepto et al., 2019). Here, we identify 99 genes not previously shown to be altered in PCOS muscle, two of them consistently upregulated in all three studies. We show that immunometabolism can be added to the list of dysfunctional pathways and present genes with large effect size that warrant further investigation. However, gene expression does not always translate to alterations in biological function and data on total and phosphorylated protein levels in skeletal muscle is limited. Therefore, functional studies and proteomics analysis is needed to validate and advance our understanding of insulin resistance in skeletal muscle in these women.

DATA AVAILABILITY STATEMENT
All datasets presented in this study are included in the article/Supplementary Material.

AUTHOR CONTRIBUTIONS
MM, ES-V, and AB designed the study and contributed to the interpretation of data. MM and AB performed data analyses and statistical analyses, and wrote the manuscript. All authors critically revised and approved the manuscript. All authors contributed to the article and approved the submitted version. The funders had no role in the study design, data collection and analysis, decision to publish, or preparation of the manuscript.