Novel Interactions Between Circulating microRNAs and Gut Microbiota Composition in Human Obesity.

DOI: https://doi.org/10.21203/rs.3.rs-66883/v1

Abstract

Background: Unbalances in microRNAs (miRNA) and gut microbiota patterns have been proposed as putative factors concerning onset and development of obesity and other metabolic diseases. However, the determinants that mediate the interactions between miRNAs and the gut microbiome impacting on obesity are scarcely understood. Thus, the aim of this article was to investigate possible interactions between circulating miRNAs and gut microbiota composition in obesity.

Method: The analyzed sample comprised 78 subjects with obesity [cases, body mass index (BMI): 30 – 40 kg/m2] and 25 eutrophic individuals (controls, BMI £ 25 kg/m2). The expression of 96 miRNAs was investigated in plasma of all individuals using miRCURY LNA miRNA Custom PCR Panels (Exiqon). Bacterial DNA sequencing was performed following the Illumina 16S protocol. The FDR (Benjamini-Hochberg test, q-value) correction was used for multiple comparison analyses.

Results: A total of 26 circulating miRNAs and 12 bacterial species were found differentially expressed between cases and controls. Interestingly, an interaction among three miRNAs (miR-130b-3p, miR-185-5p, and miR-21-5p) with Bacteroides eggerthi, and BMI levels was evidenced (r2= 0.148, P= 0.004). Those miRNAs that correlated with obesity-associated gut bacteria abundance are known to regulate target genes that participate in metabolism-related pathways, such as fatty acid degradation, carbohydrate digestion and absorption, insulin signaling, and glycerolipid metabolism.

Conclusion: This study characterized an interaction between the abundance of 4 bacterial species and 14 circulating miRNAs in relation to body adiposity. Moreover, the current study also suggests that miRNAs may serve as a communication mechanism between the gut microbiome and human hosts.

Clinical trial registration: clinicaltrials.gov (reg. no. NCT02737267).

Background

Obesity is a worldwide epidemic that arises as a chronic long-term imbalance between calorie intake and energy expenditure (Campión et al., 2009). Despite nutritional interventions and physical education programs, the prevalence of obesity is still increasing and ∼600 million people worldwide are expected to be obese by 2025 (Organization., 2018). Obesity is the result of complex and not completely understood pathological processes deriving from a crosstalk among environmental factors, genetic susceptibility, and epigenetic mechanisms (Guyenet and Schwartz, 2012; van Dijk et al., 2015).

Among the epigenetic mechanisms, microRNAs (miRNAs) are a class of small non-coding RNAs that regulate gene expression (Bartel, 2004; Esteller, 2011; Butz et al., 2016). These molecules have recognized roles in the regulation of several biological processes, including cell cycle, as well as such as cellular differentiation, proliferation, metabolism, ageing, and apoptosis (Bartel, 2004). Moreover, it is estimated that miRNAs regulate the expression of more than 60% of protein-coding genes (Esteller, 2011); and, consequently, changes in their expressions and functions have been linked to many diseases, including metabolic disorders and obesity (Maurizi et al., 2018; Lorente-Cebrián et al., 2019).

Recent findings indicate that host miRNAs contribute to the regulation of the gut microbiome, specially involving at least two main processes: (i) host-secreted miRNAs regulate the gut microbiota; and (ii) the gut microbiota affects the host via inducing special miRNA expression (Belcheva, 2017). Indeed, evidences suggest that miRNAs produced by the host’s intestinal epithelial cells (IECs) participate in shaping the gut microbiota and affect bacterial growth (Liu and Weiner, 2016). These miRNAs target bacterial mRNA, and then the host controls the gut microbiota via bacterial mRNA degradation or translational inhibition (Liu and Weiner, 2016). On the other hand, it was demonstrated, using Dicer1 knock-out mice, that miRNAs were essential for epithelial cell proliferation, differentiation, nutrient absorption, and that defective miRNA biogenesis was also responsible for impaired intestinal barrier function (McKenna et al., 2010).

Additionally, the gut microbiota regulates miRNA expression in IECs subtypes, and this regulation may alter intestinal homeostasis (Nakata et al., 2017). In this sense, it was demonstrated that the expression of some miRNAs is different among IEC subtypes and the difference depends on microbial patterns (Peck et al., 2017). Similarly, the expression of 16 miRNAs was found to be altered in the caecum of conventionally raised versus germ-free mice (Singh et al., 2012). Recently, it has been reported that the gut microbiota specifically controlled the miR-181 family expression in white adipocytes during homeostasis to modulate key pathways related to adiposity, insulin sensitivity, and white adipose tissue (WAT) inflammation in mice (Virtue et al., 2019). Furthermore, high-fat diet (HFD) feeding altered the composition of the gut microbiota, leading to aberrant overexpression of miR-181 in WAT adipocytes (Virtue et al., 2019). Altogether, these studies provide clues that gut microbiota regulates host gene expression through modulation of the host miRNA signature, and that host metabolism could be influenced by this interaction.

According to these findings, miRNAs appear to play an important role in host-to-microbe interaction and may be considered as molecular targets for novel anti-microbial therapies to be developed. However, very little is known about the interactions between miRNAs and the host microbiome in the context of obesity. Therefore, the aim of this study was to investigate interactions between circulating miRNA patterns and gut microbiota composition in obesity.

Methods

Study population

This study was designed in accordance with STROBE guidelines for reporting association studies (von Elm et al., 2008; Vandenbroucke et al., 2014). The sample comprised 78 subjects with obesity [cases, body mass index (BMI): 30–40 kg/m2] and 25 eutrophic individuals (controls, BMI ≤ 25 kg/m2). Obesity was classified following the World Health Organization (WHO) guidelines (Organization., 2018). All volunteers were enrolled between October 2015 and February 2016 in the Metabolic Unit of the Centre for Nutrition Research of the University of Navarra, Spain. Major exclusion criteria included a history of diabetes mellitus (DM), cardiovascular disease or hypertension, pregnant or lactating women, current use of lipid-lowering drugs or medications that affect body weight, and weight change ≥3 kg within three months before the recruitment.

This study followed the ethical principles for medical research in humans from the Helsinki Declaration (Association, 2013). Moreover, the research protocol was properly approved by the Research Ethics Committee of the University of Navarra (ref. 132/2015) and it is registered at clinicaltrials.gov (reg. no. NCT02737267). A written informed consent of each participant was obtained prior to enrollment in the study.

All patients underwent anthropometric and laboratory evaluations, as previously described (Lopez-Legarrea et al., 2013; Ramos-Lopez et al., 2018). The measurements of height (cm), body weight (kg), and waist circumference (WC, cm) were collected in the fasting state by trained nutritionists following validated procedures (Lopez-Legarrea et al., 2013). BMI was calculated as the ratio between weight and squared height (kg/m2). Body composition was quantified by dual-energy X-ray absorptiometry according to instructions provided by the supplier (Lunar Prodigy, software version 6.0, Madison, WI, USA). Biochemical measurements including fasting plasma glucose (FPG, mg/dl), total cholesterol (TC, mg/dl), high-density lipoprotein cholesterol (HDL-c, mg/dl), and triglycerides (TG, mg/dl) were determined in an automatic analyzer (Pentra C200, HORIBA Medical), following standardized procedures. Endocrine markers such as insulin, adiponectin, and leptin were quantified with commercial ELISA kits (Mercodia Insulin ELISA, Biovendor Human adiponectin ELISA, and Mercodia Leptin).

