During an outbreak of severe diarrhea in 3–50-day-old piglets in Yunnan province, China, between January and May 2020, 202 fecal samples were collected from 53 different pig farms across 15 prefectures/cities of Yunnan province, representing ~5–10% of the total number of farms in Yunnan, in format of rectal swabs using a directed sampling method to minimize cross-contamination, followed by shipment in sterile 15-ml falcon tubes on ice to laboratory for storage at −80°C till isolation of viral nucleic acid and virus detection.

First, the 202 clinical fecal samples were processed for extraction of viral RNA according to the user instructions of TRIpure total RNA extraction reagent (cat. no.: RP001, BioTeke, Beijing, China). Then, cDNA synthesis of the 202 RNA samples was performed by incubating the reaction at 37°C for 15 min and 85°C for 5 s using EasyScript RT/RI Emzyme Mix (cat. no.: AE311-02, Transgen, Beijing, China) and the downstream primer (5′-CGGTACGCGTAACCAGGGAAAGA-3′ for GIII or 5′-AGTTGTTCATTTYTGGCCATCC-3′ for GV, oligonucleotides purchased from Tsingke Biotech, Kunming, China) for detection of PoSaV-GIII and PoSaV-GV. The subsequent PCR was performed using 2× Phanta Max Master Mix (cat. no.: P505-01, Vazyme, Nanjing, China) and the corresponding prime pairs (Forward 5′-CCCTCATTGGACCAAGTGGGA-3′ and Reverse 5′-CGGTACGCGTAACCAGGGAAAGA-3′ for GIII, Forward 5′-ATCCCAGAGAACATGATGGC-3′ and Reverse 5′-AGTTGTTCATTTYTGGCCATCC-3′ for GV, which were designed according to the reference sequences FJ387164 and KX000383, respectively). The amplification conditions were denaturing at 95°C for 3 min, followed by 35 cycles of 95°C for 15 s, 60°C for 30 s, and 72°C for 60 s, with an additional extension step at 72°C for 5 min. PoSaV-positive samples were determined by gel electrophoresis to visualize the specific amplification products (591 bp and 370 bp in length for the GIII amplicons and the GV amplicons, respectively).

In addition, some most common diarrhea-causing bacterium and parasites, such as Escherichia coli, Salmonella, Shigella, and Coccidia, were also screened for the PoSaV-positive samples according to the standard protocol. Furthermore, co-infections with other nine diarrhea-related porcine viruses were determined for the PoSaV-positive samples using the primers validated in our recent report (42), namely, porcine sapovirus (PoSaV), porcine epidemic diarrhea virus (PEDV), transmissible gastroenteritis coronavirus (TGEV), porcine rotavirus (PoRV), porcine bocavirus (PBoV), porcine astrovirus (PAstV), porcine deltacoronavirus (PDCoV), classical swine fever virus (CSFV), pseudorabies (PRV), and porcine circovirus 2 (PCV2).

To determine the sequences of capsid protein genes (VP1 and VP2) or genome of PoSaV from the positive samples, the PCR products from the amplifications using the prime pairs listed in Supplementary Table S1 (for PoSaV-GIII strains) and Supplementary Table S2 (for PoSaV-GV strains) were gel purified, cloned into pMD18-T vector (cat. no.: 6011, Takara, Dalian, China), purified, and transformed into E. coli DH5α for bulk culture. The plasmid was then extracted from the bacterial culture according to the user instructions of EZ-10 Spin Column Plasmid DNA Minipreps Kit (cat. no.: B610413, Sangon Biotech, Shanghai, China) for subsequent Sanger sequencing confirmation at Sangon Biotech (Shanghai, China). A total of three colonies were sequenced to determine the sequence of VP1, VP2 or the sequenced genome of PoSaV to minimize the potential sequencing error.

VP1 is a highly divergent capsid protein and closely related to the viral immunogenicity, and therefore, its encoding gene has been commonly used for genetic characterization of SaV (3, 4). Some selected SaV reference sequences (listed in Supplementary Table S3) were retrieved from the GenBank database for sequence alignments and phylogenetic analyses with the sequences of capsid genes (VP1 and VP2) determined in this study (accession numbers in Supplementary Table S4 and sequences in Supplementary File). DNAStar 6.0 software with the default parameters was used to assemble the sequences and to perform the similarity analysis of the nucleic acid sequences between the selected reference sequences representative of 16 SaV genogroups (Supplementary Table S3) and the 21 strains obtained in this study (Supplementary Table S4 and Supplementary File). Since the alignments contain gaps due to indels, which are phylogenetically informative, the gaps were included into the calculations and further analyses.

