Changes in Splicing Machinery Components Influence, Precede, and Early Predict the Development of Type 2 Diabetes: From the CORDIOPREV Study.

BACKGROUND: Type-2 diabetes mellitus (T2DM) is a major health problem with increasing incidence, which severely impacts cardiovascular disease. Because T2DM is associated with altered gene expression and aberrant splicing, we hypothesized that dysregulations in splicing machinery could precede, contribute to, and predict T2DM development.
METHODS: A cohort of patients with cardiovascular disease (CORDIOPREV study) and without T2DM at baseline (at the inclusion of the study) was used (n = 215). We determined the expression of selected splicing machinery components in fasting and 4 h-postprandial peripheral blood mononuclear cells (PBMCs, obtained at baseline) from all the patients who developed T2DM during 5-years of follow-up (n = 107 incident-T2DM cases) and 108 randomly selected non-T2DM patients (controls). Serum from incident-T2DM and control patients was used to analyze in vitro the modulation of splicing machinery expression in control PBMCs from an independent cohort of healthy subjects.
FINDINGS: Expression of key splicing machinery components (e.g. RNU2, RNU4 or RNU12) from fasting and 4 h-postprandial PBMCs of incident-T2DM patients was markedly altered compared to non-T2DM controls. Moreover, in vitro treatment of healthy individuals PBMCs with serum from incident-T2DM patients (compared to non-T2DM controls) reduced the expression of splicing machinery elements found down-regulated in incident-T2DM patients PBMCs. Finally, fasting/postprandial levels of several splicing machinery components in the PBMCs of CORDIOPREV patients were associated to higher risk of T2DM (Odds Ratio > 4) and could accurately predict (AUC > 0.85) T2DM development.
INTERPRETATION: Our results reveal the existence of splicing machinery alterations that precede and predict T2DM development in patients with cardiovascular disease. FUND: ISCIII, MINECO, CIBERObn.

RCT Entities:   Population Interventions Outcomes

Evidence before this study

Early identification of individuals at high risk of type 2 diabetes mellitus (T2DM) development, especially among patients with cardiovascular disease (CVD), is critical for prevention. Traditionally, the identification strategy has been based on biomarkers such as glycated hemoglobin (HbA1c) or predictive scores like FINDRISK. However, these approaches have limitations and cannot precisely predict an individual's risk of developing T2DM. In this context, the identification of key modifiers of phenotypic plasticity that define individual susceptibility to develop T2DM may hold predictive potential. Particularly, there is emerging evidence that under adverse metabolic conditions, including obesity and T2DM, the splicing process and, consequently, the generation of splicing variants is markedly dysregulated in most tissues. However, until the present study, the association between potential changes in the regulatory machinery governing the splicing process and alterations in phenotypic plasticity —which is essential to maintain whole body homeostasis, and therefore, to prevent the development of metabolic pathologies— was still unknown. In this context, we took advantage of the emerging evidence showing that gene expression pattern of peripheral blood mononuclear cells (PBMCs) commonly reflects and accompany disease-characteristic expression patterns, and may thus serve as a general sentinel, biosensor and early indicator of the instauration of metabolic disease. Until now, it was not known whether the expression of the splicing machinery components is dysregulated in PBMCs under conditions of loss of phenotypic plasticity and, therefore, if it could be associated with the subsequent development of T2DM.

Added value of this study

We took advantage of the CORDIOPREV study (CORonary Diet Intervention with Olive oil and cardiovascular PREVention), an ongoing, prospective, randomized trial that includes the 5-year follow-up of patients with CVD, but without T2DM at the inclusion of the study (baseline), to show, for the first time, that striking changes in the expression of key splicing machinery components of PBMCs precede and can early predict the development of T2DM. Moreover, our study demonstrates that expression pattern of the splicing machinery components can be dynamically modulated during the post-prandial phase. These rapid changes likely reflect the phenotypic flexibility of these patients in response to a metabolic challenge, and, most importantly, were clearly altered at the inclusion of the study in patients that developed T2DM during the 5-year follow-up. In particular, the fasting/postprandial expression levels of a set of bioinformatically-selected splicing machinery components in the PBMCs of CORDIOPREV patients were associated to higher risk of T2DM (Odds Ratio > 4) and could accurately predict (AUC > 0.85) diabetes development in this cohort, which is significantly higher than the current standard procedures (HbA1c and FINDRISK exhibited AUCs < 0.66 in our cohort). Finally, incubation of PBMCs with serum from patients that developed, or not, T2DM revealed a potential causal relationship in that the latter was able to reproduce some of the changes observed in the expression of key splicing machinery components.

Implications of all the available evidence

Altered expression of splicing machinery components may be associated with the development of T2DM, preceding the clinical instauration of this pathology, and could serve as a sensor and early predictor for T2DM development in CVD patients. Our data demonstrate the existence of a splicing machinery-associated molecular fingerprint capable to predict the future development of T2DM in individual patients with high precision, which even outperforms the capacity of classical predictors of T2DM development. Therefore, this splicing machinery-associated molecular fingerprint could become a valuable, non-invasive, new tool for early risk assessment of T2DM in clinical practice to prevent disease development.

Alt-text: Unlabelled Box

Introduction

Type 2 diabetes mellitus (T2DM) is a major health problem [1] closely associated to cardiovascular disease (CVD), wherein patients with myocardial infarction and T2DM have a higher risk of developing a new cardiovascular event than those without T2DM [2]. Thus, early identification of individuals at high risk for T2DM development, especially among patients with CVD [3,4], is critical for prevention [5,6]. Traditionally, this strategy has been based on biomarkers [glycated hemoglobin (HbA1c)] or predictive scores (FINDRISK); however, these approaches have limitations and cannot precisely predict an individual's risk of developing T2DM [7,8].

