LncRNA Proling of Exosomes Derived From Staphylococcus Aureus-induced Bovine Mammary Epithelial Cell Lines MAC-T

Background: Exosomes (carrying proteins, miRNA, lncRNA, etc.) play a vital role in participating in the occurrence and development of mastitis. LncRNA can play a variety of regulatory roles by combining with protein, RNA, and DNA. The expression of mRNA and lncRNA in exosomes derived from bovine mammary epithelial cells infected by S. aureus is rarely understood. Methods: We used ultracentrifugation to separate exosomes derived from S. aureus- stimulated MAC-T cells. After extracting exosomal RNA, RNA sequencing was performed. We screened differentially expressed genes, mRNA and lncRNA, and predicted their functions through KEGG and GO analysis. 6 mRNA and lncRNA were randomly selected for RT-qPCR verication. Results: Analysis of the sequencing results showed that there were 186 differentially expressed genes, 431 differentially expressed mRNAs and 19 differentially expressed lncRNAs in the exosomes derived from S. aureus-infected and non-infected MAC-T acells. By predicting lncRNA target genes, it was found that 19 differentially expressed lncRNAs all acted on multiple mRNAs in cis and trans. GO analysis revealed that differentially expressed genes and lncRNA target genes played signicant roles in such metabolism, transmembrane transport, cellular response to DNA damage stimulus, response to cytokines. KEGG enrichment indicated that lncRNA target genes gathered in the TNF pathway, Notch pathway, phosphatidylinositol signaling system and p53 pathway. Conclusion: The data suggested that S. aureus induced changes in genes, mRNA, and lncRNA in exosomes secreted by MAC-T cells, which furtherly triggered other cells damage through exosomal delivery. This study would provide valuable resource for understanding the lncRNA information in exosomes derived from dairy cow mammary epithelial cells, and conduced to the study of S. aureus infection in dairy cow mammary glands.


Introduction
Mastitis is a common and widespread infectious disease in dairy farms all over the world, which characterized by in ammation of the mammary glands, resulting in a decrease in milk production and quality [1][2][3]. There are many kinds of pathogens of mastitis, and Staphylococcus aureus (S. aureus) is one of the important pathogens in subclinical mastitis [4]. S. aureus often causes clinical or subclinical mastitis, and the intramammary infections are usually characterized as mild, chronic or persistent [5].
Bovine mammary epithelial cells (BMEC) are not only essential for the synthesis of milk components, but also the main line of defense against pathogen invasion [3,6]. Milk contains a large number of exosomes, which are mainly derived from BMEC. The cargo carried by milk exosomes has been suggested to play a role in cell growth, development, immune regulation and regulation [7]. Exosomes are tiny extracellular vesicles with a diameter ranging from 40-150 nm, which originate from the internal vesicles of multivesicles in almost all cell types [8].
Typical exosomes are surrounded by phospholipid membranes, which contain proteins, lipids and nucleic acid molecules (including miRNA, cicrRNA, lncRNA, etc) [9]. With the deepening of the understanding of exosomes, exosomes are used as carriers to deliver drugs for targeted therapy. Because exosomes have been proven to tolerate the harsh environment in the intestine, are absorbed by various cell types, cross biological barriers and reach surrounding tissues [7]. Exosomes can also initiate a signal cascade of neighboring cells through the paracrine pathway. Therefore, it is very valuable to explore the RNA sequence in exosomes derived from MAC-T cells infected by S. aureus.
Long non-coding RNA (lncRNA) is a type of RNA that does not have the ability to encode proteins with more than 200 nucleotides in length [10]. LncRNA can bind to DNA, RNA, and proteins to play a role by mediating DNA methylation, mRNA degradation, chromatin remodeling, histone modi cation and protein modi cation [10,11].
There have been research evidences that lncRNA plays a key role in mastitis caused by S. aureus, Escherichia coli, Mycoplasma bovis and other pathogens though regulating the expression of speci c genes. LRRC75A-AS1 in bovine mammary epithelial cells and breast tissue was down-regulated under in ammatory conditions, and knocking out LRRC75A-AS1 in MAC-T cell line inhibited the adhesion and invasion of S. aureus and attenuated the activation of the NF-κB pathway [12]. LncRNA-TUB was higher expression in the mammary epithelial cells of dairy cows that received pro-in ammatory stimulation, and knockout of lncRNA-TUB indicated that it played a key role in the morphology, proliferation, migration and β-casein secretion of mammary epithelial cells [13].
Due to the spongy effect of exosomes, we speculated that the expression pro le of lncRNA had changed in the exosomes derived from MAC-T cells infected with S. aureus. The purpose of this study was to determine the expression pro le of lncRNA in exosomes derived from S. aureus-infected and non-infected MAC-T cells to provide new targets for the link between bovine mastitis caused by S. aureus and lncRNA for future research.

