Changes in the spike and nucleocapsid protein of porcine epidemic diarrhea virus strain in Vietnam—a molecular potential for the vaccine development?

Introduction

Porcine epidemic diarrhea (PED) characterized by watery diarrhea, vomiting, and severe dehydration in suckling piglets and led to 50%–90% mortality among susceptible piglets (Madson et al., 2014). The disease is caused by the porcine epidemic diarrhea virus (PEDV), an enveloped, positive-sense single-stranded RNA virus that belongs to the Coronaviridae family in Nidovirales order (genus Alphacoronavirus) (Lin et al., 2016). The genome of PEDV is approximately 28 kb in length. It comprises a 5′ untranslated region (5′ UTR), a 3′ UTR with a polyadenylated tail, and seven open reading frames (ORF1a, ORF1b, and ORF2-6) encoding four structural proteins (spike, S; envelope, E; membrane, M; nucleocapsid, N), two nonstructural proteins, and one accessory protein ORF3 (Song & Park, 2012).

The S protein (a glycoprotein) contains a specific receptor binding site that is an antigenic target for neutralizing antibodies and relates to its pathogenicity and immunogenicity (Pospischil, Stuedli & Kiupel, 2002). This protein plays a critical role in viral entry through the viral-cellular fusion activity and induces an immune response in the natural host during replication (Lee, 2015; Teenavechyan et al., 2016). Thus, it is used to study the genetic match between vaccine strains and circulating PEDV strains (Lee et al., 2010; Lee & Lee, 2014; Oh et al., 2014). The M protein is a surface protein and plays an important role in the process of virus-assembly and the induction of protective antibodies with neutralizing activity (Lee, 2015; Teenavechyan et al., 2016). The N protein is highly conserved and binds to virion RNA to provide a structural basis for viral transcription, replication, and assembly (Curtis, Yount & Baric, 2002; Chen et al., 2013a; Sun et al., 2015). The N protein was commonly used to diagnose infection with PEDV (Song, Moon & Kang, 2015a) and was known to protect the viral genome during the coronavirus assembly. It also affects other anti-virus responses through host immune evasion strategies (Lee, 2015; Teenavechyan et al., 2016). The epitopes on the N protein are considered to possibly cause for induction of cell-mediated immunity (CMI) (Saif, 1993). However, protein S is more antigenic than other PEDV proteins and anti-S antibodies detected in PEDV-infected pigs last longer than anti-N antibodies (Knuchel et al., 1992). For non-structural proteins, ORF1a and ORF1b are multifunctional related to viral genome replication (Brian & Baric, 2005), and the accessory ORF3 protein is known to be associated with the virus virulence (Park et al., 2007; Park et al., 2008; Wang et al., 2012; Chen et al., 2013b; Song et al., 2015b). Therefore, the ORF3 protein has been studied extensively in the molecular epidemiology of PEDV (Shirato et al., 2011; Chen et al., 2013a; Li et al., 2013; Temeeyasen et al., 2014; Song et al., 2015b). These proteins have been the targets of many studies to understand the causes of the outbreaks and develop more effective vaccines (Lee & Lee, 2014; Temeeyasen et al., 2014; Lee, 2015; Song et al., 2015b; Lin et al., 2016; Su et al., 2016).

Won et al. (2020) reviewed 299 published articles and showed that current vaccines are produced mainly based on PEDV strains: CV777 (belonged to G1a group) and SM98, DR13 (belonged to G2a group). The G2b whole-virus killed vaccines have been developed and used in the USA and the G2b live oral vaccine (based on the KNU-141113 strain) has been developed and used in Korea from 2020. However, the virus evoluted and accumulated mutations as the time passed that may lead to sub-optimally match to the actual pandemic virus of the vaccines. Therefore, complete genome sequencing of PEDV strains circulating in each country is of critical in order to develop effective vaccines (Vlasova et al., 2014; Vui et al., 2015; Jarvis et al., 2016; Fan et al., 2017; Lee & Lee, 2017; Rasmussen et al., 2018; Liang et al., 2020; Garcia-Hernandez et al., 2021; Wen et al., 2021).

