Circulation of SARS-CoV-2 Variants among Children from November 2020 to January 2022 in Trieste (Italy)

Introduction: The ongoing coronavirus disease 19 (COVID-19) outbreak involves the pediatric population, but to date, few reports have investigated the circulation of variants among children. Material and Methods: In this retrospective study, non-hospitalized pediatric patients with SARS-CoV-2-positive nasopharyngeal swabs (NPS) were enrolled at the Institute for Maternal and Child Health-IRCCS Burlo Garofolo, Trieste (Italy), from November 2020 to January 2022. SARS-CoV-2 variants were identified by in vitro viral isolation, amplification, automatic sequencing of the receptor binding domain (RBD) of the SARS-CoV-2 spike coding gene, and subsequent next-generation sequencing. The growth curves of the isolated strains were defined in vitro by infecting Vero-E6 cells and quantifying the viral load in the supernatants up to 72 h post-infection by qRT–PCR. The neutralization activity of sera obtained from a COVID-19 vaccinated subject, recovered (2020) patient, vaccinated and recovered (2021) patient, and seronegative subject was assessed by microneutralization assay against the different variants. Results: In total, 32 SARS-CoV-2-positive children, 16 (50%) females, with a median age of 1.4 years (range: 1 day–13 years), were enrolled. The D614G amino acid substitution was detected in all isolated and amplified viral strains. Of the 32 isolates, 4 (12.5%) carried a nonsynonymous nucleotide mutation leading to the N439K (3/4), lineage B.1.258 (∆H69/∆V70), and S477N (1/4) substitution. In 7/32 (21.8%) isolates, amino acid substitutions allowed the identification of a delta variant, lineage B.1.617.2-AY.43, and in 1/32 (3.1%), the Omicron strain (B.1.1.529.BA1) was identified. The growth curves of the B.1, B.1.258 (∆H69/∆V70), B.1.617.2-AY.43, and B.1.1.529.BA1 variants did not show any significant differences. A reduction in the serum neutralizing activity against B.1.258 (∆H69/∆V70) only in a vaccinated subject (1.7-fold difference), against B.1.617.2-AY.43 in a vaccinated subject and in recovered patients (12.7 and ≥2.5-fold differences, respectively), and against B.1.1.529.BA1 variant (57.6- and 1.4-fold differences in vaccinated and in vaccinated and recovered patients) were observed compared to the B.1 variant. Conclusions: SARS-CoV-2 variants carrying the B.1.258 (∆H69/∆V70) and S477N substitutions were reported here in a pediatric population for the first time. Although the growth rates of the isolated strains (B.1.258, B.1.617.2-AY.43, B.1.1.529.BA1) did not differ from the B.1 variant, neutralizing activity of the sera from vaccinated subjects significantly decreased against these variants. Attention should be devoted to the pediatric population to prevent the spread of new SARS-CoV-2 variants in an unvaccinated and predominantly naive population.


