Genome-wide gene expression array identifies novel genes related to disease severity and excessive daytime sleepiness in patients with obstructive sleep apnea

We aimed to identify novel molecular associations between chronic intermittent hypoxia with re-oxygenation and adverse consequences in obstructive sleep apnea (OSA). We analyzed gene expression profiles of peripheral blood mononuclear cells from 48 patients with sleep-disordered breathing stratified into four groups: primary snoring (PS), moderate to severe OSA (MSO), very severe OSA (VSO), and very severe OSA patients on long-term continuous positive airway pressure treatment (VSOC). Comparisons of the microarray gene expression data identified eight genes up-regulated with OSA and down-regulated with CPAP treatment, and five genes down-regulated with OSA and up-regulated with CPAP treatment. Protein expression levels of two genes related to endothelial tight junction (AMOT P130, and PLEKHH3), and three genes related to anti-or pro-apoptosis (BIRC3, ADAR1 P150, and LGALS3) were all increased in the VSO group, while AMOT P130 was further increased, and PLEKHH3, BIRC3, and ADAR1 P150 were all decreased in the VSOC group. Subgroup analyses revealed that AMOT P130 protein expression was increased in OSA patients with excessive daytime sleepiness, BIRC3 protein expression was decreased in OSA patients with hypertension, and LGALS3 protein expression was increased in OSA patients with chronic kidney disease. In vitro short-term intermittent hypoxia with re-oxygenation experiment showed immediate over-expression of ADAR1 P150. In conclusion, we identified a novel association between AMOT/PLEKHH3/BIRC3/ADAR1/LGALS3 over-expressions and high severity index in OSA patients. AMOT and GALIG may constitute an important determinant for the development of hypersomnia and kidney injury, respectively, while BIRC3 may play a protective role in the development of hypertension.

LGALS3 over-expressions and high severity index in OSA patients. AMOT and GALIG may constitute an important determinant for the development of hypersomnia and kidney injury, respectively, while BIRC3 may play a protective role in the development of hypertension.

Background
Obstructive sleep apnea (OSA) has been independently associated with endothelial dysfunction which may explain the increased risk for cardiovascular events, hypertension, and all-cause mortality in this population through chronic intermittent hypoxia with re-oxygenation (IHR) injury [1,2]. A recent meta-analysis has demonstrated that continuous positive airway pressure (CPAP) treatment significantly improved endothelial function as assessed by flow-mediated dilation [3]. It has been observed that circulating leukocytes of OSA patients exhibited markedly enhanced in vitro release of superoxide radical anions and delayed apoptosis, which may in turn lead to endothelial dysfunction [4,5]. However, underlying mechanisms by which IHR leads to endothelial dysfunction and other adverse consequences, such as hypertension, excessive daytime sleepiness (EDS), and chronic kidney disease (CKD), in OSA are largely unknown.
The aim of this study was to explore which expression signatures in peripheral blood could be representative of, or associated with OSA, by investigating the expressions of genes in peripheral blood mononuclear cells (PBMCs) that may be involved in these effects. We hypothesized that the gene expressions of PBMCs involved in the responses to chronic IHR in OSA patients would be markedly different from those in subjects with primary snoring, and that additional differences would be seen between OSA patients with and without hypertension or between OSA patients with and without EDS. Furthermore, we aimed to improve the understanding of the molecular signatures that can correlate endothelial dysfunction with long-term and sufficient CPAP treatment, with the hope that novel genes may be found to be over-or under-expressed after treatment.
A few studies have applied DNA microarray technology to investigate gene expressions in patients with OSA [6][7][8][9]. In one of these studies, which focused on gene expressions in the blood leukocytes from adult OSA patients, genes involved in modulation of reactive oxygen species, cell growth, proliferation, and the cell cycle were found to be altered after one night of sleep in four OSA patients free of any co-morbidity [6]. The extent to which the leukocyte genes play a role in OSA patients with adverse consequences, and the effect of nasal CPAP on gene signatures is unclear. Therefore, we extended our investigation into OSA patients with long-term CPAP treatment, hypertension, or EDS by analyzing whole-genome gene expression profiles of PBMC in three comparisons: (1) treatment-naïve moderate to very severe OSA patients versus subjects with primary snoring; (2) moderate to very severe OSA patients with hypertension or EDS versus those without hypertension or EDS, respectively; (3) treatmentnaïve very severe OSA patients versus those receiving at least one year of adequate CPAP treatment. We identified a novel association between several genes related to apoptosis or endothelial function and OSA or its clinical phenotypes.

