Higher entropy observed in SARS-CoV-2 genomes from the first COVID-19 wave in Pakistan

Background We investigated the genome diversity of SARS-CoV-2 associated with the early COVID-19 period to investigate evolution of the virus in Pakistan. Materials and methods We studied ninety SARS-CoV-2 strains isolated between March and October 2020. Whole genome sequences from our laboratory and available genomes were used to investigate phylogeny, genetic variantion and mutation rates of SARS-CoV-2 strains in Pakistan. Site specific entropy analysis compared mutation rates between strains isolated before and after June 2020. Results In March, strains belonging to L, S, V and GH clades were observed but by October, only L and GH strains were present. The highest diversity of clades was present in Sindh and Islamabad Capital Territory and the least in Punjab province. Initial introductions of SARS-CoV-2 GH (B.1.255, B.1) and S (A) clades were associated with overseas travelers. Additionally, GH (B.1.255, B.1, B.1.160, B.1.36), L (B, B.6, B.4), V (B.4) and S (A) clades were transmitted locally. SARS-CoV-2 genomes clustered with global strains except for ten which matched Pakistani isolates. RNA substitution rates were estimated at 5.86 x10−4. The most frequent mutations were 5’ UTR 241C > T, Spike glycoprotein D614G, RNA dependent RNA polymerase (RdRp) P4715L and Orf3a Q57H. Strains up until June 2020 exhibited an overall higher mean and site-specific entropy as compared with sequences after June. Relative entropy was higher across GH as compared with GR and L clades. More sites were under selection pressure in GH strains but this was not significant for any particular site. Conclusions The higher entropy and diversity observed in early pandemic as compared with later strains suggests increasing stability of the genomes in subsequent COVID-19 waves. This would likely lead to the selection of site-specific changes that are advantageous to the virus, as has been currently observed through the pandemic.


Results
In March, strains belonging to L, S, V and GH clades were observed but by October, only L and GH strains were present. The highest diversity of clades was present in Sindh and Islamabad Capital Territory and the least in Punjab province. Initial introductions of SARS- and S (A) clades were transmitted locally. SARS-CoV-2 genomes clustered with global strains except for ten which matched Pakistani isolates. RNA substitution rates were estimated at 5.86 x10 −4 . The most frequent mutations were 5' UTR 241C > T, Spike glycoprotein D614G, RNA dependent RNA polymerase (RdRp) P4715L and Orf3a Q57H. Strains up until June 2020 exhibited an overall higher mean and site-specific entropy as compared with sequences after June. Relative entropy was higher across GH as compared with GR and L clades. More sites were under selection pressure in GH strains but this was not significant for any particular site.

