Characterization of microRNA expression in B cells derived from Japanese black cattle naturally infected with bovine leukemia virus by deep sequencing

Bovine leukemia virus (BLV) is the causative agent of enzootic bovine leukosis (EBL), a malignant B cell lymphoma. However, the mechanisms of BLV-associated lymphomagenesis remain poorly understood. Here, after deep sequencing, we performed comparative analyses of B cell microRNAs (miRNAs) in cattle infected with BLV and those without BLV. In BLV-infected cattle, BLV-derived miRNAs (blv-miRNAs) accounted for 38% of all miRNAs in B cells. Four of these blv-miRNAs (blv-miR-B1-5p, blv-miR-B2-5p, blv-miR-B4-3p, and blv-miR-B5-5p) had highly significant positive correlations with BLV proviral load (PVL). The read counts of 90 host-derived miRNAs (bta-miRNAs) were significantly down-regulated in BLV-infected cattle compared to those in uninfected cattle. Only bta-miR-375 had a positive correlation with PVL in BLV-infected cattle and was highly expressed in the B cell lymphoma tissue of EBL cattle. There were a few bta-miRNAs that correlated with BLV tax/rex gene expression; however, BLV AS1 expression had a significant negative correlation with many of the down-regulated bta-miRNAs that are important for tumor development and/or tumor suppression. These results suggest that BLV promotes lymphomagenesis via AS1 and blv-miRNAs, rather than tax/rex, by down-regulating the expression of bta-miRNAs that have a tumor-suppressing function, and this downregulation is linked to increased PVL.


Introduction
In this study, we performed comparative analyses of B cell miRNA expression in cattle uninfected and naturally infected with BLV. In these cattle, the relationships between the miRNA expression and BLV proviral load, tax/rex gene expression, and AS1 gene expression were investigated.

Materials and methods
Blood and tissue sample collection, serum isolation, and DNA/RNA extraction Blood was collected from the jugular vein of 16 BLV-naturally infected and 6 BLV-uninfected Japanese Black cattle, bred at the Iwate University Field Science Center. The BLV provirus genomes of all cattle were examined by quantitative PCR (qPCR), as described below, and ELISA using anti-BLV antibodies according to the manufacturer's instructions (JNC Inc., Tokyo, Japan). Lymphoma tissues were also obtained from five cattle diagnosed with EBL at the Iwate University Field Science Center. Details of animals used in this study are shown in Table 1. All procedures and animals used in this study were approved by the Iwate University Animal Care and Use Committee (no. A201704).
Genomic DNA was extracted from EDTA-treated whole blood with a magLEAD 1 12gC instrument (Precision System Science, Chiba, Japan) immediately after the blood collection. RNA was extracted from whole blood collected in PAXgene Blood RNA tubes (PreAnalytix, Hombrechtikon, Switzerland) and stored at -70˚C for several months after the blood collection. RNA was also extracted from the bovine B cell leukemia cell line KU-17 [27] with TRIzol Reagent (Invitrogen, Carlsbad, CA, USA). These DNA and RNA extraction procedures were performed according to manufacturer's instructions. Extracted DNA and RNA were stored at -20˚C and -70˚C, respectively, until analyzed.

Isolation of B cells from peripheral blood mononuclear cells
EDTA-treated whole blood was layered over 60% percoll (GE Healthcare, Tokyo, Japan) in Leucosep tubes (Greiner Bio-One, Kremsmunster, Austria) and the peripheral blood mononuclear cells (PBMCs) were isolated via density gradient centrifugation for 20 min at 1,000 g. The isolated cells (10 8 cells) were incubated with 1,000 μL of anti-bovine IgM mouse monoclonal antibody (diluted 1:100) (BIG73A; VMRD, Pullman, WA, USA), diluted with MACS buffer [2 mM EDTA, 0.5% BSA in PBS (pH 7.2)], at 4˚C for 15 min. The cells were then incubated with anti-mouse IgG microbeads (Miltenyi Biotec, Bergisch Gladbach, Gemany) at 4˚C for 15 min. The cells were passed through a cell strainer (EASY strainer; pore size 40 μm, Greiner Bio-One) and applied to a MACS LS column (Miltenyi Biotec) in the magnetic field of a MACS separator (Miltenyi Biotec). After washing three times with MACs buffer, the column was removed from the MACS separator and the magnetically labeled cells were flushed into a collection tube. Approximately 3 x 10 7 -7 x 10 7 PBMCs were recovered.