T2DM development is critically affected by the loss of phenotypic flexibility (i.e. the difficulty to cope with stressors to maintain metabolic homeostasis), wherein the identification of key modifiers of phenotypic plasticity that define individual susceptibility to develop T2DM may hold predictive potential [9]. Particularly, there is emerging evidence that under adverse metabolic conditions (obesity, insulin resistance, etc.) the splicing machinery is markedly dysregulated in most tissues [[10], [11], [12]], including peripheral blood mononuclear cells (PBMCs); and that dysregulations in splicing processes are associated with development of several pathologies [10,13,14]. Specifically, splicing process is catalyzed by the minor and major spliceosomes, which act on different types of introns [15]. The spliceosome is an intricate macromolecular complex, whose functional core is comprised by several small nuclear ribonucleoproteins (snRNPs) subunits, which dynamically interact to regulate the splicing process. In addition, the activity of the spliceosome is modulated by >300 splicing factors (SFs) that specifically recognize certain sequences in exons and introns. Consequently, the dysregulation of the expression and/or function of certain spliceosomal components may drive the aberrant alteration of the normal splicing process [16]. Indeed, correct function of the splicing machinery, i.e. spliceosome components and SFs, is an essential mechanism to maintain whole body homeostasis [11,12,17,18]. Therefore, perturbed splicing may play a major role as a pathogenic factor for and may serve as a predictive marker of the development of T2DM [12,19].

In this sense, since gene expression patterns in PBMCs commonly reflect and accompany disease-characteristic expression patterns [20], the hypothesis of this study is that the pattern of expression of certain splicing machinery elements in the PBMCs, especially during the post-prandial response, when changes in phenotypic flexibility are more evident, could be associated with the risk and could early predict the development of T2DM. To test this idea, we took advantage of the CORDIOPREV study (CORonary Diet Intervention with Olive oil and cardiovascular PREVention), an ongoing, prospective and randomized trial that includes the follow-up of patients with CVD and, therefore, with high risk to develop T2DM.

Material and methods

Study cohort

This work has been conducted in the context of the CORDIOPREV study (Clinical Trials Registry NCT0092493741) [21]. The study protocol was approved by the Human Investigation Review Committee of the Reina Sofia University Hospital (HURS, Cordoba, Spain), according to institutional and Good Clinical Practice guidelines. The CORDIOPREV study is a prospective, randomized, controlled trial that includes 1002 coronary heart disease (CHD) patients, who had their last coronary event over six months before joining the study. Details of the study, including inclusion and exclusion criteria and cardiovascular risk factors of the patients have been published elsewhere [21]. Briefly, patients aged 20–75 years, with established CHD but without clinical events in the last six months with no other serious illnesses and a life expectancy of at least five years were eligible, with the intention of following a long-term monitoring study. In addition to conventional treatment for CHD, subjects were randomized in two different dietary models (Mediterranean and low-fat diets). The intervention phase is still in progress, and will have a follow up of seven years. Specifically, 462 of those patients were non-T2DM at the inclusion of the study (baseline) [22] and, from those, 107 patients developed T2DM (Incident-T2DM cases) after a mean follow-up of 60 months, according to all the American Diabetes Association (ADA) diagnosis criteria, evaluated on the basis of glucose tolerance tests performed each year (Supplementary Fig. 1). T2DM was diagnosed if one or more of following criteria were present in the study subjects: fasting plasma glucose (FPG) concentration ≥ 126 mg/dL, impaired fasting glucose (IFG); FPG ≥200 mg/dL after 2 h of oral glucose test (OGTT), impaired glucose tolerance (IGT); glycated hemoglobin (HbA1c) ≥6.5%. Specifically, among the 107 incident-T2DM patients, 43 subjects were diagnosed during the first year of follow-up, 24 during the second year, 11 during the third year, 19 during the fourth year and 10 during the fifth year (Supplementary Fig. 1). The remaining 355 subjects did not develop T2DM during the study period. In the present study, all the incident-T2DM cases (n = 107) were included together with 108 matched, randomly selected controls, who did not develop T2DM during the follow-up (non-T2DM subjects) (Supplementary Table 1). The random selection of the non-T2DM subjects was performed using computational stratified sampling from the 355 non-T2DM subjects of the CORDIOPREV study according to the following clinical, anthropometric and biochemical variables: diet, age, gender, fasting plasma glucose, body mass index, LDL-cholesterol and HDL-cholesterol. To implement this type of sampling, the target population was first divided into separate strata and then, samples were randomly selected within each stratum through simple automatic sampling by using the R software.

Metabolic study design

Oral glucose tolerance test (OGTT) and fat-rich meal tests were implemented in all patients to dynamically determine the metabolic status of the patient, as previously reported [21]. Briefly, OGTT (75 g dextrose monohydrate in 250 ml water, NUTER. TEC GLUCOSA, Subra, Toulouse, France) was started at 8:00 am, and plasma samples were collected at 0, 30, 60, 90 and 120 min to determine plasma glucose and insulin levels. Fat-rich meal test (0.7 g fat and 5 mg cholesterol per kg body weight with 12% saturated fatty acids, 10% polyunsaturated fatty acids, 43% monounsaturated fatty acids, 10% protein, and 25% carbohydrates) was performed at 8:00 am. Blood samples were collected before the meal and after 4 h. Biochemical determination of metabolic parameters and calculation of insulin resistance and sensitivity indexes were performed as previously reported [21]. Further details are provided in Supplementary Material and Methods.

Blood sampling and processing to isolate PBMCs

