Identification and characterization of differentially expressed exosomal microRNAs in bovine milk infected with Staphylococcus aureus

MicroRNAs (miRNAs) in milk-derived exosomes may reflect pathophysiological changes caused by mastitis. This study profiled miRNAs in exosomes from both normal milk and mastitic milk infected by Staphylococcus aureus (S. aureus). The potential targets for differentially expressed (DE) miRNAs were predicted and the target genes for bta-miR-378 and bta-miR-185 were also validated. Total RNA from milk exosomes was collected from healthy cows (n = 3, the control group) and S. aureus infected cows (n = 6, the SA group). Two hundred ninety miRNAs (221 known and 69 novel ones) were identified. Among them, 22 known and 15 novel miRNAs were differentially expressed. Target genes of DE miRNAs were significantly enriched in intracellular protein transport, endoplasmic reticulum and identical protein binding. The expression of two miRNAs (bta-miR-378 and bta-miR-185) with high read counts and log2 fold changes (> 3.5) was significantly higher in mastitic milk infected with S. aureus. One target gene (VAT1L) of bta-miR-378 and five target genes (DYRK1B, MLLT3, HP1BP3, NPR2 and PGM1) of bta-miR-185 were validated. DE miRNAs in exosomes from normal and S. aureus infected milk were identified. The predicted targets for two DE miRNAs (bta-miR-378 and bta-miR-185) were further validated. The linkage between the validated target genes and diseases suggested that we should pay particular attention to exosome miRNAs from mastitic milk in terms of milk safety.


Background
MicroRNAs (miRNAs) are short noncoding (~22 nucleotide in length), regulatory RNAs that modulate gene expression at the post-transcriptional level, mostly via binding to perfectly/partially complementary sites at the 3′-UTR of target mRNAs [1]. Among different body fluids, milk contains the highest amount of miRNAs [2]. Milk is an essential source of nutrients to all mammalian offspring. Bovine milk and dairy products have long traditions in human nutrition. In addition to providing nutrition, milk has long been known to protect the infant from infections and to play developmental functions integral to the infant, in which miRNAs are likely to be highly involved [3].
The majority of milk's miRNAs are transported and protected by the lipid bilayer of extracellular vesicles, predominantly exosomes of about 100 nm in diameter secreted by mammary epithelial cells [4]. Exosomes are cell-derived vesicles that are present in all biological fluids including blood, saliva, urine, amniotic fluid, bronchoalveolar lavage fluid and milk [5,6]. Milk exosomes have been reported in cows [7], buffalos [8], goats [9], pigs [10], marsupial tammar wallabies [11] and humans [12]. Exosomes protect miRNA molecules from effects of RNase digestion and low pH [13]. Thus, miRNAs in milk exosomes may be transferred into the gastrointestinal tract of infants and likely contribute to infant development and protection against infections [14].
Cells can take up exosomes through a variety of endocytic pathways, including clathrin-dependent endocytosis, clathrin-independent pathways such as caveolin-mediated uptake, macropinocytosis and phagocytosis [15]. Uptake of milk exosomes including their miRNAs has been demonstrated in human colon carcinoma Caco-2 cells and rat intestinal epithelial cell (IEC) IEC-6 cells [16]. Further, orally administered exosomes escaped re-packaging in the intestinal mucosa, and accumulated in liver and spleen. The same group later reported that labelled RNA derived from the milk exosomes was detected in mouse brain, kidney and lung [17]. Porcine milk exosomes promoted IEC proliferation in mice and increased mouse villus height, crypt depth and ratio of villus length to crypt depth of intestinal tissues were associated with miRNA-mediated gene regulatory changes in IECs [18]. In another study, oral delivery of bovine milk exosomes ameliorated experimentally induced arthritis [19]. Correctively, these data suggest that miRNA in milk exosomes can get into the body.
Accumulating evidence suggests that exosomal miRNAs play crucial roles in numerous diseases such as hepatocellular carcinoma [20], breast cancer [21] and Alzheimer's disease [22]. Secretion of milk exosomes is affected by bacterial infections in mammary glands. Staphylococcus aureus (S. aureus) is one of the most important etiologic agents for chronic bovine mastitis. Our previous in vitro study showed that 5 miRNAs (miR-2339, miR-21-3p, miR-92a, miR-23a and miR-365-3p) were up-regulated in bovine mammary epithelial cells when challenged with S. aureus [23]. In bovine mammary gland infected with S. aureus, a total of 77 miRNAs showed significant differences compared to the control group [24]. Previous studies have also investigated milk exosomal miRNAs following S. aureus induced bovine mastitis [25,26]. However, no study has focused on miR-NAs in exosomes derived from milk naturally infected with S. aureus. More importantly, previous studies focused on the profiling of miRNAs in milk exosomes, without experimental confirmation of the predicted target genes by bioinformatics. In addition, how miRNAs in exosomes affect milk safety has not been considered.
The objective of this study, therefore, was to characterize the miRNA expression profiles comprehensively in exosomes from normal and uninfected milk (the control group) and S. aureus infected milk (the SA group) and to predict potential targets for DE miRNAs and explore their possible functions.