MicroRNA library preparation
Total RNA, containing miRNA, was extracted from B cells (10 7 cells per animal) and lymphoma tissues using miRNeasy Mini Kits (Qiagen K.K., Tokyo, Japan). The RNA integrity number (RIN) was determined on an Agilent RNA 6000 Nano Bioanalyzer (Agilent Technologies, Santa Clara, CA, USA). The 3' and 5' adaptors were ligated to the total RNA extracted from the isolated B cells with TruSeq Small RNA Library Preparation Kits (Illumina, San Diego, CA, USA). For 3' adaptor ligation, total RNA was incubated at 70˚C for 2 min and then transferred to ice. Subsequently, the following reagents were added to the mixture: 5 μL of 1 μg total RNA, 1 μL of RNA 3' adaptor, 2 μL of Ligation Buffer, 1 μL of RNase Inhibitor, 1 μL of 10× T4 RNA Ligase 2, Deletion Mutant (Epicentre, Madison, WI, USA). The reaction was incubated at 28˚C for 60 min, 1 μL of stop solution (Stop oligo) was added, and the reaction mixture was then incubated at 28˚C for 15 min. For 5' adaptor ligation, a 5' RNA adaptor was denatured by heating at 70˚C for 2 min and was then transferred on ice. The following reagents were added to the 3' adaptor ligation mixture: 1 μL of RNA 5' adaptor, 1 μL of ATP (10 mM), 1 μL of T4 RNA ligase, and 11 μL  Table 2.

Reverse transcription of adapter ligation products
RNA RT Primer (1 μL) was added to 6 μL of the adaptor ligation mixture, described in the previous section, heated at 70˚C for 2 min, and then immediately placed on ice. The following reagents were then added to the ligation mixture: 2 μL of 5× First Strand Buffer, 0.5 μL of 12.5 mM dNTP mix, 1 μL of 100 mM DTT, 1 μL of RNase Inhibitor, and 1 μL of SuperScript II Reverse Transcriptase (Thermo Fisher Scientific, Waltham, MA, USA). The reaction mixture was then incubated at 50˚C for 60 min. The sequence of the RNA RT Primer is shown in Table 2.

PCR amplification and purification of PCR products
The following reagents were added to 12.5 μL of reverse transcription reaction mixture described in the previous section: 25 μL of PCR mix, 2 μL of miRNA PCR primer, 2 μL of miRNA PCR primer Index, and 8.5 μL of nuclease-free purified water to make the total reaction mixture up to 50 μL. The PCRs were performed under the following conditions: initial denaturation at 98˚C for 30 s; followed by 15 cycles of heat denaturation at 98˚C for 10 s, annealing at 60˚C for 30 s, and extension at 72˚C for 15 s; then a final extension at 72˚C for 1 min. The sequences of the primers used are shown in Table 2. The PCR products (145 bp to 160 bp) were purified on 6% Novex TBE gels (Life Technologies, Waltham, MA, USA), following the manufacturer's instructions. The PCR products were evaluated with a microchip based capillary electrophoresis system (MultiNA, Shimadzu, Tokyo, Japan), and the concentrations were measured on a Qubit fluorometer (Thermo Fisher Scientific).

