Plasma Proteomics in Healthy Subjects with Differences in Tissue Glucocorticoid Sensitivity Identifies A Novel Proteomic Signature

Significant inter-individual variation in terms of susceptibility to several stress-related disorders, such as myocardial infarction and Alzheimer’s disease, and therapeutic response has been observed among healthy subjects. The molecular features responsible for this phenomenon have not been fully elucidated. Proteomics, in association with bioinformatics analysis, offer a comprehensive description of molecular phenotypes with clear links to human disease pathophysiology. The aim of this study was to conduct a comparative plasma proteomics analysis of glucocorticoid resistant and glucocorticoid sensitive healthy subjects and provide clues of the underlying physiological differences. For this purpose, 101 healthy volunteers were given a very low dose (0.25 mg) of dexamethasone at midnight, and were stratified into the 10% most glucocorticoid sensitive (S) (n = 11) and 10% most glucocorticoid resistant (R) (n = 11) according to the 08:00 h serum cortisol concentrations determined the following morning. One month following the very-low dose dexamethasone suppression test, DNA and plasma samples were collected from the 22 selected individuals. Sequencing analysis did not reveal any genetic defects in the human glucocorticoid receptor (NR3C1) gene. To investigate the proteomic profile of plasma samples, we used Liquid Chromatography–Mass Spectrometry (LC-MS/MS) and found 110 up-regulated and 66 down-regulated proteins in the S compared to the R group. The majority of the up-regulated proteins in the S group were implicated in platelet activation. To predict response to cortisol prior to administration, a random forest classifier was developed by using the proteomics data in order to distinguish S from R individuals. Apolipoprotein A4 (APOA4) and gelsolin (GSN) were the most important variables in the classification, and warrant further investigation. Our results indicate that a proteomics signature may differentiate the S from the R healthy subjects, and may be useful in clinical practice. In addition, it may provide clues of the underlying molecular mechanisms of the chronic stress-related diseases, including myocardial infarction and Alzheimer’s disease.


Introduction
Glucocorticoids exert numerous actions in target tissues through binding to their cognate receptor, the glucocorticoid receptor (GR), which functions as a ligand-activated transcription factor influencing the transcription rate of a large number of glucocorticoid responsive genes [1][2][3]. Moreover, the ligand-bound GR may influence gene expression in a positive or negative fashion through protein-protein interactions with other transcription factors, including the nuclear factor-κB (NF-κB) [4], signal transducers and activators of transcription (STATs) [5], and the activator protein-1 (AP-1) [4,6]. In addition to their genomic actions, glucocorticoids may exert rapid nongenomic actions mediated by a membrane-anchored GR, which triggers the activation of the mitogen-activated protein kinases (MAPK) [7,8] or the phosphatidylinositol 3-kinase (PI 3 K) [8,9]. Finally, accumulating evidence suggests that glucocorticoids may signal through mitochondrial GRs, thereby regulating fundamental cellular function, including apoptosis, energy production, and stress response [10].
In humans, tissue glucocorticoid sensitivity is complex and multifactorial, depending on polymorphisms [11] or genetic defects in the NR3C1 gene that encodes the human GR (hGR) [12], several tissue-specific hGR protein isoforms [13,14], or the hGR "interactome" influenced by an ever-increasing number of interacting proteins, miRNAs [15] and non-coding RNAs [15,16]. All these regulatory factors influence substantially tissue responsiveness to glucocorticoids in clinical practice. Indeed, a large number of patients with inflammatory, allergic, hematologic, and lymphoproliferative diseases are treated with synthetic glucocorticoids due to their anti-inflammatory and immunosuppressive effects [17]. However, the therapeutic outcome might be compromised by the adverse effects of hypercortisolism [18]. Therefore, the identification of novel protein biomarkers that could stratify healthy adults as most glucocorticoid-sensitive or most glucocorticoid-resistant is currently an unmet clinical need.
Proteomics approaches allow the identification and quantification of proteins from biological samples and could lead to the identification of biomarkers associated with the pathogenesis, diagnosis, and treatment of many diseases [19]. Of note, untargeted proteomics has not yet been used to investigate differences among healthy individuals with variations in tissue glucocorticoid sensitivity. The aim of our study was to study any differences at the proteomics level among healthy participants with increased or decreased sensitivity to glucocorticoids, and to identify a proteomic signature that could be used to distinguish most glucocorticoid-sensitive (S) from most glucocorticoid-resistant (R) healthy adult subjects, who do not harbor any genetic defects in the human glucocorticoid receptor (NR3C1) gene, using untargeted plasma proteomics.