Identification of S. aureus in bovine milk
Based on colony counting and PCR results for nuc and bacterial 16S rRNA genes, 13.95% (42/301) milk samples were infected with S. aureus. At the cow level, the infection rate was 31.58% (24/76) (Additional file 8: Table S1).

Isolation and detection of exosome from bovine milk
Bovine milk exosomes with approximately 100 nm in diameter were observed (Additional file 1: Figure S1a). Exosomes with a particle diameter in the range of 20 nm to 200 nm amounted to 84.1% of the total (Additional file 1: Figure S1b). The expression of CD63 and CD81 on the surface of exosomes was at a positive rate of 72.0 and 77.9%, respectively (Additional file 2: Figure S2).

Characterization of bovine milk exosomal miRNAs
Average RNA contents of the exosomes from 40 mL of the control or S. aureus infected milk samples were 1301 ± 38.7 ng (n = 3) and 1223 ± 56.6 ng (n = 6), respectively. Bovine milk exosomal RNA contained little or no 28S and 18S ribosomal RNA (data not shown).
The total raw read count from the sequencing of 9 libraries was 101,392,712 with an average of 11,265,857 reads per sample. After removing linker reads, reads containing N and poly A/T structure, length-anomalous reads, lowquality reads and reads greater than 35 nt or less than 17 nt, the resulting high quality clean data accounted for 83 to 96% of the original raw read counts. The majority of retained reads were 22 nt in length (Fig. 1a).
A total of 221 known and 69 novel miRNAs satisfied the conditions of having at least 1 transcript per million clean tags and were present in a minimum of four libraries. These 290 miRNAs were used for differently expressed (DE) analysis (Additional file 11: Table S4).

Target genes for bta-miR-378 and bta-miR-185 were validated
Read counts of DE miRNAs showed that bta-miR-378, bta-miR-185 and bta-miR-146b were 3 top DE miRNAs in the SA group as compared to the control group. Considering the potential importance based on read counts and log2foldchange values, the potential targets for bta-miR-378 and bta-miR-185 were further validated (Fig. 4).