MicroRNA deep sequencing and analysis
MicroRNA analysis was performed using a MiniSeq Sequencing System (Ilumina). The libraries were diluted to 1 nM with 10 mM Tris HCl (pH 8.5) and made up to 5 μL each, to which, 5 μL of 2-fold diluted sodium hydroxide solution (Fluka Analytical, St. Gallen, Switzerland) was added. The mixture was incubated for 5 min at room temperature, and then 5 μL of 200 mM Tris HCl (pH 7) was added and the reaction mixture kept on ice. The mixed library reaction was diluted to 1.8 pM with hybridization buffer (Ilumina), and then 500 μL of the reaction mixture was applied to a MiniSeq High Output Reagent Cartridge (Ilumina). Deep sequence analysis was performed according to the manufacturer's recommended protocol for small RNA sequencing. Afterwards, sequencing reads were processed with CLC Genomics Workbench software (Ver. 9.5.5; Qiagen KK) to obtain the final miRNA counts for each sample (see Qiagen tutorial manual for small RNA Analysis using Illumina Data for detail; https:// resources.qiagenbioinformatics.com/tutorials/Small_RNA_analysis_Illumina.pdf). Briefly, adapter sequences were removed from the partial adapter sequences of the FASTQ file. The adapter trimming parameters were set to default values; i.e., mismatch cost and gap cost were 2 and 3, respectively; match threshold was selected to "allow end matches"; and the minimum score at the end was set to 6. Subsequently, for sequence filtering, the minimum and maximum length values were used as default values; i.e., reads below length of 15 and above length of 55 were discarded, and the sample threshold for the minimum sampling count was set to 1. The number of copies of each of the resulting small RNAs was counted. To annotate the small RNA sample, the bovine miRNA database in miRBase 22 [28] (http://www.mirbase.org/blog/ 2018/03/mirbase-22-release/) was downloaded. The trimmed sequences were compared VIC-TCTGTCCATTGTCCC-MGB For mRNA and miRNA expression

Quantification of BLV provirus
We performed duplex quantitative PCR (qPCR) that targeted the BLV tax/rex gene region, and bovine RPPH1 gene as an internal control. The qPCR was performed under the following conditions: initial denaturation at 95˚C for 20 s, followed by 40 amplification cycles of denaturation at 95ºC for 1 s and annealing/extension at 60ºC for 20 s. Reaction mixtures consisted of 5 μL of template genomic DNA derived from whole blood, 10 μL of Premix Ex Taq Table 2. The qPCRs were performed on a QuantStudio™ 3 Real-Time PCR System (Applied Biosystems, Life Technologies, Foster City, CA, USA). Standard curves were generated by creating 10-fold serial dilutions of standard plasmids that contained the relevant BLV tax/rex or bovine RPPH1 genes, amplified by the appropriate PCR primers. The standards for calibration ranged from 10 0 to 10 5 copies/reaction and were run in duplicate. The number of BLV copies was indicated as proviral load per 10 ng DNA. The percent of BLV-infected cells was calculated by the following equation (as there were two copies of the RPPH1 gene per cell): ½% of BLV-infected cells ¼ BLV tax=rex copy number � ðRPPH1 copy number � 2Þ � 100�:

Quantification of mRNA and miRNA expression by quantitative RT-PCR (qRT-PCR)
SYBR Prime Script RT-PCR Kits (Takara Bio) were used for tax/rex, AS1, and bovine HMG box-containing protein 1 (HBP1) mRNA and miScript II RT Kits (Qiagen KK) were used for bta-miR-375. For tax/rex, AS1, and bovine HBP1 mRNA, reverse transcription reaction mixtures consisted of 400 ng/μL of template RNA, derived from whole blood; 4 μL of 5× prime-Script Buffer; 1 μL of 50 μM Oligo-dT primer; 1 μL of 100 μM random 6-mer primer; 1 μL of PrimeScript RT Enzyme Mix 1; and RNase-free water to make the total reaction volume up to 10 μL. The reaction was incubated at 37˚C for 15 min, and then heated at 85˚C for 5 s for enzyme inactivation and placed on ice. For bta-miR-375, reverse transcription reaction mixtures consisted of 300 ng of template RNA, derived from isolated B cells and lymphoma tissues; 4 μL of 5x miScript HiFlex Buffer; 2 μL of 10x miScript Nucleic Mix; 2 μL of miScript Reverse Transcriptase Mix, and RNase-free water up to a total reaction volume of 20 μL. The reaction mixture was incubated at 37˚C for 60 min, and then heated at 95˚C for 5 min for enzyme inactivation and placed on ice. The concentrations of cDNAs obtained were calculated by absorbance at 260 nm on a NanoDrop One (Thermo Fisher Scientific K.K).
Quantitative RT-PCRs that targeted tax/rex and AS1 mRNAs were performed with GAPDH mRNA as the internal control. Bovine HBP1 mRNA was targeted with beta actin (ACTB) mRNA as an internal control [29], and bta-miR-375 was targeted with bta-miR-16a-5p as an internal control. For tax/rex, AS1, and bovine HBP1 mRNA, reaction mixtures consisted of 40 ng/5 μL of cDNA, 0.8 μL each of 10 μM forward and reverse primers, 10 μL of SYBR Premix Ex Taq (Takara Bio), 0.4 μL of Rox Reference Dye (Takara Bio), and 3 μL of sterilized ultrapure water. PCRs were performed under the following conditions: initial denaturation at 95˚C for 30 s, 40 cycles of denaturation at 95˚C for 5 s and annealing/extension at 60˚C for 30 s. The gene copy number was calculated via the standard curve method. For bta-miR-375, miScript SYBR green PCR Kits (Qiagen KK) were used. Reaction mixtures consisted of 3 ng of cDNA, 2.5 μL of 10x miScript Universal Primer, 2.5 μL of microRNA-specific primer, 12.5 μL of 2x QuantiTect SYBR Green PCR Master Mix, and sterilized ultrapure water up to a total reaction volume of 25 μL. PCRs were performed under the following conditions: initial denaturation at 95˚C for 15 min; 40 cycles of denaturation at 94˚C for 15 s, annealing at 55˚C for 30 s, and extension at 70˚C for 30 s. The relative miRNA expression levels were calculated using the ΔΔCT comparative method by Quantstudio™ design and analysis software (Version 2.4, Thermo Fisher Scientific). The sequences of the primers used are shown in Table 2.