Venous blood from the participants (12 h overnight fast) was collected in tubes containing EDTA during the fat-rich meal test, at 0 and 4 h. PBMCs were isolated as previously described [21,23].

RNA extraction, quantification and reverse transcribed

Total RNA from PBMCs was isolated using Direct-zol RNA kit (Zymo Research, Irvine, CA, USA) following manufacturer's instructions. The amount of RNA recovered was determined and its quality assessed by the NanoDrop2000 spectrophotometer (Thermo Fisher). Specifically, all the RNA samples passed the quality controls, being the 260/280 and 230/260 absorbance ratios among 1.8–2.0. One μg of RNA was reverse transcribed (RT) to cDNA using random hexamer primers with the First Strand Synthesis Kit (Thermo Fisher).

In vitro culturing and treatment of PBMCs

PBMCs from n = 7 healthy subjects (Supplementary Table 2) were extracted as described above, and then cultured in serum-free RPMI medium (Lonza). 500,000 cells/well from each subject were seeded on ultra-low attachment multi 12-well plates (Corning Costar, Sigma) and treated per duplicate with 10% baseline fasting and postprandial serum derived from control and incident-T2DM patients (specifically, we used serum form individuals that developed T2DM during the first two years of follow-up). After the incubation periods (4 and 24 h), PBMCs were centrifuged and RNA was extracted and isolated using TRI-reagent (Sigma) [23].

Analysis of splicing machinery components by qPCR dynamic array based on microfluidic technology

A 48.48 Dynamic Array based on microfluidic technology (Fluidigm, San Francisco, CA, USA) was implemented to determine the expression of 48 transcripts in 48 samples, simultaneously. Specific primers for human transcripts including components of the major (n = 13) and minor spliceosome (n = 4), associated SFs (n = 28) and three housekeeping genes were specifically designed (Supplementary Table 3). The panel of splicing machinery components was selected on the basis of two main criteria: 1) the relevance of the given spliceosome components in the splicing process (such as the components of the spliceosome core) and 2) the demonstrated implication in the regulation of splicing variants implicated in the pathophysiology of T2DM (as is the case of the 28 splicing factors selected in this study).

Primers were selected using Primer3 software with selection parameters set to identify primer pairs that: 1) span an intron (when possible), 2) differ by no more that 1 °C in annealing temperature, 3) are at least 20 bp in length, 4) have a GC content between 45 and 55%, but 5) exclude primers that may form primer-dimers. Sequences of selected primers were used in BLAST (NCBI) searches to check for potential homology to sequences other than the designated target. Initial screening of primer efficiency using real-time detection was performed by amplifying 2-fold dilutions of RT products, where optimal efficiency was demonstrated by a difference of one cycle threshold between dilutions and a clear melting peak followed by a graded temperature-dependent dissociation to verify that only one product was amplified. The thermocycling profile consisted of: 1) 95 °C for 1 min; 2) 35 cycles of denaturing (95 °C for 5 s) and annealing/extension (60 °C for 20s); and 3) a last cycle where final PCR products were subjected to graded temperature-dependent dissociation (60 °C to 95 °C, increasing 1 °C/3 s). PCR products were then column-purified (FAVORGEN Biotech, Vienna, Austria) and sequenced to confirm target specificity. After confirmation of primer efficiency and specificity, the concentration of purified products was determined, and PCR products were serial diluted to obtain standards containing 1, 101,102, 103, 104, 105, and 106 copies of the synthetic template. Standards were then amplified by qPCR, and standard curves were generated using Stratagene Mx3000p software. The slope of a standard curve for each template examined was approximately −3.33 (R2 ≈ 1), indicating that the efficiency of amplification of our primers was 100%, meaning that all templates in each cycle were copied.

Preamplification, exonuclease treatment and qPCR dynamic array based on microfluidic technology were implemented following manufacturer's instructions using the Biomark System and the Real-Time PCR Analysis Software (Fluidigm). Additional details are provided in Supplementary Material and Methods. Finally, in the case of PBMC samples from non-T2DM and Incident-T2DM patients, the expression level of each transcript was adjusted by a normalization factor (NF) obtained from the expression levels of two different housekeeping genes [Beta-actin (ACTB) and Glyceraldehyde-3-Phosphate Dehydrogenase (GAPDH)] using Genorm 3.3 [24], while ACTB was the housekeeping gene used to normalize the expression level of each transcript in the in vitro experiment. This selection was based on the stability of these housekeeping genes among the experimental groups to be compared, wherein the expression of these housekeeping genes was not significantly different among groups.

Statistical and bioinformatical analysis

Data were evaluated for heterogeneity of variance using the Kolmogorov–Smirnov test and are expressed as mean ± SEM. Statistical analysis was carried out using unpaired t-test, Mann Whitney U test or one-way ANOVA followed by Dunnett's test. Significant correlations were studied using bivariate Spearman correlation methods. The odds ratio (OR) was calculated by logistic regression analysis by comparing T1 vs. T3. The subjects were categorized into tertiles of expression levels as low (T1), intermediate (T2) and high (T3). Predictive models were constructed by logistic regression (first with the SPSS software and later validated with the R language), Random Forest and C4.5 (both with R language) algorithms and followed by cross-validation studies by using the full-dataset of variables or a selection of them (obtained by feature ranking processes) as described in Supplementary Material and Methods. AUCs from ROC curves were compared by DeLong test [25]. P-values smaller than 0.05 were considered statistically significant. When appropriate, analyses were adjusted for age, diet, gender, BMI, glycated hemoglobin, HDL and triglycerides. Statistical analyses were carried out with GraphPad Prism 6 (La Jolla, CA, USA) and SPSS 17.0 (IBM).