Insulin resistance was estimated by the homeostatic model assessment-insulin resistance (HOMA-IR) index according to the following formula: (fasting insulin (mU/L) × plasma glucose (mmol/L)/22.5), while the triglyceride-glucose (TyG) index was calculated as: ln [fasting triglycerides (mg/dl) × fasting plasma glucose (mg/dl)/2], as described elsewhere (Navarro-González et al., 2016).

A validated semiquantitative food frequency questionnaire was used to evaluate habitual consumption (daily, weekly, monthly, or never) of 137 foods during the previous year (de la Fuente-Arrillaga et al., 2010). Energy and nutrient intakes were further calculated with an ad hoc computer program based on the standard Spanish food composition tables (Moreiras, 2018). The physical activity level was estimated using a validated questionnaire (Martínez-González et al., 2005). The volume of activity was expressed in metabolic equivalents (METs, kcal/kg/h), as described elsewhere (Basterra-Gortari et al., 2009).

MicroRNA expression analysis

miRNA isolation and reverse transcription-quantitative PCR

Total RNA was extracted from 200µ l EDTA-plasma using the miRNeasy Serum/Plasma Advanced kit (Qiagen. Hilden, Germany) according to the manufacturer’s recommendations. RNA spike-in was added to each sample (RNA Spike-In Kit, Qiagen. Hilden, Germany). The purity and concentration of RNA samples were measured using the NanoDrop ND-1000 Spectrophotometer (Thermo Fisher Scientific. Massachusetts, USA) (Bustin et al., 2009).

Relative expression of the 86 miRNAs was analyzed in plasma from all subjects using the Custom Pick-&-Mix microRNA PCR Panel v5 (Qiagen. Hilden, Germany). Moreover, 9 controls (reference genes + Spike-in controls) and a blank were also included in each plate, as shown in Supplementary Table 1. The selection of these miRNAs was based on the available literature (Vienberg et al., 2017; Lorente-Cebrián et al., 2019) and by searching on the miRWalk 2.0 database (Dweep et al., 2011) those miRNAs potentially associated with obesity in humans.

Total RNA (4µ l) was reverse-transcribed in 10µ l reaction using the miRCURY LNA Universal RT microRNA PCR, Polyadenylation, and cDNA synthesis kit II (Qiagen, Vedbaek, Denmark). RT-qPCR experiments were performed using a CFX384 Real-time system (Bio-Rad. USA). The cycling conditions were used: 95 °C for 2 min, followed by 40 cycles at 95 °C for 10 s and 56 °C for 1 min. Relative expressions were calculated using the 2ΔΔCq method (Rao et al., 2013).

Quality control and normalization

Quality control was carried out using synthetic spike-in RNAs to analyze the robustness of the RNA isolation process and quality of extracted miRNA. The RNA isolation controls (UniSp2, UniSp4, and UniSp5; Qiagen. Denmark) were added to the thawed plasma before the isolation process, aiming to detect differences in the extraction efficiency. The cDNA synthesis control (UniSp6, Qiagen) and cel-miR-39-3p were added to the reverse transcription reaction to determine the effectiveness of this process. Furthermore, UniSp3 was included in all plates and used as an inter-plate calibrator and PCR amplification control.

Hemolysis was assessed by the ratio between hsa-miR-451a (which is expressed in erythrocytes) and hsa-miR-23a-3p (which is relatively stable in plasma and not affected by hemolysis) as described elsewhere (Blondal et al., 2013). The difference in expression values between these 2 miRNAs provides a good measure of the extent of hemolysis, with values > 5 suggesting erythrocyte miRNA contamination. Only samples without hemolysis (values < 5) were included in the study. (Blondal et al., 2013). The assay cut-off was 35 cycles, and miRNAs expressed in at least 20% of the total sample (Gevaert et al., 2018). All individual samples were run on a predefined assay panel of 96 specific human miRNAs (Supplementary table 1). The miRNAs with complete data were used for the global mean method for normalization of the data, since this approach was found to be the most stable normalizer (Mestdagh et al., 2009).

miRNA target prediction and pathway enrichment analysis