Introduction
Severe Acute Respiratory Syndrome coronavirus 2 (SARS-CoV-2), the causative agent of COVID-19 was first reported in Wuhan, China [1] and has been to date reported in 190 million cases globally. SARS-CoV-2 RNA belongs to the Sarbecovirus subgenus of Betacoronavirus [2]. The first Wuhan strain identified in January 2020 was an L clade isolate and by the middle of 2020, SARS-CoV-2 was seen to diversify into S, V and G and its sub clades [3,4]. SARS-CoV-2 has 8 coding and 6 non-coding genes and earlier, the greatest genetic variation observed in Nucleocapsid (N) and Orf1ab regions [5] however, since the end of 2020, evolution in Spike glycoprotein (S) regions has resulted in a number of variants of concern which are spreading rapidly across the globe [6]. Variants observed in Orf1ab, Orf3a, N and S genes are associated with evolutionary changes and thirteen signature single nucleotide variations (SNVs) divide strains into L, S, V, I, and G and its subclades as defined by GISAID [7]. Additional mutations have further characterize SARS-CoV2 genomes into GH, GR and GV clades [8].
In Pakistan, the first wave of COVID-19 occurred between March and July 2021, with a peak in mid-June where 4000-6000 positive COVID-19 cases were diagnosed each day [9]. There was a decrease in cases with a low in early September followed by a resurgence of COVID-19 cases with a second wave from October 2020 until January 2021. There was a brief reprieve in February with a third wave of cases between March and May 2021 [10]. COVID-19 cases started to rise again in June and currently in July 2021, Pakistan is well into the fourth wave [11]. Up to 17 July 2021, approximately 986,668 COVID-19 cases have been diagnosed with 22,760 deaths [12]. Of these, 354,103 (36%) COVID-19 cases were from Sindh province with 21% were from the city of Karachi. Provincial distribution of cases was found to be 349,890 (35%) from Punjab, 140,293 (14%) from Khyber Pakhtunkhwa Province (KPK), 84,399 (9%) from Islamabad Capital Territory (ICT), 28,884 (3%) from Baluchistan, 21,811 (2%) from Azad Jammu Kashmir and 7288 (1%) from Gilgit Baltistan. To date, the case fatality rate (CFR) for SARS-CoV-2 in Pakistan has been 2% with some regional variations [13]. Morbidity due to COVID-19 has been relatively lower as compared with many other countries [14,15]. This could be due to both pathogen and host-related factors.
Study of genomic variation of SARS-CoV-2 strains helps understand the epidemiology of COVID-19. Global SARS-CoV-2 genomic sequencing efforts have contributed data of thousands (2,503,415) of strains into public database such as, GISAID and Nextstrain, NCBI SARS-CoV-2 Resources (https://www.ncbi.nlm.nih.gov/sars-cov-2/). These have allowed the interrogation of viral diversity with associated disease transmission in different countries [16].
There is limited genomic epidemiological data available for SARS-CoV-2 strains in Pakistan. The first introduction of SARS-CoV-2 was made by a traveller from Iran in March 2020. The pandemic was initially associated with travelers but local transmission was identified within two weeks of the first known COVID-19. G and S clade strains were shown to be present in the first wave of COVID-19 in Pakistan [17]. Genetic studies of SARS-CoV-2 from the second wave in Pakistan have identified B.1 and its sub-clades together with B.6 lineage strains to be present [18].
Our clinical laboratory at the Aga Khan University Hospital, Karachi, Pakistan tested 196,588 respiratory samples for SARS-CoV-2 by PCR and reported 40,709 (21%) PCR positive COVID-19 cases between March and December 2020. Here we investigated the genomic diversity of SARS-CoV-2 studying the phylogeny of isolates from the first COVID-19 wave and compared them with later strains. Further, we studied diversity, site-specific mutation and entropy across the genome in strains isolated before and after June 2020, representing earlier and later periods.

Sample selection
This study was approved by the Ethical Review Committee at the Aga Khan University (AKU), Karachi, Pakistan as a waiver of informed consent for the study. Samples used were those archived as de-identified samples previously reported positive for SARS-CoV-2 by the Clinical Laboratories, Aga Khan University Hospital (AKUH), Karachi, Pakistan. Laboratory data (including age, gender) was utilized where available.
Nasopharyngeal swab specimens were confirmed positive for SARS-CoV-2 by reverse transcription (RT) polymerase chain reaction (PCR) using the SARS-CoV-2 Cobas 6800 Roche assay at the Section of Molecular Pathology, AKUH, Karachi, Pakistan. A random, convenience sampling was done to include samples from March until October 2020.
Inclusion criteria were; samples amplified at CT value of 30 and below were selected, specimens collected between March and October 2020. Exclusion criteria were: samples with a crossing threshold (CT) value greater than 30.

Sequencing of SARS-CoV-2 strains
A total of 70 samples were selected for SARS-CoV-2 genome sequencing. RNA was extracted from selected specimens using the QiaAmp RNA minikit (Qiagen, USA) and used as input for NGS library preparation.
In total, we successfully sequenced thirty-two SARS-CoV-2 specimens (twenty-one whole genome sequences and eleven partial genome sequences). Whole genomes of eight isolates were available through sequencing using the Nextera XT DNA Library Preparation kit (Illumina) was used for. As described previously, first-strand cDNA was synthesized with Super-Script III Reverse Transcriptase (SSIII), Thermo Fisher Scientific, USA, followed by Second strand cDNA synthesized with DNA polymerase I, Large Fragment, Klenow (Invitrogen, USA), (11). We obtained sequence data for twenty-four SARS-CoV-2 isolates using the Tru-Seq1 Stranded Total RNA Library Preparation kit (Illumina). All normalized libraries were pooled and spiked with PhiX control prior to sequencing was performed on the Illumina Miniseq platform using a 300 cycle Miniseq Reagent Kit v2 (Illumina). In total we obtained twenty-one full and eleven partial genome sequences.  Table).