Results

Description of the cohort

The group of patients who developed T2DM during the study (Incident-T2DM; n = 107) presented higher weight, BMI, HbA1c levels, fasting insulin and TGs at baseline compared to non-T2DM controls (n = 108) (Table 1). HOMA-IR, HIRI, DI and ISI were also significantly different at baseline between both groups (Table 1). Remarkably, the cohort of n = 108 randomly selected non-T2DM patients showed comparable levels of all the parameters determined except for age and DI compared to the total population (n = 355) of non-T2DM individuals (Supplementary Table 1).

Table 1

Baseline characteristics of subjects who did not develop T2DM (Non-T2DM) vs subjects who developed T2DM (Incident-T2DM) after a median follow-up of 5 years.

Variables	Non-T2DM	Incident-T2DM	p-value
n	108	107
Sex (male; female)	93; 15	89; 18	0.550
Diabetes family history n (%)	35 (32.4)	41 (38.3)	0.220
Age (years)	60.30 ±0.806	58.75 ±0.873	0.191
Waist circumference (cm)	102.48 ±0.958	105.28 ±1.08	0.053
Weight (kg)	81.92 ±1.194	85.69 ±1.47	0.047
BMI (kg/m2)	30.16 ±0.362	31.39 ±0.466	0.038
Glucose (mg/dl)	94.4 ±0.952	96.18 ±1.403	0.208
HbA1c (%)	5.88 ±0.032	6.03 ±0.033	0.001
Insulin (mU/l)	8.07 ±0.514	10.5 ±0.656	0.004
TG (mg/dl)	109.24 ±4.699	132.60 ±6.608	0.004
Total cholesterol (mg/dl)	159.55 ±3.027	164.97 ±3.409	0.235
c-LDL (mg/dl)	91.20 ±2.38	93.4 ±2.657	0.538
c-HDL (mg/dl)	44.58 ±0.899	43.52 ±1.039	0.440
NEFA (mmol/L)	0.286 ±0.015	0.317 ±0.016	0.174
Apo A1 (mg/dl)	133.5 ±2.093	135.15 ±2.312	0.596
Apo B (mg/dl)	71.57 ±1.934	76.22 ±1.835	0.083
hs-CRP (mg/dl)	2.428 ±0.32	2.878 ±0.292	0.300
HOMA-IR	2.5424 ±0.126	3.3734 ±0.302	0.012
HIRI	1024.55 ±50.85	1370.214 ±120.93	0.009
MISI (x102)	0.021 ±0.002	0.019 ±0.002	0.402
DI	0.8948 ±0.041	0.7685 ±0.041	0.030
ISI	4.0815 ±0.256	3.2758 ±0.186	0.012
IGI	1.0646 ±0.103	0.6633 ±0.294	0.200

Values expressed as mean ± SEM. BMI: Body mass index; HbA1c: Glycated hemoglobin; TG: Triglycerides; c- LDL: Low density lipoprotein; c-HDL: High density lipoprotein; NEFA: Non-esterified fatty acids; Apo A1: Apolipoprotein A1; Apo B: Apolipoprotein B; hs-CRP: High sensitivity C-reactive protein; HOMA-IR: Homeostasis model assessment- insulin resistance; HIRI: Hepatic insulin resistance index; MISI: Muscle insulin sensitivity index; DI: Disposition index; ISI: Insulin sensitivity index; IGI: Insulinogenic index.

Expression of splicing machinery components is different between incident-T2DM and non-T2DM patients

The expression pattern of several spliceosome components and SFs was dysregulated in the PBMCs of incident-T2DM compared to non-T2DM patients at baseline (Fig. 1a). Specifically, PBMCs of incident-T2DM patients exhibited significantly lower levels of RNU2, RNU4, RNU6ATAC, SNRNP200, ESRP1, SRSF1 and SRSF5 (Fig. 1b).

Fig. 1

Baseline expression pattern of the selected spliceosome components and splicing factors in the PBMCs of Incident-T2DM and control non-T2DM patients. a) Fold-change levels between Incident-T2DM and non-T2DM subjects, represented in red (increase) or blue (decrease). Specific spliceosome components or splicing factors significantly altered are highlighted in bold. b) mRNA expression levels [adjusted by a normalization factor (NF) calculated from the expression level of GAPDH and ACTB] of specific spliceosome components (first row) and splicing factors (second row) in the PBMCs from non-T2DM and Incident-T2DM subjects. Values represent the mean ± SEM. Asterisks indicate values that significantly differ from non-T2DM patients (t-test: *, p < .05).

Most notably, expression of many of these spliceosome components and SFs was dynamically, and differentially, regulated in the PBMCs of these individuals during a fat-rich meal test (a stress condition that challenges phenotypic flexibility, enabling its analysis) (Supplementary Fig. 2). Specifically, some spliceosome components and SFs were selectively altered during the post-prandial phase compared to the baseline state in a specific group of patients (e.g. RNU12 was increased in non-T2DM patients, whereas SNRNP70 was increased in incident-T2DM patients), leading to the appearance of significant differences between both groups of patients (Supplementary Fig. 2a). Interestingly, several of the changes observed at baseline (RNU2, RNU4, RNU6ATAC, ESRP1 and SRSF1) were also evident in the postprandial state, while others did not (SNRNP200 and SRSF5) (Supplementary Fig. 2B and C, respectively).