Selection of the Healthy Subject Cohort
Fifty males and fifty-one females (n = 101) of a mean age of (±SD) 26.5 (±5.0) years, were recruited to participate in our study prospectively. The participants had unremarkable medical history and were not receiving any medications, including oral contraceptives for females. To study the inter-individual variation in tissue glucocorticoid sensitivity among the participants, the very low dose dexamethasone suppression test was applied, as previously described in studies with healthy subjects [20,21]. According to this protocol, all participants were given per os 0.25 mg of dexamethasone at midnight and were advised to wake up at 07.00 h. The following morning (08.00 h), serum cortisol (261.8 ± 206.9 nmol/L) and plasma ACTH (16.1 ± 12.2 pg/mL) were measured. Depending on their serum cortisol concentrations, all participants were rank ordered. The 11 subjects who had the lowest serum cortisol concentrations were classified as the most glucocorticoid sensitive (S) group, while the 11 subjects with the highest serum cortisol concentrations were classified as the most glucocorticoid-resistant (R) group ((mean serum cortisol concentrations ± SD: 34.4 ± 15.0 nmol/L in the S participants, mean serum cortisol concentrations ± SD: 622.4 ± 93.7 nmol/L in the R participants, p < 0.001); (mean plasma ACTH concentrations ± SD: 2.8 ± 2.4 pg/mL in the S participants, mean plasma ACTH concentrations ± SD: 31.6 ± 10.6 pg/mL in the R participants, p < 0.001)) ( Table 1). Thirty days following the very low dose dexamethasone suppression test, plasma, DNA and RNA samples were collected from the 22 individuals for further analyses.

Approval of the Study
The study was approved by the Aghia Sophia Children's Hospital Committee on the Ethics of Human Research EB-PASCH-MoM: 13/02/2013, Re: 1490-21/01/2013). Written informed consent was obtained by all subjects prior to their participation in the study.

Sample Collection
Whole blood samples (total volume: 5 mL) were obtained from the study participants at 08:00 h in EDTA-containing tubes, centrifuged immediately following collection, and stored at −80 • C until further analysis. Two aliquots of 150 µL per blood sample were shipped in dry ice to the proteomics laboratory.

Assays
The standard hematological, biochemical, and endocrinological investigations were determined, as previously described [22]. Taking into account the reproductive cycle of the females, gonadotropin and sex steroid concentrations were determined between the 3rd and 5th day of the menstrual cycle.

Sequencing of the NR3C1 Gene
Genomic DNA was isolated from whole blood samples from the 22 selected participants using the Maxwell 16 instrument for automated DNA extraction (Promega Corp., Madison, WI, USA). The protein-coding sequences of the NR3C1 gene, including the junctions of introns and exons, were PCR-amplified and sequenced using the Big Dye Terminator cycle sequencing kit (Applied Biosystems, Carlsbad, CA, USA) on an ABI 3100 sequencer (ABI 3100, Applied Biosystems, Carlsbad, CA, USA), as previously described [22,23].

Sample Preparation for Proteomics Analysis
Plasma samples (200 µg total protein per sample) were processed with the filter aided sample preparation (FASP) method as described previously [24], with minor modifications [25]. Briefly, proteins were reduced with DTE (0.1 M), alkylated with iodoacetamide (0.05 M), and digested overnight by trypsin in 50 mM NH 4 HCO 3 pH 8.5. The peptides were lyophilized and kept at −80 • C.

LC-MS/MS Analysis
Samples were resuspended in 200 µL mobile phase A (0.1% FA). A 5 µL volume was injected into a Dionex Ultimate 3000 RSLS nano flow system (Dionex, Camberly, UK) configured with a Dionex 0.1 × 20 mm, 5 µm, 100 Å C18 nano trap column with a flow rate of 5 µL/min. The analytical column was an Acclaim PepMap C18 nano column 75 µm × 50 cm, 2 µm 100 Å with a flow rate of 300 nL/min. The trap and analytical column were maintained at 35 • C. The chromatography method used had a duration of 4 h and peptide elution was achieved by the following gradient: 2% to 33% mobile phase B (ACN 100%, Formic acid 0.1%). MS/MS analysis was performed with a Q Exactive operated with Data Dependent Acquisition mode.