In Vietnam, PEDV caused an outbreak for the first time in 2009 and then occurred again in 2013 (Diep et al., 2018). Several PEDV vaccines have been developed and being used in Vietnam based on vaccine strains CV777/CN/KT323979, SM98/Korea/GU937797, and DR13/Korea/JQ023162 (belonged to G1a and G2a groups). Although a vaccine campaign for piglets has been conducted since 2011 for PEDV control, diarrhea disease still exists in many provinces and causes serious damage to the livestock industry in recent years in Vietnam. The low effectiveness of the vaccines could be due to the genetic differences between vaccines and field epidemic strains, emphasizing the necessity of novel vaccines against new viral strains. Vui et al. (2015) has sequenced the first complete genome of PEDV from outbreaks in Vietnam in 2013 and suggested that the Vietnamese PEDV isolates were new variants. Other studies mainly focused on the spike gene sequences indicated that there have been remarkable changes in the Vietnamese PEDV strains collected from 2012 to 2016 (Kim et al., 2015; Diep et al., 2018; Than et al., 2020) and these may decrease the vaccine efficiencies. Although PEDV strains remain the main cause of piglet losses in large litters in Vietnam; however, very limited molecular data of their genotypic and genealogy of PEDV are available from Vietnam. Since the first complete genome sequence of Vietnamese PEDV strains announced in 2015, no more complete sequence of Vietnamese PEDV strain was reported. In this study, we report the molecular characteristics and the changes in the spike and nucleocapsid protein of a PEDV strain isolated from pigs causing severe diarrhea in 2018 in Hung Yen province, Northern Vietnam. Our study emphasizes the importance and urgency of studying the genetic diversity of current PEDV strains and their evolution in order to develop a suitable vaccine strategy for each geographic region.

Materials and Methods

Results

Phylogenetic analysis

The phylogenetic tree constructed based on the full-length genome sequence of 37 PEDV references indicated that the PEDV strains could be divided into two major groups G1 and G2 (Fig. 1). Group 1 contained five strains including CV777 (Belgium/2001/AF353511), FR001 (France/2014/KR011756), Minnesota52 (USA/2013/KJ645704), and Ohio126 (USA/2014/KJ645702). Group 2 was divided into two subgroups G2a and G2b. The G2a group comprised of six strains from countries including Japan (JPN/2014/LC063813, JPN/2013/LC063814), Korea (SM98/2010/GU937797), China (SXYl/2016/MF462814), Thailand (CBR1/2014/KR610993), and USA (PC22A/2013/KU893861). The G2b group, which is considered highly virulent, comprised of strains from Vietnam (IBT/VN/2018, KCHY, VAP, and JFP), China, Korea, USA, Germany, and Belgium. Remarkably, the PEDV strain isolated from North of Vietnam and CH Hubei/CN/2016/KY928065 strain—a highly pathogenic strain from China in 2016, showed the highest genetic identity with each other compared to other strains.

The phylogenetic tree based on the S protein sequences also showed that all the strains could be divided into two groups (Fig. 2). Group 1 (G1) consisted of two subgroups: G1a and G1b. G1a subgroup included the classical non S-INDEL strains (Korea/2010/GU937797, JPN/2011/AB548622, and JPN/2011/AB548623). G1b subgroup included the classical S-INDEL strains (Belgium/2015/KR003452, CV777/Belgium/AF353511, CV777/CN/KT323979, SM98/Korea/KJ857455, DR13-V/Korea/DQ462404, DR13/Korea/DQ862099, HUA PED67/VN/KP455319, and some strains from France, Germany, USA, China, Japan, and Korea). Group 2 (G2) contained two subgroups: G2a and G2b. Among them, the Vietnamese strains (IBT/VN/2018, HUA PED45, HUA PED47, KCHY, VAP, and JFP) belonged to the G2b subgroup including Asian S-INDEL strains such as the vaccine strain AJ1102/CN/JX188454 and China/2012/JX112709. G2a subgroup including 24 strains from USA, Mexico, Korea, China, Taiwan, and Thailand belonged to Northern American and Asian S-INDEL subgroup.

Multiple sequence alignments of S, N and ORF3 genes