Comparison between the expression profiles of PBMCs from non-T2DM and incident-T2DM patients at post-prandial state revealed a drastic dysregulation of many spliceosome components and SFs after this meal challenge (Fig. 2a). Specifically, expression levels of RNU2, RNU4, RNU6, RNU4ATAC, RNU6ATAC, RNU12, NOVA1, ESRP1 and SRSF1 were lower, whereas those of SNRNP70 were higher in PBMCs from incident-T2DM patients vs. non-T2DM during the post-prandial phase (Fig. 2b). Remarkably, these fasting and postprandial changes observed between non-T2DM and incident-T2DM patients were also observed when only considering the male population, which were age and BMI matched (data not shown).

Fig. 2

Baseline expression pattern of the selected spliceosome components and splicing factors in the PBMCs of incident-T2DM and control non-T2DM patients during the postprandial state. a) Fold-change levels between Incident-T2DM and non-T2DM subjects, represented in red (increase) or blue (decrease). Specific spliceosome components or splicing factors significantly altered are highlighted in bold. b) mRNA expression levels [adjusted by a normalization factor (NF) calculated from the expression level of GAPDH and ACTB] of specific spliceosome components and SFs in the PBMCs of non-T2DM and Incident-T2DM subjects. Values represent the mean ± SEM. Asterisks indicate values that significantly differ from non-T2DM patients (t-test: *, p < .05; **, p < .01; ***, p < .001).

Expression of splicing machinery components is different according to the year of T2DM diagnosis

Subjects diagnosed of T2DM in the first 2-years of follow-up exhibited more pronounced changes in the expression pattern of spliceosome components and SFs as compared to non-T2DM controls and subjects diagnosed in subsequent years, at both baseline fasting and post-prandial states (Supplementary Figs. 3 and 4). Specifically, fasting RNU4, and post-prandial RNU4, RNU4ATAC and RNU6ATAC expression levels were altered in patients who developed T2DM during the first and second years of follow-up compared to non-T2DM patients (Supplementary Figs. 3 and 4). In contrast, fasting RNU2, SNRNP200, RNU6ATAC and SRSF5 and post-prandial RNU2, RNU6, NOVA1, ESRP1 and SRSF1 expression levels were lower only in patients that developed T2DM during the first year of follow-up compared to non-T2DM patients (Supplementary Figs. 3 and 4). Also, post-prandial expression levels of RNU1 and RNU12 were significantly altered during the second year of follow-up in patients that developed T2DM compared to non-T2DM patients (Supplementary Fig. 4).

Baseline expression of splicing machinery components was correlated with key clinical parameters

Baseline expression of several spliceosome components and SFs was correlated with clinically relevant parameters, thus suggesting a potential pathophysiological role. Although the baseline fasting or postprandial expression levels of the splicing components were not associated with BMI, baseline fasting levels of ESRP1, SRSF1 and RNU4 was directly correlated with fasting c-HDL and NEFA levels. Furthermore, fasting levels of SRSF1 and RNU4 were directly correlated with postprandial NEFA levels and fasting levels of SNRNP200 were directly correlated with postprandial glucose. Fasting SRSF5 expression levels were directly correlated with fasting and postprandial c-HDL, postprandial Apo A1 and ISI values, and inversely correlated with fasting HbA1c and postprandial triglycerides and C-reactive protein. Additionally, baseline postprandial expression levels of ESRP1, RNU4, RNU6 and NOVA1 were inversely correlated with HbA1c levels, while postprandial levels of RNU2 were directly correlated with DI and postprandial RNU12 levels inversely correlated with ISI. Finally, postprandial SRSF1 levels were directly correlated with postprandial NEFA levels (Supplementary Table 4).

Expression of some splicing machinery components was associated to the risk of developing T2DM

The possible association between the expression of the components of the splicing machinery and the risk of developing T2DM was assessed by logistic regression analysis (odds ratio, OR). This analysis revealed that low baseline fasting expression levels of RNU4 and RNU2, as well as low baseline postprandial levels of RNU4ATAC, NOVA1, RNU6ATAC, RNU6, RNU12, SRSF1 and RNU4 were strongly associated to the relative risk of T2DM development when adjusting by age, gender, diet, BMI, HbA1c, c-HDL and TG levels (Table 2).

Table 2

Association between the PBMC expression of the components of the splicing machinery and the relative risk of developing T2DM by logistic regression analysis of relative risk (odds ratio, OR).

	OR	95% C.I.		p-value
		Inferior	Superior
Fasting RNU4	2.521	1.117	5.688	0.026
Fasting RNU2	2.283	1.012	5.153	0.047
Postprandial RNU4ATAC	4.456	1.821	10.903	0.001
Postprandial NOVA1	4.099	1.836	9.154	0.001
Postprandial RNU6ATAC	3.762	1.706	8.298	0.001
Postprandial RNU6	3.762	1.706	8.298	0.001
Postprandial RNU12	2.274	1.041	4.967	0.039
Postprandial SRSF1	2.204	1.007	4.825	0.048
Postprandial RNU4	2.109	0.963	4.619	0.062

Subjects were categorized in tertiles according to the expression level of each spliceosome component or splicing factor as follows: low expression levels (T1), medium expression levels (T2) and high expression levels (T3), and the OR estimated between T1 and T3 for each element of interest. OR: Odds Ratio; C.I.: Confidence intervals.

Alterations in the expression of the splicing machinery components predict the development of T2DM