Introduction
The ongoing coronavirus disease  outbreak is caused by the severe acute respiratory syndrome coronavirus-2 (SARS-CoV-2). The first cases of atypical pneumonia were reported on 31 December 2019, in Wuhan, China, and as of 26 December 2021, 281 million cases and 5.41 million deaths have been reported worldwide [1,2]. Regarding the pediatric population, few studies have indicated the infection rate among children, and 1-10% of total COVID-19 cases have been observed in this population [3][4][5]. The low rate of SARS-CoV-2 detection via nasopharyngeal swabs (NPS) may be due to mild/asymptomatic infection that escapes detection, the lower susceptibility to infection among children than adults (as children express fewer ACE2 receptors), or the prioritization of testing adults [6][7][8]. A multinational, multicenter cohort study conducted among 25 different countries in Europe provided an overview of the clinical impact of SARS-CoV-2 infection. The study found that COVID-19 is milder and has lower hospitalization rates among children than adults [9]. Regarding neonates, SARS-CoV-2 infection seems to be mainly transmitted from mothers, even if vertical transmission is still unclear [10,11]. Family clusters might be the main source of infection for children [12], while educational settings increase the circulation of the virus among adolescents [13]. In addition to adults, children are a reservoir for SARS-CoV-2 transmission, even if they are asymptomatic or affected by mild illness, thereby promoting the evolution of the virus [14][15][16]. To date, few studies have examined the circulation of SARS-CoV-2 variants of concern (VOCs) and variants of interest (VOIs) in the pediatric population. In Israel, when the alpha variant circulated from December 2020 to February 2021, the virus transmission rate increased among children under 10 years of age [13]. In contrast, in a study conducted by Howard-Jones, the SARS-CoV-2 Delta strain was shown to have a lower level of transmission among children than among adults [17]. To obtain more data on the circulation of SARS-CoV-2 variants in children, the current retrospective study examined 32 pediatric SARS-CoV-2-positive patients with nasal-pharyngeal swabs (NSPs) collected between 18 November 2020, and 7 January 2022, at the Institute for Maternal and Child Health-IRCCS Burlo Garofolo, Trieste, Friuli-Venezia Giulia, Italy. The in vitro growth curves of the isolated strains were assessed, and the microneutralization activity of antibodies elicited by the COVID-19 vaccine and natural antibodies after infection against these variants was evaluated.

Study Design
From 18 November 2020 to 7 January 2022, 32 NPSs were positive for SARS-CoV-2 molecular detection and collected from 32 non-hospitalized children at the Institute for Maternal and Child Health-IRCCS Burlo Garofolo, Trieste, Friuli-Venezia Giulia, Italy. Ethical approval was obtained from the local committee (IRB-Burlo Garofolo IRCCS Trieste, Approval Code: 18.11.20, Approval Date: September 2020), and informed consent was signed by the patients' parents.

RNA Isolation from NPS and SARS-CoV-2 Detection
Automated extraction of viral RNA was performed from 200 µL of NPS using the commercial Maxwell ® RSC Viral Total Nucleic Acid Purification Kit (Promega, Madison, WI, USA), according to the manufacturer's instructions. The presence of the SARS-CoV-2 genome was evaluated by real-time PCR using the Neoplex™ COVID-19 Detection Kit (Genematrix, Korea), which simultaneously targeted the RdRp and N viral genes. The test was performed according to the manufacturer's instructions. The reaction was made up in 20 µL, containing 5 µL of 1-step master mix, 5 µL of COVID-19 PPM, 5 µL of RNase-free water, and 5 µL of RNA. The thermal cycles were as follows: 30 min at 50 • C, 15 min at 95 • C, and 45 cycles at 95 • C for 15 s and at 61 • C for 60 s.

Quantitative Reverse Transcription Real-Time PCR (RT-qPCR)
The SARS-CoV-2 N1 gene amount was assessed by means of a quantitative reverse transcription real-time PCR (qRT-PCR) assay using the Applied Biosystem 7500 Real-Time PCR system instrument (ThermoFisher, Waltham, MA, USA) and the CDC and WHO 2019-Novel Coronavirus (2019-nCoV) Real-Time RT-PCR Diagnostic Panel for the detection of nucleic acids from 2019-nCoV (N1 gene) [18,19] (Table 1). The reaction was made up in 25 µL, containing 0.4 µM forward primer, 0.4 µM reverse primer, 0.1 µM probe, 12.5 µL 2x RT-PCR Buffer, 1 µL 25x RT-PCR Enzyme Mix, and 5 µL of RNA using the AgPath-ID One-Step kit (ThermoFisher, Waltham, MA, USA). The thermal cycles were as follows: 30 min at 45 • C, 10 min at 95 • C, and 45 cycles at 95 • C for 15 s and at 55 • C for 30 s. The standard curve was constructed using serial dilution of the pEX-A128-nCoV_all plasmid (Eurofins, Luxemburg) containing part of the SARS-CoV-2 genome. Samples were analyzed in duplicate, and a negative control consisting of a mix supplemented with water was added. The limit of viral genome detection was evaluated as 3 copies/reaction.

