Genomic evidence of multiple SARS-CoV-2 introductions into Morocco

Data set

The Morocco related-sequences used in this study have been deposited in the GISAID database by the Laboratory of the Royal Gendarmerie of Morocco (LRAM) (one sequence), the Pasteur Institute of Morocco (IPM) (17 sequences), the Anoual Laboratory in Casablanca (one sequence), and the SeqCovid Spanish Project (two sequences). The sequences were collected from the GISAID database and analyzed in comparison with 500 selected genomes, that were collected between the 15^th of February 2020 (2 weeks before the first detected case in Morocco) and the 30^th of March 2020 (2 weeks after the international borders of Morocco were closed; lockdown). Only high-quality sequences have been included in this study. All sequences and metadata are provided as Extended data⁵ and can be found at the GitHub repository https://github.com/covidmor/covid19_morocco.

Variant calling and network analysis

The Moroccan sequences were mapped to the Wuhan reference genome (NC_045512) using Bowtie 2.3.5.1⁶. Variants were called using mpileup 1.7 and bcftools 1.9⁷. The variants were annotated using the China National Center for Bioinformation Annotator. Variant network analysis was performed using Gephi 0.9.2.

Phylogenetic analysis

A total of 500 sequences from different countries were retrieved from the GISAID database (www.gisaid.org) and aligned using Muscle 3.8.1551 multiple sequence alignment⁸; see Extended data, metadata⁵ for details on each. A maximum likelihood model was created using RaxML 8.2.12⁹ with a bootstrap of 100 using the Wuhan reference genome (NC_045512) as outgroup. The phylogenetic tree was generated using Figtree 1.4.4¹⁰.

Phylogenetic tree analysis

To put the complete genomes of SARS-CoV-2 that were isolated from 19 Moroccan patients, a patient from Melilla and a Moroccan patient from Cadiz (Andalusia, Spain), into the context of the global pandemic, they were aligned together with the dataset of 500 SARS-CoV-2 complete genomes from different countries. Extended data, Figure 1⁵ shows the estimated maximum likelihood phylogeny. To have a good presentation of the phylogeny, we created a circular presentation subset of 120 sequences from our dataset (collected from 18^th of February 2020 (15 days before the first Moroccan case) to 15^th of March 2020 (Date of closing international borders in Morocco)) (Figure 1b).

Figure 1. Phylogenetic tree analysis.

(a) Maximum likelihood phylogeny, constructed using complete SARS-CoV-2 genomes from the GISAID database. The blue sequences represent the Moroccan isolates (19), a strain isolated from Melilla and a strain isolated from a Moroccan patient in Cadiz, Andalusia. The phylogenetic tree is divided into four clades: the red part for the GR strains, the green one for the GH strains, the yellow one for the G strains, and the pink one for the S strain. This Figure represents selected sub-tree extracts from the extended phylogenetic tree in Extended data, Figure 1⁵. A1–3: represent clusters of Moroccan strains originated from Italy and belonging to the GR clade; B: represents a cluster of Moroccan strains from a French origin and belonging to the GH clade; C: represents the Moroccan strain isolated from Ouarzazate, it originated from Portugal and belongs to the G clade; D: represents a cluster of Moroccan strains originated from France and belonging to the G clade; E: represents a Moroccan strain that originated from Italy and it belongs to the G clade; F: represents a Moroccan strain from a Spanish origin and belonging to the G clade, G: represents the sequence isolated from Melilla belonging to the G clade with a Spanish origin, H: represents the sequence isolated from Cadiz; it belongs to the S clade and it’s close to strains from Spain. (b) Circular representation of the phylogenetic tree of a subset of 120 sequences collected before closing Moroccan international borders, constructed using complete SARS-CoV-2 genomes from the GISAID database. The blue sequences represent the Moroccan isolates (19). The phylogenetic tree is divided into three clades: the red part for the GR strains, the green one for the GH strains, and the yellow one for the G strains. The clade S containing the Cadiz sequence is absent as we removed this sequence from this phylogenetic analysis.

The extended phylogenetic tree is subdivided into seven clades that correspond to the main seven SARS-CoV-2 strain types GR (red), GH (green), G (yellow), S (pink), V (magenta), L (light gray) and O (others; dark gray). The four strain types S, G, GR, and GH are scattered all over the world (Extended data, Figure 2)⁵. According to GISAID (update of the 12^th of June 2020), they made their first appearance in February 2020, being G the strain type that gave origin to the GH and GR ones, while the S strain was mainly found in Spain¹¹.