The S gene sequence of the IBT/VN/2018 strain was 4,158 nt in length and encoded a protein consisting of 1,385 aa. The genetic identity of the nucleotide sequence and amino acid sequence of the S gene of the Vietnamese strain was 92.2%–98.1% and 91.1%–98.4% compared to the strains from other countries, respectively (Table S2). The highest identity was found between IBT/VN/2018 strain and CN/2012/JX112709 strain, while the lowest identity was found with Korea/2010/GU937797 strain. The S protein of the IBT/VN/2018 strain contained a signal peptide (aa 1–19), four neutralizing epitopes (aa 501–640, 751–758, 766–774, and 1,370–1,376), a transmembrane domain (aa 1,328–1,350), and a short cytoplasmic domain. Four epitopes were determined at aa positions 501–640 (COE), ⁷⁵¹YSNIGVCK⁷⁵⁸ (SS2), ⁷⁶⁶LSQSGQVKI⁷⁷⁴(SS6), and ¹³⁷⁰GPRLQPY¹³⁷⁶ (2C10). In the epitope SS6, substitution ⁷⁶⁶P > L⁷⁶⁶ was identified in IBT/VN/2018, HUA PED45 and HUA PED47 strains isolated in Vietnam only while it was not detected in other strains investigated.

The S protein can also be divided into two domains: S1 (aa 1–789), S2 (aa 790–1,385), or subdomains as described by Li et al. (2016) (Fig. 3). In the S⁰ subdomain of the IBT/VN/2018 strain, two substitutions (¹³⁵DN¹³⁶ > ¹³⁵SI¹³⁶ and ¹⁴⁴N > D¹⁴⁴) were found when compared to the vaccine strains (CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162) (Fig. S1). In the S^A subdomain, there were eight substitutions between the IBT/VN/2018 strain and the vaccine strains (CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162) at aa positions: ²⁹⁴I > M²⁹⁴, ³¹⁸A > S³¹⁸, ³³⁵V > I³³⁵, ³⁶¹A > T³⁶¹, ⁴⁹⁷R > T⁴⁹⁷, ⁵⁰¹SH⁵⁰² > ⁵⁰¹IY⁵⁰², ⁵⁰⁶I > T⁵⁰⁶, ⁶⁸²V > I⁶⁸², and ⁷⁷⁷P > L⁷⁷⁷. The IBT/VN/2018 strain also had two insertions (at aa position ⁵⁹NQGV⁶² and ¹⁴⁵N) and one deletion (¹⁶⁸DI¹⁶⁹) in S protein when compared with the CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162 strains (Fig. 4).

The S2 domain presents the typical structural features of class I fusion proteins, including a hydrophobic fusion peptide (FP, residues 891–908), two heptad repeat regions (HR1, residues 978–1,117 and HR2, residues 1,274–1,313), and a C-terminal transmembrane domain (TM, residues 1,328–1,350). In the S2 domain of IBT/VN/2018 strain, eight substitutions at aa positions ⁸⁵⁷V > A⁸⁵⁷, ¹⁰⁰⁹L > M¹⁰⁰⁹, ¹⁰⁸⁹S > L¹⁰⁸⁹, ¹²⁰⁷T > D¹²⁰⁷, ¹²²¹F > Y¹²²¹, ¹²²⁹S > G¹²²⁹, ¹²⁵¹D > E¹²⁵¹, and ¹²⁷⁹P > S¹²⁷⁹ were detected in comparison to the CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162 strains. Among them, two substitutions were found in subdomain S^HR1 at aa position ¹⁰⁰⁹L > M¹⁰⁰⁹ and ¹⁰⁸⁹S > L¹⁰⁸⁹, and one substitution at aa ¹²⁷⁹P > S¹²⁷⁹ in subdomain S^HR2.

The N protein of the IBT/VN/2018 strain had a genetic identity at the nt and the aa sequence with other strains ranging from 92.5%–99.0% and 93.4%–98.8%, respectively (Table S3). The IBT/VN/2018 strain had the highest genetic identity with CH Hubei/CN/2016/KY928065 strain and the lowest identity with the vaccine strain CV777/CN/KT323979. Notably, two substitutions at aa ³⁶⁴V > I³⁶⁴ and ³⁷⁸N > S³⁷⁸ were only detected in the IBT/VN/2018 strain (Fig. S2).