Discussion
Milk provides important nutrients that are of benefit to most people throughout life. Due to the direct effects of the protein, fat, lipid, vitamin, and mineral fractions, milk has a specific growth promoting effect in children in both developing and developed countries [27].
Pasteurization is widely used in commercial milk production and it destroys all known pathogens and most of the spoilage bacteria in raw milk. Nowadays, there is compelling evidence that milk exosomes are retained in pasteurized commercial milk [17] and reach the systemic circulation and tissues of the human milk consumer [28]. Furthermore, pasteurization did not affect the profile expression of miRNA in bovine milk [29]. Bovine milk exosomes miRNAs which resist the harsh conditions in the gastrointestinal tract [30] are taken up via receptor mediated endocytosis by intestinal epithelial cells [16] and vascular endothelial cells [31]. More importantly, in vivo studies confirmed that milk exosomes miRNAs could reach distant tissues [19] and human plasma [32]. In our study, two miRNAs (bta-miR-378 and bta-miR-185) with high read counts were up-regulated significantly in exosomes from bovine milk infected by S. aureus. These two miRNAs have been reported to be associated with health. MiRNA-378 facilitates the development of hepatic inflammation and fibrosis [33]. In addition, expression of the miR-378 was reported to promote tumor growth [34]. MiR-185-5p may inhibit ameloblast and osteoblast differentiations and result in cleidocranial dysplasia [35], and promotes lung epithelial cell apoptosis [36]. How these miRNAs affect health parameters is not clear and it is plausible that their target genes are involved.
Consistent with two previous studies [25,26], the expression level of bta-miR-148a was the highest among all milk-derived exosomal miRNAs in our study. In these two previous studies, the miRNAs with the greatest differences in expression in milk-derived exosomes after S. aureus infection were bta-miR-142-5p [25] and bta-miR-223 [26], respectively. While the expression level of bta-miR-142-5p was also up-regulated significantly in our study, it was not the most differentially expressed one. In addition, expression of bta-miR-223 was not significantly changed in our study. These discrepancies between our study and the other two studies could be due to the fact that exosomes were isolated from mastitic milk naturally infected with S. aureus in this study, while the other two studies used milk samples from the mammary gland challenged with S. aureus.
We have functionally validated VAT1L as a target gene of bat-miR-378. Through a network-based analysis of three independent schizophrenia genome-wide association studies, Chang et al. reported that VAT1L may be one of the genes associated with schizophrenia [37]. In addition, DYRK1B, HP1BP3, MLLT3, NPR2 and PGM1 were identified as the target genes of bta-miR-185 in our study. Surprisingly, deficiency of these target genes also leads to a variety of diseases. DYRK1B belongs to the Dyrk family of proteins, a group of evolutionarily conserved protein kinases that are involved in cell differentiation, survival, and proliferation [38]. Mutations in DYRK1B were associated with a clinical phenotype that is characterized by central obesity, hypertension, type II diabetes and early-onset coronary artery disease [38]. HP1BP3 was identified as a novel modulator of cognitive aging and HP1BP3 protein levels were significantly reduced in the hippocampi of cognitively impaired elderly humans relative to cognitively intact controls [39]. Targeted knockdown of HP1BP3 in the hippocampus induced cognitive deficits [40]. MLLT3 gene is required for normal embryogenesis in mice, and an MLLT3 null mutation caused perinatal lethality [41]. A loss-offunction mutation of the AF9/MLLT3 gene was hypothesized to relate to neuromotor development delay, a Values represent the amount of expression of the miRNA in the SA group divided by the amount in the control group, followed by a base-2 logarithmic value $ miRNA only expressed in the SA group $$ miRNA only expressed in the control group cerebellar ataxia and epilepsy [42]. Homozygous inactivating mutations of NPR2 caused a severe skeletal dysplasia, acromesomelic dysplasia and Maroteaux type [43]. PGM1 deficiency has been described in a patient with myopathy and exercise induced hypoglycemia [44,45]. PGM1 deficiency causes a non-neurological disorder of glycosylation as well as a rare muscular glycolytic defect [46]. In addition to bta-miR-378 and bta-miR-185, several other miRNAs were also differentially expressed, including miR-1, miR-122, miR-1246, miR-142-5p, miR-146a,  miR-154, miR-184, miR-196 and miR-205. They were also associated with a variety of human diseases. For example, circulating miR-122 is strongly associated with the risk of developing metabolic syndrome and type II diabetes [47]. MiR-196 was overexpressed in the inflammatory intestinal epithelia of individuals with Crohn's disease [48]. While it is beyond the scope of this study to confirm the linkage between discussed miRNAs and health parameters, the above discussion about increased expression of certain miRNAs in exosomes from S. aureus infected milk argues strongly to be vigilant about the safety of mastitis milk, even after pasteurization.