Figure 2. Variant network analysis of the 19 Moroccan genome sequences, one sequence from Melilla and another one from a Moroccan host in Cadiz, Andalusia.

Every dot represents a sequence, and each edge on the network connects, using a shared variant, a sequence to one or more other sequences. The red circle group Moroccan sequences originated from Italy, the green one group sequences originated from France, the purple one assemble sequences originated from Spain, and the blue one represents the sequence originated from Portugal; we used different color gradients to separate sequences that belong to different clades in the same cluster. A: Strain isolated from Cadiz, it belongs to the S clade and it shares no variants with Moroccan strains, including the Spanish ones. M: represents the strain isolated from Melilla, it belongs to the G clade in a cluster of Spanish sequences. IPM3: is a Moroccan strain isolated by Pasteur Institut of Morocco (IPM), it belongs to the G clade and it’s coming from Spain. O: represents the strain isolated by LRAM from Ouarzazate, it's the only sequence with a Portuguese origin and it belongs to the G clade. IPM9: represents a sequence collected by IPM, it originated from Italy and belongs to G clade. C; IPM1; IPM2; IPM7; IPM10; IPM13; IPM15-17: represent a strain isolated by Anoual laboratory in Casablanca (C) and sequences collected by IPM, they belong to the GR clade and they all are from Italian origin. IPM4; IPM11-12: represent strains isolated by IPM, they originated from France and they belong to the GH clade. IPM5-6; IPM8; IPM14: represent strains isolated by IPM they are from a French origin and they belong to the G clade.

Zooming into the trees (Extended data, Figure 1)⁵, Figure 1a shows how the Moroccan isolates group into four independent clusters. In the GR clade, the sequences GISAID IDs: EPI_ISL_459973, EPI_ISL_459966 and EPI_ISL_459965 (the first infected patient in Morocco) cluster together with sequences from Italy (GISAID ID: EPI_ISL_417922), Thailand (GISAID ID: EPI_ISL_430837), Mexico (GISAID ID: EPI_ISL_452139), and Iceland (GISAID ID: EPI_ISL_417829) (Figure 1a: A1–3). Interestingly, the hosts of these four sequences were all travelling in Italy during the same period (Extended data, Metadata)⁵. This suggests that the first Moroccan patient (GISAID: EPI_ISP_459965), that returned to Morocco from recent travel to Bergamo, Lombardy (the most affected region in Italy), and the three other individuals from Mexico, Iceland and Thailand (GISAID: EPI_ISL_452139, EPI_ISL_417829 and EPI_ISL_430837, respectively) were probably infected by the same person in Italy (GISAID: EPI_ISL_417922). From this cluster, we infer, based on the sample collection dates (Extended data, Metadata)⁵, that the first Moroccan patient (GISAID: EPI_ISP_459965) probably infected the two other Moroccan patients whose infecting virus sequences are found in the same cluster (GISAID: EPI_ISL_459973 and EPI_ISL_459966). The Moroccan sequences (GISAID IDs: EPI_ISL_458150, EPI_ISL_459980 and EPI_ISL_459984) that cluster together in the GR clade, are close to sequences from Iceland (GISAID IDs: EPI_ISL_417861, EPI_ISL_417863 and EPI_ISL_417852) whose carriers have travel history to Italy. Furthermore, all these three Icelandic sequences are related to one Italian sequence (GISAID ID: EPI_ISL_460081). For their part, the Moroccan sequences EPI_ISL_459983, EPI_ISL_459977 and EPI_ISL_459978 cluster with another sequence from Italy (GISAID ID: EPI_ISL_452186). Surprisingly, many sequences from the GR clade were isolated from patients with known travel history to Italy (Figure 1a: A1–3). Therefore, we can conclude that all the viral sequences in this clade, where 9 Moroccan viral sequences clustered, were of Italian origin.

The four Moroccan sequences that cluster together within the G type strain (GISAID IDs: EPI_ISL_459968, EPI_ISL_459974, EPI_ISL_459972 and EPI_ISL_459981) were close to French sequences (GISAID ID: EPI_ISL_417333, EPI_ISL_443270 and EPI_ISL_416752). These four Moroccan sequences were collected the 17^th of March, the 20^th of March, the 20^th of March and the 19^th of April 2020, respectively. Interestingly, the third Covid-19 case in Morocco was a French male tourist, in his 30s that arrived in Morocco on the 7^th of March 2020 as reported by the Moroccan Health authorities. He, therefore, was most likely the source of the Moroccan viral strains that belong to the G type. The case reported on the 19^th of April 2020 (sequence with the GISAID ID: EPI_ISL_416752) seems to be the result of a local community transmission.