Statistical analysis
Differences in expression of BLV miRNAs (blv-miRNAs) and bovine miRNAs (bta-miRNAs) between BLV-infected and BLV-uninfected cattle were assessed using Mann-Whitney test. Correlations between parameters in BLV-infected cattle were assessed by the Spearman's correlation coefficients. Differences in bta-miR-375 expression among BLV negative, BLV-positive, and EBL cattle were assessed by Kruskal-Wallis test with Steel-Dwass post-hoc test. Differences in expression of BLV miRNAs (blv-miRNAs) and bovine miRNAs (bta-miRNAs) between BLV AS1 high expression cattle and low expression cattle were assessed using a Mann-Whitney test. These data analyses were performed by R, a language and environment for statistical computing (R Core Team, 2020. URL https://www.R-project.org/). Statistical significance was determined as p < 0.05.

MicroRNA sequencing reads in B cells of BLV-infected and uninfected cattle
B cells were isolated from 16 BLV-infected and 6 uninfected cattle at purity levels between 82% and 97%. The RINs of RNA samples derived from B cells were 6.9 to 10. The numbers of miR-NAs that were read in these RNA samples were between 1.33 × 10 6 and 4.12 × 10 6 . Among these miRNAs, 614 bovine-derived miRNAs (bta-miRNAs) were detected out of 1,064 currently registered in the database (miRBase) (S1 Table). In addition, the 10 BLV provirusderived miRNAs (blv-miR), which were previously reported [22], were also detected (S2 Table).

Correlation between miRNA copies and BLV tax/rex and AS1 expression in BLV-infected cattle
The expression of BLV tax/rex and AS1 genes in PBMCs were quantified by qRT-PCR. BLV tax/rex and AS1 mRNA copy numbers were correlated against blv-miRNAs and bta-miRNAs copies. BVL-infected cattle had between 1.6 and 91 BLV tax/rex mRNA copies per 10 4 B cells and 2.1 to 1,388 AS1 mRNA copies per 10 4 B cells. There was no significant correlation between tax/rex and AS1 mRNA expression (S1A and S1B Fig).