Sanger Sequencing of the SARS-CoV-2 Receptor Binding Domain (RBD) Genomic Region
The receptor-binding domain (RBD) of the SARS-CoV-2 spike-encoding region was amplified and sequenced by PCR and subsequent Sanger automatic sequencing. A one-step RT-PCR kit (Thermo Fisher, Waltham, MA, USA) was employed for the amplification step using the set of primers summarized in Table 2. The 50 µL comprised 25 µL of 2× reaction mix, 2 µL of SuperScript III RT/Platinum Taq Mix (Thermo Fisher, Waltham, MA, USA), 1 µL of forward primer (10 µM), 1 µL of reverse primer (10 µM), and 5 µL of RNA sample. The PCR steps were as follows: 1 cycle at 50 • C for 15 min and at 95 • C for 2 min, followed by 40 cycles at 95 • C for 15 s, 54 • C for 30 s, and 68 • C for 40 s, followed by a last step at 68 • C for 5 min. The amplification length product (1000 bp) was checked on 1.5% agarose gels in 0.5× TBE buffer by electrophoresis. PCR products were purified using the QIAquick Gel Extraction Kit (QIAgen, Hilden, Germany) in accordance with the manufacturer's instructions. Automatic sequencing of the fragments was performed by an external facility (Eurofins Genomics, Munich, Germany) using both forward and reverse primers with the Sanger method. The sequence similarity was searched through the NCBI site (USA) using BLAST [20], comparing the amplified fragments to the reference Wuhan-Hu-1 strain (NC_045512.2). The predicted amino acid sequences were obtained by blastx [20] and compared to the Wuhan-Hu-1 RBD protein sequence with blastp [20].

Next Generation Sequencing
RNA samples that underwent to Next-Generation Sequencing (NGS) were quantitated using a Qubit 3.0 Fluorometer (ThermoFisher Scientific, Waltham, MA, USA) and Qubit RNA HS Assay Kit (Invitrogen Q32855). A RIN for each sample was determined by running 2uL on an Agilent 4200 TapeStation using High-Sensitivity RNA Screen Tapes (Agilent 5067-5579). Then, 5-7 ng of each sample was used to performed libraries following the SMART-Seq Stranded Kit User Manual (Takara Bio 634444). Briefly, cDNAs were first generated from all total RNA fragments after RNA fragmentation. Then, Illumina adapters, indexes, and library purification followed. Final library amplification and purification were performed after enzymatically removing ribosomal fragments originating from rRNA molecules using specific probes. Libraries were quantitated using the Qubit dsDNA HS Assay Kit (Invitrogen, Q32854). The average size of each library was determined by running each sample on Agilent 4200 TapeStation using High-Sensitivity D5000 Screen Tapes (Agilent 5067-5593). Each library was diluted to 4 nM in Resuspension Buffer (Illumina). All libraries were pooled at a concentration of 4 nM before denaturing with 0.2 N NaOH for 5 min at RT. The libraries were neutralized using 200 mM Tris-HCl, pH 7 before being diluted to 20 pM in HT1 buffer (Illumina). Finally, denatured libraries were loaded into an Illumina NextSeq 550 v2 High-Output reagent cartridge at a final concentration of 1.4 pM in HT1 buffer. Denatured PhiX control library was also included at 5% concentration. Dual indexes and paired-end 150 bp sequencing was performed on a High-Output flow cell, yielding about 20M read pairs per sample.

Phylogenetic Analysis
The S gene sequences included in this study were aligned using the ClustalW program, implemented in the BioEdit software [21]. Intragroup nucleotide identity was calculated using the Sequence Identity Matrix tool of BioEdit software [21]. The mean intragroup nucleotide identities were expressed as crude percentages with the corresponding minimum and maximum range values.
The phylogenetic analysis was conducted by means of MEGA software, version 6.0 [22]. The best substitution model was selected by analyzing the sequence dataset with the Models Tool available in MEGA software, version 6.0 [22]. The phylogenetic tree was constructed by means of the Neighbor-Joining method and Tamura 3-parameter model, including the S gene sequences of SARS-CoV-2 study strains and reference strains retrieved from the Global Initiative on Sharing All Influenza Data (GISAID) EpiFlu database [23]. The robustness of the phylogenetic tree was estimated by means of a bootstrap analyses with 1000 replicates [24].