Exosomes isolation
Cell line exosomes were separation by ultracentrifugation as previously described [14]. The cell supernatant was sequentially subjected to 300 × g, 10 min; 2000 × g, 10 min; 10,000 × g, 30 min in a large-capacity high-speed centrifuge (AVANTI J-26XPN, Beckman, USA). After collecting the supernatant, ultra-high-speed centrifuge with SW32Ti rotor (optima XE-90, Beckman, USA) was hired to ultracentrifuge at 100,000 × g for 70 min. Pellet was collected and resuspended in PBS, and then perform ultracentrifugation again at 100,000 × g for 70 min. Finally, the exosomal pellet was resuspended in PBS and collected for subsequent analysis.

Transmission electron microscopy
Exosomes were mixed with an equal volume of 4% paraformaldehyde. 100 µl mixed liquid was sucked off and added dropwise to the copper net, then phosphotungstic acid was added for negative staining. After drying, transmission electron microscope was used to observe exosomal morphology at 80 kV.
Particle size distribution assay Exosomes obtained by ultracentrifugation was diluted by PBS to the optimal concentration for analysis. ZETASIZER Nanoparticle Tracking Analysis (NTA) instrument (Malvern, Worcestershire, England) was performed to analyze particle size of exosome. According to the Stokes-Einstein principle, dynamic light scattering detected Brownian motion of particles to generate trajectories and displacements to determine particle size distribution.