Subjects
The study was approved by the Institutional Review Board of Chung Gung Memorial Hospital, Taiwan. The study participants were recruited from the pulmonary clinics and health examination center of Kaohsiung Chung Gung Memorial Hospital January 2012 through December 2014. Written informed consent was obtained from each subject participating in the study, who was aged 20 years or older. All study participants underwent full-night in-laboratory polysomnography examinations as described previously, and OSA was diagnosed according to the AASM guideline [10]. The exclusion criteria included ongoing infections, autoimmune disease, use of immunosuppressive agent in the past 6 months, narcolepsy, severe obesity (body mass index, BMI, ≧35 kg/m 2 ), old age (>65 year-old), and those with a BMI < 21 kg/m 2 . Study cohorts used for the whole-genome microarray gene expression and protein expression experiments included 48 (cohort 1) and 68 subjects (cohort 2: cohort 1 with expanded sample size), respectively. All these participants were classified into the following four groups based on apnea hypopnea index (AHI) and long-term use of CPAP: subjects with primary snoring (PS; AHI<5), treatment-naïve patients with moderate to severe OSA (MSO; 15<AHI≦50), treatment-naïve patients with very severe OSA (VSO; AHI>50), and very severe OSA patients on long-term CPAP treatment (VSOC; AHI>50 and regular use of CPAP: >4 hours/night, > one year). Nocturnal hypoxemia was evaluated in terms of the percentage of total minutes of recording time with oxyhemoglobin saturation <90% (%time <90% SaO2), and the number of dips >4% of basal SaO2%//h (oxygen desaturation index, ODI). The Epworth Sleepiness Scale (ESS) recorded at the PSG exam was used to measure sleep propensity, and EDS was defined as ESS>10. Hypertension was defined as baseline blood pressure>140/90 mmHg. Heart disease included ischemic heart disease, cardiac arrythmia, and congestive heart failure. CKD was defined as estimated glomerular filtration rate <60 mL/min/1.73 m 2 for ! 3 months.

