METHODS
Sample collection
The questing ticks were collected by the flagging method, and blood-sucking ticks were collected from cattle in Shulan. These ticks were identified to species as described elsewhere (Deng, 1978; Shao, Zhang, Li, Huang, & Yan, 2020). Every 10 ticks were pooled into a 1.5 mL Eppendorf tube according to the collection sites and species and stored at -80°C until use.
RNA library construction and sequencing
After washing with 75% ethanol and RNA/DNA-free water, pooled ticks in tubes were added with 800 μL Dulbecco’s modified Eagle’s minimum essential medium (DMEM) and two stainless steel beads (3 mm diameter), and crushed using the Tissuelyser (Jingxin, Shanghai, China) at 70 Hz for 2 min. The lysates were centrifuged at 12000 rpm for 10 min at 4 °C, and the supernatant was further pooled for library construction according to the collection sites and species (Supplementary Table S1). After digested with micrococcal nuclease (NEB, USA) in 37 °C for 2 h, the pooled samples were used for viral RNA extraction with the TIANamp Virus RNA kit (TIANGEN, Beijing, China).
The extracted RNA was subjected to metagenomic sequencing at Tanpu Biological Technology Co., LTD (Shanghai, China). Briefly, the RNA from each pool was used for library preparation with the NEBNext® UltraTM RNA Library Prep Kit for Illumina® (NEB, USA) according to the manufacturer’s instructions. After adapter ligation, 10 cycles of PCR amplification were performed for target enrichment. The libraries were pooled at equal molar ratio, denatured and diluted to optimal concentration, and sequenced with an Illumina NovaSeq 6000 System.
Transcriptome analysis
Transcriptome analysis was conducted as described elsewhere (Gordon et al., 2020). Briefly, after trimming and removing low quality reads, the paired-raw reads were purified by removing ribosomal RNA, host contamination, and bacteria sequences using BBMap program (https://github.com/BioInfoTools/bbmap), and assembled into contigs with SPAdes v3.14.1 (https://github.com/abl ab/spades) and SOAPdenovo v2.04 (https://github.com/aquaskyline/SOAPdenovo-Trans) (Prjibelski, Antipov, Meleshko, Lapidus, & Korobeynikov, 2020; Xie et al., 2014). After being compared with the nonredundant nucleotide (nt) and protein (nr) database downloaded from GenBank using BLAST+ v2.10.0, the assembled contigs were filtered to remove the host and bacterial sequences. The relative abundance of the identified viruses was determined by mapping the reads back to the assembled contigs using Bowtie2 v2.3.3.1.
Viral genome confirmation and annotation
The assembled contigs were compared with NCBI nucleotide and viral refseq database using BLAST (V2.10.0+), and used as a reference for designing specific primers to confirm and analyze the sequences of terminal ends, using the nested reverse transcription-polymerase chain reaction (RT-PCR) and the rapid amplification of cDNA ends (RACE) as described elsewhere (Ma et al., 2021; Y. C. Wang et al., 2021). Detailed information on the primers for the detection or whole genome amplification are shown in Supplementary Table S2. Potential open reading frames (ORFs) in the viral sequences were predicted using ORFfinder (https://www.ncbi.nlm.nih.gov/ orffinder/).
Virus classification
All the viruses identified in this study were classified according to the latest International Committee on Taxonomy of Viruses (ICTV) report of virus taxonomy (https://talk.ictvonline.org/ictvreports/ictv_online_report/). A novel viral species should be satisfied with one of the following conditions as described before (Xu et al., 2021). namely, (i) <80% nucleotide (nt) identity across the complete genome; or (ii) <90% amino acid (aa) identity of the RNA-dependent RNA polymerase (RdRp) domain with the known viruses. All novel viruses were named as the collection sites that the virus was first identified, followed by common viral names according to their taxonomy. All the viral strains would be marked with ‘Northeatern (NE)’ to distinguish them from the virus strains identified in other studies.
Phylogenetic analyses
To confirm the phylogenetic relationships of the viruses discovered in this study, representative reference viral sequences were retrieved from the GenBank database (Supplementary Table S3), and aligned using ClustalW available within MEGA 7.0. Phylogenetic analyses were conducted with the aligned sequences using the maximum-likelihood method in MEGA version 7.0 with the best-fit substitution model for each alignment (Kumar, Stecher, & Tamura, 2016). A bootstrapping analysis of 1000 replicates were conducted in the analysis, and the bootstrap values more than 70 were shown in the trees.