The Moroccan strain isolated from Ouarzazate (South of Morocco) (GISAID ID: EPI_ISL_451400) is close to the Portuguese sequences (GISAID ID: EPI_ISL_413648 and EPI_ISL_453947) in the G type suggesting a third introduction path to Morocco from Portugal. For their part, the Moroccan sequence (GISAID ID: EPI_ISL_459967), that also belongs to the G type, is close to the Spanish sequences (GISAID IDs: EPI_ISL_455349 and EPI_ISL_455332) making Spain another point of introduction to Morocco. Furthermore, the Moroccan sequence (GISAID ID: EPI_ISL_459975) is closely related to the Italian G strain (GISAID: EPI_ISL_417921), which represents another Italian introduction to Morocco. It’s also similar to a sequence from Iceland (GISAID ID: EPI_ISL_417730), where the host patient has travel history to Italy (Figure 1a: E) and could have been either infected in Italy or by another asymptomatic carrier coming from there.

Three other viral sequences isolated from Moroccan patients (GISAID IDs: EPI_ISL_459976, EPI_ISL_459982 and EPI_ISL_459979) appear in the clade of the GH strain types where they cluster with French sequences (GISAID IDs: EPI_ISL_416748 and EPI_ISL_416501) (Figure 1a: B) suggesting a further introduction to Morocco, once more from France.

The viral sequence isolated from a patient in Melilla (GISAID ID: EPI_ISL_455344) is clustered within the G type, with Spanish sequences (GISAID IDs: EPI_ISL_419235, EPI_ISL_450337 and EPI_ISL_419236). While the sequence isolated from a Moroccan patient in Cadiz, Andalusia (GISAID ID: EPI_ISL_452463) (Figure 1a: G) that belongs to the S type, also clusters with Spanish sequences; however, this sequence is quite far from the other Spanish related sequences. Thus, supporting the apparent neutrality of the host’s genetics for the transmission and potential selection on the type of strain.

Finally, similarity analysis of the Moroccan sequence (GISAID ID: EPI_ISL_467299) collected on May 21^st, 2020 and sequenced by the LRAM Laboratory, is close to the first identified French cluster, suggesting that the virus strain circulating in Morocco did not experience any major mutations and, once again, confirms the local community-based transmission of the virus. This sequence was not included in the overall analysis, as it was deposited when the writing of this work was finalizing.

Thus, the phylogenetic analysis shows that the virus outbreak in Morocco was likely the result of multiple introductions. We can highlight four independent introductions to Morocco from Italy, France, Spain, and Portugal. This finding has obvious implications for the epidemiological tracing of the pathogenic agent’s initial introduction that caused the current outbreak in Morocco. Generally, the viral sequences isolated from Moroccan patients showed close relationships primarily to European strains. The geographic nearness and tourist and migratory connections, therefore, play a key role in the spread of the virus. We caution that further analyses are needed to evaluate the statistical robustness of the interference suggested herein.

Network analysis

To further understand the evolution of the SARS-CoV-2 virus within the Moroccan population and trace the infection pathways, we performed a Variant Network Analysis using Gephi. We used the complete sequences of the 21 SARS-CoV-2 genomes, 19 from Morocco, and two sequences from Melilla and Andalucia (Figure 2). The variant network in Figure 2 shows three main clusters, one isolated sequence from Portugal and another one from Cadiz sharing no variants with the other sequences. The first cluster comprising the following viral strains, all closely related to each other and diagnosed in Casablanca (west Morocco) (GISAID ID: EPI_ISL_458150), IPM_1 (GISAID ID: EPI_ISL_459965), IPM_2 (GISAID ID: EPI_ISL_459966), IPM_13 (GISAID ID: EPI_ISL_459980), IPM_15 (GISAID ID: EPI_ISL_459978), IPM_7 (GISAID ID: EPI_ISL_459973), IPM_16 (GISAID ID: EPI_ISL_459977), IPM_10 (GISAID ID: EPI_ISL_459984) and IPM_17 (GISAID ID: EPI_ISL_459983) were closely related to each other. This is in agreement with the phylogenetic tree where all these sequences belong to the GR clade and appear to have been originated from Italy (Figure 1 and Extended data, Figure 1)⁵. The second cluster (Blue) made out of the sequences labeled IPM_5 (GISAID ID: EPI_ISL_459974), IPM_6 (GISAID ID: EPI_ISL_459968), IPM_8 (GISAID ID: EPI_ISL_459972) and IPM_14 (GISAID ID: EPI_ISL_459981) belong to the clade G and has a French sequence (GISAID ID: EPI_ISL_418222) as a root in the phylogenetic tree. The sequences IPM_4 (GISAID ID: EPI_ISL_459976), IPM_11 (GISAID ID: EPI_ISL_459982) and IPM_12 (GISAID ID: EPI_ISL_459979) belong to the GH clade and seem to have originated from another French sequence (GISAID: EPS_ISL_418219). Interestingly, this sequence belongs to the G clade and is close to a strain that was isolated from an Icelandic patient that was most probably infected during his stay in Italy (GISAID ID: EPI_ISL_417730). Another separated Moroccan sequence, labelled IPM_3 (GISAID ID: EPI_ISL_459967), is closer to Spanish sequences (GISAID ID: EPI_ISL_455349, EPI_ISL_455332). The sequence of the viral strain from Ouarzazate (GISAID ID: EPI_ISL_451400) is closely related to a Portuguese sequence (GISAID ID: EPI_ISL_413648), just as shown in the phylogenetic tree. Yet, all these sequences belong to the G clade, which has an Italian sequence as a root. Hence, the Variant Network Analysis is in accordance with the phylogenetic tree and further supports our deduction of multiple SARS-CoV-2 introductions to Morocco from at least four European countries. Thus, while geography marks the variation that the virus undergoes in a way that we can identify geographically separated strains, population interconnection, via travel and migratory flows, seems more important than the geographical nearness for the worldwide spread of the virus.