Processes of RNA isolation and cRNA synthesis
Peripheral whole blood (20 ml) was collected at AM 6:00 to 8:00 after written informed consent was obtained from all the recruited participants. The PBMCs were isolated by Ficoll-Hypaque gradient centrifugation (HISTOPAQUE 1 -119, Sigma-Aldrich, Inc., St. Louis, MO USA) within 90 min of drawing blood, washed in PBS, and then stored in RNAlater (Ambion Inc., Austin, TX, USA) at -80˚C until RNA isolation. An RNeasy 1 Plus Mini Kit (Qiagen, Hilden, Germany) was used for isolation of high quality total RNA, and treated with DNase according to the manufacture protocol. RNA samples were run on a RNA 6000 Nano Gel System (Agilent Technologies Inc., Palo Alto, CA, USA) using an Agilent 2100 Bioanalyzer (Agilent) to determine the quality of RNA. Only samples with A260/A280 ratios of 1.9 to 2.1 were used for further analysis. A total of 300 ng RNA was used for synthesis of first strand cDNA and transcription of cRNA using an Illumina Totalprep RNA Amplication kit (Ambion, Inc.).
Gene expression profiling and microarray data analysis in the study cohort 1 Illumina (San Diego, CA) HumanRef-12-version 2 bead microarrays were used with 750 ng labeled cRNA for each sample according to the manufacturer's protocol. Human Ref-12 version 2 arrays consist of 27,455 probes representing 21,910 unique human genes on an eightstrip format array. All expression dataset has been deposited in the NCBI Gene Expression Omnibus with accession number of GSE75097.The preprocessing, quality control, background subtraction, quantile normalization and log2 transformation of array data for cross comparison were performed using BeadStudio software. Statistical analysis of the microarray data was further performed, using GeneSpring software version 11 (Agilent Technologies Inc., Santa Clara, CA, USA) as previously described [11]. Because the sample sizes were relatively small and statistics could be influenced by outliers, a non-parametric U test for unpaired comparisons of the two independent groups was applied. The Benjamini-Hochberg false-discovery rate correction method was used for controlling false positives and a corrected p-value cutoff of 0.05 was used to select the sets of significantly up-and down-regulated genes. The gene sets were then grouped into functional categories according to the Gene Ontology Biological Processes Classification. The pathway enrichment process was used to find direct relationships between genes of interest. This was performed in GeneSpring software with the "simple and direct interaction" algorithm.
Measurement of protein expression levels of five selected genes from PBMC samples in the study cohort 2 PBMCs were lysed in radio-immuno-precipitation assay-buffer containing a protease inhibitor cocktail (Sigma-Aldrich). Protein lysate normalized to 20 ng total protein was used to measure angiomotin (AMOT) P130, pleckstrin homology domain containing, family H member 3 (PLEKHH3), baculoviral IAP repeat-containing 3 (BIRC3; LGALS3), adenosine deaminase RNA-specific (ADAR1) P150, and galectin-3 internal gene (GALIG) levels by a commercial enzyme-linked immunosorbent assay (ELISA) kit, where AMOT kit was from USCN Business Co. (China), PLEKHH3 kit was from Sunlong (China), BIRC3/ADAR kits were from Cusabio (China), and LGALS3 kit was from R&D Systems (Minneapolis, MN). Briefly, 40 μl dye (Bio-Rad Protein Assay Dye Reagent Concentrate #500-0006) was added to 10 μl bovine serum albumin with serial dilution to generate a standard curve for total protein concentration by measuring values at OD 595nm.
In vitro blood cell culture under IHR conditions PMBC from six healthy subjects (600 μL per well, and adjusted to 1×10 6 cells per ml) were exposed to normoxia (NOX) or IHR in a custom-designed, incubation chambers which are attached to an external O2-CO2-N2 computer-driven controller, as described previously [12]. Air-phase set point consisted of a 35-min hypoxic period (0% O2 and 5% CO2), followed by 25 min of re-oxygenation (21% O2 and 5% CO2), using the BioSpherix OxyCycler C42 system (BioSpherix, Redfield, NY), 7 hours each day for 4 days. Control cells were maintained in NOX conditions for the same durations. Previous studies have shown that a 30-40% decreases in blood SaO2 could be achieved in the conditioned media by 25 min of continuous exposure of cells to 0%O2 and 5%CO2 [13,14]. Protein expression levels of the five selected genes in the PBMC samples were determined using ELISA method, as described above.

Statistical analysis
Continuous values were presented as the mean ± standard deviation (SD). ANOVA test followed by post hoc analysis with Bonferroni test was used for comparing mean values of more than two experimental groups in case of homogeneous data, while Brown-Forsythe test followed by post hoc analysis with Games-Howell test was used in case of non-homogeneous data. Chi-square tests were used to assess the differences of category values between different groups. In subgroup analyses, multivariate linear regression model was used to adjust for confounding factors, including age, gender, BMI, co-morbidities (hypertension, diabetes mellitus, stroke, heart disease, and CKD), smoking and alcoholism history, and to obtain adjusted p values. Pearson correlation test was used to assess the correlation between two continuous variables. All tests were two tailed and the null hypothesis was rejected at p < 0.05. A statistical software package (SPSS, version 15.0, Chicago, IL) was used for all analyses.

Demographic data of the participants
The baseline, sleep, and biochemistry data of the study cohort 1 and 2 are listed in Table 1.
The study population was all residents in Taiwan. Age, BMI, male gender ratio, smoking history, alcoholism history, and co-morbidities were all matched among the four subgroups, except that more patients in the VSOC group had diabetes mellitus than those in the other 3 groups in the study cohort 2. Differentially expressed genes (DEG) identified in the whole-genome gene expression experiment in the validation cohort 1. Microarrays from the 48 subjects that passed the quality-control filters were included in this study. From the 22,150 probes, 285 transcripts were differentially expressed between the treatment-naïve OSA patients (MSO and VSO groups) and PS group using a Mann-Whitney unpaired test (comparison I). To assess the effect of long-term CPAP treatment on PBMC, gene expression profiles in the VSOC group was compared with that in the VSO group, and differential expression of 1583transcripts was identified (comparison II). The intersection of comparison I and II result in eight DEG up-regulated with OSA and down-regulated with CPAP treatment, and five DEG down-regulated with OSA and up-regulated with CPAP treatment (Fig 1 and Table 2). Gene Ontology and gene interaction analyses revealed that AMOT, BIRC3, and ADAR are all involved in anti-apoptosis signaling, and LGALS3 serves as pro-apoptosis molecule, while AMOT, BIRC3, and PLEKHH3 are all involved in angiogenesis or tight junction (Fig 2A and 2B). To further clarify the effects of hypertension or EDS on the gene signature of PBMCs, gene expression profiles in the treatment-naïve OSA patients with hypertension or EDS was compared with those without hypertension or EDS, respectively. Six genes, including AMOT and PLEKHH3, were upregulated ( Table 3) in the treatment-naïve OSA patients versus PS subjects and further up-regulated in OSA patients with EDS as compared to that in those without EDS.LGALS3(fold change 2.23, p = 0.003) and metallopeptidase inhibitor 2 (TIMP2; fold change 4.44, p = 0.043) were up-regulated in the treatment-naïve OSA patients versus PS subjects, and further up-regulated in OSA patients with hypertension(LGALS3: fold change 4.37, p = 0.044; TIMP2; fold change 5.42, p = 0.018) as compared to that in those without hypertension ( Table 3).

Effects of in vitro IHR on protein expressions of the five selected genes
To determine whether IHR per se can affect expressions of the five selected genes, PBMC from six healthy subjects were exposed in vitro to either 7 cycles of IHR per day for 4 days or 4 days of continuous normoxic condition. Short-term IHR treatment in vitro resulted in a significant increase of ADAR1 protein expression (0.057±0.013 versus 0.025±0.008, p value<0.001, Fig  4H) as compared with normoxic condition, whereas no significant change was found for the other 4 protein expressions.

Discussion
In this study, we identified a novel association between several genes related to apoptosis or endothelial function and OSA or its clinical phenotypes through whole-genome gene expression microarray analyses. Moreover, we verified AMOT, PLEKHH3, ADAR, BIRC3, and LGALS3 protein over-expressions in the treatment-naïve OSA patients, and discovered a correlation between AMOT/BIRC3/LGALS3 protein expression and the presence of EDS/hypertension/CKD, respectively. Furthermore, we found an immediate effect of in vitro IHR on ADAR1 P150 protein over-expression in a cell culture model.  Endothelial tight junctions form a seal between polar cells, isolating the lumen of the blood vessel from the surrounding tissue and restricting the diffusion of solutes from the blood to the surrounding cells. AMOT is a trans-membrane receptor for the angiostatic factor angiostatin. Alternative splicing of AMOT mRNA results in two protein isoforms (p80 and p130), which exert very distinct roles during angiogenesis [15]. AMOT p80 stimulates endothelial cell migration and angiogenesis in response to angiostatin, whereas AMOT p130 localizes to actin and tight junctions [16,17]. AMOT is also involved in the maturation of epithelial tight junctions as part of the zonular signalosome formation [17,18]. In this study, AMOT gene and P130 protein expressions were both up-regulated in the treatment-naïve OSA patients, whereas AMOT gene expression was down-regulated and P130 protein expression was further up-regulated in the OSA patients with long-term CPAP treatment. We speculate that post-transcriptional modifications, maturation, alternative splicing, and degradation processes, may contribute to this discrepancy. The further increase of AMOT P130 protein expression in the VSOC group may be attributed to long-term air flow stimulation of upper airway epithelial cells with CPAP use. In line with our findings, AMOT expression has been shown to be a good indicator of plasticity of the vascular network in skeletal muscle, and P130 expression can be reduced with exercise training in an obese rat model [19].
Pleckstrin homology (PH) domains family function as a versatile protein-protein interaction platform and are integrated in an increasing number of available multidomain structures [20]. For example, PLEKHA7 is a recently identified protein of the epithelial zonula adhaerens, and stabilizes it by modulating the dynamics of assembly and disassembly of the tight junction barrier, through E-cadherin protein complex-and microtubule-dependent mechanisms [21]. Zonula adhaerens are topologically associated with tight junctions in the apical junctional complex at the apicolateral border of both epithelial and endothelial polar cells. In this study, we found a positive correlation between AHI and PLEKHH3 expression level, which was reversed with CPAP treatment. We think that PLEKHH3, another PH family member, may play a pivotal role in chronic IHR-related endothelial dysfunction in OSA patients through regulating tight junction formation of endothelial cells specifically.  Anti-apoptotic marker, BIRC3, has been shown to mediate the pro-survival and inflammatory responses induced by the docosahexaenoic acid / neuroprotectin D1 pathway under oxidative stress in an ischemia-reperfusion stroke model [22].BIRC3 has been demonstrated to be a downstream effector of HIF-1 signaling involved in the survival response of endothelial cells to hypoxia [23]. BIRC3 is also a negative regulator of the non-canonical NF-κB signaling pathway and mutated primarily in patients with aggressive chronic lymphocytic leukemia [24]. In this study, BIRC3 was up-regulated in the treatment-naive OSA patients, and down-regulated with CPAP treatment. BIRC3 under-expression was associated with the occurrence of hypertension in OSA patients. We speculate that BIRC3 may play a role in protecting from hypertension in OSA patients through inhibiting endothelial cell apoptosis and NF-κB signaling.
Adenosine deaminase acting on RNA1 (ADAR1) catalyzes the C6 deamination of adenosine (A) to produce inosine (I) in regions of double-stranded RNA, known as A-to-I RNA editing. Alternative splicing gives rise to transcripts that encode twoADAR1 protein size isoforms. ADAR1 p150 is an interferon-inducible dsRNA adenosine deaminase found in the cytoplasm and nucleus, mounting pro-viral and anti-apoptotic responses, whereas ADAR1 p110 is constitutively expressed in the nucleus [25]. Among the biologically relevant substrates of ADAR1 that result in amino acid coding changes following editing are transcripts for the 5-HT2c-R neurotransmitter receptor for serotonin [26]. In recent years, this modification has been discovered to occur not only in coding RNAs but also in non-coding RNAs, such as microRNAs, small interfering RNAs, transfer RNAs, and long non-coding RNAs [27].The malfunction of this editing machinery is associated with various human diseases, such as neurodegenerative, cardiovascular, and carcinogenic diseases. For the first time, we found an association between ADAR1 up-regulation and AHI, and demonstrated a direct link between IHR and ADAR1 P150 protein over-expression in vitro.
LGALS3 (GALIG) is a novel cell death gene encoding mitogaligin, a protein promoting cytochrome c release upon direct interaction with the mitochondria or nucleus [28,29].
LGALS3 pro-apoptotic gene is up-regulated during neutrophils apoptosis and underexpressed in acute myeloid leukemia cells [30]. In this study, we found that LGALS3 was up regulated in moderate to very severe OSA patients, especially in those with CKD. Impaired renal function has been observed in OSA patients with metabolic syndrome, hypertension or heart failure, possibly through reduced endothelial nitricoxide synthase expression [31,32]. Our findings suggest that LGALS3 may play a crucial role in mediating renal dysfunction in OSA patients.
There are several limitations to the present study. First, the cause and effect relationship could not be determined in this cross-sectional clinical study design, and further studies are required to elucidate underlying mechanisms for these five novel biomarkers. However, the in vitro experiment demonstrated immediate over-expression of ADAR1 P150 under short-term IHR stimuli. Second, gene expression levels of the five selected genes were not examined in the PBMC samples because of inadequate RNA samples. However, their protein expressions showed corresponding changes, indicating a functional role of these novel genes in mediating pathogenesis of OSA and its adverse consequences. Third, the sample size of each subgroup is relatively small, and many confounding factors may affect the expression levels. Further sleepiness. (F) PLEKHH3 protein expression was increased in the VSO group and decreased in the VSOC group. PLEKHH3 protein expression was positively correlated with (G) AHI and (H) AMOT P130 protein expression. *p<0.001 compared between the VSO and PS groups. **P<0.001, compared between the VSO and MSO groups. ***P<0.001, compared between the VSO and VSOC groups. #p<0.001, compared between the MSO and PS groups.
https://doi.org/10.1371/journal.pone.0176575.g003 studies with sufficiently large sample sizes are required for the internal and external validity and the reliability of the results.
In summary, we reported a novel association of increased AHI in OSA patients of Asian origin with over-expressions of the AMOT, PLEKHH3, ADAR1, BIRC3, and LGALS3 genes in blood immune cells. The findings extend reports linking AMOT with EDS in OSA patients, and provide direct evidence that perturbation of BIRC3 and LGALS3 signaling may play an important role in the mediation of hypertension and CKD in OSA patients, respectively.