Virus Isolation
SARS-CoV-2 isolation and expansion were conducted on the Vero E6 cell line, starting from NPS. Briefly, a total of 4 × 10 5 cells/well were seeded in a 6-well culture multi-well plate and left to adhere overnight in complete medium at 37 • C with 5% CO 2 until 80% confluence was reached. Cell infection was obtained by adsorption of 500 µL of NPS for 120 min at 37 • C with 5% CO 2 . After 2 h of adsorption, the inoculum was removed, the cells were washed twice with Dulbecco's phosphate-buffered saline (DPBS) (EuroClone, Italy), and 2 mL of fresh complete medium was added. Vero E6 cells were maintained for 3 days or more, checking for cytopathic effect (CPE) daily. The procedure was carried out in a biosafety level 3 (BSL3) facility. Infected cell medium was then recovered and stored at −80 • C until use.

Plaque Reduction Assay
The virus titer of supernatant from the infected cells was assessed by means of a plaque reduction assay. Briefly, 4 × 10 5 cells/well were seeded in a 6-well culture multi-well plate and left overnight in complete medium at 37 • C with 5% CO 2 until 80% confluence was reached. Then, 500 µL of 10-fold serially diluted supernatants was inoculated onto Vero E6 monolayer cells in duplicate for 2 h. The viral inoculum was removed, and cells were overlaid with complete medium with 3% agarose and maintained for 48 h or 72 h at 37 • C with 5% CO 2 . The cell monolayer was then fixed with paraformaldehyde 3.7% and stained with methylene blue. The titer of the virus was calculated as plaque-forming units (PFU) per milliliter (PFU/mL) and used to define the Multiplicity of Infection (MOI) for the subsequent experiments.

Variant Strain Growth Curves: Kinetics of Infection
A total of 4.5 × 10 4 Vero E6 cells/well were seeded at 80% confluence in a 48-well culture multi-well plate and left overnight in complete medium at 37 • C with 5% CO 2 . Cell infection was performed by inoculating SARS-CoV-2 isolates at aMOI of 0.01 for 2 h at 37 • C with 5% CO 2 . Growth curves were evaluated at 24, 48, and 72 h post-infection by means of RT-qPCR, as described above. Two replicates for each time points were tested.

Tissue Culture Infectious Dose Assay
The 50% Tissue Culture Infectious Dose (TCID 50 ) assay was used to measure the infectivity of the isolated SARS-CoV-2. Briefly, virus serial 1:10 dilutions were performed from 1:10 to 1:10 8 . The prepared dilutions were then transferred on Vero E6 cells, seeded the day before, and incubated for 2 h at 37 • C with 5% CO 2. Then, complete medium was added, and the plates were kept at 37 • C with 5% CO 2 for 5-7 days. The CPE end-point method was used for the viral titration, and CPE was counted visually. TCID 50 values were obtained applying the Reed and Muench formula [25,26].

Microneutralization Assay
For the microneutralization assay, 1.5 × 10 4 cells/well were seeded in a 96-well culture multi-well plate and left overnight in complete medium at 37 • C with 5% CO 2. A microneutralization assay was carried out using four different serum samples that were collected: (a) 6 months after the infection from a COVID-19 adult patient infected by the B.1 strain and recovered in May 2020; (b) 2 weeks after the administration of the second dose of Comirnaty (Pfizer-BioNTech) vaccine from a uninfected, vaccinated adult subject; (c) from a recently recovered (December 2021) patient (2 weeks after the infection with B.1.1.529) who became infected after receiving the third dose of the Comirnaty (Pfizer-BioNTech) vaccine; (d) from a SARS-CoV-2-seronegative subject. Serum samples were heat-inactivated at 56 • C for 40 min and serially diluted two-fold in complete medium, from 1:2 to 1:1024, in quintuplicate. In total, 200 TCID 50 of virus was mixed with serum dilutions and incubated at 37 • C with 5% CO 2 for 1 h. Then, 30 µL of virus-serum mix was seeded on a Vero E6 cell monolayer, grown on 96-well microplates, and incubated at 37 • C with 5% CO 2 for 2 h. The complete medium was then added to the cells, which were inspected for CPE with an inverted optical microscope 6 days post-infection. A neutralizing antibody titer that reduces the number of infected wells by 50% (NT 50 ) was defined by the means of the Reed and Muench method [26,27].