Variant calling and phylogenetic analysis
We performed whole genome analysis on thirty-two sequences from this study in addition to fifty-eight SARS-CoV-2 genomes from Pakistan available from GenBank (S1 Table). These comprised 79 full-length and 11 partial genomes. Fifty-four full length genomes were those isolated in the period March until June, 2020. Twenty-five full length genomes were from between July and October, 2020.
Phylogenetic analysis was performed using the 79 available full length SARS-CoV-2 genomes (S1 Table) along with the 449 full-length SARS-CoV-2 reference sequences (S1 File) from different pandemic countries obtained from the NCBI SARS-CoV-2 Resources (https:// www.ncbi.nlm.nih.gov/sars-cov-2/) were subjected to Multiple Sequence Alignment (MSA) along using MAFTT online server [21]. The MSA was subsequently used to generate a Maximum Likelihood (ML) phylogenetic tree using PhyML 3.0 (http://www.atgc-montpellier.fr/ phyml/) with a GTR-based nucleotide substitution model and aLRT SH-Like branch support. The root of the tree and branch length variance was determined using the TreeRate tool [22] by applying a generalized midpoint rooting strategy. The tree was visualized and edited in Figtree software (http://tree.bio.ed.ac.uk/software/figtree/). Mean and individual pairwise distance between SARS-CoV-2 sequences from our study and previously deposited Pakistani SARS-CoV-2 sequences was calculated using MEGA 7 [23].

Entropy and site selection (dN-dS) analysis
Genome-wide Shannon entropy and site selection (dN-dS) analyses were performed to evaluate genomic variability between SARS-CoV-2 genomes from Pakistan isolated up until June 2020 (representing strains from the first COVID-19 wave) and those isolated later. Further, we also performed dN-dS entropy analysis for each clade of strains separately. However, as both entropy and site selection analysis require at least three taxa to give a meaningful output, these analyses could be performed for L, GH, and GR clade strains only and not the V and S clade strains in our study. Shannon entropy for each clade was carried out using Bioedit entropy function. The statistical significance in mean entropy value between 1 st and 2 nd wave sequences was evaluated using the paired T-test, while statistical significance in mean entropy value between different clades was evaluated using the One-way Anova test. Both tests were performed using 95% confidence interval and p<0.05 as significant value. The calculations were performed using GraphPad Prism tool.
For site selection analysis, codon alignment was performed using MEGA 7, using Muscle algorithm, and the codon aligned file was subsequently used for site selection using SNAP tool available at the Los Alamos Database (www.hiv.lanl.gov) and SLAC tool available at DataMonkey [25]. It is important to note that SLAC was pre-cited as the best model for our data set based on algorithm selection criteria available on DataMonkey website. Normalized dN-dS data was plotted and analyzed for the presence of significant positively or negatively selected sites.

Description of SARS-CoV-2 isolates
We determined the genomic epidemiology of ninety SARS-CoV-2 strains isolated between March and October 2020. Overall, the strains were found to be belong to clades L (n = 13, 14.5%), S (n = 5, 5.5%), V (n = 2, 2.2%), O (n = 1, 1.1%), G (n = 1, 1.1%) clade with sub-types GH (n = 59, 65%) and GR (n = 9, 10%), S1 Table. A month-wise distribution of the strains showed eleven isolates sequences from March, nineteen from April until May, forty-two from June until July and eighteen between August and October 2020. G and L clade strains were present throughout the period studied (S1 Fig). Strains belonging to the GH clade became predominant from April onwards until October. The frequencies of L and S clades reduced over the study period. GR and V clades were found in March and then in June-July. Notably, the diversity of clades was reduced over the study period in that by October only L, GR and GH (as the predominant clade) were present.
Four SARS-CoV-2 GH clades strains identified in March were from travelers from Iran and Turkey. Apart from these, the remaining eight-six strains were from across the four different provinces of the country without any known travel history. These cases of likely local transmission were from the province of Sindh (n = 64: clades G, n = 1; GH, n = 46; GR, n = 1, L n = 10; S n = 3; V n = 2), Punjab (n = 5; clade GR, n = 5), ICT (n = 17; S n = 2; GH n = 12; GR n = 3) and KPK (n = 4, GH n = 1; L n = 3).

SARS-CoV-2 lineage analysis
Phylogenetic analysis. Phylogenetic analysis of 79 full length SARS-CoV-2 genomes from Pakistan and 449 global isolates was conducted (Fig 1). Of the 21 genomes sequenced in this study (AKU), seven clustered with sequences from Saudi Arabia and India (Fig 1, orange and purple branches, respectively). Four clustered with sequences from US (Fig 1, pink branches), and ten (AKU-2, -3, -21, -24, -25, -26, -33, -46, -47, -56) clustered with previously deposited sequence from Pakistan (Fig 1, green color branches). These ten AKU strains were all from Karachi (S1 Table), suggesting that the viral strains circulating in the city were predominantly similar to those circulating in other parts of Pakistan. The mean pairwise genetic distance between our sequences; sequences previously deposited from Pakistan, and sequences from India, Saudi Arabia, USA, and Australia was found to be 0.00, indicating phylogenetic relatedness between the genomes.
Nextstrain analysis. Expanded phylogenetic analysis to examine the genetic divergence of strains was conducted against a representative subset of 825 Asian SARS-CoV-2 genomes present in Nextstrain database using data for the period January 2020 to January 2021. Pakistani SARS-CoV-2 genomes aligned throughout the phylogenetic tree indicating multiple introductions of the SARS-CoV-2 in the country (Fig 2). 71 sequences from Pakistan are displayed in the global phylogenetic tree of SARS-CoV-2, primarily classified on the basis of clade (19A, 19B, 20A, 20B, 20C and 20D). Out of 71 sequences from Pakistani population, 48 sequences were clade 20A, seven were clade 20B, one was clade 20C and seven belonged to clade 20D. Further, there are 10 and 3 Pakistani sequences, of clades 19A and 19B respectively (S1 Table).  Table) along with 58 other Pakistani and 449 full-length SARS-CoV-2 reference sequences (S1 File). AKU study sequences are indicated in blue, while other Pakistan sequences are shown in green. AKU sequences clustered with those from India (purple), Saudi Arabia (orange), United States (pink), China (light blue), Australia (turquoise), Bangladesh (yellow), France (light green) and also with other sequences from Pakistan (green). The root of the tree was determined using TreeRate tool by applying generalized midpoint rooting strategy. Nodes with significant (>0.90) aLRT-SH like support values are colored maroon. The tree was visualized and edited in Figtree software.