Flow cytometry analysis
Exosomes acquired by ultracentrifugation separation were incubated with CD63 (#557288, BD Bioscience) or CD81 (#555676, BD Bioscience), and exosomes without antibody incubation was used as a negative control. BD accuri C6 ow cytomenter was used for sample detection according to the operating instructions, and then analyzed with FlowJo software.

RNA extraction and lncRNA sequencing
Total RNA was extracted from exosomal samples using Trizol reagent (Invitrogen, Carlsbad, CA). Total RNA is used as input material for constructing cDNA library. Mainly, exosomal trace RNA puri cation, cDNA synthesis through reverse transcription, end repair and adaptor, PCR ampli cation and other steps to build a library, and then go through the library quality inspection and sequence on an Illumina HiSeq 4000 platform (Illumina Inc., San Diego, CA).

Sequence data quality control
The original image le (BCL) obtained by sequencing was converted into raw data in FASTQ format after base recognition. Quality analysis on the original data was performed to assess whether it was suitable for bioinformatics analysis. Quality analysis mainly included sequencing quality and base composition analysis. Clean data was obtained by removing linker sequence and low quality sequence from raw data. The quality of clean data was then evaluated using FASTQC v0.10.1.
Reference genome alignment and genome-wide reads distribution map Clean reads were obtained through ltering sequencing data by Trimmomatic (Version 3) and removing ribosomal RNA by Bowtie2 (v2.2.5). Bowtie2 (v2.2.5) was used to generate the index of the reference genome, and alignment of paired-end clean reads to the reference genome was performed with HISAT2 (v2.1.0). The assembly of mapped transcripts for each sample was performed with Scripture (beta2) and Cu inks (v2.1.1) in a reference-based approach [15].

Differential expression analysis
Quantitative analysis of genes and transcripts was performed via htseq v0.11.2. DESeq2 (1.18.1) is used to determine differential expression in transcripts or gene expression data. It is determined that the transcripts or genes with P value < 0.05 are differentially expressed in biological replication. Q value < 0.05 and |log2 (fold change)|<1 are considered as thresholds for signi cant differential expression in abiotic replication.

Prediction of lncRNA-gene interactions
In order to predict the targets of lncRNAs in MAC-T cells and thus understand potential functions of lncRNAs. Method as in previous research [15]. Firstly, we searched the 10 k/100 k coding genes upstream and downstream of lncRNA to nd cis role lncRNA. Then construction logarithm was modeled as log (P / (1-P)) = a + bx + cy, where a, b and c were model parameters, y represented the percentage of GC content in the sequence, and P represents lncRNA binding to a 10 kb long sequence the odds ratio. Predict statistical signi cance by using Wald test on Z statistics (p value < 0.05). Only lncRNAs with p value < 0.05 and z score > 0 were annotated as trans lncRNA.

Kyoto encyclopedia of genes and genomes (KEGG) biological pathway enrichment and Gene ontology (GO) functional enrichment analysis
The GOseq R package was used to analyze GO enrichment for differentially expressed genes or lncRNA target genes [15]. GO terms with Q values < 0.05 were considered signi cantly enriched by differentially expressed genes.
The summary differential expression genes and lncRNA target genes in KEGG pathways enrichment was analyzed by KOBAS (v3.0) software.

RT-qPCR validation
Total exosomal RNA samples from S. aureus-infected and non-infected MAC-T cells for the lncRNA-seq were analyzed by RT-qPCR. Then cDNA was obtained by reverse transcriptase reagent according to the instruction manual. The qPCR was performed using the SYBR Green Plus Reagent Kit with the Light Cycler 96 instrument (Roche, Basel, Switzerland) following the instructions of the manufacturer. GAPDH was used as an internal reference gene. The fold change of mRNA transcripts and lncRNAs was calculated using 2 −∆∆Cq methods.

Statistical analysis
IBM SPSS 20 was used to perform statistical analyses. One-Way Analysis of Variance (ANOVA) was performed to detect statistically differences in lncRNAs and mRNAs between S. aureus infected and non-infected groups. p < 0.05 and p < 0.01 were considered as signi cant differences.

Identi cation of exosomes
In order to directly observe the morphology of exosomes, we employed transmission electron microscopy to observe negatively stained exosomes samples. As shown in Fig. 1A, it had a very obvious single-layer membrane structure under the electron microscope, which appeared as a saucer-shaped or disc-shaped vesicle with one side recessed. NTA analysis was performed to further determine the size of exosomes. The main peak of the particle size of the exosomes we separated was 70.55 nm, and the 80.3% particle size was between 40-150 nm (Fig. 1B). Flow cytometry was used to detect the surface characteristic markers of exosomes. The results revealed that the positive rates of CD81 and CD63 were 92.7% and 85.9%, respectively (Fig. 1C).

Sequencing data statistics of exosomes derived from MAC-T cells
In the current study, a total of 4 cDNA libraries were constructed by isolating total RNA from exosomes derived from S. aureus-infected and non-infected MAC-T cells. The raw reads in the cDNA library of different samples were cont1, 67572136; cont2, 66474707; infect1, 69033293; infect2, 70292999. GC content (%) percentages were found to be cont1, 48.61; cont2, 47.63; infect1, 49.75; infect2, 49.98 (Table 1). We conducted quality control analysis on raw data, and the results are as follows. The quality of the bases in raw reads was counted in each sequencing cycle. The base composition of each cycle is shown in Fig. 2A. As shown in Fig. 2B, the blue line represented the average base quality of the cycle in the base quality distribution box plot. Base error rate is limited to 0.11-0.13 in raw reads (Fig. 2C). We ltered raw data to get clean data, and clean reads, GC%, Q20, Q30, clean bases were counted in Table 2, and the quality control results of clean reads were shown in Fig. 2D, E and F. In addition, when the clean reads were aligned with the bovine reference genome using HISAT2 software by an improved BWT algorithm (FM index), it was determined that the Total Mapped Reads or Fragments belonging to all the samples was above 80% and was mapped in the reference genome (Table 3).

Distribution of reads in the chromosomes and known types of RNAs
After mapped in the reference genome, we counted the position information of the genomes corresponding to all reads to evaluate the depth of sequencing data coverage. The distribution of reads on each chromosome was displayed in Fig. 3A. Subsequently, the reads aligned on the chromosome were annotated to exonic, intronic, and intergenic. Reads compared to the genome was counted the distribution (Fig. 3B). And the reads distributions in the known RNA types were shown in Fig. 3C.

Differentially expressed genes analysis in exosomes derived from S. aureus-infected and non-infected MAC-T cells
DESeq2 was used to screen differentially expressed genes. We selected the differentially expressed genes between samples through the two levels of multiple of difference (|log2(Fold Change)|>1) and signi cance level (Q value < 0.05) [16][17][18]. The statistics of the number of differentially expressed genes in exosomes derived from S. aureusinfected and non-infected MAC-T cells were shown in Fig. 4A. There were a total of 186 genes expression disorders, of which 31 genes were signi cantly up-regulated and 155 genes were signi cantly down-regulated. We used the RPKM/COUNT of the differential gene as the expression level. Hierarchical Clustering analysis was performed to cluster genes with the same or similar expression patterns into clusters, and used different colored regions to represent different clustering grouping information to determine clustering mode of different sample control modes [15,19]. In this study, the hierarchical clustering analysis of differentially expressed genes is shown in Fig. 4B.
GO enrichment analyses were conducted to search for the biological processes, cellular components, and molecular functions of differentially expressed genes in more detail [15,20]. GO analysis was performed on the up-regulated and down-regulated genes. According to the enrichment point, all disordered genes are classi ed according to GO terms ( Fig. 4C and D). Moreover, KEGG pathway analysis results were calculated according to the enrichment points, and the best gene pathways related to the up-regulated and down-regulated genes were listed in Fig. 4E and F.

Differentially expressed mRNAs and lncRNAs analysis in exosomes derived from S. aureus-infected and noninfected MAC-T cells
The expression level of the transcript was calculated using the RPKM value. We found that there were 78 upregulated and 353 down-regulated known mRNAs, and 8 up-regulated and 11 down-regulated lncRNAs ( Fig. 5A  differentially expressed mRNAs and lncRNAs were also clustered in the hierarchical cluster analysis map (Fig. 5C, D).

Prediction of lncRNA target genes and mRNAs
The mode of action of lncRNA regulating target genes could be divided into two categories: cis regulation (In Cis) and trans regulation (In Trans). We made predictions about how the differentially expressed lncRNAs regulated the target genes (Additional le 1 & 2). Differentially expressed lncRNAs were selected to draw a network diagram of the interaction between these lncRNAs and their target genes (Fig. 6), so as to provide reference and help for the overall analysis of the functions of lncRNAs in samples.
GO enrichment and KEGG enrichment analysis of lncRNA target genes and mRNAs GO analysis was performed on the target gene, and take the P < 0.05 of GO annotation enrichment as the signi cance threshold. GO terms of biological processes, cellular components and molecular functions were shown 10 GO terms, respectively (Fig. 7A). The Additional le for detailed biological processes, cell components and molecular functions were seen in Additional le 3, 4 and 5. KEGG analysis was performed on the target gene, and use the P < 0.05 of the enrichment degree of KEGG pathway as the signi cance threshold to obtain the analysis result. KEGG enriched target gene pathways were classi ed according to environmental information processing, human diseases, metabolism, and biological systems. It was found that lncRNA target genes were most closely related to human diseases (Fig. 7B). KEGG pathway analysis results are shown in the bubble diagram of the KEGG pathway (Fig. 7C). Multiple target genes are involved in the TNF signaling pathway, Notch signaling pathway, Phosphatidylinositol signaling system and p53 signaling pathway in the KEGG target gene enrichment pathway.
Results of RT-qPCR validation was in accordance with RNA-Seq The results of qPCR were found to be in accordance with RNA-seq (Fig. 8A, B), suggesting that RNA-seq results were reliable.

Discussion
Mastitis is an in ammation of the udders of dairy cows that is almost always caused by pathogenic microorganisms [21]. It seriously harms dairy farming and the dairy industry. Cow mammary epithelial cells can produce exosomes and excrete with milk during lactation, and exosomes can also communicate between cells through paracrine pathways. The cargo carried in the exosomes is the switch that initiates communication. The aim of this study was to investigate the expression of mRNA and lncRNA in exosomes derived from bovine mammary epithelial cells infected by S. aureus (Fig. 9).
In this study, we used ultracentrifugation to separate the exosomes in the MAC-T cell culture medium. We directly observed saucer-shaped vesicles through transmission electron microscopy. The NTA instrument detected 80% of the extracellular vesicles with a size of 40-150 nm. In addition, the positive rates of CD81 and CD63, which are characteristic markers on the surface of exosomes, were more than 80% by ow cytometry. All the above data implied that exosomes were successfully isolated.
After the total exosomal RNA extracted, high-throughput sequencing was performed on the Illumina platform to reveal the lncRNA and gene expression pro les in the exosomes derived from S. aureus infected and non-infected MAC-T cells. It was found that186 differentially expressed genes, 431 differentially expressed mRNAs and 19 differentially expressed lncRNAs were screened out. The function of most of these lncRNAs was unclear. Usually, we predicted the function of lncRNA based on its closely related coding genes [22]. In this study, we rst estimated the target genes of differentially expressed lncRNA using cis and trans, then GO and KEGG analysis were performed on the target gene to nd the potential role of lncRNA [15]. We found that these differentially expressed target genes were mainly related to in ammation or apoptosis signaling pathways through enrichment analysis of GO and KEGG pathways and coding-non-coding co-expression network analysis.
Differentially expressed lncRNA XR_001500758.2 and XR_003029422.1 were enriched in TNF signaling pathway, Notch signaling pathway, Phosphatidylinositol signaling system p53 signaling pathway, RIG-I-like receptor signaling pathway, MAPK signaling pathway and PI3K-Akt signaling pathway. These signaling pathway are closely related to in ammation or apoptosis. The target genes predicted for XR_001500758.2 included ATM, FAS, SERPINB5, PPIP5K2, SNW1, MLKL, IL1A, PLA1A. The target genes predicted for XR_003029422.1 included MAML2, GORAB, DGKH, DGKI, SELE, ITGA1. ATM protein kinase is a product encoded by the telangiectasia ataxia mutated gene ATM, which senses DNA damage, transmits DNA damage signals to downstream target genes, initiates stress systems and produces cycle arrest, cell repair and apoptosis [23,24]. FAS, adapter proteins associated with the death domain (FADD and TRADD) belongs to TNFR family. FAS could activate NF-kappa B, MAPK and PI3K-Akt signaling pathway. Based on predictions, we found that lncRNA had multiple targets and was involved in the complex process of regulating the immune response.

Conclusion
In conclusion, the molecular mechanism of S. aureus infecting dairy cow mammary epithelial cells involves the paracrine mode of exosomal cargo lncRNA playing a regulatory role in the TNF pathway, Notch pathway, Exosomal particle size analysis was measured by NTA instrument. (C) Exosomal characteristic markers CD81 and CD63 were detected by ow cytometry. NC was negative control, P1 was gate.    Differentially expressed mRNAs and lncRNAs. (A) Volcano plot of differential expression mRNA transcripts. (B) Volcano plot of differential expression lncRNA transcripts. Red color represents genes with higher expression, green color represents genes with lower expression. (C) Clustering of differential expressed mRNA. (D) Clustering of differential expressed lncRNA. Hierarchical clustering based on FPKMs, where log10(FPKM+1) is used for clustering. Red color represents genes with higher expression, while blue color represents genes with lower expression.
Page 18/21 Figure 6 Network diagram of the interaction between lncRNA and target genes.