In our cohort, incident-T2DM and non-T2DM patients exhibited different baseline characteristics (Table 1); however, HbA1c and/or FINDRISK, the classic predictors of T2DM development, exhibited low predictive potential (HbA1c: AUC = 0.643, p < .001, TPR = 0.669, TNR = 0.565; FINDRISK: AUC = 0.548, p = .231, TPR = 0.622, TNR = 0.393; HbA1c and FINDRISK: AUC = 0.643, p < .001, TPR = 0.651, TNR = 0.606) (Table 3). In this sense, although the expression levels of spliceosome-associated elements altered in PBMCs of incident-T2DM compared to non-T2DM patients generated significant but low (<0.65) AUCs in ROC curves (Supplementary Table 5), the combination of all these elements led to more significant and accurate ROC curves. Specifically, the fingerprints comprised by the baseline expression of spliceosome components and SFs during fasting, post-prandial or their combination generated significant ROC curves with AUC > 0.85 using different algorithms (Table 3 and Supplementary Table 6), which was validated (AUCs ranging 0.65–0.81) by cross-validation analysis (Supplementary Table 6). In particular, the fingerprint comprised by the fasting and postprandial expression of spliceosome components and SFs presented an AUC = 0.813 (TPR = 0.802, TNR = 0.689) in the cross-validation analysis (Fig. 3a – green line), which was significantly higher than the capacity of the clinically-relevant HbA1c and FINDRISK to predict T2DM development (Fig. 3a), as demonstrated by DeLong test comparing the AUCs (p < .05 vs. HbA1c or FINDRISK).

Table 3

Capacity of the molecular fingerprint comprised by baseline fasting and/or postprandial levels of spliceosome components and splicing factors as T2DM predictive models by logistic regression and ROC curve analysis.

Model	AUC	p-value
Splicing machinery components baseline expression during fasting	0.894	0.000
Splicing machinery components baseline expression during postprandial phase	0.853	0.000
Splicing machinery components baseline expression during fasting and postprandial phase	1	0.000
HbA1c	0.643	0.000
FINDRISK	0.548	0.231
HbA1c + FINDRISK	0.643	0.000
Splicing machinery components baseline expression during fasting + HbA1c	0.898	0.000
Splicing machinery components baseline expression during postprandial phase + HbA1c	0.867	0.000
Splicing machinery components baseline expression during fasting and postprandial phase + HbA1c	1	0.000
Splicing machinery components baseline expression during fasting + FINDRISK	0.895	0.000
Splicing machinery components baseline expression during postprandial phase + FINDRISK	0.856	0.000
Splicing machinery components baseline expression during fasting and postprandial phase + FINDRISK	1	0.000

Logistic regression models considering the fasting and/or postprandial baseline levels of the measured elements, alone or in commination with the classic predictors of T2DM (HbA1c and FINDRISK). AUC: Area under curve; HbA1c: Glycated hemoglobin.

Fig. 3

Spliceosome components and splicing factors-based predictive models generated by Random Forest computational algorithm and ROC curve analysis. ROC curves parameters were calculated for the predictive models generated by Random Forest algorithm considering the expression of all the splicing machinery elements determined at fasting (orange), postprandial (blue) or their combination (green) in non-T2DM and Incident-T2DM patients (a) or using a selection of the most relevant and discriminatory splicing machinery components (b). Specifically, the subset of specific splicing machinery components were RNU4ATAC, TIA1, KHDRBS1, SRSF10, PTBP1, RAVER1, RNU2, RNU5, SRSF9, U2AF2, RBM45, SRSF4, RBM3 for baseline fasting, RNU4ATAC, RNU6ATAC, RAVER1, SF3B1, SRSF3, NOVA1, SRM160, SRSF6, ESRP1, U2AF1 for baseline postprandial and fasting SRSF3, SRSF10, SRSF6, TRA2B, PTBP1, SF3B1 and postprandial RNU4ATAC, RAVER1, RNU4ATAC, NOVA1, RNU4, RNU6, PRPF8 for combined analysis. The same ROC curves were calculated considering the patients that developed T2DM during the first two years of follow-up (c and d). In this case, the subset of specific splicing machinery components were RNUU4ATAC, PTBP1, TRA2A, RM17, RNU12, TIA1, SRSF5, RNU2 for baseline fasting, RNU4ATA, PTBP1, MAGOH, SRSF9, RAVER1, PRP8, PRPF40A, SRRM1, SRSF6, SNRNP200, TIA1, RNU2 for baseline postprandial, and baseline fasting RNU4ATAC, TRA2A, SRSF5, RBM17, SRSF10, SRSF3, RNU2 and postprandial RNU4ATAC, RNU4, RNU6ATAC, MAGOH, RAVER1, PRPF40A, RBM4, U2AF2, SRSF10, RNU11 TRA2B, SND1 for combined analysis. ROC curves of HbA1c and FINDRISK were also estimated. HbA1c: glycosylated hemoglobin; AUC: Area under the curve; TPR: True positive ratio; TNR: True negative ratio; p: p value.

More remarkably, a subset of splicing machinery elements selected by computational approaches (feature ranking analysis) exhibited even higher predictive capacity compared to the whole dataset. Indeed, a predictive model comprised by fasting expression levels of RNU4ATAC, SRSF3, SRSF6, SRSF10, TRA2B, PTBP1, SF3B1 and post-prandial levels of RNU4, RNU6, RNU4ATAC, RAVER1, NOVA1, and PRPF8 exhibited AUC = 0.881 (TPR = 0.801, TNR = 0.700) in cross-validation analysis (Fig. 3b – green line), which was validated by different modeling methods, clearly outperformed the capacity of the classic predictors of T2DM development (p < .001 by DeLong test vs. HbA1c or FINDRISK AUCs) (Fig. 3 b).