Statistical Analysis
Pearson correlation was performed to evaluate the difference in the variants' growth curves. NT 50 values against the isolated variants were analyzed using the Mann-Whitney nonparametric test. p values < 0.05 were considered significant.
The phylogenetic relationship between study strains was assessed by calculating the pairwise p-distance implemented in the MEGA software, version 6.0 [22]. The mean nucleotide genetic distances were expressed as crude rates with the corresponding standard deviation (sd).

Sequencing of the SARS-CoV-2 RBD Spike Region
SARS-CoV-2 RBD sequences, obtained by means of amplification and subsequent Sanger automatic sequencing, were aligned with the reference SARS-CoV-2 Wuhan-Hu-1 strain (NC_045512.2). The alignment results are reported in Table 3. Briefly, in all sequences (32/32, 100%), the D614G substitution was detected, lineage B.1. A total of 6 out of 32 (18.8%) sequences showed a nucleotide point mutation leading to a synonymous amino acid substitution, while 4/32 (12.5%) sequences carried a total of 5 nucleotide mutations, 3 of which (C22879A) led to the N439K amino acid substitution and 1 of which (G22992A) led to the S477N amino acid substitution. In total, 7/32 (21.8%) strains showed the presence The nucleotide mutations and amino acid mutations are shown in Table 4. The nucleotide sequences were deposited in the GISAID EpiCov TM online database.   We decided to perform the full genome sequencing of the viral strains that were more represented in the sampled population: the N439K, present in 3/32 patients, and the Delta-like, present in 7/32 patients. Full genome sequences obtained by the NGS of a single SARS-CoV-2 isolated strain, N439K, and a single Delta-like strain, which were deposited at GISAID (N439K: EPI_ISL_7208675; Delta: EPI_ISL_7698448).

Phylogenetic Analysis
The phylogenetic analysis of 31/32 study SARS-CoV-2 nucleotide sequences revealed an intragroup nucleotide identity of 99.6% (range: 97.5-100%) and a mean genetic distance (p-distance) of 0.004 ± 0.010. This analysis demonstrates that the study sequences were genetically related.
However, as shown in the phylogenetic tree in Figure 1, the study sequences were distributed in three major clades. In total, 23/31 sequences were "no-VOC," showing an intragroup nucleotide identity of 99.9% (range: 99.5-100%) and a p-distance of 0.001 ± 0.000. These sequences were strictly related each other and were detected from November 2020 to February 2021. In addition, 7/31 sequences belonged to the lineage B.1.617.2, also known as Delta VOC. These sequences showed an intragroup nucleotide identity of 100%, and were detected in July 2021 (1/7) and in November 2021 (6/7). Finally, 1/31 sequences belonged to the lineage B.1.1.529, also known as Omicron VOC. This sequence was detected in January 2022.

SARS-CoV-2 Isolates Growth Curves
Isolation and growth curves of SARS-CoV-2 strains carrying nonsynonymous amino acid mutations were performed. A comparison was performed with the SARS-CoV-2 B.1 strain (EPI_ISL_584051) previously isolated in our laboratory [29]. No significant differences were observed between the four growth curves (p value > 0.05) (Figure 2). The viral loads increased starting from 24 h post-infection for the four strains.

Microneutralization Assay
The results obtained from the microneutralization assay were evaluated by inspecting the inhibition of the CPE for each serum dilution by an inverted optical microscope. The NT 50 for the SARS-CoV-2 tested strains, mixed with the sera from a Comirnaty (Pfizer-BioNTech) COVID-19 vaccinated subject, COVID-19 patient recovered in 2020, patient who was infected after the third dose of Comirnaty (Pfizer-BioNTech) and recovered recently (December 2021), and seronegative patient were determined by means of the Reed and Muench method (Figure 3). NT 50