Conclusions
In conclusion, we characterized the miRNA profiles in exosomes derived from control and S. aureus infected bovine milk, and 37 miRNAs (22 known and 15 novel) were significantly differentially expressed between the control group and SA group. This is the first report of functional validation of VAT1L, and DYRK1B, MLLT3, HP1BP3, NPR2 and PGM1 as target genes for bta-miR-378 and bta-miR-185, respectively. Finally, the potential safety hazards of mastitic milk were discussed, in the context of miRNAs within milk exosomes.

Milk sample collection and bacteria identification
Milk samples from seventy-six 3-to 4-year-old Holstein cows in the middle stage of lactation from 4 dairy farms (the Shaanxi Academy of Agricultural Sciences Farm, the farm of Delikang Dairy Co., Ltd., the Duzhai dairy farm and the Cuidonggou dairy farm) in Shaanxi Province were collected for this study with the approval of the Animal Use and Care Committee of Northwest A&F University (NWAFAC3751). Milk samples from all four quarters of each cow were aseptically collected and stored at − 80°C.
To select the samples for the control and the SA groups, 100 μL of each milk sample were plated on the Plate Count Agar (BD Diagnostics, Sparks, MD, USA) and incubated at 32°C for 48 h. The control group samples (n = 3) were randomly selected from the samples among which the colony count was zero. The milk samples with colony counts more than 1000 were marked as samples with bacterial infection for further detection. To rule out the interference caused by Escherichia coli (E. coli) infection for subsequent experiments, samples were cultured on the BactiCard™ E. coli (Thermo Oxoid Remel, Lenexa, USA). The milk samples without E. coli infection were selected for the identification of S. aureus by the Baird-Parker agar (Oxoid, Basingstoke, Hampshire, UK) as described previously [49]. Briefly, aliquots of individual milk samples were added to an equal volume of a double-strength enrichment broth (a trypticase soy broth supplemented with 10% NaCl and 1% sodium pyruvate) (Oxoid, Basingstoke, Hampshire, UK). After 24 h incubation at 35°C, the enrichment broth was streaked onto the Baird-Parker (Oxoid) agar containing 30% egg yolk with 1% tellurite (Oxoid) and onto the phenol red mannitol salt agar plates. Following 48 h incubation at 35°C, the colonies on the plates were counted, and one to three presumptive staphylococcal colonies from each plate were transferred to trypticase soy agar plates. Yellow colored colonies from the phenol red mannitol salt agar plates were assumed to be S. aureus. Further identification of these presumptive staphylococcal colonies was first based on conventional methods including Gram stain staining, colony morphology, a catalase test and a coagulase test with rabbit plasma. The culture result was further confirmed by using a polymerase chain reaction (PCR) assay targeting S. aureus-specific region of the thermonuclease gene (nuc) [49] and bacterial 16S rRNA genes. The samples with nuc positive S. aureus, which was also confirmed by sequencing the PCR products of 16S rRNA genes, were selected for the SA group (n = 6).

Preparation and purification of milk exosomes
Milk exosomes from the SA group (n = 6) and the control group (n = 3) were isolated by differential centrifugation as described previously [50]. Briefly, milk samples were centrifuged at 5000×g for 60 min at 4°C to remove milk fat and milk somatic cell. For the removal of casein and other cell debris, skimmed milk samples were subjected to three successive centrifugations at 4°C for 1 h each at 12,000×g, 35,000×g and 70,000×g (Beckman Coulter, USA). The whey was collected and centrifuged at 135,000×g at 4°C for 90 min (Beckman Coulter) to remove large particles and micro-vesicles. The supernatant was carefully collected and filtered through a 0.22 μm syringe-driven filter unit (Merck KGaA, Darmstadt, Germany). The percolate was collected and centrifuged at 150,000×g for 90 min at 4°C (Beckman Coulter). The exosome pellet was re-suspended in 1 mL sterile PBS and filtered through a 0.22 μm syringe-driven filter unit (Merck KGaA). Finally, exosomes were stored in aliquots of 200 μL at − 80°C until being used.