Finally, the capacity of splicing machinery alterations to predict T2DM development was even superior when considering only the patients that developed T2DM during the first 2 years of follow-up. Indeed, the fasting and/or postprandial baseline fingerprints comprised by spliceosome components and SFs generated ROC curves with AUC ranging 0.630–0.851 by cross-validation analysis (Fig. 3c and Supplementary Table 6), which were again clearly higher than those generated by HbA1c or FINDRISK (p < .05 by DeLong test in all cases). Similarly, models generated after selection of the most relevant elements exhibited AUCs>0.8 when predicting T2DM development in the first two years of follow-up (Fig. 3d), which were clearly higher than those generated by HbA1c or FINDRISK (p < .001 by DeLong test in all cases).

Incident-T2DM patient serum altered spliceosome components expression in PBMCs of healthy patients

To test whether changes in the expression of spliceosome components and SFs in PBMCs could be caused by the metabolic dysregulation occurring in incident-T2DM patients, PBMCs from healthy subjects (Supplementary Table 1) were incubated with serum from non-T2DM and incident-T2DM patients. Remarkably, 24 h incubation with baseline fasting serum from Incident-T2DM patients markedly reduced RNU4 expression compared to the serum from non-T2DM patients (Fig. 4a); whereas expression of other spliceosome components and SFs altered at fasting baseline was not significantly modulated (Supplementary Fig. 5). More strikingly, 4 h incubation with baseline postprandial serum from Incident-T2DM patients induced a significant reduction of RNU4 and RNU12 expression and a non-significant reduction trend of RNU2 (p = .09) compared to non-T2DM treated samples (Fig. 4b). It should be noted that PBMCs survival was evaluated to assess the possible effect of human sera and that after 24 h of culture with baseline fasting and postprandial serum derived from control and incident-T2DM patients, the survival of the PBMCs was minimally affected (>95% of cell survival in all cases) by the different sera used and that there were not significant changes among the different types of serum (data not shown).

Fig. 4

Expression levels of key spliceosome components in PBMCs from healthy subjects treated with non-T2DM and Incident-T2DM derived baseline fasting during 24 h (a) or postprandial serum during 4 h (b). mRNA expression levels, adjusted by ACTB expression levels, of specific spliceosome components. Values represent the mean ± SEM. Asterisks indicate values that significantly differ from PBMCs treated with non-T2DM derived serum (t-test: *, p < .05; ***, p < .001).

Discussion

In this study, we analyzed the expression pattern of selected splicing machinery elements in PBMCs from non-T2DM patients with CVD and, therefore, with high risk to develop T2DM, who were followed-up during 5-years (individuals that suffered a cardiovascular event and were included in the CORDIOPREV study [21]). This revealed, for the first time, that PMBCs expression levels of certain splicing machinery components were significantly altered in Incident-T2DM at the inclusion of the study (before the development of T2DM) and associated with T2DM development. Indeed, the molecular fingerprints comprised by the fasting and, especially, by the postprandial levels of certain of these splicing machinery components were able to differentiate between patients who subsequently developed T2DM (Incident-T2DM) from those patients who did not develop the disease with an AUC > 0.8 after cross-validation analysis, which is significantly higher than the current standard procedures (HbA1c and FINDRISK exhibited AUCs<0.66 in our cohort). In addition, patients with low PBMCs expression levels of specific splicing machinery components, including RNU2, RNU4 or RNU12, which we have shown to be modulated by the serum of Incident-T2DM patients, were at higher risk (OR > 4) of T2DM development as compared with those patients with high PBMC levels of these components.

In line with our findings, previous reports have indicated that under adverse metabolic conditions (e.g. obesity, insulin resistance, etc.) the splicing machinery is markedly altered in most tissues [[10], [11], [12]] and associated with the development of several pathologies [10,13,14]. Actually, alternative splicing seems to reside at the crossroad between hyperinsulinemia, insulin resistance, obesity and T2DM [11,12,17], and, consequently, the correct function of the splicing machinery (spliceosome components and SFs) is essential to maintain whole body homeostasis [18]. However, we provide herein novel, primary evidence demonstrating that the alteration in the expression of certain spliceosome components and SFs precedes the development of T2DM in CVD patients and, most notably, our results provide the first indication that the expression of certain spliceosome components and SFs in PBMCs can be dynamically modulated, according to the metabolic status of the individual, and specially under stress conditions as those observed during the post-prandial response, thus suggesting an adaptive capacity of the splicing machinery in order to respond to metabolic disturbances.