Discussion
In this study, we performed deep sequencing analysis to comprehensively compare miRNAs expressed in B cells derived from BLV-infected healthy cattle and those derived from BLVuninfected cattle and determined the correlations between B cell miRNAs and the pathogenesis of BLV. Furthermore, the correlations between B cell miRNAs and BLV proviral load (PVL), BLV tax/rex and AS1 mRNA expression were also investigated.
Ten BLV provirus-derived microRNAs (blv-miRNAs) were detected in B cells derived from BLV-infected cattle, and these blv-miRNAs accounted for 38% of all detected miRNAs. This is in agreement with a study that reported that approximately 40% of total miRNAs were  blv-miRNAs in B cell lymphoma derived from sheep experimentally infected with BLV [22]. These results suggest that blv-miRNAs are constantly expressed at a high rate in B cells in healthy BLV-infected cattle, from the asymptomatic stage to the onset of EBL. The risk of EBL onset in BLV-infected cattle harboring higher PVLs is higher than that in BLV-infected cattle harboring lower PVLs [98]. The read counts of four blv-miRNAs (blv-miR-B1-5p, blv-miR-B2-5p, blv-miR-B4-3p, and blv-miR-B5-5p) had strong positive correlations with PVL. Among these blv-miRNAs, blv-miR-B4-3p had the highest sequencing reads, which is in agreement with a previous study that showed blv-miR-B4-3p had the highest number of read counts in B-cell lymphoma derived from sheep experimentally infected with BLV [22]. The blv-miR-B4-3p has seven bases in common with the 5' flanking region of the host genome-derived miR-29 family (miR-29a, miR-29b, and miR-29c) [19], suggesting that blv-miR-B4-3p and the miR-29 family have similar functions. The miR-29 family is involved in cell proliferation, apoptosis, angiogenesis, and metastasis in a variety of human tumor cells [46,99]. The blv-miR-B4-3p also promotes cell proliferation by down-regulating the expression of a transcription repressor HMG box-containing protein 1 (HBP1), which suppresses the cell cycle of ovine malignant B cell lymphoma in vitro [19,22,100]. However, the results of this study show that HBP1 gene expression was not decreased in B cells derived from healthy BLV-infected cattle. These results suggest that BLV provirus-derived miRNAs, including blv-miR-B4-3p, modulate proliferation and apoptosis of BLV-infected B cells in an HBP1-independent manner and contribute to the increased PVL seen prior to the onset of EBL.
Of the 49 bta-miRNAs that had significant differences in their read counts between BLVinfected and uninfected cattle, 32 bta-miRNAs significantly correlated with PVL (Table 4); 31 bta-miRNAs had a negative correlation with PVL and 1 (bta-miR-375-3p) had a positive correlation with PVL. Of the 31 miRNAs that had negative correlations with PVL, 3 bta-miRNAs (bta-miR-28-5p, bta-miR-150-5p and bta-miR-197-3p) had a strong negative correlation (r < -0.7, p < 0.01). MiR-28 controls cell proliferation, is down-regulated in B-cell lymphomas [75], and reduces HTLV replication and infection [101]. The role of miR-150 in human cancer is context-dependent as this miRNA can have either oncogenic or tumor suppressor activity in cells that originate from different tissues. This is highlighted by the upregulated expression of miR-150 in B cells from chronic lymphocytic leukemia (CLL) [102,103] but downregulated expression in chronic myeloid leukemia [104,105] and mantle cell lymphoma [106]. MiR-197 functions as a tumor suppressor in multiple myeloma and hepatocellular carcinoma and as a key repressor of the p53-dependent apoptotic cascade in lung cancer [67,68,107]. Moreover, 12 of the 13 bta-miRNAs with relatively weaker negative correlation coefficients (-0.7 < r < 0.6) function as tumor suppressors and/or oncogenes. In particular, miR-146a has been deregulated in HTLV-1-transformed T-cells [108]. Taken together, these results suggest that increased PVL down-regulates the expression of bta-miRNAs, the majority of which have functions involved in suppressing cell proliferation and viral replication.
Little is known about the function of the BLV AS1 gene, which is encoded by the minus strand of BLV provirus. AS1 transcripts are not present in the cytoplasm and AS1 protein has not been identified [20], suggesting that the AS1 gene functions as transcripts (RNA), but not as protein. In this study, there was a positive correlation between AS1 transcripts and the expression levels of two of the five blv-miRNAs (blv-miR-B1-5p and blv-miR-B2-5p), whereas there was a positive correlation between tax/rex transcripts and the read counts of the five blv-miRNAs. Furthermore, the correlation coefficient with AS1 was lower than that with tax/rex. Although the reason is unknown, the interaction with bta-miRNAs might be different between tax/rex transcribed from the 5' flanking region and AS1 transcribed from the 3' flanking region. In addition, blv-miR-B1-5p and blv-miR-B2-5p might more strongly influence AS1 transcription than do blv-miR-B4-3p, blv-miR-B2-3p, and blv-miR-B5-5p.
Of the 31 bta-miRNAs that had a negative correlation with PVL, 24 bta-miRNAs also had a negative correlation with AS1 mRNA expression, and the majority of the bta-miRNAs function as tumor suppressors or oncogenes (Table 4). In particular, bta-miR-150-5p, an important tumor suppresser of leukemia/lymphoma that targets Nanog (a homeobox transcription regulatory factor involved in stem cell pluripotency) [36,112], had a strong negative correlation with AS1 mRNA expression. MiR-150 is expressed at high levels in mature T and B cells, is downregulated in regulatory T cells (Tregs) through the action of Foxp3 [113], is downregulated in HTLV 1-infected cells, and is upregulated in adult T cell leukemia/lymphoma (ATLL) cells [114,115]. Our data for BLV-infected B cells is consistent with the down-regulation of miR-150 in HTLV 1-infected cells. In addition, the 3' UTR of HIV-1 mRNA is targeted by miR-150 and miR-28, and these interactions influence the ability of the virus to effectively infect CD4 + T cells, monocytes, and macrophages [116,117]. MiR-150 specifically targets the signal transducer and activator of transcription 1 (STAT1) 3' UTR, reducing STAT1 expression and dampening STAT1-dependent signaling in human T cells [118]. HTLV-I-transformed and ATL-derived cells have reduced levels of miR-150 expression, which coincides with increased STAT1 expression and STAT1-dependent signaling. STAT1 plays a role in immune modulatory functions, anti-viral responses, apoptosis, and anti-proliferative responses [119]. In addition, STAT1 can act as a potent tumor promoter of leukemia development [120]. Interestingly, HBZ interacts with STAT1 and enhances its transcriptional activities [121]. Assuming that AS1 has the same function as HBZ, AS1 might activate STAT1 and promote lymphomagenesis. MiR-484 and miR-1307 also function as tumor suppressors or oncogenes in several cancers [96,97,122,123]. In particular, miR-484 is down-regulated in malignant B cell lymphoma [124]. Therefore, two of the three miRNAs that have strong negative correlations with AS1 mRNA expression were associated with lymphomagenesis.
Moreover, 12 of the 15 bta-miRNAs with relatively weaker negative correlations to AS1 mRNA expression (-0.7 < r < -0.6) also function as tumor suppressors, oncogenes, or both. Our study shows that both miR-532-5p and miR-106b-5p are down-regulated and have a negative correlation with AS1 mRNA expression, which is consistent with the results of another study that showed that miR-532-5p and miR-106a-5p are significantly down-regulated in HTLV-1 asymptomatic carriers [125]. MiR-106b targets the cell cycle regulatory gene p21 (CDKN1A) and is also specifically downregulated in HIV-1 infected CD4 + T cells [126]. MiR-197 induces apoptosis and suppresses tumorigenicity in multiple myeloma [67]. MiR-221 inhibits erythroleukemic cell growth [127]. MiR-425 inhibits proliferation of CLL cells [128]. MiR-342 suppresses the proliferation and invasion of acute myeloid leukemia [51]. MiR-27a functions as a tumor suppressor gene in acute leukemia [66]. MiR-191 displays tumor-type specific roles in tumorigenesis, as miR-191 inhibits cyclin-dependent kinase 6 (CDK6) expression in thyroid follicular cancer [129]. Our results suggest that AS1 may function to down-regulate the expression of bta-miRNAs that suppress cell proliferation, BLV replication, or both. In contrast, miR-146a is an NF-κB-dependent gene and is important in the control of Toll-like receptor and cytokine signaling [130]. In addition, miR-146a is highly expressed in HTLV-1-infected T-cell lines and is directly induced by Tax protein through the activation of NF-κB signaling [108,131]. However, our results were inconsistent with those found in HTLV 1-infected cells. BLV tax and HTLV tax may have different functions for miR-146a. Therefore, these cellular miRNAs may also be pivotal in BLV latency and tumorigenesis.
AS1 mRNA copy number was negatively correlated with six bta-miRNAs (bta-miR-92b-3p, bta-miR-423-3p, bta-miR-155-5p, bta-miR-26a-5p, bta-miR-532-5p, and bta-miR-107-3p), which were not associated with PVL. Five miRNAs (miR-92b, miR-423-3p, miR-155-5p, miR-26a-5p and miR-523-5p) function as both oncogene and tumor suppressor genes. In particular, miR-155-upregulation has been reported in HTLV-1 cell lines and adult T-cell leukemia (ATL) patients [114,132]; however, our results are inconsistent with this and showed that expression levels of miR-155-5p were significantly decreased. MiR-26a-5p is frequently downregulated in various types of cancer, suggesting that these miRNAs function as tumor suppressors by targeting multiple oncogenes, whereas there are some reports that miR-26a promotes tumorigenesis [34, 133,134]. Since sequencing reads of bta-miR-26a, as well as that of bta-miR-191-5p, were very high in B cells from both BLV-uninfected and infected cattle compared to that of other miRNAs, miR-26a seems to be necessary for B cell proliferation, survival, or both. Interestingly, the target sequence of miR26a/b exists in the 3'-UTR of Cell-Division Cycle (CDC)6, and CDC6 gene expression is suppressed by miR-26a/b [135]. CDC6 protein is a key factor for loading the helicase mini-chromosome maintenance (MCM) proteins onto replication origins for the assembly of the pre-replicative complex (pre-RC) at the M-to-G1 phase transition, which is required to establish replication licensing [136,137]. Overexpression of CDC6 gene has been shown to contribute to oncogenesis [138]. Therefore, it is possible that these five bta-miRNAs are affected by AS1 specifically, as there are no associations between these bta-miRNAs and PVL. MiR-26a may be an important miRNA for BLV induced lymphomagenesis.
Three bta-miRNAs (bta-miR-363-3p, bta-miR-874-3p, and bta-miR-339b) were negatively correlated with PVL; however, these were not associated with either tax/rex or AS1 mRNA copy number. Although the reason is unknown, these bta-miRNAs might be affected by PVL via other accessory genes, such as G4 or R3, rather than tax/rex and AS1.
In addition to bta-miR-375 expression significantly correlating with PVL in healthy BLVinfected cattle, at the onset of EBL, bta-miR-375 expression increased to significantly higher levels than those in healthy BLV-infected and uninfected cattle. Several organs express miR-375, which is significantly down-regulated in multiple types of cancer, although it has been found to be up-regulated in prostate and breast cancers [41]. This particular miRNA is a crucial regulator of phagocyte infiltration and the subsequent development of a tumor-promoting microenvironment [139]. In EBL, miR-375 up-regulation may be important for tumor development. Furthermore, our result has confirmed that bta-miR-375 expression levels can be used to distinguish between healthy BLV-infected and EBL cattle. This indicates that bta-miR-375 may be used as a diagnostic biomarker of EBL onset.
The deletion of the miR-15/16 cluster accelerates the proliferation of both human and mouse B cells by modulating the expression of genes so as to control cell cycle progression. In addition, the miR-15/16 cluster has been shown to be deleted or its expression down-regulated in two-thirds of B cell chronic lymphocytic leukemia (B-CLL) cases, which is characterized by the clonal expansion of CD5 + B cells and is similar to that seen in EBL [2,140,141]. In this study, however, the expression of bta-miR-16a was the most stable in B cells, among all 22 cattle, both BLV-infected and BLV-uninfected cattle, and was used as the internal control to normalize the read counts of other bta-miRNAs. Although bta-miR-16b was significantly downregulated in BLV-infected cattle, its down-regulation was not affected by PVL or AS1 expression. Furthermore, the expression levels of bta-miR-15a and -15b in the B cells of BLV-infected cattle did not different significantly from those in BLV-uninfected cattle. This suggests that the bta-miR-15/16 cluster may not be involved in B cell lymphoma caused by BLV.
In conclusion, our deep sequencing analysis demonstrated that BLV provirus-derived blv-miRNAs are preferentially expressed in B cells and correlate with PVL in healthy BLV-infected cattle. In contrast, the expression of some bovine-derived bta-miRNAs, which are believed to be involved in tumor and/or tumor suppression, were significantly down-regulated. These results suggest that BLV promotes lymphomagenesis by down-regulating the expression of bta-miRNAs that have tumor-suppressing functions. However, this lymphomagenesis promotion involves AS1 and blv-miRNAs rather than the tax/rex genes and is associated with increased PVL. Further studies are needed to investigate the molecular function of blv-miR-NAs and bta-miRNAs in the pathogenesis of EBL induced by BLV.