Estimation of mutation rates
We investigated the stability of the SARS-CoV-2 genomes by determining their mutation rates. Pairwise time tree distances and substitution rates were estimated 5.68 x 10 −4 substitution per site per year (16.98 substitution per year) for the Pakistan tree compared to 8 × 10 −4 substitution per site per year (23.92 substitution per year) globally, (Fig 3). This demonstrates the substantial variation in SARS-CoV-2 phylogenies between March and October, 2020.
Seventy-nine different non-synonymous SNVs were observed in the Orf1ab gene encoding the non-structural proteins with the greatest number observed in the nsp3 (n = 21), nsp2  (n = 16) and RdRp (n = 12) regions, respectively. Within structural proteins, 59 different nonsynonymous SNVs were observed, with greatest numbers seen in N gene (n = 21), S gene (n = 15) and ORF3a (n = 7). In N gene, the most frequent ns-SNV were R209I and S194L where in the S-gene the most frequent ns-SNV was D614G (S2 Table).

Entropy and site selection pressure analysis
In the next step, we investigated the correlation between genetic diversity, entropy and site selection pressures across SARS-CoV-2 genomes. Entropy gives an estimate of the probability of acquiring mutations in a given set of genomic sequences [26].
The first COVID-19 wave in Pakistan peaked in June 2020, hence, we compared SARS--CoV-2 strains within our study period as those isolated before and after June 2020. There were fifty-seven genomes from March until the thirtieth of June. Thirty-three genomes were isolated from between first of July and until the end of October.
We first investigated the difference in genome-wide entropy between Pakistani SARS-CoV-2 sequences from before and after June, representing samples from the first COVID-19 wave in Pakistan and later isolates. Our results revealed that sequences collected before June 2020 exhibited an overall higher mean (Fig 5A, p<0.001) and site-specific entropy (Fig 5B) across genome as compared with sequences collected after June, 2020. Indicating, greater overall stability in the genomic pressures of SARS-CoV-2 strains in the second wave. Analysis of entropy prevalent clades L, GR and GH revealed clade GH to have an overall higher mean entropy (p<0.001) across the genome as compared with GR and L clade genomes (Fig 6A). In the GH clade genomes we observed multiple sites exhibiting entropy above 0.5 and up to 0.9 (Fig 6B). This was followed by second highest entropy across genomes of the L clade (p<0.001), while GR clade (p<0.05) genomes had fewer high-entropy sites and an overall lower pan-genomic entropy as compared to clade GH genomes.
A subsequent site-selection analysis revealed several genomic locations under positive or negative selection in sequences before and after June 2020 (Fig 7A and 7B), however, none of the sites were under statistically significant selection pressure. In agreement with the entropy analysis, clade GH exhibited higher number of sites under selection pressures, followed by clade L and GR (S3-S5 Figs). However, none of the sites in any clade were found to be under statistically significant selection pressure.