To further test the multiple introductions to Morocco from the aforementioned four European countries, we performed a Variant Network analysis using all the Moroccan sequences as one group and compared them with randomly selected sequences from different countries (Italy, France, Spain and Portugal, in addition to Austria, Australia, and Brazil). These results show that the sequences detected during the early stages of the epidemic in Moroccan share variants with Italian (seven variants), French (four), Portuguese (three), and Spanish (three) strains. They also share one variant with strains from Austria, one with Brazil and another one with Australia (Figure 3). Interestingly, there was one shared variant between all countries located in the Spike gene (23403A>G) which induces the amino acid change D614G; it began spreading in Europe in early February and, when introduced to new regions, it rapidly became the dominant form. This mutation is suspected to increase the transmissibility of the virus¹².

Figure 3. Variant network analysis elucidating the shared variants between sequences from different countries.

Each dot represents a group of 19 sequences randomly selected from each country. The Moroccan node represents all 19 sequences isolated in Morocco. Morocco shares more variants with Italy (seven), France (four), Portugal (three) and Spain (three) than with Brazil (one), Austria (one) and Australia (one).

Variant calling

The genome-wide SNPs and the corresponding amino-acid positions and variations of the virus proteins are described in Table 1. Our results showed that all 20 sequences including the sequence coming from the Melilla patient, shared four mutations in common and 13 novel mutations exclusively present in the 19 sequences isolated in Morocco (Table 2). The four common mutations happened to be the same recurrent mutations described in several reports¹³. The mutation in the leader sequence (241C>T) is one of the most common mutations in the SARS-CoV-2 genome; it affects an important genomic site for discontinuous sub-genomic replication. The mutation in the 5’UTR genomic region co-evolved with two other mutations, 14408C>T and 23403A>G. They both affect critical RNA replication proteins (241C>T, 14408C>T) and the ACE2 receptor binding protein S (23403A>G). We noticed that these four mutations are prevalent in the virus isolates from Europe, where the infections seem to be more severe. In fact, Yin et al.¹⁴ suggested that these mutations could increase the transmissibility of the virus. The sharing of these mutations between European and Moroccan isolates further indicates that strains circulating in Morocco came mainly from Europe. Moreover, at least one variant of the Moroccan sequences of SARS-CoV-2 is shared with one of the other countries (Table 3). The mutations found in the sequences identified in Morocco, and the identification of amino acid change, should be further investigated in order to understand whether they affect the transmission and clinical characteristics of the SARS-CoV-2 virus circulating in Morocco.

Source data

Sequences analyzed in this study were downloaded from the GISAID database. The identity of each sequence is given in the Extended data, Metadata file⁵.

Extended data

Zenodo: Genomic evidence of multiple SARS-CoV-2 introductions into Morocco. http://doi.org/10.5281/zenodo.3909463⁵.

This project contains the following extended data:

All supplementary material and figures can be found at GitHub: https://github.com/covidmor/covid19_morocco.

Extended data are available under the terms of the Creative Commons Zero "No rights reserved" data waiver (CC0 1.0 Public domain dedication).