Potential targets of selected miRNAs were searched using miRWalk 3.0 (http://zmf.umm.uni-heidelberg.de/apps/zmf/mirwalk2/, accessed 4th August 2020). To better understand the biological relevance of the miRNAs target genes, a network analysis was executed using PathDIP (accessed 4th August 2020, (Rahmati et al., 2017)). A hypergeometric test was used to calculate the statistical significance of the enriched pathways, and P-values were corrected for multiple tests using the Benjamini-Hochberg procedure, which provides a False Discovery Rate (FDR) adjusted-P-value (q-value). Pathways associated with a q-value < 0.05 were considered significantly enriched.

Gut Microbiota Analysis

Fecal sample collection and DNA extraction

Volunteers self-collected fecal samples using OMNIgene•GUT kits from DNA Genotek (Ottawa, Canada), according to the standardized instructions provided by the fabricant. The isolation of DNA from fecal samples was performed with the QIAamp® DNA kit (Qiagen. Hilden, Germany) according to the manufacturer’s protocol.

16S rRNA sequencing and sequence analysis

Bacterial DNA sequencing was performed by the Servei de Genòmica i Bioinformàtica from the Autonomous University of Barcelona (Barcelona, Spain). The Illumina 16S protocol, based on the amplification of the V3-V4 variable regions of the 16S rRNA gene, was followed. Paired-end sequencing was performed on the MiSeq System (Illumina. California, USA). In the process, two PCR reactions were carried out. In the first one, 12.5 ng of genomic DNA and the 16S-F and 16S-R primers were used. After the cleaning process, 5 µl of the first PCR was used in the second PCR. The primers used in this PCR were part of the Nextera XT DNA Index Kit (96 indexes, 384 samples) FC-131-1002 (Illumina. California, USA). After each PCR, the quality of the process was checked in a Labchip Bioanalyzer (Agilent Technologies Inc. Delaware, USA). Once all the samples were obtained, up to 40 samples were multiplexed in each run of 2 × 300 cycles for obtaining around 500,000 reads per sample. For this purpose, equimolar concentrations of each of the samples were mixed and the pool diluted up to 20 pM. A total of 3 runs were performed on the MiSeq sequencer with the MiSeq Reagent Kit v3 (600 cycle) MS-102-3003. The maximum of reads obtained was 1,867,496 and the minimum 5,019; the mean was 164,387.

The 16S rRNA sequences were filtered following quality criteria of the OTU processing pipeline LotuS (release 1.58) (Hildebrand et al., 2014). This pipeline includes UPARSE de novo sequence clustering (Edgar, 2013) and removal of chimeric sequences and phix contaminants for the identification of Operational Taxonomic Units (OTUs) and their abundance matrix generation (Rideout et al., 2014; Pichler et al., 2018). OTU refers to organisms clustered by similarities in DNA sequence (Hao et al., 2011). Finally, taxonomy was assigned using BLAST (Altschul et al., 1990) and HITdb (Ritari et al., 2015) achieving up to species-level sensitivity. The abundance matrices were first filtered and then normalized in R/Bioconductor (Lawrence et al., 2013; Pasolli et al., 2017) at each classification level: OTUs, phylum, genus, family, order, class, and species. Briefly, taxa less than 10% of frequency in our population were removed for the analysis and a global normalization was performed using the library size as a correcting factor and log2 data transformation.

To evaluate alpha diversity, the Shannon index was calculated (Shannon, 1997). To assess beta diversity, permutational multivariate analysis of variance (PERMANOVA) was used to analyze whether the structures of gut microbiota were significantly different among groups based on the Jaccard and Bray–Curtis distance matrices (Anderson et al., 2006).

Statistical analysis

Normalized data (RQ expression levels) were initially analyzed, with an estimation and comparison of expression levels between groups. Normal distribution of data was assessed using the Kolmogorov-Smirnov and Shapiro-Wilk tests. Variables with normal distribution are presented as mean ± standard deviation (SD). Variables with skewed distribution were log-transformed prior to analysis and are presented as median (25th – 75th percentiles). Categorical data are shown as percentages. Clinical and laboratory characteristics, miRNA expressions, and gut microbiota abundance were compared among groups using Student’s t-test or χ2 tests, as appropriate. Correlations between quantitative variables were assessed using Pearson’s correlation tests.

All classical statistical analyses were performed using the SPSS statistical package (v.20.0) for Windows (SPSS Inc, Chicago, IL, USA) and PAST v3.24 (University of Oslo, Norway) for statistical analyses of biodiversity. FDR correction was used to account for multiple comparisons using the Benjamini–Hochberg method (q-value < 0.05). The network visualization of miRNA-microbe was generated using Cytoscape v.3.7.1 (Shannon et al., 2003). One heatmap plot of the correlation values were produced using MORPHEUS web tool (Morpheus, https://software.broadinstitute.org/morpheus).

Results

Clinical and laboratory characteristics of individuals included in the study.

Clinical, laboratorial, and nutritional characteristics of cases with obesity and eutrophic controls are shown in Table 1. There were no differences between cases and controls regarding age, gender, and energy intake. Moreover, both groups had a comparable dietary macronutrient composition. As expected, subjects with obesity presented higher waist circumference, glucose, total cholesterol, and triglyceride levels compared to normal weight individuals. Additionally, cases also presented elevated levels of metabolic markers such as insulin, leptin, TyG, and HOMA-IR indexes and, lower levels of METs compared to controls.

Table 1

Clinical, dietary, and laboratory characteristics of the sample included in the study.

Characteristics

Subjects with Obesity (cases, n = 78)

Eutrophic individuals (controls, n = 25)

P-value*

Age (years)

46.6 ± 9.4

44.7 ± 9.1

0.106

Gender (% male)

36.1

40.0

0.443

Anthropometric and clinical data

   

BMI (kg/m2)

32.9 ± 2.4

18.6 ± 2.1

-

WC (cm)

104.9 ± 10.2

75.2 ± 7.6

0.0001

SBP (mmHg)

131 ± 16

111 ± 10

0.0001

DBP (mmHg)

86 ± 9

70 ± 8

0.0001

Metabolic profile

     

FPG (mg/dl)

97.4 ± 11.9

85.3 ± 6.8

0.0001

TC (mg/dl)

222.5 ± 40.1

192.6 ± 37.1

0.001

HDL-c (mg/dl)

54.2 ± 14.0

61.6 ± 12.7

0.022

TG (mg/dl)

101.6 ± 54.1

65.6 ± 25.0

0.002

TyG index

4.6 ± 0.3

4.2 ± 0.2

0.0001

HOMA-IR index

1.6 (1.1–2.8)

0.6 (0.4–1.0)

0.0001

Adiponectin (ng/ml)

10.9 (7.9–13.5)

12.2 (9.3–15.7)

0.067

Insulin (mU/L)

6.8 (4.7–11.5)

3.2 (2.8–4.8)

0.0001

Leptin (ng/ml)

33.1 (17.2–46.8)

4.9 (2.1–11.7)

0.0001

Body composition

     

Fat mass (%)

34.7 ± 6.5

13.6 ± 5.7

0.0001

Lean mass (%)

57.0 ± 11.7

47.6 ± 12.2

0.001

Dietary intake and energy expenditure

   

Energy (Kcal)

2,961 ± 1051

2,588 ± 701

0.101

Carbohydrates (%)

41.4 ± 7.1

44.8 ± 6.4

0.034

Protein (%)

16.7 ± 2.9

15.9 ± 3.4

0.245

Fat (%)

40.1 ± 6.4

37.7 ± 4.9

0.100

Fiber (g/day)

27.9 ± 11.4

32.7 ± 11.5

0.070

METs (kcal/kg/h)

17.0 (7.5–27.0)

33.2 (20.0–44.4)

0.001

Variables are shown as mean ± SD, median (25th–75th percentiles) or %, as appropriate. *P-values were computed using χ2 or Student´s t-test, as appropriated.
BMI: body mass index; DBP: diastolic blood pressure; FPG: fasting plasma glucose; HDL-c: high-density lipoprotein cholesterol; HOMA-IR index: homeostatic model assessment-insulin resistance index; METs: metabolic equivalents; SBP: systolic blood pressure; TC: total cholesterol; TG: triglycerides; TyG index: triglyceride glucose index; WC: waist circumference.

Quality control of miRNA expression

The RNA spike-in expressions presented low variation in Cq among samples in RNA isolation and cDNA synthesis, demonstrating that extraction, reverse transcription, and qPCR were effective, and none of the samples contained inhibitors. As expected, the expression of UniSp2, UniSp4, UniSp5, and UniSp6 did not differ between groups (cases vs. controls). UniSp5 was expressed in all analyzed samples, demonstrating that miRNAs expressed in low levels was not lost during isolation. The ratio between miR-451a and miR-23a-3p ranged between 5 and − 1, indicating that the samples were not affected by hemolysis. Generally, these results showed a good and similar level of sample quality and reproducibility of the miRNA profiling processes.

MicroRNAs differentially expressed in plasma of patients with obesity.

Expression of 86 target miRNAs was evaluated in plasma of subjects with obesity and in normal weight individuals. Of these 86 miRNAs, 61 were expressed in at least 20% of the sample with Cq values ≥ 35. Of these 61 miRNAs, 26 were differentially expressed between cases and controls after FDR correction (Table 2).

Table 2

Relation of 26 microRNAs whose expression profile in plasma is significantly different between cases with obesity and eutrophic controls.

miRNA

Subjects with obesity (cases, n = 78)

Eutrophic individuals (controls, n = 25)

P-value*

q-value**

miR-103a-3p

0.169 (0.057–0.443)

0.567 (0.263–1.426)

0.006

0.020

miR-107

0.201 (0.067–0.520)

0.614 (0.309–1.386)

0.014

0.038

miR-130a-3p

19.321 (6.029–32.947)

47.904 (31.421–79.317)

0.005

0.004

miR-130b-3p

0.208 (0.096–0.433)

0.442 (0.221–1.426)

0.003

0.015

miR-140-3p

0.237 (0.076–0.578)

0.872 (0.548–1.561)

0.0001

0.0012

miR-142-5p

0.118 (0.036–0.214)

0.279 (0.162–0.941)

0.002

0.007

miR-144-3p

1.308 (0.276–2.872)

6.903 (2.036–11.542)

0.0001

0.0012

miR-148a-3p

0.259 (0.096–0.479)

0.647 (0.212–1.135)

0.006

0.019

miR-181a-5p

0.549 (0.229–1.107)

1.633 (0.369–3.021)

0.006

0.008

miR-183-5p

0.374 (0.244–0.649)

0.775 (0.476–1.928)

0.001

0.009

miR-185-5p

0.198 (0.083–0.480)

0.665 (0.288–1.272)

0.034

0.040

miR-200c-3p

0.540 (0.257–1.118)

1.001 (0.454–2.161)

0.037

0.044

miR-205-5p

0.163 (0.122–0.559)

0.711 (0.319–1.322)

0.005

0.020

miR-21-5p

0.257 (0.137–0.615)

0.559 (0.261–1.238)

0.036

0.041

miR-210-3p

0.110 (0.071–0.245)

0.503 (0.117–0.672)

0.030

0.078

miR-221-3p

0.212 (0.065–0.454)

0.555 (0.268–1.549)

0.004

0.013

miR-222-3p

0.311 (0.158–0.712)

0.859 (0.470–1.396)

0.005

0.019

miR-15a-5p

0.118 (0.060–0.322)

0.356 (0.155–0.877)

0.022

0.054

miR-22-3p

0.126 (0.056–0.286)

0.320 (0.127–1.027)

0.012

0.034

miR-29c-3p

0.137 (0.068–0.350)

0.453 (0.202–1.269)

0.007

0.020

miR-30a-5p

0.629 (0.294–1.347)

1.426 (0.515–2.062)

0.043

0.048

miR-30c-5p

0.274 (0.093–0.625)

0.694 (0.235–1.462)

0.042

0.050

miR-33a-5p

1.268 (0.573–2.425)

4.54 (1.009–23.596)

0.012

0.016

miR-375

0.228 (0.124–0.609)

0.765 (0.484–1.563)

0.0001

0.002

miR-424-3p

0.892 (0.698–1.135)

2.000 (0.834–3.127)

0.016

0.030

miR-486-3p

0.301 (0.190–0.570)

0.645 (0.238–1.326)

0.004

0.009

Data are shown as median (25th–75th percentiles) of n-fold values. *P-values were obtained using Student t test using the log-transformed variable. **P-values were corrected using false discovery rate (FDR; q-value).

All these 26 miRNAs were negatively correlated with BMI (P ≤0.05). Moreover, miR-107, miR-130a-3p, miR-140-3p, miR-142-5p, miR-144-3p, miR-181a-3p, miR-21-5p, miR-221-3p, miR-375, and miR-424-3p expressions were negatively associated with glucose levels (P < 0.05). Otherwise, miR-200c-3p and miR-375 positively correlated with HDL-c levels (r = 0.232, P = 0.047; and r = 0.295, P = 0.015, respectively). Regarding hormone levels, miR-140-3p and miR-144-3p were negatively associated with leptin levels (r= -0.222, P = 0.034; and r= -0.245, P = 0.025, respectively), miR-144-3p and miR-183-5p were inversely correlated with insulin (r= -0.245, P = 0.032; and r= -0.264, P = 0.033, respectively) and adiponectin levels were positively associated with miR-375 (r = 0.272, P = 0.025) and miR-424-3p (r = 0.405, P = 0.012).

Gut microbiota profile in subjects with obesity compared to eutrophic individuals

The effect of obesity on gut microbiota composition was investigated at the genus and species levels. The levels of eighteen bacterial genera were significantly different when comparing obese and normal weight individuals, being nine bacterial genera significantly increased in obese subjects when compared to controls (Fig. 1A). Twelve bacterial species were statistically different between obese and normal weight individuals, being ten of them more abundant in subjects with obesity compared to eutrophic individuals (Fig. 1B and Table 3).

Table 3

Bacterial species whose abundance is statistically different when comparing cases with obesity and eutrophic controls.

Bacteria

Subjects with obesity (cases, n = 78)

Eutrophic individuals (controls, n = 25)

P-value*

q-value**

Abiotrophia defectiva

0.0001 (0.0001–0.2825)

0.0001 (0.0001–0.001)

0.0001

-

Actinomyces odontolyticus

0.0001 (0.0001–1.643)

0.0001 (0.0001–0.652)

0.010

-

Allisonella histaminiformans

0.0001 (0.0001- 5.260)

0.0001 (0.0001–0.001)

0.0001

-

Bacteroides eggerthii

3.573 (1.903–10.086)

9.769 (3.32–11.892)

0.039

0.042

Barnesiella intestinihominis

2.885 (0.730–5.809)

1.889 (0.001–3.361)

0.019

0.036

Dorea longicatena

9.959 (9.295–10.557)

9.547 (8.887–10.268)

0.039

0.045

Haemophilus parainfluenzae

5.834 (2.990–7.775)

8.254 (6.181–10.708)

0.001

0.002

Howardella ureilytica

2.439 (0.0001–6.303)

0.0001 (0.0001–2.476)

0.010

0.023

Lactobacillus curvatus

0.0001 (0.0001–1.626)

0.0001 (0.0001–0.001)

0.0001

-

Megamonas funiformis

0.0001 (0.0001–0.974)

0.0001 (0.0001–0.001)

0.0001

-

Mitsuokella jaladudinii

0.0001 (0.0001–2.521)

0.0001 (0.0001–0.001)

0.002

-

Odoribacter laneus

0.0001 (0.0001–2.707)

0.0001 (0.0001–1.015)

0.027

-

Data are shown as median (25th − 75th percentiles). *P-values were obtained using Student t test using the log-transformed variable. **P-values were corrected for false discovery rate (FDR; q-value) using the Benjamini-Hochberg test.

Shannon index, which reflects the alpha diversity, was not different between obese and normal weight groups (Supplementary Fig. 1). However, the beta diversity values of gut microbiota, based on Jaccard index (PERMANOVA, P = 0.025; Supplementary Fig. 2A) and Bray-Curtis dissimilarity (PERMANOVA, P = 0.015; Supplementary Fig. 2B), was significantly different between groups.

Crosstalk between host miRNAs and gut microbiota

To further investigate the relationships between circulating miRNAs and the gut microbiota composition, interactions between bacteria and miRNAs differentially expressed in obesity were analyzed. At the genus level, of the 18 genera differently expressed in obesity, 9 were significantly correlated with the expression of 10 miRNAs out of 26 miRNAs differently expressed in subjects with obesity (Fig. 2A). Fourteen of these miRNAs were significantly associated with 4 bacterial species (Dorea longicatena, Banesiela intestinihominis, Bacteroides eggerthii, and Haemophillus parainfluenzae), as illustrated in Fig. 2B and Fig. 2C.

A diagram was built to visualize the relationships between the miRNAs and their significantly correlated bacteria (Fig. 2C). The correlation network shows a highly interconnected relationship between these miRNAs and bacterial species. Interestingly, B. eggerthii negatively correlated with miR-103a-3p, miR-21-5p, miR-130a-3p, miR-185-5p, miR-144-3p, miR-210-3p, miR-33a-5p, miR-15a-5p, miR-130b-3p, miR-183-5p, miR-221-3p, miR-222-3p, and miR-142-5p. Moreover, an interaction among miRNAs, B. eggerthi and BMI levels was found. Individually, the expression of miR-103a-3p (r2 = 0.1229, P = 0.051), miR-130b-3p (r2 = 0.0933, P = 0.021), miR-185-5p (r2 = 0.0894, P = 0.035), miR-21-5p (r2 = 0.1124, P = 0.008), and miR-210 (r2 = 0.0866, P = 0.052) interacted with these bacterial species and BMI levels. Furthermore, an interaction among three miRNAs levels (miR-130b-3p, miR-185-5p, and miR-21-5p), B. eggerthi and BMI levels was also evidenced (r2 = 0.148, P = 0.004). Interestingly, there was also an interaction among B. eggerthi, adiponectin levels, and miR-183-5p (r2 = 0.1294, P = 0.009).

In the same way, B. intestinihominis abundance was negatively associated with miR-107, miR-103a-3p, miR-222-3p, and miR-142-5p expressions. The expression of miR-15a-5p was inversely associated with the abundance of H. parainfluenzae, and an interaction with insulin levels (r2 = 0.0592, P = 0.027) was found. In contrast, D. longicatena was positively associated with miR-21-5p, miR-130a-3p, miR-185-5p, and miR-144-3p. However, interactions among the bacterial abundance, miRNA expression, and BMI levels were not found for these three bacterial species. No association among the bacterial species, miRNAs, and leptin was found.

Predicted functions of miRNAs correlated with obesity-associated bacteria

Target gene prediction of the 14 miRNAs that correlated with the 4 bacterial species associated with obesity were investigated (Supplementary Table 2). Of the total 9,584 genes identified as potential targets of these miRNAs, 5,381 were found to be regulated by two or more miRNAs; however only 719 were experimentally validated (Supplementary Table 2). After that, functional enrichment analysis of miRNA targets was carried out to explore biological pathways possibly regulated by this set of miRNAs. A total of 248 pathways were significantly enriched (q-value < 0.05) for these miRNAs (Supplementary Table 3). However, considering only the experimentally validated target genes, 98 pathways were significantly enriched (Supplementary Table 3).

As shown in Fig. 3, H. parainfluenzae, D. longicatena, B. intestinihominis, and B. eggerthii correlated with miRNAs associated with pathways related to obesity and metabolic processes, including carbohydrate and lipid turnover, endocrine and inflammatory signaling pathways. More specifically, the target genes of miRNAs associated with the four bacterial species related to obesity participate in the fatty acid degradation, mineral absorption, carbohydrate digestion and absorption, insulin signaling pathway, and glycerolipid metabolism.

Discussion

Over the past decade, there has been an increasing attention about the crucial roles played by miRNAs in a wide variety of cellular processes (Lorente-Cebrián et al., 2019). In the present study, 26 miRNAs were found as differentially expressed in plasma of subjects with obesity compared to normal weight individuals. Furthermore, the expression of 14 miRNAs (miR-107, miR-103a-3p, miR-142-5p, miR-222-3p, miR-221-3p, miR-183-5p, miR-183-5p, miR-130b-3p, miR-15a-5p, miR-33a-5p, miR-210-3p, miR-144-3p, miR-185-5p, miR-130a-3p, and miR-21-5p) was linked with the relative abundance of 4 bacterial species that also significantly differed between cases and controls (D. longicatena, B. intestinihominis, B. eggerthii, and H. parainfluenzae).

These miRNAs that interacted with obesity-associated bacteria regulate the expression of genes that participate in several metabolism and obesity-related pathways, such as carbohydrate and lipid metabolism, endocrine and inflammatory signaling pathways. Indeed, evidence suggests that the majority of miRNAs do not regulate a specific or individual target gene, but rather they modulate the expression of large number of genes in networks, demonstrating their importance in the regulation of several metabolic processes (Bartel, 2004; Virtue et al., 2019).

Additionally, an interaction between BMI levels, B. eggerthii abundance, and the expression of three miRNAs (miR-130b-3p, miR-185-5p, and miR-21-5p) was also evidenced. Interestingly, B. eggerthii is one of the intestinal bacteria that metabolize phenolic acids, which are regarded as beneficial for human health (Russell et al., 2013). In a recent study, B. eggerthii abundance was significantly higher in children with obesity and positively correlated with body fat percentage, but negatively with insoluble fiber intake in Mexican children (López-Contreras et al., 2018). On the other hand, this bacterium was found to be underrepresented after sleeve gastrectomy surgery (Medina et al., 2017).

Of the three miRNAs associated with the abundance of B. eggerthii and BMI levels, miR-185-5p and miR-21-5p were also correlated with D. longicatena abundance. Furthermore, miR-185-5p was described as involved in oxidative stress, obesity, and DM in many studies [reviewed at (Matoušková et al., 2018)]. Moreover, miR-185-5p was identified as a regulator of de novo cholesterol biosynthesis and low density lipoprotein uptake (Yang et al., 2014). However, we could not find in the literature evidences of association between this miRNA and gut microbiota.

Regarding miR-21-5p, the 16S rRNA sequencing revealed significant differences in the composition of WT and miR-21−/− intestinal microbiota in a dextran sodium sulphate (DSS)-induced colitis mouse model (Johnston et al., 2018). Moreover, commensal bacteria induced the expression of miR-21 in IECs, with implications in the regulation of intestinal epithelial permeability (Nakata et al., 2017). Otherwise, miR-130b-3p was only correlated with B. eggerthi abundance and there is evidence that the expression of this miRNA was influenced by microbial status in intestinal epithelial stem cell of conventionalized mice compared to germ-free mice (Peck et al., 2017), showing that the miRNAs expression could be modulated by gut microbiota.

Moreover, an association among B. eggerthi abundance, miR-183-5p expression, and adiponectin levels was also found. Previous findings demonstrated that miR-183 may contribute to adipocyte differentiation, adipogenesis, and development of fat cells (John et al., 2012; Sedgeman et al., 2019). Additionally, miR-183 was identified as a novel positive regulator during 3T3-L1 adipogenesis. Both gain-of-function and loss-of-function assays showed that miR-183 promoted 3T3-L1 adipocyte differentiation, lipid accumulation, and adipogenesis by enhancing the expressions of peroxisome proliferator activated receptor gamma (PPARγ), CCAAT enhancer binding protein alpha (C/EBPα), adiponectin, and fatty acid synthase (FAS) (Chen et al., 2014).

The expression of miR-15a-5p was found associated with H. parainfluenzae abundance and insulin levels in our study. miR-15a positively regulates insulin biosynthesis by inhibiting endogenous uncoupling protein 2 (UCP2) expression, leading to higher ATP levels in islets and improving glucose-stimulated insulin secretion. Moreover, circulating levels of miR-15a were found downregulated before the onset of type 2 DM (T2DM) (Zampetaki et al., 2010), and also in incident-T2DM subjects compared to controls, with intermediate values in the pre-DM and incident pre-DM patients (Jiménez-Lucena et al., 2018).

Regarding gut microbiota composition, our results evidenced that obesity had no significant impact in alpha diversity, indicating that microbial species diversity is relatively stable in response to obesity. However, obesity influenced the beta diversity of human gut microbiota compared to the control group, suggesting that this disease is accompanied by species replacement (changes in species taxa) and species sorting (changes in abundance).

According to a meta-analysis of metagenomic datasets obtained from fecal samples of healthy human adults living in different world regions, Bacteroides and Barnesiella genera are markers of Western populations (Mancabelli et al., 2017). Barnesiella spp (represented mainly by the specie Barnesiella intestinihominis) were identified only in populations living in developed countries, suggesting that their presence was promoted by the urbanization/ industrialization process and Western-type diet (Mancabelli et al., 2017).

In agreement with our results, the levels of Dorea genera were previously reported to be higher in overweight children compared to normal weight counterparts (Karvonen et al., 2019). Moreover, this association was stronger for non-white children than for white children, and also stronger for boys than for girls (Karvonen et al., 2019). Interestingly, a recent study in an early-life HFD mouse model found that the this diet increased the relative abundances of Dorea genus (Villamil et al., 2018).

Our investigation has strengths and limitation. The strengths include study and data analyses of a very-well characterized cohort of subjects with obesity and eutrophic subjects was analyzed. Moreover, several quality controls for miRNA extraction, cDNA synthesis, and PCR process were implemented. Additionally, robust bioinformatic analyses were performed to explore the pathways where these miRNAs target genes are participating, explaining the association with obesity. Likewise, we highlighted candidates for potentially linking host miRNAs and gut microbiota, which can be directly validated and explored in model systems.

Even though these methods are powerful, this evaluation has some limitations. First, it is important to note that our study uses 16S rRNA gene sequencing to characterize microbiome taxonomic composition. Second, the results from bioinformatics are predictions and may not represent the real biological system. Third, our approach identifies correlations and not causal relationships. Even though a hypothesis-driven approach was performed, selection of only miRNAs previously associated with obesity or metabolism makes possible type I or type II errors due to multiple comparisons. These limitations should be considered when interpreting the results. Although limitations exist in the current data, the patterns uncovered here are important for understanding the contribution of miRNAs and gut microbiota in obesity.

Conclusion

This current research characterized a global relationship between microbial community composition and miRNA expression in plasma of subjects with obesity compared to normal weight individuals. Indeed, our study featured an interaction between B. eggerthi abundance and circulating miRNA expression in the control of body adiposity. The current study also adds to the growing body of literature that suggests that miRNAs may serve as a communication mechanism between the gut microbiota and human hosts.

Abbreviations

Body mass index

BMI

Diabetes mellitus

DM

False discovery rate

FDR

Fasting plasma glucose

FPG

High-density lipoprotein

HDL

High-fat diet

HFD

Homeostatic model assessment insulin resistance

HOMA-IR

Intestinal epithelial cells

IECs

Metabolic equivalents

METs

microRNA

miRNA

Operational taxonomic units

OTUs

Total cholesterol

TC

Triglyceride-glucose

TyG

Triglycerides

TG

White adipose tissue

WAT

Declarations

Acknowledgements. Authors thank Marta Cuervo, Leticia Goni, Laura Olazaran, and Iosune Zubieta for the recruitment and nutritional assistance, Ana Lorente for laboratory support, and Anna Barceló for the contribution for the 16S rRNA sequencing.

Authors' contributions. T.S.A. designed the study, analyzed and interpreted the data, and drafted the manuscript. A.C. collected and analyzed the microbiota data. J.I.R. designed the study, interpreted the data, and critically reviewed the manuscript. F.I.M. designed and supervised the study, interpreted the data, and critically reviewed the manuscript. J.AM. designed and supervised the study, interpreted the data, and critically reviewed the manuscript.

Funding. This research was supported by CIBERobn (grant number: CB12/03/30002), Government of Navarra (Obekit-PT024 and Microbiota-PI035 projects), and the Spanish Ministerio de Ciencia, Innovación y Universidades (reference RTI2018-102205-B-I00). T.S.A is recipient of scholarships from Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (PDE/CAPES, grant number: 88881.170123/2018-01).

Ethics approval and consent to participate. The research was conducted in accordance with the rules of the Helsinki Declaration. The research protocol was properly approved by the Research Ethics Committee of the University of Navarra (ref. 132/2015) and it is registered at clinicaltrials.gov (reg. no. NCT02737267). A written informed consent of each participant was obtained prior to enrollment in the study.

Consent for publication. All authors approved the final version and agreed to be accountable for all aspects of the work regarding accuracy and integrity aspects. All authors agree to publish this article in the journal of Molecular Medicine.

Competing interests. All authors declare no competing interests. The datasets used during the current study are available from the corresponding author on reasonable request.

References

  1. Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403–10.
  2. Anderson MJ, Ellingsen KE, McArdle BH. Multivariate dispersion as a measure of beta diversity. Ecol Lett. 2006;9:683–93.
  3. Association WM. World Medical Association Declaration of Helsinki: ethical principles for medical research involving human subjects. JAMA. 2013;310:2191–4.
  4. Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116:281–97.
  5. Basterra-Gortari FJ, Bes-Rastrollo M, Pardo-Fernández M, Forga L, Martinez JA, Martínez-González MA. Changes in weight and physical activity over two years in Spanish alumni. Med Sci Sports Exerc. 2009;41:516–22.
  6. Belcheva A. 2017. MicroRNAs at the epicenter of intestinal homeostasis. Bioessays 39.
  7. Blondal T, Nielsen J, Baker S, Andreasen A, Mouritzen D, Wrang Teilum P, M. and Dahlsveen IK. Assessing sample and miRNA profile quality in serum and plasma or other biofluids. Methods. 2013;59:1–6.
  8. Bustin SA, Benes V, Garson JA, Hellemans J, Huggett J, Kubista M, Mueller R, Nolan T, Pfaffl MW, Shipley GL, Vandesompele J, Wittwer CT. The MIQE guidelines: minimum information for publication of quantitative real-time PCR experiments. Clin Chem. 2009;55:611–22.
  9. Butz H, Kinga N, Racz K, Patocs A. Circulating miRNAs as biomarkers for endocrine disorders. J Endocrinol Invest. 2016;39:1–10.
  10. Campión J, Milagro FI, Martínez JA. Individuality and epigenetics in obesity. Obes Rev. 2009;10:383–92.
  11. Chen C, Xiang H, Peng YL, Peng J, Jiang SW. Mature miR-183, negatively regulated by transcription factor GATA3, promotes 3T3-L1 adipogenesis through inhibition of the canonical Wnt/β-catenin signaling pathway by targeting LRP6. Cell Signal. 2014;26:1155–65.
  12. de la Fuente-Arrillaga C, Ruiz ZV, Bes-Rastrollo M, Sampson L, Martinez-González MA. Reproducibility of an FFQ validated in Spain. Public Health Nutr. 2010;13:1364–72.
  13. Dweep H, Sticht C, Pandey P, Gretz N. miRWalk–database: prediction of possible miRNA binding sites by "walking" the genes of three genomes. J Biomed Inform. 2011;44:839–47.
  14. Edgar RC. UPARSE: highly accurate OTU sequences from microbial amplicon reads. Nat Methods. 2013;10:996–8.
  15. Esteller M. Non-coding RNAs in human disease. Nat Rev Genet. 2011;12:861–74.
  16. Gevaert AB, Witvrouwen I, Vrints CJ, Heidbuchel H, Van Craenenbroeck EM, Van Laere SJ, Van Craenenbroeck AH. MicroRNA profiling in plasma samples using qPCR arrays: Recommendations for correct analysis and interpretation. PLoS One. 2018;13:e0193173.
  17. Guyenet SJ, Schwartz MW. Clinical review: Regulation of food intake, energy balance, and body fat mass: implications for the pathogenesis and treatment of obesity. J Clin Endocrinol Metab. 2012;97:745–55.
  18. Hao X, Jiang R, Chen T. Clustering 16S rRNA for OTU prediction: a method of unsupervised Bayesian clustering. Bioinformatics. 2011;27:611–8.
  19. Hildebrand F, Tadeo R, Voigt AY, Bork P, Raes J. LotuS: an efficient and user-friendly OTU processing pipeline. Microbiome. 2014;2:30.
  20. Jiménez-Lucena R, Camargo A, Alcalá-Diaz JF, Romero-Baldonado C, Luque RM, van Ommen B, Delgado-Lista J, Ordovás JM, Pérez-Martínez P, Rangel-Zúñiga OA, López-Miranda J. 2018. A plasma circulating miRNAs profile predicts type 2 diabetes mellitus and prediabetes: from the CORDIOPREV study. Exp Mol Med 50, 168.
  21. John E, Wienecke-Baldacchino A, Liivrand M, Heinäniemi M, Carlberg C, Sinkkonen L. Dataset integration identifies transcriptional regulation of microRNA genes by PPARγ in differentiating mouse 3T3-L1 adipocytes. Nucleic Acids Res. 2012;40:4446–60.
  22. Johnston DGW, Williams MA, Thaiss CA, Cabrera-Rubio R, Raverdeau M, McEntee C, Cotter PD, Elinav E, O'Neill LAJ, Corr SC. Loss of MicroRNA-21 Influences the Gut Microbiota, Causing Reduced Susceptibility in a Murine Model of Colitis. J Crohns Colitis. 2018;12:835–48.
  23. Karvonen AM, Sordillo JE, Gold DR, Bacharier LB, O'Connor GT, Zeiger RS, Beigelman A, Weiss ST, Litonjua AA. Gut microbiota and overweight in 3-year old children. Int J Obes (Lond). 2019;43:713–23.
  24. Lawrence M, Huber W, Pagès H, Aboyoun P, Carlson M, Gentleman R, Morgan MT, Carey VJ. Software for computing and annotating genomic ranges. PLoS Comput Biol. 2013;9:e1003118.
  25. Liu S, Weiner HL. Control of the gut microbiome by fecal microRNA. Microb Cell. 2016;3:176–7.
  26. Lopez-Legarrea P, de la Iglesia R, Abete I, Bondia-Pons I, Navas-Carretero S, Forga L, Martinez JA, Zulet MA. Short-term role of the dietary total antioxidant capacity in two hypocaloric regimes on obese with metabolic syndrome symptoms: the RESMENA randomized controlled trial. Nutr Metab (Lond). 2013;10:22.
  27. Lorente-Cebrián S, González-Muniesa P, Milagro FI, Martínez JA. MicroRNAs and other non-coding RNAs in adipose tissue and obesity: emerging roles as biomarkers and therapeutic targets. Clin Sci (Lond). 2019;133:23–40.
  28. López-Contreras BE, Morán-Ramos S, Villarruel-Vázquez R, Macías-Kauffer L, Villamil-Ramírez H, León-Mimila P, Vega-Badillo J, Sánchez-Muñoz F, Llanos-Moreno LE, Canizalez-Román A, Del Río-Navarro B, Ibarra-González I, Vela-Amieva M, Villarreal-Molina T, Ochoa-Leyva A, Aguilar-Salinas CA, Canizales-Quinteros S. Composition of gut microbiota in obese and normal-weight Mexican school-age children and its association with metabolic traits. Pediatr Obes. 2018;13:381–8.
  29. Mancabelli L, Milani C, Lugli GA, Turroni F, Ferrario C, van Sinderen D, Ventura M. Meta-analysis of the human gut microbiome from urbanized and pre-agricultural populations. Environ Microbiol. 2017;19:1379–90.
  30. Martínez-González MA, López-Fontana C, Varo JJ, Sánchez-Villegas A, Martinez JA. Validation of the Spanish version of the physical activity questionnaire used in the Nurses' Health Study and the Health Professionals' Follow-up Study. Public Health Nutr. 2005;8:920–7.
  31. Matoušková P, Hanousková B, Skálová L. 2018. MicroRNAs as Potential Regulators of Glutathione Peroxidases Expression and Their Role in Obesity and Related Pathologies. Int J Mol Sci 19.
  32. Maurizi G, Babini L, Della Guardia L. 2018. Potential role of microRNAs in the regulation of adipocytes liposecretion and adipose tissue physiology. J Cell Physiol.
  33. McKenna LB, Schug J, Vourekas A, McKenna JB, Bramswig NC, Friedman JR, Kaestner KH. MicroRNAs control intestinal epithelial differentiation, architecture, and barrier function. Gastroenterology. 2010;139(64):1654. ., -, 1664.e1.
  34. Medina DA, Pedreros JP, Turiel D, Quezada N, Pimentel F, Escalona A, Garrido D. Distinct patterns in the gut microbiota after surgical or medical therapy in obese patients. PeerJ. 2017;5:e3443.
  35. Mestdagh P, Van Vlierberghe P, De Weer A, Muth D, Westermann F, Speleman F, Vandesompele J. A novel and universal method for microRNA RT-qPCR data normalization. Genome Biol. 2009;10:R64.
  36. Moreiras O, Cabrera A, LuisaCuadrado C. 2018. Tablas de Composición de Alimentos. Guía de prácticas., 19 ed.
  37. Nakata K, Sugi Y, Narabayashi H, Kobayakawa T, Nakanishi Y, Tsuda M, Hosono A, Kaminogawa S, Hanazawa S, Takahashi K. Commensal microbiota-induced microRNA modulates intestinal epithelial permeability through the small GTPase ARF4. J Biol Chem. 2017;292:15426–33.
  38. Navarro-González D, Sánchez-Íñigo L, Pastrana-Delgado J, Fernández-Montero A, Martinez JA. Triglyceride-glucose index (TyG index) in comparison with fasting plasma glucose improved diabetes prediction in patients with normal fasting glucose: The Vascular-Metabolic CUN cohort. Prev Med. 2016;86:99–105.
  39. Pasolli E, Schiffer L, Manghi P, Renson A, Obenchain V, Truong DT, Beghini F, Malik F, Ramos M, Dowd JB, Huttenhower C, Morgan M, Segata N, Waldron L. Accessible, curated metagenomic data through ExperimentHub. Nat Methods. 2017;14:1023–4.
  40. Peck BC, Mah AT, Pitman WA, Ding S, Lund PK, Sethupathy P. Functional Transcriptomics in Diverse Intestinal Epithelial Cell Types Reveals Robust MicroRNA Sensitivity in Intestinal Stem Cells to Microbial Status. J Biol Chem. 2017;292:2586–600.
  41. Pichler M, Coskun Ö, Ortega-Arbulú AS, Conci N, Wörheide G, Vargas S, Orsi WD. A 16S rRNA gene sequencing and analysis protocol for the Illumina MiniSeq platform. Microbiologyopen. 2018;7:e00611.
  42. Rahmati S, Abovsky M, Pastrello C, Jurisica I. pathDIP: an annotated resource for known and predicted human gene-pathway associations and pathway enrichment analysis. Nucleic Acids Res. 2017;45:D419–26.
  43. Ramos-Lopez O, Riezu-Boj JI, Milagro FI, Cuervo M, Goni L, Martinez JA. 2018. Prediction of Blood Lipid Phenotypes Using Obesity-Related Genetic Polymorphisms and Lifestyle Data in Subjects with Excessive Body Weight. Int J Genomics 2018, 4283078.
  44. Rao X, Huang X, Zhou Z, Lin X. An improvement of the 2ˆ(-delta delta CT) method for quantitative real-time polymerase chain reaction data analysis. Biostat Bioinforma Biomath. 2013;3:71–85.
  45. Rideout JR, He Y, Navas-Molina JA, Walters WA, Ursell LK, Gibbons SM, Chase J, McDonald D, Gonzalez A, Robbins-Pianka A, Clemente JC, Gilbert JA, Huse SM, Zhou HW, Knight R, Caporaso JG. Subsampled open-reference clustering creates consistent, comprehensive OTU definitions and scales to billions of sequences. PeerJ. 2014;2:e545.
  46. Ritari J, Salojärvi J, Lahti L, de Vos WM. Improved taxonomic assignment of human intestinal 16S rRNA sequences by a dedicated reference database. BMC Genom. 2015;16:1056.
  47. Russell WR, Duncan SH, Scobbie L, Duncan G, Cantlay L, Calder AG, Anderson SE, Flint HJ. Major phenylpropanoid-derived metabolites in the human gut can arise from microbial fermentation of protein. Mol Nutr Food Res. 2013;57:523–35.
  48. Sedgeman LR, Michell DL, Vickers KC. Integrative roles of microRNAs in lipid metabolism and dyslipidemia. Curr Opin Lipidol. 2019;30:165–71.
  49. Shannon CE. 1997. The mathematical theory of communication. 1963. MD Comput 14, 306 – 17.
  50. Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, Amin N, Schwikowski B, Ideker T. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13:2498–504.
  51. Singh N, Shirdel EA, Waldron L, Zhang RH, Jurisica I, Comelli EM. The murine caecal microRNA signature depends on the presence of the endogenous microbiota. Int J Biol Sci. 2012;8:171–86.
  52. van Dijk SJ, Tellam RL, Morrison JL, Muhlhausler BS, Molloy PL. Recent developments on the role of epigenetics in obesity and metabolic disease. Clin Epigenetics. 2015;7:66.
  53. Vandenbroucke JP, von Elm E, Altman DG, Gøtzsche PC, Mulrow CD, Pocock SJ, Poole C, Schlesselman JJ, Egger M, Initiative S. Strengthening the Reporting of Observational Studies in Epidemiology (STROBE): explanation and elaboration. Int J Surg. 2014;12:1500–24.
  54. Vienberg S, Geiger J, Madsen S, Dalgaard LT. MicroRNAs in metabolism. Acta Physiol (Oxf). 2017;219:346–61.
  55. Villamil SI, Huerlimann R, Morianos C, Sarnyai Z, Maes GE. Adverse effect of early-life high-fat/high-carbohydrate ("Western") diet on bacterial community in the distal bowel of mice. Nutr Res. 2018;50:25–36.
  56. Virtue AT, McCright SJ, Wright JM, Jimenez MT, Mowel WK, Kotzin JJ, Joannas L, Basavappa MG, Spencer SP, Clark ML, Eisennagel SH, Williams A, Levy M, Manne S, Henrickson SE, Wherry EJ, Thaiss CA, Elinav E, Henao-Mejia J. 2019. The gut microbiota regulates white adipose tissue inflammation and obesity via a family of microRNAs. Sci Transl Med 11.
  57. von Elm E, Altman DG, Egger M, Pocock SJ, Gøtzsche PC, Vandenbroucke JP, Initiative S. The Strengthening the Reporting of Observational Studies in Epidemiology (STROBE) statement: guidelines for reporting observational studies. J Clin Epidemiol. 2008;61:344–9.
  58. Yang M, Liu W, Pellicane C, Sahyoun C, Joseph BK, Gallo-Ebert C, Donigan M, Pandya D, Giordano C, Bata A, Nickels JT. Identification of miR-185 as a regulator of de novo cholesterol biosynthesis and low density lipoprotein uptake. J Lipid Res. 2014;55:226–38.
  59. Wordl Health Organization. WHO. 2018. Obesity and overweight.
  60. Zampetaki A, Kiechl S, Drozdov I, Willeit P, Mayr U, Prokopi M, Mayr A, Weger S, Oberhollenzer F, Bonora E, Shah A, Willeit J, Mayr M. Plasma microRNA profiling reveals loss of endothelial miR-126 and other microRNAs in type 2 diabetes. Circ Res. 2010;107:810–7.