SARS-CoV-2 Isolates Growth Curves
Isolation and growth curves of SARS-CoV-2 strains carrying nonsynonymous amino acid mutations were performed. A comparison was performed with the SARS-CoV-2 B.1 strain (EPI_ISL_584051) previously isolated in our laboratory [29]. No significant differences were observed between the four growth curves (p value >0.05) (Figure 2). The viral loads increased starting from 24 h post-infection for the four strains.

Microneutralization Assay
The results obtained from the microneutralization assay were evaluated by inspecting the inhibition of the CPE for each serum dilution by an inverted optical microscope. The NT50 for the SARS-CoV-2 tested strains, mixed with the sera from a Comirnaty (Pfizer-BioNTech) COVID-19 vaccinated subject, COVID-19 patient recovered in 2020, patient who was infected after the third dose of Comirnaty (Pfizer-BioNTech) and recovered recently (December 2021), and seronegative patient were determined by means of the Reed and Muench method ( In summary, a reduction in serum neutralizing activity against variants was observed, especially against B.1.617.2+AY.43 (GISAID ID: EPI_ISL_7698448) and B.1.1.529.BA1 (EPI_ISL_10156328), but the differences were not significant (p value > 0.05).

Discussion
The continuing spread of SARS-CoV-2 in the human population induces the emergence of new variants that could impact the epidemiological situation of geographical areas where these variants become predominant. The transmissibility, severity of these strains, and immune response to the infection may be altered [30]. Children are a poorly studied population, and very few studies have investigated the circulation of SARS-CoV-2 and its variants in this population. Children are reservoirs of the virus, contributing to the spread of SARS-CoV-2 variants [16]. Currently, COVID-19 vaccines are available and recommended for children aged 5 and above [31], but no vaccines have been authorized for children aged 0-4 years. In summary, a reduction in serum neutralizing activity against variants was observed, especially against B.1.617.2+AY.43 (GISAID ID: EPI_ISL_7698448) and B.1.1.529.BA1 (EPI_ISL_10156328), but the differences were not significant (p value > 0.05).

Discussion
The continuing spread of SARS-CoV-2 in the human population induces the emergence of new variants that could impact the epidemiological situation of geographical areas where these variants become predominant. The transmissibility, severity of these strains, and immune response to the infection may be altered [30]. Children are a poorly studied population, and very few studies have investigated the circulation of SARS-CoV-2 and its variants in this population. Children are reservoirs of the virus, contributing to the spread of SARS-CoV-2 variants [16]. Currently, COVID-19 vaccines are available and recommended for children aged 5 and above [31], but no vaccines have been authorized for children aged 0-4 years.
In this regard, the focus of this paper was to analyze the sequences of SARS-CoV-2 isolated from SARS-CoV-2-positive NPSs collected from 32 children to investigate the circulation of variants in a geographic area of northeastern Italy.
The patients were all infected with the SARS-CoV-2 B.1 variant, which appeared in North Italy in April 2020 and was the predominant variant, together with its derivatives, until February 2021 [32]. The B.1 strain shows the D614G substitution, able to confer moderate advantages in viral infectiveness and transmissibility, compared to the Wuhan-Hu-1 strain [33,34].
Among the isolates, the N439K substitution was found in three samples. Based on the Pangolin nomenclature, the SARS-CoV-2 strain presenting this mutation was assigned to the B.1.258 lineage. The strain was detected for the first time in Scotland in March 2020 [35], and spread quickly in many European countries [36]. The N439K substitution occurs in the RBD, the highly variable region of the SARS-CoV-2 spike protein, which interacts with ACE2 receptors for viral entry and is the main target of neutralizing antibodies [36][37][38]. This mutation seems to increase the affinity of virus S protein to the ACE2 receptor, and reduce the activity of designed monoclonal antibodies and serum polyclonal antibodies from recovered COVID-19 patients [35,36]. Consequently, N439K might evade the immune system either in recovered COVID-19 patients or in vaccinated subjects. Our results are in line with the published literature, indicating that the neutralizing activity of sera from the COVID-19 vaccinated subject decreased against the B.1.258 strain up to two-fold compared to the B.1 strain.
NGS also revealed the presence of the ∆69/70 mutation in the B.1.258 strain, a mutation that emerged independently in the other five strains (namely B.1.1.7, B.1.1.298, B.1.160, B.1.177, B.1.375), always after the acquisition of receptor binding motif replacements, such as N439, in several European countries in the beginning of August 2020. During the fall of 2020, the B.1.258 variant gained significant prevalence in multiple countries, and Brejova and colleagues speculated that its diffusion was associated with increasing hospitalization trends, starting before the appearance of the B.1.1.7 variant [39]. Different in vitro studies have reported contrasting results regarding the increased replication of the strains carrying the deletion [40,41]. In this regard, we did not observe a significant difference in the replication trend of the isolated SARS-CoV-2 strains.
Of note, NGS analysis also revealed the presence of a stop codon in the ORF7a coding region (G58) due to the G/T mutation at nt 27505. ORF7a, together with ORF3, ORF6, and ORF8, is an accessory viral protein involved in SARS-CoV-2 immune escape mechanisms since they are able to antagonize type I IFN production and signaling [42][43][44].
Mutations are not infrequent in the ORF7 region. Nemudriy and colleagues reported that C-terminal mutations are often reported in patients worldwide. In particular, mutations leading to the truncation of the ORF7a products have been previously observed, and they have led to a loss of protein function that abrogates its interaction with IFN signaling [45].
To the best of our knowledge, this is the first study to report the detection of the N439K substitution in a pediatric population. No significant differences were observed among the clinical outcomes of the three children positive for N439K and the others, confirming the findings reported by Thomson and colleagues, who found no association between this substitution and increases in COVID-19 symptomatology and severity [36]. In contrast, we did not find any differences in the viral load between the N439K and D614G strains [36].
The variant carrying the S477N substitution spread in some European countries, with a higher prevalence in France during the summer of 2020 [46] and a frequency of detection ranging from 4-7% [47]. In our study, the S477N substitution was detected in one sample, with a frequency of 3.1%, resembling its frequency rate in Europe, which ranges from 4 to 7% [47]. As for N439K, S477N occurs in the RBD region and seems to confer enhanced binding ability to the ACE2 receptor [35], but it does not seem to be associated with the loss of neutralizing activities of convalescent serum antibodies [35,48,49].
In our cohort, and compared to the D614G strain, the S477N viral load was slightly lower, while no difference in the clinical presentation was observed.
Seven isolates in our cohort of children were Delta isolates collected from July 2021 to November 2021. Starting in May 2021, the Delta variant became the predominant variant in many European countries, drastically enhancing the spread of the virus and causing a high number of infections [50]. The clinical symptoms of SARS-CoV-2 Delta-positive children were similar to those affecting other positive children, suggesting that the variant does not affect the clinical course in the pediatric population independently by the variant [51].
As reported in another study [52], we observed a marked reduction (up to almost 13-fold) in the effectiveness of the Comirnaty vaccine against the delta-like variant in vitro. A decrease in the neutralization ability against this variant was also observed testing the sera from both the recovered and the vaccinated and recovered patients (up to 5-fold and up to 2.5-fold, respectively).
Since January 2022, the Omicron variant has become the predominant one in Italy, as well as the rest of Europe. We collected the NPS from one asymptomatic, Omicron-infected child, and we observed that the B.1.1.529.BA1 variant could be neutralized only after the booster dose of vaccine, as already reported in the literature [53,54]. No or very weak neutralization was indeed observed in the recovered patient, and in the subject with two doses of vaccine, respectively.

Conclusions
Although we are aware that the conclusions arose from a limited number of patients, it should be underlined that studying the pediatric population is not simple, since children have a smaller number, lower incidence, and milder symptoms than adult patients [55,56].
Further investigations on a large group of children are needed to assess the circulation of SARS-CoV-2, especially because children contribute to the spread of the virus like adults. The vaccination of young people will be an important step to stop the spread of the virus.