The identity of nucleotide sequence and the amino acid sequence of the ORF3 region of the IBT/VN/2018 strain was 89.1%–99.4% and 90.1%–100% compared to other strains, respectively (Table S4). Based on the nucleotide sequence of the ORF3 region, the IBT/VN/2018 strain was most closely related to the CH Hubei/CN/2016/KY928065 strain, sharing a sequence identity of 99.4%. Remarkably, it shared the lowest nucleotide identity (89.1%) with the vaccine strain CV777/CN/KT323979 and JS2008/CN/KC109141 strain. In the ORF3 region, four substitutions were detected at aa positions ²⁵L > S²⁵, ⁷⁰I > V⁷⁰, ¹⁰⁷C > F¹⁰⁷, and ¹⁶⁸D > N¹⁶⁸ in IBT/VN/2018, KCHY, VAP, and JFP strains and some other strains such as GD1/CN/2011/JX647847, GDA/CN/2012/JX112709, CH Hubei/CN/2016/KY928065, and SBR1/Thailand/2014/KR610993 (Fig. S3).

Discussion

In this study, the complete genome of the PEDV strain (IBT/VN/2018) isolated in a severe outbreak in North Vietnam 2018 has been sequenced and analyzed. The genome sequence analysis indicated that the IBT/VN/2018 strain carried specific structural characteristics of members of the Coronaviridae family, with genetic identities ranging from 95.7% to 99.3% compared to other PEDV strains (Table S1). The PEDV, based on the S gene sequence, was classified into two genogroups (G1 and G2) and each can sub-divided into groups (G1a, G1b, G2a, and G2b). Classical strains are designated G1a, and the new variant strains (with insertion and deletion in the S gene) belong to G1b. The G2a and G2b include the highly virulent strains in Asia and North America. The phylogenetic tree constructed from the whole genome sequences (Fig. 1) and the S gene sequences (Fig. 2) showed that the IBT/VN/2018 strain belonged to the G2b subgroup including Northern American and Asian S-INDEL strains. According to previous studies, strains that belonged to the G2b and Northern American and Asian S-INDEL subgroup were highly virulent (Lee, 2015; Lin et al., 2016). This result is consistent with previous studies on PEDV strains isolated in Northern provinces of Vietnam which also showed that Vietnamese PEDV strains belong to the G2b subgroup (Kim et al., 2015; Diep et al., 2018).

The S protein composes of 1,383 to 1,386 amino acids consisting of the S1 subunit (aa 1–789) and the S2 subunit (aa 790–1,383) (Fig. 3), depending on the strain (Sun et al., 2007). S protein is known to mediate viral entry and inducing neutralizing antibodies in the natural host. The S1 subunit is the extracellular domain and can bind to target cell receptors (Deng et al., 2016). It is important for cell membrane fusion and virus entry and it is the antigenic target of neutralizing antibodies (Sun et al., 2008; Lee et al., 2011; Deng et al., 2016). It contains a putative signal peptide (aa 1–24), a large extracellular region contains two subdomains: NTD (aa 21–324) and CTD (aa 253–638) (Li, 2015), a single transmembrane domain (aa 1,334–1,356), and a short cytoplasmic tail (Oh et al., 2014). Four neutralizing epitopes (aa 499–638, 748–755, 764–771 designated as COE, SS2, SS6 in domain S1, and 2C10 1,368–1,374 in domain S2) have been determined on the surface of S protein (Sun et al., 2008; Li et al., 2016; Okda et al., 2017). The two regions including SS2 and 2C10 are highly conserved. In contrast, there are two positions at aa 499 and 520 in the COE subdomain and one position at aa 766 in the SS6 subdomain are highly variable (Chang et al., 2002; Cruz, Kim & Shin, 2008; Sun et al., 2008; Lara-Romero et al., 2018). These findings are the basis for the development of new effective vaccines in the future (Yu et al., 2020). Further analysis of sequences in functional protein regions of the IBT/VN/2018 strain showed that at the N-terminal of S protein there were two aa insertions and one aa deletion compared to the CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162 strains (Fig. 4). The substitution at aa ⁷⁶⁶P > L⁷⁶⁶ was also found in the IBT/VN/2018 strain in the epitope SS6. These changes may lead to the ineffectiveness of the drugs and vaccines.