The construction of phylogenetic trees from the aligned nucleotide sequences using bootstrap method with 1,000 replicates was achieved using the ClustalW alignment program included in MEGA 7.0 software package (43) with p-distance as the substitution model (22) and using neighbor-joining method. Of note, p distance may not be the best method to estimate genetic distances as the model fails to consider some biologically relevant factors; however, this model is still used in this study, for direct comparison with some previously published similar studies (5, 22). Recombination Detection Program (RDP, version 4) (44) was used to perform recombination analyses using the unweighted pair-group method (UPGMA) (45). All the seven methods (i.e., RDP, GENECONV, BootScan, Maxchi, Chimaera, SiScan, and 3Seq) in the RDP package were included in the recombination analysis, and only when all of them detected a recombination event, can it be accepted as a recombination. The p-values reported in the study correspond to the ones from RDP method. Phylogenetic trees and RDP plots are present for each recombinant strain.

A considerable number of pig farms in Yunnan province reported severe diarrhea in piglets between January and May 2020, which spread in different batches and across adjacent pens, causing huge economic losses to the farms. Over the period, the pig farms collected 202 stool samples in total from diarrheic piglets (3-−6 samples per pig farm) for pathogen diagnoses at College of Animal Veterinary Medicine, Yunnan Agricultural University.

First, we checked the samples for some most common diarrhea-causing bacterium and parasites, such as Escherichia coli, Salmonella, Shigella, and Coccidia. It turned out that Salmonella, Shigella, and Coccidia were not detected from the samples. Only Escherichia coli were isolated from some of the samples; however, the Escherichia coli-killing antibiotics (such as gentamicin and enrofloxacin) from drug sensitivity test did not show noticeable therapeutic effects after being applied to the diarrheic piglets for treatment, indicating that the direct link between Escherichia coli infection and diarrhea in piglets was not supported.

Then, RT-PCR detection method showed that the PoSaV-positive rate in diarrheic piglets varied from 22.2 to 64.7% depending on the prefectures or cities, with a total positive rate of 35.2% (71/202, Table 1). In spite of co-infection of PoSaV with some other diarrhea-related porcine viruses, such as PEDV, PoRV, and PAstV, 15 of the 71 PoSaV-positive samples (21.1%) were not co-infected with any of the other nine diarrhea-related porcine viruses, suggesting that PoSaV is a contributing factor for diarrhea in swine herds, in agreement with a very recent report (46). As there is no report yet about PoSaV infection in Yunnan, in this report, we focus on the PoSaV epidemiology and the sequence characterization.

Not surprisingly, the three prefectures or cities (Chuxiong, Yuxi and Kunming) which are the three biggest swine producing areas in Yunnan province showed highest PoSaV detection rates, when excluding Pu'er city from discussion due to the very limited sample number (only 4 samples in total from Pu'er). Interestingly, besides 70 PoSaV-GIII-positive samples, a stool sample turned to be PoSaV-GV-positive, which is the first case for detection of SaV-GV in pigs from China.

In this study, 2–4 PoSaV-positive samples from each prefecture/city were randomly selected for nucleic acid sequence determination of the structural capsid protein genes VP1 and VP2, or of the sequenced PoSaV genome. Subsequently, nucleic acid sequences of capsid protein genes of 17 GIII PoSaVs, whole genomic sequences of three PoSaV-GIII strains (YNLH, YNAN, and YNJD) and one PoSaV-GV strain (YNQB) were determined using the primer pairs listed in Supplementary Tables S1, S2, respectively. The 21 Yunnan PoSaV strains obtained in this study are documented in GenBank under the accession numbers MW285639–MW285642 and MW296248–MW296266 (Supplementary Table S4 and Supplementary File).

Next, phylogenetic tree analysis using MEGA7.0 software was performed for the 21 Yunnan PoSaV strains and 60 selected PoSaV reference sequences representing 16 PoSaV genogroups (Supplementary Table S3). In both phylogenic trees constructed on either only VP1 gene sequences (Figure 1A) or longer sequences (VP1 and VP2) (Figure 1B), 20 identified Yunnan PoSaV strains (indicated by solid blue squares) were clustered together with the selected PoSaV-GIII reference strains (Figure 1, block highlighted in light green). On the other hand, one Yunnan PoSaV strain (YNQB, indicated by a red triangle) was clustered together with the 9 selected PoSaV-GV strains (Figure 1, highlighted in dark blue).