MS Data Processing
Proteome Discoverer 1.4 was used for processing the raw mass spectrometry data with the Sequest search algorithm and the Uniprot FASTA file for Homo sapiens (20,243 entries downloaded on 15 December 2017). Cysteine cabamidomethylation was considered as a fixed modification, whereas methionine oxidation was considered as variable modification. The mass tolerance for precursors was set to 10 ppm, whereas for fragment ions it was 0.05 Da. The False discovery rate (FDR) was set to 0.01.

Data Processing (Methods)
Data processing, visualizations, and statistical analysis were conducted in the R environment (version 3.6) for Windows. Raw spectral areas per sample were merged with an in-house script and were subjected for sample normalization according to X = [X/sum(Xi)] × 10 6 , with X being the raw area of a given protein ID for a given sample, sum(Xi) the sum of all raw protein areas for the same sample and X the normalized protein area. To increase the power of the downstream analysis, proteins showing null abundance in more than 65% of samples in both biological groups ("hypersensitive" and "resistant"; described in the results section) were excluded. The non-parametric Mann-Whitney test was utilized for assessing statistical significance of the continuous variables. Heatmap was constructed after row z-normalization with the package ComplexHeatmap and the volcano plot with the package EnhancedVolcano. Pathway analysis was conducted in the Cytoscape plug-in ClueGO+CluePedia [26] and significantly enriched pathways were defined with a twosided hypergeometric test corrected with Benjamini-Hochberg (p < 0.05). Grouping of the enriched terms was based on 50% gene similarity. A random forest-based classifier was trained to distinguish between the two biological groups using the package randomForest. The complete study population was utilized for training and feature selection involved proteins with no missing values across samples and specifically the union of statistically significant proteins (defined by Mann-Whitney p < 0.05; n = 11) and of proteins with fold change >1.5 (n = 3 extra proteins). Tuning of classifier parameters was conducted iteratively and was based on minimizing the out of bag error [27], which is an estimate of the correct classification rate of random partitions of the training set. For optimal tuning, number of trees and variables used at each split was set to 48 and 3, respectively, with a stepfactor of 0.5 and improve of 0.05. Variable importance was assessed in terms of the increase in the misclassification rate after permutation and comparison of the accuracy to the standard model (Mean Decrease in Accuracy and Mean Decrease in Gini index; [28]) and also in terms of the mean minimal depth achieved over all trees, which indicates the impact of a variable in the classification of observations (the lower the minimal depth of a variable the highest the number of observations that were classified on the basis of that variable), with the packages caret and randomForestExplainer.
Summaries and detailed data on protein intensities and statistics can be found in the Supplementary Tables.

Clinical Characteristics, Hematological, Biochemical, and Endocrinological Parameters in the Most Glucocorticoid Sensitive (S) and Most Glucocorticoid Resistant (R) Healthy Subjects
The clinical characteristics (sex, weight, height, and body mass index (BMI)), serum cortisol and plasma ACTH concentrations at the time following the very low dose dexamethasone suppression test are presented in Table 1. All the hematological, biochemical and endocrinological parameters of the S and R healthy subjects measured one month after the very-low dexamethasone suppression test are presented in Table 2. No statistically significant differences were found between the S and R groups (p-value > 0.05).

The S and R Healthy Participants Did Not Harbor Any Genetic Defects in the NR3C1 Gene
To study whether the observed variation in tissue glucocorticoid sensitivity among the S and R healthy participants might be caused by genetic defects in the NR3C1 gene, the protein-coding exons and the intron-exon junctions were first PCR amplified and then sequenced. We have not identified any point mutations, deletions, insertions, or polymorphisms in the NR3C1 gene of the S and R participants.