In particular, this study demonstrated that the expression of some splicing machinery components is altered in control individuals during the postprandial phase, which is consistent with previous results demonstrating a regulatory response of specific splicing variants to the postprandial environment [26,27]. Therefore, these data suggest that the changes observed in the splicing machinery during the postprandial phase may be responsible for the regulation of the expression of particular splicing variants under these conditions [[12], [19], [28]], which could be essential to the appropriate response of the organism to metabolic challenges and disturbances. Moreover, our results also demonstrate for the first time that the response of key splicing machinery components to metabolic insults is altered in individuals who will develop T2DM (Incident-T2DM patients), but especially in those developing T2DM during the first two years of follow-up. Therefore, since postprandial alterations are closely related to the phenotypic flexibility, which is strongly linked to T2DM development [9], our data primarily demonstrate that the alteration in the splicing machinery precedes the instauration of T2DM, thereby suggesting its putative implication as a driving force in the development of this pathology. Based on all the information mentioned above, it is tempting to propose that the splicing machinery could be acting as a biosensor of the whole body metabolism to adapt cell gene expression to the pathophysiological conditions, and that its dysregulation could lead to an unbalance in the landscape of splicing variants present in a given cell at a given moment [12,28,29], which may be associated to the instauration of T2DM [10,30]. This idea is further supported by two pieces of evidence presented herein. First, we have demonstrated that low fasting or postprandial expression levels of certain splicing machinery components drastically increase the relative risk of T2DM development in these patients, suggesting that a dysregulated expression of certain splicing machinery components could augment the risk of developing T2DM. Secondly, our in vitro studies demonstrate that the expression of relevant spliceosome components, specially RNU2, RNU4 or RNU12, which are key elements responsible for the appropriate function of the spliceosome [31], can be modulated by baseline fasting and/or postprandial serum from Incident-T2DM patients. This observation might suggest the existence of specific factors in the serum of these patients capable to modulate the expression of relevant spliceosome components and, therefore, the function of the splicing machinery. In this sense, previous studies have found a relationship between circulating factors and the modulation of SFs in different tissues. For example, it has been described that insulin signaling can up-regulate the expression of the splicing factor SRSF1 in pancreatic beta cells, inducing the splicing of the insulin receptor to generate the INSR-B isoform [32]. The same study also found a regulation of the protein levels of the splicing factor MBNL1 by high glucose levels. Furthermore, the splicing of the Fatty acid desaturase 3 has been observed to be modulated in the liver of baboons in response to different diets and in human liver HepG2 cells after treatment with polyunsaturated fatty acids [33]. Obviously, further studies will be required to attain a more comprehensive understanding of the changes in the splicing process and the contribution of the dysregulation in the splicing machinery to the generation of alternative spliced isoforms bearing pathological implications; as well, additional, and more detailed and refined experiments will be needed to assess the specific contribution of each blood cell subset to the effects observed herein. Ultimately, all of this information would help to elucidate the nature of the elements causing the changes in the splicing machinery.

In conclusion, although this study was implemented in CVD patients from the CORDIOPREV study [21], which limits our findings to people with these characteristics, and precludes generalization to healthy people, the data presented herein strongly support the notion that altered expression of spliceosome components and SFs variants may be associated with the development of T2DM, preceding the clinical instauration of this pathology and, therefore, could likely serve as a sensor and early predictor for T2DM development in CVD patients. Certainly, our data demonstrate the existence of a spliceosome-associated molecular fingerprint capable to predict the future development of T2DM in individual patients with high precision (AUC = 0.881, TPR = 0.801, TNR = 0.700), which even outperforms the capacity of classical predictors of T2DM development, such as HbA1c or FINDRISK. Therefore, this splicing machinery-associated molecular fingerprint could become a valuable, non-invasive, new tool for early risk assessment of T2DM in clinical practice to prevent disease development.

Role of the funding sources

This work was funded by Instituto de Salud Carlos III, co-funded by European Union (ERDF/ESF, “Investing in your future”) [PIE14/00005, PIE 14/00031, PI16/00264, CP15/00156, PI17/002287], MINECO (BFU2016-80360-R, TIN2017–83445-P, PI13/00023, AGL2012/39615, AGL2015-67896-P), Junta de Andalucía (BIO-0139, CTS-1406, CTS-525, PI-0541-2013, CVI-7450) and CIBERobn. CIBER is an initiative of Instituto de Salud Carlos III, Ministerio de Sanidad, Servicios Sociales e Igualdad, Spain. The CORDIOPREV study is supported by the Fundación Patrimonio Comunal Olivarero, Junta de Andalucía (Consejería de Salud, Consejería de Agricultura y Pesca, Consejería de Innovación, Ciencia y Empresa), Diputaciones de Jaén y Córdoba, Centro de Excelencia en Investigación sobre Aceite de Oliva y Salud and Ministerio de Medio Ambiente, Medio Rural y Marino, Gobierno de España.

The following funding sources were essential for data collection, analysis, interpretation, and patient recruitment: PIE14/00005, PIE 14/00031, TIN2017-83445-P, PI13/00023, AGL2012/39615, AGL2015–67896-P, CTS-525, CVI-7450.

The following funding sources were essential for data collection, analysis, and interpretation: CP15/00156, PI16/00264, BFU2016–80360-R, BIO-0139, CTS-1406, PI17/002287, PI-0541-2013.

The corresponding authors declare that they had full access to all the data in the study and had the final responsibility for the decision to submit for publication.

Conflicts of interest

The authors declare that they have no conflict of interest.

Authors contributions

Manuel D. Gahete contributed to the conception and design of the work, literature search, acquisition, analysis, and interpretation of the data, drafting and revising the work and final approval of the manuscript.

Mercedes del Rio-Moreno contributed to the design of the work, literature search, acquisition, analysis, and interpretation of the data, drafting and revising the work and final approval of the manuscript.

Antonio Camargo contributed to the conception and design of the work, analysis and interpretation of the data, drafting and revising the work and final approval of the manuscript.

Juan F Alcala contributed to the conception and design of the work, acquisition and analysis of the data and final approval of the manuscript.

Emilia Alors-Perez contributed to the acquisition and analysis of the data, drafting the work, and final approval of the manuscript.

Javier Delgado-Lista contributed to the conception and design of the work, acquisition of the data and final approval of the manuscript.

Oscar Reyes contributed to the analysis and interpretation of the data, drafting the work and final approval of the manuscript.

Sebastian Ventura contributed to the analysis and interpretation of the data, revising the work, and final approval of the manuscript.

Pablo Perez-Martinez contributed to the conception and design of the work, literature search and interpretation of the data; revising the work and final approval of the manuscript.

Justo P. Castaño contributed to the conception and design of the work, literature search and interpretation of the data, revising the work and final approval of the manuscript.

Jose Lopez-Miranda contributed to the conception and design of the work, literature search and interpretation of the data, revising the work and final approval of the manuscript.

Raúl M. Luque contributed to the conception and design of the work, literature search and interpretation of the data, drafting and revising the work and final approval of the manuscript.