Discussion
Our study provides insights into the SARS-CoV-2 strains circulating between the first wave of COVID-19 in Pakistan and those after this period. It associates genetic diversity of the mostly commonly found GH clade strains with greater relative entropy without a site-specific selection bias.
SARS-CoV-2 genomic surveillance data available at GISAID provide key information regarding the viral diversity of SARS-CoV-2 genomes across the globe. New mutations are continually being identified but the biological effects of most of them remain unclear. Recently reported spectrum of mutations mainly spike 69/70 deletion, E484K and N501Y variants in B.1 lineage strains, have revealed a role for mutations in driving virulence of the virus by impacting host infection, transmission, diagnostics, and vaccine escape [27,28].
Coronavirus replication is error-prone, with a high mutation rate and a nucleotide mutation rate estimated at 4x10 -4 substations/site/year [29]. For SARS-CoV-2, the substitution rate is considerably higher i.e., 9.90 × 10−4 substitutions/site/year (6.29 ×10 −4 to 1.35 × 10 −3 ) based on one study, whilst it is thought to be 5.3504 × 10 −3 and 5.35 × 10 −3 based on other studies [30]. We found the SARS-CoV-2 genome substitution rate to be 5.68 x 10-4, which lies within the range previously calculated. The variability in substitution rates in different studies might be due to estimation model bias or rapid evolution of virus [30].
Pakistan SARS-CoV-2 genomes were genetically diverse and clustered with those from Saudi Arabia, India, USA, Australia, Italy and China. Comparison with a sub-sample of Asian isolates revealed that there was clustering with virus genomes from Japan, which might indicate that strains from these countries might be genetically similar and evolving at a similar rate. Ten sequences from this study clustered with previously deposited sequenced from Pakistan suggesting limited diversity between strains circulating in Pakistan. However, this data may also be impacted by the limited genomic information on SARS-CoV-2 strains currently available for the country.
We show that in Karachi there was a shift from S, L, V and G clades in March to predominantly the GH clade strains from May onwards [17]. Traveler-associated strains belonged to S and GH clades. GH (B.1.255) and S (A) clade strains were introduced through travelers from Iran and Turkey, corroborating with previous reports [31]. An L lineage isolate (B.6) was identified from a religious pilgrim who attended a super spreader event in Punjab Province which had visitors from abroad including China. L lineage isolates exhibiting L3606F mutations were first reported on January 2020 in China and associated with super spreader events in the USA, Singapore, Japan, and Europe [3]. GH (B.1., B.1.36, B.1.160, B.1.255, B.1.275) and L (B, B.6) strains associated with local transmission persisted between March and October. Data from Karachi concur with previous reports where B.1 sub-lineages have been reported [32]. Variant analysis revealed Orf1ab, nsp3 and N region to have the highest number of mutations, as shown previously [5,33]. N gene P13L was present in both S clade and A lineage isolates in this study, it has been previously reported from the UK and Australia [34]. N gene variants at codon 202, 203, 204 were in 9% of the isolates whereas and 209 were found in 20% of isolates further reinforcing the increased mutability observed in this protein; these variants are associated with a split of G into GH and GR clades.
The most common non-synonymous variants observed were S D614G, Orf1ab RdRp gene P4715L and, Orf3a Q57H which have previously been reported from Europe and the United States [33]. Five isolates with the ORF8 L84S mutation also had N gene S202N. The L84S mutation is common amongst S clade isolates found in Europe where it was found to co-evolve with mutations such as P323L [34]. The L lineage initially split equally into G and V subclades, with G reaching 50% of the viruses isolated in March 2020 and then splitting further into GR and GH subclades [35]. We observed a GH isolate from May to have Orf1ab L3606F mutation co-occurring with P323L. L3606F is not typical of GH the clade, and is a transitory mutation defining L to S clade divergence, possibly identifying further evolutionary transitions in this isolate [3].
The D614G mutation in S gene was present in the majority (76%) of the genomes studied here and reported globally [17]. The Spike glycoprotein assists viral entry into host cells by binding to the ACE2 receptor [36,37]. The substitution of Glycine at the 614 mutations may have introduced structural instability into the spike protein [38]. It has been associated with increased virulence and transmission of SARS-CoV2, most probably due to a higher viral load [39,40]. G clade isolates comprise the predominant proportion of isolates in Europe and North America. The D614G mutation has also been associated with greater mortality observed in Belgium, Spain, Italy, France, Netherlands and Switzerland [41]. Importantly, eight different singly occurring S gene ns-SNV were observed in the study isolates. One strain had two S gene mutations, D614G and L5F. This double mutation has been shown to increase the infectivity of human cell lines as compared to the reference Wuhan strain [6].
The RdRp variant P4715L was present in 72% of all isolates. P4715L and D614G containing isolates have shown significant positive correlations with fatality rates in many countries [42]. RdRp is the target for polymerase inhibitors and mutations in RdRp can potentially decrease drug-RdRp complex binding affinity leading to resistance and a differential effect of antiviral treatments [43]. The frequency of P4715L has been shown to differ between SARS-CoV-2 genomes isolated from the USA (63.0%) and China (11.2%) and have been attributed to the differential efficacy observed in the clinical trials of the antiviral agent Remdesivir [44].
Variants in the non-coding 5'UTR and 3' UTR regions of the SARS CoV2 virus have been reported to affect viral replication and transcription. The 241C > T mutation has been reported to result in the strong binding of TARDBP (RNA/DNA-binding protein) to the 5'UTR region of the SARS CoV2 virus [45]. This variant has been implicated in facilitating the translation of viral proteins resulting in its effective propagation within the human host. Interestingly, variant (241C > T) of the 5'UTR region often coexists with spike glycoprotein variant (S protein, D614G) [46]. This coexistence is also evident in our study as we found 4 GH clade strains and a V clade strain to have the 5' UTR variant +241 C>T. Two S clade strains had the 3' UTR + 29742 G>A. This mutation is thought to affect binding of the miR-1307 which regulates the stability of the RNA at a post-transcriptional level and therefore, it potentially weakens the host immune response against the virus [47].
SARS-CoV-2 genomic surveillance data into GISAID provides key information regarding the viral diversity of genomes across the globe. New mutations are continually being identified but the biological effects of most of them remain unclear. Recently reported spectrum of mutations mainly spike 69/70 deletion, E484K and N501Y variants in B.1 lineage strains, have revealed a role for mutations in driving virulence of the virus by impacting host infection, transmission, diagnostics and vaccine escape [27,28].
Coronavirus replication is error-prone, with a high mutation rate estimated at 4x10 -4 substations/site/year [29]. For SARS-CoV-2, the substitution rate is considerably higher i.e., 9.90 × 10 −4 substitutions/site/year (6.29 ×10 −4 to 1.35 × 10 −3 ) based on one study, whilst it is thought to be 5.3504 × 10 −3 and 5.35 × 10 −3 based on other studies [30]. We found the SARS-CoV-2 genome substitution rate to be 5.86 x 10 −4 , which lies within the range previously calculated. Such a high substitution rate can have a significant impact on the evolution of viral clades and emergence of new SARS-CoV-2 variants. The variability in substitution rates in different studies might be due to estimation model bias or rapid evolution of virus [30].
We investigated the correlation between viral diversity and sequence variability using entropy and site selection analysis. Genomes from clade GH were found to have an overall higher pan-genomic entropy as well as higher evidence of pervasive selection as compared with clades GR and L. Increase site specific pressure and entropy, especially across Orf1ab and N regions in GH clade is consistent with the mutation data, and might suggest that the high selection pressure of GH clade is attributed to high transmission opportunities and/or dissemination into genetically diverse population observed in the transient phase of viral evolution before the virus adapts to the population/host selection pressures [48].
A limitation of this study is the relatively small sample size. However, as the Aga Khan University Hospital initiated SARS-CoV-2 PCR based diagnostics in Karachi, Sindh at the start of the pandemic in February 2020, a large proportion of test samples from Sindh province were referred to our laboratory. Therefore, the genomes described here are considered representative of the SARS-CoV-2 strains circulating in the province during and after the first COVID-19 wave.

Conclusions
Our results highlight the value of genomic surveillance in understanding the evolution of SARS-CoV-2 strains in Pakistan. Our data shows that the predominant GH clade strains were more genetically diverse but did not display site specific pressures. Importantly, there was a higher mean and site-specific entropy in SARS-CoV-2 genomes isolated before as compared with after June 2020. This suggests that the viral genome achieved stability after the initial early period of the pandemic. It will be important to continue genomic epidemiological studies surveillance of SARS-CoV-2 to understand transmission patterns in the country.