Finally, in agreement with the data from specific PCR amplification and sequence determination, phylogenetic tree analyses revealed that PoSaV infection is prevalent in Yunnan pig farms, with PoSaV-GIII to be the predominant genogroup. Furthermore, a PoSaV-GV strain was identified from Yunnan, which is the first report in China.

At present, the 19 genogroups of SaVs can be classified into at least 52 genotypes based on complete sequences of VP1 genes using a pairwise distance cutoff value of ≤0.169 to distinguish different genotypes or clusters (1, 37). To further scrutinize the genetic diversity of the 20 Yunnan PoSaV-GIII strains, the complete sequences of VP1 genes of these 20 strains were subjected to phylogenetic analysis together with 38 selected PoSaV-GIII reference sequences retrieved from NCBI (Supplementary Table S3). Phylogenetic analysis revealed that the 58 PoSaV-GIII VP1 gene sequences could be divided into 11 genotypes (Figure 2), with the 20 Yunnan PoSaV-GIII strains dispersed into 5 genotypes (genotypes 3, 5, 9, 10, and 11) according to the previous classification (29). Sequence dissimilarity comparison (Supplementary Figure S1) confirmed that the sequence discrepancy within a cluster was ≤0.169 and that the sequence disparity between different clusters was mostly >0.169, in agreement with the widely accepted genotype definition (1, 37). Furthermore, when correlating the identification years and places of the PoSaV-GIII strains (Supplementary Table S3) with the genotype assignment (Figure 2), for most of the clusters, we did not find clear clustering profiles based on the years or places of strain identification. Besides, none of the 20 Yunnan PoSaV-GIII strains grouped together with the Cowden strain (KT922087 in genotype 1), and the 4 Yunnan PoSaV-GIII strains in genotype 3 grouped together with 12 PoSaV-GIII strains identified from USA, Spain, Korea, Japan, and China between 2012 and 2019. Taken together, the in-depth phylogenetic analysis revealed high genetic diversity of Yunnan PoSaV-GIII strains.

For the 20 Yunnan PoSaV-GIII strains identified in this study, the sequence assembly revealed that the length of VP1 gene is 1635 nt, whereas the length of VP2 gene is 516 (10 strains with 9 nt deletions, highlighted in the red square in Figure 3A) or 525 nt (another 10 strains without deletions, highlighted in the black square in Figure 3A). In addition, the genome sequences were determined for 3 PoSaV-GIII strains with different lengths (7341 nt for YNJD and YNAN, whereas 7350 nt for YNLH). Analysis using MegAlign function incorporated in software DNASTAR 6.0 revealed that the nucleic acid sequence similarity between the 20 Yunnan PoSaV-GIII capsid protein genes (VP1 and VP2) ranged from 75.3% (YNSJ1 vs. YNYM2) to 99.9% (YNDY1 vs. YNDY2) as highlighted in red squares (Supplementary Figure S1), on the other hand, when cross-comparing the 20 Yunnan PoSaV-GIII sequences with the 38 selected PoSaV-GIII reference sequences, as highlighted in black squares (Supplementary Figure S1), the sequence disparity of capsid protein genes (VP1 and VP2) is between 3.5 (YNLH-MK965898) and 25.1% (YNSJ1 vs. MK965899, YNLH vs. LC215876, and YNLH vs. LC965875). Overall, high sequence discrepancy with deletions and insertions again implied high genetic diversity of Yunnan PoSaV-GIII strains.