In this study, many substitutions were found on the S protein of the IBT/VN/2018 strain in comparison to the CV777/CN/KT323979, AJ1102/CN/JX188454, SM98/Korea/KJ857455, and DR13/Korea/JQ023162 strains. Among 19 aa substitutions in S protein, the IBT/VN/2018 strain had eleven changes that were similarly in other Vietnamese PEDV strains (KCHY, VAP, JFP and/or HUA PED45, HUA PED47, HUA PED67, and HUA PED96) and some strains from China (CN/2011/JQ638920, CN/2012/JX112709, and CN/2013/KF761675), Korea (Korea/2010/GU937797 and Korea/2011/KC879280), and Thailand (Thailand/2011/KC764959). Remarkably, we found seven substitutions that were only found in the IBT/VN/2018 strain in comparison to other PEDV strains including ¹⁴⁴N > D¹⁴⁴, ³¹⁸A > S³¹⁸, ³³⁵V > I³³⁵, ⁵⁰¹SH⁵⁰² > ⁵⁰¹IY⁵⁰², ⁶⁸²V > I⁶⁸², ¹⁰⁰⁹L > M¹⁰⁰⁹, and ¹⁰⁸⁹S > L¹⁰⁸⁹ (Fig. S1). The results of phylogenetic analysis showed that the IBT/VN/2018 strain was closely related to PEDV strains in Asia but differ from US and European strains. This can be speculated that the Vietnamese strains have been genetically changes to adapt to environmental conditions and it leads to the reduction of the effectiveness of the vaccine. Our hypothesis was further supported by a recent study which revealed that the mutations in the neutralizing epitope regions in the S gene cause inefficiencies in vaccination (Li et al., 2013). These regions were detected at aa 7–146 and 271–278 in the neutralizing epitopes (Li et al., 2013). However, more studies need to be carried out to confirm this assumption.

S⁰ also known as the N-terminal region is a functional receptor for the porcine epidemic diarrhea virus. PEDV uses the N-terminal region as the major receptor for cell entry (Li et al., 2017). The N-terminal region binds to sugar which acts as its co-receptor. This process is the first important step to help the virus penetrate cells (Deng et al., 2016). Recent reports suggest that any amino acid mutation can change the virulence of Coronavirus (Peng et al., 2011; Wu et al., 2012; Teenavechyan et al., 2016; Qin et al., 2019; Wabalo et al., 2021). Consequently, the N-terminal domain can be used as a vaccination strategy to prevent PEDV infections. The S⁰ subdomain of the IBT/VN/2018 strain had two aa substitutions at positions ¹³⁵DN¹³⁶ > ¹³⁵SI¹³⁶ and N¹⁴⁴ > D¹⁴⁴ which may make virus entry into cells easier. In other words, these mutations can enhance the pathogenicity of the Vietnamese strain. However, this needs to be experimentally confirmed.

In the IBT/VN/2018 strain, we also identified two aa substitutions in the S^HR1 subdomain at aa positions ¹⁰⁰⁹L > M¹⁰⁰⁹ and ¹⁰⁸⁹S > L¹⁰⁸⁹, and one substitution at aa ¹²⁷⁹P > S¹²⁷⁹ in subdomain S^HR2. In Coronavirus, membrane fusion is initiated by the invention and intervention of the fusion protein (FP) into the target cell membrane, then the fusion protein combined with HR1 and HR2 regions form a stable structure. This process took the virus’s transmembrane domain accessed and fused it with the membrane of the host cell (Eckert & Kim, 2001; Li et al., 2016). Therefore, the substitutions in the functional areas may lead to change this process and interferes with the viral entry.