Proteomics Data Analysis
In total, there were 2737 proteins identified and quantified in at least one of the analyzed samples (listed in the Supplementary Tables). Based on their response to cortisol, samples were labeled either as "most sensitive" (mean serum cortisol concentrations ± SD: 30 ± 20 nmol/L; n = 11 samples) or "most resistant" (mean serum cortisol concentrations ± SD: 620 ± 90 nmol/L; n = 11 samples), and statistical analysis was conducted for the comparison of resistant vs sensitive. After selecting those proteins with presence in at least 35% of the samples in one of the two groups (n = 466 features), differentially expressed proteins were defined as the subset of statistically significant changes (Mann-Whitney p-value < 0.05) with a fold change greater than 1.5 (or less than 0.67). This counted for 66 proteins with higher abundance in the resistant and 110 proteins with higher abundance in the sensitive group ( Figure 1). Among them, there were 21 proteins being present exclusively in only one of the two groups (Table 3). Over-representation test (ClueGO plug-in) for proteins being overexpressed in the sensitive group (FC < 0.67; n = 110 proteins) yielded some significantly enriched Reactome pathways (Table 4), including erythrocyte gas exchange and platelet activation and aggregation. However, there were no significant pathway enrichments for proteins being at higher abundance in the resistant group (FC > 1.5; n = 66 proteins). Detailed protein data for this analysis can be found in the "35% threshold" spreadsheet of the Supplementary Tables.
In order to predict the response to cortisol prior to administration, a random forest classifier was developed out of the proteomics data, to distinguish between sensitive and resistant samples. As a training set, we used the complete set of samples (n = 22) aiming to validate the classifier in new independent cohorts. For optimal training, feature selection involved proteins with no missing values across samples and specifically those passing Mann-Whitney p < 0.05, as well as those with a fold change >1.5 (and <0.67; n = 14). After tuning for optimal parameters, the classifier showed promising results in correctly assigning random partitions of the training data to the studied groups, achieving an overall accuracy score of 0.86. In detail, there were only three misclassified cases, one for the sensitive and two for the resistant groups. The individual importance of each protein in the model was evaluated in terms of the Mean Decrease in Accuracy, the Mean Decrease in Gini index and the mean minimal depth. Out of the 14 proteins utilized for training, apolipoprotein A4 (APOA4) and gelsolin (GSN) were the most important variables in the classification (Figure 2), and warrant further investigation to determine their prognostic capacity.

Discussion and Conclusions
In this study, 22 individuals were selected as the most glucocorticoid sensitive (S) (n = 11) and most glucocorticoid resistant (R) (n = 11) from a large sample of 101 healthy adults, according to their 08.00 h serum cortisol concentrations following the very low

