Emergence of Avian coronavirus Escape Mutants Under Suboptimal Antibody Titers.

To perform a quasispecies assessment of the effect of vaccine combinations and antibody titers on the emergence of Avian coronavirus (AvCoV) escape mutants, 5-week-old males from a commercial chicken breeder lineage were vaccinated intramuscularly with one dose of a monovalent (genotype GI-1) or a bivalent (genotypes GI-1 and GI-11 (n = 40 birds/group) AvCoV vaccine. Seven birds were kept as controls. Six weeks later, pools of sera of each group were prepared and incubated at virus neutralization doses of 10 and 10-1 with the Beaudette strain (GI-1) of AvCoV in VERO cells. Rescued viruses were then submitted to genome-wide deep sequencing for subconsensus variant detection. After treatment with serum from birds vaccinated with the bivalent vaccine at a titer of 10-1, an F307I variant was detected in the spike glycoprotein that mapped to an important neutralizing region, which indicated an escape mutant derived from natural selection. Further variants were detected in nonstructural proteins and non-coding regions that are not targets of neutralizing antibodies and might be indicators of genetic drift. These results indicate that the evolution of AvCoV escape mutants after vaccination depends on the type of vaccine strain and the antibody titer and must be assessed based on quasispecies rather than consensus dominant sequences only because quasispecies may be otherwise undetected.

Introduction

Avian coronavirus AvCoV (Nidovirales: Coronaviridae: Coronavirinae: Gammacoronavirus), the causative agent of avian infectious bronchitis (IB), is a pantropic coronavirus that occurs as a host of serotypes classified into 6 genotypes and 32 lineages that use chickens as the main hosts (Jones 2010; Valastro et al. 2016).

The AvCoV genome is a 27.6 kb positive-sense single stranded RNA that codes for nonstructural proteins (NPSs) 2 to 16 of the replicase, with NPS12 representing the RNA-dependent RNA-polymerase and NPS14 representing the 3’–5’ exonuclease, as well as the envelope protein spike glycoprotein S, envelope protein E, membrane protein M and nucleocapsid N protein, which is found in the helical nucleocapsid (Denison et al. 2011; Masters 2006; Minskaia et al. 2006).

Accessory proteins 3a, 3b, 5a and 5b may accelerate growth kinetics and act as type I IFN antagonists (Laconi et al., 2018), while the 5′UTR (untranslated region), which contains the leader and the body TRS (transcription regulatory sequence), and 3'UTR secondary structures regulate the synthesis of sub-genomic mRNAs (Masters 2006). ORFX has been described between genes M and 5a (Gomaa et al. 2008), but its function is still unknown.

For RNA viruses, such as AvCoV, a population of mutants thrives under the consensus sequence level and present single and multiple nucleotide substitutions and indels at different frequencies, thus producing a spectrum of mutants called the quasispecies that evolved by natural selection and genetic drift (Lauring and Andino 2010).

Although AvCoV variants have already been described for vaccine strains and in challenged birds (McKinley et al. 2008; Toro et al. 2012), these reports were based only on partial spike gene sequences by Sanger sequencing and could not provide a more robust mechanistic explanation for AvCoV quasispecies evolution.

AvCoV serotypes show a 20–25% aa difference in the S1 ectodomain of the spike glycoprotein, but as few as 10 to 15 aa changes are sufficient to produce new serotypes; moreover, several antigenic domains are mapped along this subunit (Cavanagh 2007), thus allowing for the emergence of escape mutants due to a low antigenic identity between vaccines and field viruses.

Considering the lack of knowledge on the molecular evolution of AvCoV that may lead to escape mutants, the present study was performed to assess in vitro the mutant spectra evolution of AvCoV under neutralization with different doses of sera of chickens vaccinated with homologous and heterologous vaccine strains.

Materials and Methods

Virus and Vaccines

The Beaudette strain of AvCoV of the GI-1 lineage at the 15th passage in VERO cells with a titer of 102 TCID50/100 µl (reference virus) was used for all in vitro assays. Inactivated vaccines containing the M41 strain of the GI-1 lineage Nobilis® REO + IB + G + ND (MSD Animal Health) (AvCoV monovalent vaccine) and the M41 strain plus the AvCoV BR-I strain of the GI-11 lineage CEVAC® MAXIMUNE PRO (Ceva) (AvCoV bivalent vaccine) were used for the production of anti-AvCoV sera.

Production of Sera Against AvCoV Vaccines in Chickens

This assay was approved by the Ethics Committee on the Use of Animals of the School of Veterinary Medicine, University of São Paulo, under protocol #1306080221. Eighty-seven-day-old male chicks from a commercial breeder lineage were used for this study, housed in an experimental poultry facility with multiple pens and provided water and feed ad libitum. Because live virus was not involved, isolators were not used for this step.

At 5 weeks of age, birds were vaccinated with either the mono (GI-1) or the bivalent (GI-1/GI-11) vaccine (n = 40 birds/group) at a dose of 0.5 ml/bird via chest muscle injection. The remaining seven birds were maintained as controls. Each group was kept in a separate pen, and 6 weeks after vaccination, sera from all groups were collected and tested for AvCoV antibodies with an ELISA Infectious Bronchitis Virus Antibody test kit (Biochek) and then pooled as GI-1 sera, GI-1/GI-11 sera, and control sera, with pooled titers of 9167, 7687 and 0, respectively. All birds were negative for AvCoV antibodies before vaccination.

For the titration of virus neutralizing antibodies, the three pools of sera were inactivated at 56 °C for 30 min, and tenfold dilutions of each pool in FCS-free MEM (100 µl/dilution, eight replicates/dilution) were added to 96-well plates and mixed with 100 µl of the reference virus at 102 TCID50/100 µl.

After incubation at 37 °C/5% CO2 for one hour, the total volume of each well was transferred to a 96-well plate with 48-h-old VERO cells after growth medium removal and rinsing with FCS-free MEM. The plate was incubated at 37 °C/5% CO2 for 48 h, and Reed-Muench titers were calculated based on the absence (neutralization) of syncytial CPE (cytopathic effect).

The GI-1 and GI-1/GI-11 serum pool titers were 102.71 and 102.11 and normalized to the VN titer of 10 by dilution in FCS-free MEM to titers similar to that reported in chickens (Finney et al. 1990). The control serum pool showed no detectable neutralizing antibody (titer = 0).

Escape Mutants Assay

In a 12-well plate, 500 µl of the reference virus at 102 TCID50/100 µl was mixed with 500 µl of GI-1 or GI-1/GI-11 sera at VN titers of 10 and 100.1 in wells A1-B1 and A2-B2, respectively. In addition, 500 µl of the control sera pool diluted to simulate the 10 and 100.1 dilutions of the hyperimmune sera was also mixed with 500 µl of VERO P15 at 102 TCID50/100 µl in wells C1 and C2, respectively. As controls, 500 µl MEM plus 500 µl of the reference virus at 102 TCID50/100 µl and 1 ml of the reference virus at 102 TCID50/100 µl (reference virus only) were added to wells A3 and C3, respectively.

The plate was incubated at 37 °C/5% CO2 for one hour, and the total (1 ml) of each well was transferred to the corresponding well in a 12-well plate with 48-h-old VERO cells after removing the growth medium and rinsing with FCS-free MEM. The plate was incubated at 37 °C/5% CO2 for 48 h and then frozen at -80 °C.

Determination of Virus Loads

Total RNA was extracted from the supernatants of the reference virus and the wells from the escape mutant assay plate after clarification at 1000×g/4 °C/10 min using RNeasy® Mini kit (Qiagen).

RT-qPCR for AvCoV genomic RNA was carried out using a Power SYBR® Green RNA-to-Ct™ 1-Step kit (Applied Biosystems) following the manufacturer’s instructions and primers described by Callison et al. (2006) targeting the 5ʹUTR, with the assays performed in triplicate. RT-qPCR of β-actin as an endogenous control was performed using a Power SYBR® Green RNA-to-Ct™ 1-Step kit (Applied Biosystems) and IDT Integrated DNA Technologies proprietary primers 5'ACAGAGCCTCGCCTTTG3'/5'CCTTGCACATGCCGGAG3', with the reactions performed in duplicate.

RT–qPCR was carried out in a 7500 Real-Time PCR system (Applied Biosystems) (48 °C/30 min; 95 °C/10 min; 45 cycles of 95 °C/15 s and 60 °C/40 s; and melting curve analysis). Absolute quantification (number of copies/µl of sample) of AvCoV was obtained by comparison with a tenfold dilution standard curve with a plasmid containing the corresponding 5′UTR sequence of AvCoV ranging from 1 × 103 to 1 × 109 copies per reaction (slope = − 3.536 and y-intercept = 43.665) and normalized to the β-actin Cq value of each sample.

Deep Sequencing

To assess full genome variants for the reference virus and the wells from the escape mutant assay plate, samples were filtered through 0.45 μm syringe filters and the filtrates were treated with Turbo™ DNase (Ambion) and RNase™ Cocktail Enzyme mix (Ambion) according to the manufacturer’s instructions. Total RNA was extracted using TRIzol® LS (Ambion) combined with RNeasy® Mini kit (Qiagen).

Ds-cDNAs were prepared using random primers, SuperScript™ III Reverse Transcriptase (Invitrogen) and Klenow DNA Polymerase (3ʹ–5ʹ exo-)™ (Ambion) according to the manufacturer’s instructions. The amount of AvCoV genomes in each ds-cDNA was estimated by qPCR for AvCoV as described in the previous section and with a Qubit® dsDNA BR Assay Kit (Thermo Fisher Scientific). Libraries were prepared with Nextera® Sample Preparation (Illumina), purified with magnetic beads on an AMPure® XP system (Beckman Coulter Life Sciences) and finally quantified by qPCR using a KAPA® Fast Universal system (Sigma-Aldrich). Reads were obtained in a NextSeq500 MID Output 300 system (Illumina) (2 × 150 bp).

Dominant genomes for each treatment (considered passage 16 of the Beaudette strain) were assembled with CLC Genomics Workbench 20 (Qiagen) and annotated using the complete genome of the reference virus and the Beaudette strain compete genome (NC_001451.1) as references. For each treatment, the reads were then mapped again using the specific consensus for that passage and a low-frequency variant analysis was run in CLC using the following parameters: 100 coverage, 10 counts, 5% frequency, 1% significance and 20 central and neighborhood qualities.

Variants generated using the control serum at both dilutions (wells C1 and C2), MEM plus reference virus (well A3), reference virus only (well C3) and the original reference virus at passage 15 in VERO cells were excluded from those found for the treatments with GI-1 and GI-1/GI-11 sera to account for nonspecific virus neutralization, the dilution factor, mutations due to passage only and pleosiomorphic states, respectively.

All reads were deposited in the GenBank SRA database under accession number PRJNA736341.

Results

Full genomes obtained for the eight wells of the escape mutant assay were all 27,602 nt in length, with coverages ranging from 15,781 to 64,084 (Table 1). For the reference virus (Beaudette strain at VERO passage 15, virus load = 5.44E07 genome copies/µl), a genome of the same size was found, with a coverage of 72,443 (GenBank accession # MZ368698). The other full genomes were not submitted to GenBank to avoid redundancy.

Table 1

Virus load (genome copies) after qPCR, coverage, and total number of variants after the escape mutant assay with either GI-1, GI-1/GI-11, control sera at different dilutions and MEM (Minimal Essential medium) using the Beaudette GI-1 strain of AvCoV at VERO passage 15 and 102 TCID50/100 µl as a reference virus

Treatment on the reference virus	Virus load/µl sample	Coverage	Number of variants
GI-1 serum 10	1.57E7	28,895	60
GI-1 serum 100.1	2.14E7	57,502	59
GI-1/GI-11 serum 10	4.01E6	15,781	54
GI-1/GI-11 serum 100.1	7.78E5	49,754	51
Control serum	9.28E5	36,972	51
Control serum 1:10	9.32E5	47,834	51
MEM	4.87E8	64,084	81
Reference virus only	2.28E6	37,523	45

A comparison of the full dominant genomes recovered after the escape mutants assay showed that despite minor SNPs in the 5 and 3′ UTRs, all sequences were 100% identical except for the treatment with the undiluted control serum (9.28E5 genome copies/µl), which produced a synonymous A351T in NPS3 and a G29A in NPS9 resulting in the R10K amino acid change.

The highest number of variants was found for the reference virus treated 1:1 with MEM (n = 81, with a virus load of 4.87E8 genome copies/µl), while the lowest number was found for the challenge virus with no treatment (n = 45, with a virus load of 2.28E6 genome copies/µl) (Table 1). Treatments with high and low doses of the GI-1 serum resulted in a higher number of variants (n = 60 and 59, respectively) than the treatments with GI-1/GI-11 (n = 54 and 51, respectively), as shown in Table 1.

The same number of variants (n = 51) was found for the two doses of the control serum (Table 1). The Beaudette strain at VERO passage 15 (reference virus) showed a total of 40 variants. A more complete description of all variants can be found in Supplementary file 1.

After all variants related to nonspecific virus neutralization, the dilution factor, mutations due to passage only and pleosiomorphic states were excluded from the variants found for the GI-1 and GI-1/GI-11 treatments, the numbers of variants remaining were 10 and 4 (virus loads of 1.57E7 and 2.14E7 genome copies/µl) for the GI-1 serum at VN titers of 10 and 100.1, respectively (Table 1), and 8 and 9 (virus loads of 4.01E6 and 7.78E5 genome copies/µl) for the GI-1/GI-11 treatment at VN titers of 10 and 100.1, respectively (Table 2). Fisher’s exact test at α = 0.05 (https://www.socscistatistics.com) showed no statistically significant differences for the numbers of these variants.

Table 2

Frequencies (f) of variants found after the escape mutant assay of the Beaudette strain of the lineage GI-1 of AvCoV using chicken sera against lineages GI-1 (monovalent) or GI-1/GI-11 (bivalent) at doses 10 and 100.1 for each serum, the respective genome region they map to and amino acid changes. Variants shared with treatment with control sera, MEM (minimum essential medium), the original reference virus at passage 15 in VERO cells and no treatment (challenge virus only) were excluded

Position	Dominant	GI-1 10		GI-1100.1		GI-1/GI-11 10		GI-1/GI-11 100.1		Region	aa change
		Variant	f	Variant	f	Variant	f	Variant	f
24	A		0		0		0	C	6.8	5' UTR	na
30	TA	AT	19.13	AC	18.85	AT	21.94		0	5' UTR	na
30	TAC	CTT	5.16		0		0		0	5' UTR	na
31	AC		0		0	TA	7.5		0	5' UTR	na
41	T	A	5.69		0		0		0	5' UTR	na
49	A		0		0	G	5.46		0	5' UTR	na
52	TT	AA	6.41		0		0		0	5' UTR	na
63	A		0		0	T	5.08		0	5' UTR	na
120	T	A	5.13		0		0		0	5' UTR	na
140	del		0		0	ATGAG	6.27		0	5' UTR	na
381	A		0		0		0	G	6.28	5' UTR	na
461	C	G	6.47	G	8.02		0		0	5' UTR	na
466	A	T	5.33	T	7.53		0	T	5.15	5' UTR	na
3830	G	T	5.1		0		0		0	nsp3	R428L
3853	G	C	5.11		0		0		0	nsp3	E436Q
9953	T		0		0	A	5.3		0	nsp6	V56D
18,048	G	T	8.97		0		0		0	nsp14	G375C
20,193	T		0		0		0	A	5.63	nsp16	Y231N
20,362	A		0		0		0	G	5.13	nsp16	D287G
20,367	GGT		0		0		0	ACC	5.45	nsp16/S	nsp16 G289T/ S V3P
21,280	T		0		0		0	A	5.56	S	F307I
23,952	A		0		0		0	G	5.88	E	I35V
23,959	T		0		0		0	del	5.27	E	*50
24,276	G		0	del	6.45		0		0	E	*31
25,440	T		0		0	A	6.06		0	NC	na
25,455	T		0		0	G	5.07		0	NC	na

NC non-coding, na non applicable

*Stop codon

Only three of these variants were shared among two or more treatments, and they were all in the 5′UTR. Early stop codons were detected at codons 31 and 50 of the envelope protein gene after treatments with GI-1 at a titer of 100.1 and GI-1/GI-11 at a titer of 100.1, respectively (Table 2 and Fig. 1).

Fig. 1

Schematic representation of Avian coronavirus genome showing the proteins in which variants with amino acids changes were detected after the escape mutant assay of the Beaudette strain of the lineage GI-1 of AvCoV (reference virus) using chicken sera against lineages GI-1 (monovalent) at doses 10 (black circles) and 100.1 (white circles) or GI-1/GI-11 (bivalent) at doses 10 (black triangles) and 100.1 (white triangles). Numbers bellow each amino acid change are frequencies of occurrence. The amino acids positions refer to the full length of each respective protein. The genes of each protein are not drawn to scale

Variants in the spike glycoprotein were found exclusively for the GI-1/GI-11 treatment at a titer of 100.1 (V3P and F307I, Table 2). The only other putatively viable variant in a structural protein was found in the envelope protein (I35V) after this same treatment (Table 2 and Fig. 1).

Other exclusive variants mapped to NPSs 3, 6, and 14 for the treatment with GI-1 serum at VN titer 10, NPS16 for GI-1/GI-11 serum at VN titer 100.1 and to the non-coding region between ORF X and gene 5 for GI-1/GI-11 serum at VN titer 10 (Table 2).

Discussion and Conclusion

Using chicken mono- and bivalent anti-AvCoV sera produced after vaccination with inactivated vaccines, subdominant variants could be detected along the genome after a single round of virus neutralization in vitro.

The F307I substitution in the S1 portion of the spike protein detected after the treatment of the challenge virus with the bivalent GI-1/GI-11 serum at a titer of 100.1 was unique to this treatment and mapped to a neutralizing region of the spike (Ignjatovic and Sapats 2005; Kant et al. 1992; Moore et al. 1997), which indicated that this variant could be an escape mutant.

A G variant at position 23952, which resulted in an I35V amino acid change in the envelope protein E exclusive to the GI-1/GI-11 serum at a titer of 100.1, might have been caused by a lack of complete virus neutralization. The E protein plays a major role in coronavirus virion assembly (Vennema et al., 1996; Westerbeck and Machamer 2015) and has been predicted to harbor B-cell epitopes (Lu et al. 2021).

Variants in the 5ʹ UTR have been reported to the porcine HKU15 Deltacoronavirus in clinical samples (Woo et al. 2017), but their consequences on virus fitness have not yet been understood. The A63T 5ʹ UTR variant found after treatment with GI-1/GI-11 serum at a titer of 10 (Table 2) mapped to the transcription regulatory sequence (TRS) (Keep et al. 2020) and could thus affect mRNA transcription on the mutant spectrum.

Though NSPs and non-coding RNA are not expected to be directly affected by neutralizing antibodies immune pressure, variants in NSPs 2, 12 and 15 have been detected in a COVID-19 patient during convalescent plasma therapy besides S variants, however no sustained change in virus load could be detected (Kemp et al. 2021).

One of the limitations of this study was that the ELISA used to titrate the anti-AvCoV antibodies in mono- and bivalent sera could not differentiate between antibodies against the GI-1 and GI-11 genotypes or measure the effects of genotype-specific antibodies on the emergence of escape mutants as it is targeted to generic rather than type-specific antigens. Moreover, considering that the virus neutralization test is genotype-specific, another limitation is that although antibodies against the GI-1 genotype were measured in sera using the GI-1 Beaudette strain, a cell culture-adapted GI-11 strain is not available to measure antibodies against this type and determine the emergence of escape mutants using it as a challenge virus because these field strains of AvCoV are difficult to grow in mammalian cells.

Another limitation is that haplotype calling to distinguish the population of individual genomes in which different variants co-occur was not carried out due to the short-reads approach used for deep sequencing, restricting the discussion to the impacts of individual variants.

The results indicate that at low antibody titers, avian coronavirus escape mutants may arise due to polymorphisms on the spike glycoprotein; moreover, further polymorphisms in nonstructural proteins and non-coding regions could have a role in virus fitness. As no statistically significant difference for the numbers variants was found when the treatments were compared, the fluctuation in virus loads was not influenced by the number of variants. However, the fitness of these escape mutants after further passage in cell culture with and without more rounds of neutralization has not been clarified. A more diverse set of vaccine strains to prepare sera would also allow for a better assessment of the relationship between serum specificity and escape mutant predictions.

Below is the link to the electronic supplementary material.

Frequencies of variants found after the escape mutant assay of the Beaudette strain of the lineage GI-1 of AvCoV (reference virus) using chicken sera against lineages GI-1 (monovalent) or GI-1/GI-11 (bivalent) at doses 10 and 100.1 for each serum, control serum, MEM and the reference virus only and the respective genome region they map to. Supplementary file1 (XLSX 42 kb)