Recent reports showed that the S protein, the N protein and the accessory ORF3 protein also play an important role in regulating the virulence of PEDV strains and can cause severe damage by evading host immune mechanisms (Zhang & Yoo, 2016; Kim et al., 2020; Li et al., 2020). The N protein is identified as a multifunctional region and participating in many stages in viral replication and regulating functions of PEDV (Zuniga et al., 2010; Liwnaree et al., 2019). Two novel epitopes at aa 18–133 and 252–262 were identified on the N protein (Wang et al., 2016a). The accessory protein ORF3 was identified as an ion channel and has many regulatory functions (Wang et al., 2012; Wang et al., 2016b; Kaewborisuth, He & Jongkaewwattana, 2018). Kaewborisuth, He & Jongkaewwattana (2018) indicated that the ORF3 protein can work together with the S protein for PEDV assembly at the viral replication step. The ORF3 protein is also related to the virulence of PEDV (Song et al., 2003; Park et al., 2007; Chen et al., 2013b) through its regulatory process in viral production (Wang et al., 2012). In general, the ORF3 protein is an important virulence gene for PEDV pathogenicity, and molecular epidemiology studies of PEDV (Pospischil, Stuedli & Kiupel, 2002; Song et al., 2003; Shirato et al., 2011). Two substitutions at aa positions ³⁶⁴N > I³⁶⁴ and ³⁷⁸N > S³⁷⁸ and four substitutions at aa positions ²⁵L > S²⁵, ⁷⁰I > V⁷⁰, ¹⁰⁷C > F¹⁰⁷, and ¹⁶⁸D > N¹⁶⁸ were found in the N and ORF3 region of PEDV-VN strains (IBT/VN/2018, KCHY, VAP, and JFP) (Figs. S2, S3) suggesting that these PEDV-VN strains have evoluted and changed their pathogenicity.

The results of homologous modeling showed that the acquired mutations found in S protein of the IBT/VN/2018 strain including mutations in the S⁰ domain does not remarkably affect its overall 3D structure (Fig. S4). This data is well consistent with result of a previous study which indicated that the deletion of S⁰ domain does not impart any macroscopic changes in spike protein conformation (Kirchdoerfer et al., 2021).

It has been known that glycosylation plays an important role receptor binding, virus entry, protein proteolysis, and protein transport, thereby altering the virulence and immune evasion of virus (Vigerust & Shepherd, 2007). The glycosylation sites have been reported in PEDV in recent years (Kim et al., 2015; Wen et al., 2021). Six glycosylation sites (25^G+, 123^N+, 62^{V +}, 144^D+, 1009^M+, and 1279^L+) were detected in the ORF3, N, and S proteins of the IBT/VN/2018 strain when compared with the vaccine strains. Among that three novel glycosylation sites (144^D+, 1009^M+, and 1279^L+) were detected in the IBT/VN/2018 strain. As a consequence, the amino acid substitutions in the S, N and ORF3 proteins may help the IBT/VN/2018 strain evade the host immune system and probably change the virulence of the virus.

Conclusion

Our study showed that the IBT/VN/2018 strain belonged to the G2b subgroup that along with the Northern American and Asian S-INDEL strains and are considered as a highly virulent group. Remarkably, eight amino acid substitutions (²⁹⁴I > M²⁹⁴, ³¹⁸A > S³¹⁸, ³³⁵V > I³³⁵, ³⁶¹A > T³⁶¹, ⁴⁹⁷R > T⁴⁹⁷, ⁵⁰¹SH⁵⁰² > ⁵⁰¹IY⁵⁰², ⁵⁰⁶I > T⁵⁰⁶, ⁶⁸²V > I⁶⁸², and ⁷⁷⁷P > L⁷⁷⁷) were found in S^A subdomain in the Vietnamese strain compared with the vaccine strains. In addition to these substitutions, three novel N- and O-glycosylation sites (144^D+, 1009^M+, and 1279^L+) were detected in the S protein of IBT/VN/2018 strain. The continual mutations in these genes may have generated a novel antigenic strain and help the virus escape from the host immune response induced by the vaccine. Our results highlight the importance of molecular characterization of PEDV strains circulating in Vietnam and provide a molecular potential for the development of an effective vaccine to control PEDV infections of pigs in Vietnam.

Supplemental Information