Discussion and Conclusions
In this study, 22 individuals were selected as the most glucocorticoid sensitive (S) (n = 11) and most glucocorticoid resistant (R) (n = 11) from a large sample of 101 healthy adults, according to their 08.00 h serum cortisol concentrations following the very low dose dexamethasone suppression test. This test has been widely used in studies aimed to identify mild differences in tissue glucocorticoid sensitivity due to interindividual variation among healthy adults [20,21]. Since serum cortisol concentrations are within normal range in these studies, the application of a cortisol cut-off was out of context; therefore, the 101 participants were rank-ordered according to their cortisol concentrations and the two normal "extreme" groups were selected for further analyses. Sequencing analysis of the NR3C1 gene did not identify any genetic defect in the 22 participants, indicating that the above-mentioned differences at the proteomics level could not be attributed to increased or decreased expression or activity of the human glucocorticoid receptor.
Untargeted plasma proteomics analysis revealed 110 up-regulated and 66 downregulated proteins in the S compared to the R group. Importantly, the majority of the up-regulated proteins in the S group were implicated in erythrocyte gas exchange (take up carbon dioxide and release oxygen), suggesting a state of the organism that is more capable to respond to stressful stimuli. This result is in line with the basic concepts of physiology of the stress system [29][30][31][32]. Upon exposure to external or internal stressors, a highly conserved neuroendocrine system, the stress system, is activated to restore the threatened or perceived as threatened internal balance or homeostasis. To achieve homeostasis, the stress system, through concurrent activation of the hypothalamic-pituitary-adrenal (HPA) axis and the locus caeruleus-autonomic nervous system, leads to behavioral and physical adaptation changes, thereby increasing the chances for survival [29][30][31][32]. Thus, both oxygen and nutrients are redirected to the central nervous system and the stressed body site(s) through increased respiratory rate, elevated cardiovascular tone and catabolism [29,30]. The S group had most of the up-regulated proteins involved in erythrocyte take up of carbon dioxide and release of oxygen, probably as part of the physical adaptation changes observed in stress-related homeostasis.
In addition to the increased gas exchange of erythrocytes, the S group was characterized by increased expression of proteins involved in platelet activation, aggregation, and degranulation. Of note, there were not any statistically significant differences between the number of platelets of the S participants compared to that of the R participants, indicating that the S group displayed a cellular "predisposition" for increased coagulation. A growing body of evidence suggests that natural and synthetic glucocorticoids contribute substantially to platelet activation and aggregation [33]. Indeed, the serum/glucocorticoid regulated kinase 1 (SGK1), which is a regulator of ORAI1 protein that mediates the storeoperated Ca 2+ entry in platelets, has been shown to be influenced by several hormones, including glucocorticoids [34]. Furthermore, endogenous hypercortisolism has been found to increase thromboxane A2 (TXA2) biosynthesis, leading to enhanced platelet aggregation [35]. The increased platelet activity and aggregation underlies the increased coagulation causing blood vessel infarctions, including myocardial and/or brain infarctions, often observed in chronically stressed individuals. Although the S participants did not have endogenous or exogenous hypercortisolism, we speculate that their increased stress response might predispose them to develop blood vessel infarction upon exposure to chronic stress.
Our study also succeeded in identifying a protein signature that could distinguish the S from the R participants. Out of the 14 proteins utilized for training in the random forest classifier, APOA4 and gelsolin GSN were the most important variables in the classification. The expression of ApoA4 was significantly higher in obese and overweight children compared with their normal counterparts. ApoA4 has an antiatherogenic function [36], is increased in obese subjects and declines with weight reduction [37]. Indeed, fasting plasma ApoA4 concentrations are significantly elevated in obesity and decrease to almost half of baseline concentrations during weight reduction [37]. In addition to APOA4, GSN was also found as one of the most important variables in the classification. GSN is one of the most abundant actin-binding proteins and plays a crucial role in cellular mechanisms and interactions [38]. Since GSN participates in several immunologic functions, and interacts with different types of immune cells, this protein has been identified as a potential therapeutic target. More specifically, previous studies have identified the anti-amyloidogenic role of GSN in Alzheimer's disease (AD) [39][40][41]. GSN can reduce amyloid burden by acting as an inhibitor of amyloid beta-protein (Aβ) fibrillization, and as an antioxidant and antiapoptotic protein. When GSN was administered or overexpressed in AD transgenic mice, the amyloid load was significantly reduced and the Aβ level was decreased, suggesting that GSN might be implicated in the treatment of AD [39]. The anti-amyloidogenic activity in the S participants is in line with the results of our previous published study, in which we used transcriptomics in these 22 participants, and found that the S participants had decreased expression of genes involved in Alzheimer's disease [23].
To the best of our knowledge, this is the first pilot study that identified differences at the proteomics level among healthy adults with variations in tissue glucocorticoid sensitivity. However, the study has some limitations. The sample of participants was relatively small to extract concrete conclusions. Undoubtedly, the proteomics analysis must be repeated in a larger cohort in order to confirm the differences reported between S and R individuals. In the future, it is planned to validate the differential expression of selected proteins and the deregulation of key pathways (erythrocyte gas exchange, and platelet activation) by independent analytical approaches (ELISA, flow cytometry, etc.). Another limitation of our study was the use of IMMULITE instead of HPLC for measuring serum cortisol concentrations. The use of IMMULITE is the standard procedure that we normally follow in our every-day clinical practice. Finally, the dexamethasone concentrations of the 101 participants were not measured following the very-low dexamethasone suppression test because we carried out the protocol performed by Donn and collaborators, who found a new glucocorticoid sensitivity-determining gene using gene expression profiling [21]. It is worth mentioning that in most published studies that use the very-low dose dexamethasone suppression test, the investigators have not routinely determined dexamethasone concentrations in their participants.
In conclusion, untargeted proteomics analysis revealed differences among healthy subjects with differences in glucocorticoid sensitivity. The S individuals displayed higher erythrocyte gas exchange and increased platelet activation and aggregation that might lead to increased risk of stress-related disorders, including myocardial and brain infarctions.

Conflicts of Interest:
The authors declare no conflict of interest.