Further sequence evaluation using MEGA 7.0 software revealed that 10 Yunnan PoSaV-GIII strains had a 9 nt deletion (7218–7226 nt) in the VP2 region (Figure 3A, highlighted in red square), whereas the other 10 did not possess this deletion at the same location (Figure 3A, highlighted in black square), resulting in deletion of three amino acids (2406–2408 aa, Figure 3B), in agreement with a previous report (5) and further complicating the PoSaV genetic diversity. Additionally, when looking into the origin of the 20 Yunnan PoSaV-GIII strains (Table 1), we noticed that the PoSaV-GIII strains determined from the three biggest swine producing prefectures or cities (Chuxiong, Yuxi, and Kunming) in Yunnan province contained both types, with or without the 9 nt deletion in gene VP2: two strains from Chuxiong (YNDY1 and YNDY2), two strains from Kunming (YNAN and YNAN1) and three strains from Yuxi (YNTH1, YNYM3, and YNYM4) showed the deletion of 7218–7226 nt in VP2 gene region, whereas one strain from Chuxiong (YNLF1), one strain from Kunming (YNDC1), and two strains from Yuxi (YNYM1 and YNYM2) showed no deletion of 7218–7226 nt. The data suggested that different sub-genotypes of PoSaV-GIII strains may be prevalent in some certain areas.

One PoSaV-GV strain (YNQB) was identified from Yunnan in this study, which is the first report for PoSaV-GV strain in China. Sanger sequencing and sequence assembly revealed that the length of the sequenced genome, VP1 gene, and VP2 gene of YNQB is 7496, 1736, and 495 nt, respectively. Alignment of the sequenced genome of YNQB strain with all the 9 SaV-GV strains we extracted from NCBI revealed that the sequenced genome of the Yunnan PoSaV-GV strain (YNQB) from this study shared nucleic acid sequence similarity ranging from 56.9% (YNQB vs. JN420370 of sea lion origin) to 85.5% (YNQB vs. KX000383 of swine origin) (Figure 4A). According to the previous classification (47), SaV-GV strains could be classified into five genotypes (from GV-1 to GV-5) and the YNQB strain fell into cluster 5 together with another two SaV-GV strains of swine origin (Figure 4A).

Meanwhile, in the VP1 gene region of the YNQB strain, we observed 3 nt insertion (position 5827–5829 nt), 9 nt deletion (position 6376–6384 nt), and 3 nt insertion (position 6481–6483 nt) (Figure 4B), resulting in a corresponding insertion of one amino acid (position 1943 aa), deletion of three amino acids (position 2126–2128 aa) and insertion of one amino acid (position 2161 aa) (Figure 4C), which are also presented in two out of the four SaV-GV strains of swine origin (LC483440 and KX000383). Of note, as shown in Figure 4A, two PoSaV-GV strains (AB521771 and AB521772, both from Japan in 2008) exhibit similar nucleic acid sequence profile (68.5–69.0% sequence similarity) with SaV-GV strains of human origin (AB775659, AY645856, MN161594, and NC027026), but shares lower sequence similarity (56.7–67.2%) with the other PoSaV-GV strains (LC4822440, KX000383, and YNQB).

Recombination, both inter- and intra-genogroup, is common for SaVs (1, 5, 22, 48). To the best of our knowledge, there was only one report from China about intra-genogroup (PoSaV-GIII) recombination, where KT922089 and KF204570 were determined as the two parental strains (29) for KX688107 identified from Shanghai in 2015. To perform the recombination analyses for the four Yunnan PoSaV strains (YNLH, YNAN, YNJD, and YNQB), whose genome sequences were determined, using software RDP v.4 (44), analyses with stringent threshold (high confidence for all the computation programs) revealed that two PoSaV-GIII strains out of the four genomes (three PoSaV-GIII strains and one PoSaV-GV strain) may be resulted from multiple recombination events. As shown in Figure 5, YNJD strain may come from the recombination between strains with similar sequences to major parent MK965898 (regions of 1–3803 and 4339–7582 nt) and minor parent LC215880 (region of 3804–4338 nt) with high confidence 2.344 × 10⁻⁴⁵. Meanwhile, as shown in Supplementary Figure S2, YNAN strain may originate from the recombination between strains with similar sequences to major parent MK965898 (regions of 1–5195 and 6352–7583 nt) and minor parent FJ387164 (region of 5196–6351 nt) with high confidence 1.303 × 10⁻³⁵. In agreement with most of other SaV recombination events, the breakpoints for Yunnan PoSaV-GIII strain YNJD and previously reported Shanghai PoSaV-GIII strain KX688107 (29) located in the RdRp-capsid junction region, which is critical for virus replication. However, we noticed that the breakpoint for Yunnan PoSaV-GIII strain YNAN located in the VP1 region, which is important for immunogenicity. In summary, the distinct recombination background with different parental strains and different breakpoint regions for the three recombinants from China further illustrated the high genetic diversity of Yunnan PoSaV-GIII strains.