Identification of bovine milk exosomes
Dynamic light scattering analysis was used for analyzing nanoparticle sizes. Aliquot of 200 μL stored exosomes was diluted to 1 mL volume with sterile PBS stored on ice. The exosome solution was slowly injected into the sample cell of the Malvern Zetasizer Nano ZS90 (Malvern Panalytical Ltd., United Kingdom) system and measurements were taken according to manufacturer's instructions.
For transmission electron microscopy (TEM), milk exosomes were fixed in 3% (w/v) glutaraldehyde and 2% paraformaldehyde in a cacodylate buffer, pH 7.3. The fixed exosomes were then applied to a continuous carbon grid and negatively stained with 2% uranyl acetate. The samples were examined with a HT7700 transmission electron microscope (HITACHI, Japan).

Extraction of total RNA from bovine milk exosomes
Total RNA was extracted from bovine milk exosomes using the Trizol reagent (TAKARA, Japan) according to manufacturer's protocol and dissolved in RNase free water. Quality and quantity of RNA was examined using a NanoDrop 2000/2000C (Thermo Fisher Scientific, Waltham, MA, USA) and integrity was detected using agarose gel electrophoresis.

MiRNA library preparation and sequencing
Deep sequencing was performed on all 9 individual samples. For each library, 1 μg of high-quality RNA per sample was used as the input material for a small RNA library construction using the NEXTflex™ Small RNA Sequencing Kit V3 (Illumina, San Diego, CA) according to the manufacturer's instruction. Small RNA libraries were gel purified and pooled together in equimolar concentrations and subjected to 50 bp single read sequencing on an Illumina HiSeq 2500 system (Illumina, San Diego, CA). Read quality (adaptor removal and size selection) was assessed using FastQC v0.11.5 (http:// www.bioinformatics.babraham.ac.uk/projects/fastqc/) and the cutadapt [51].

Known miRNAs identification and novel miRNAs discovery
Identification of known miRNAs was performed with miRBase v21 (http://www.mirbase.org/) [52], while novel miRNA discovery was achieved with miRDeep2 v2.0.0.8 (https://github.com/rajewsky-lab/mirdeep2) [53]. The core and quantifier modules of miRDeep2 were applied to discover novel miRNAs in the pooled dataset of all the libraries while the quantifier module was used to profile the detected miRNAs in each library. The amount of miRNA expression was calculated by the transcript per million (TPM) metric, which is calculated as number of reads per miRNA alignment/number of reads of the total sample alignment× 10 6 . MiRDeep2 score > 1 was used as a cuff point for identification of novel miRNAs. Subsequently, a threshold of ≥1 TPM of total reads and present in ≥4 libraries was applied to remove lowly expressed miRNAs. MiRNAs meeting these criteria were further used in downstream analyses including differential expression analyses.
Differential miRNA expression and predicted target genes DE miRNAs were detected with DeSeq2 (v1.14.1) (https://bioconductor.org/packages/release/bioc/html/ DESeq2.html) [54]. Following normalization, miRNAs read counts in the SA group were compared with corresponding values in the control group. Significant DE miRNAs between the control and the SA groups were defined as having a Benjamini and Hochberg [55] corrected p-value< 0.05.
In order to investigate the potential functions of DE miRNAs, their target genes were predicted using the mi-Randa algorithm [56]. Predicted target genes with tot scores above 150 and tot energy below − 15 were further used for pathway analyses. The Database for Annotation, Visualization and Integrated Discovery (DAVID) was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway annotations of their target genes [57].
Firefly and Renilla luminescent signals arising from transfected cells were quantified according to the manufacturer's instructions using a Dual Luciferase assay system (Promega) with a Multilabel Counter luminometer (Varioskan Flash, Thermo Fisher Scientific). Renilla luciferase activities to firefly luciferase activities in cells transfected with an empty psiCHECK-2 vector without a 3′-UTR fragment was set to 100%. The experiment was repeated 3 times.

Statistical analysis
Data were analyzed with the SPSS 17.0 software (SPSS Inc., Chicago, IL, USA). The statistical significance between experimental groups was analyzed using One-Way ANOVA. p < 0.05 and p < 0.01 were defined to be statistically significant and extremely significant, respectively.