Les estimations à l'échelle du génome de la recombinaison, de la mutation et de la sélection positive inspirent divers moteurs de Mycobacterium bovis

Merci de votre visite sur Nature. La version de votre navigateur ne prend pas en charge les styles CSS. Pour une expérience optimale, nous vous recommandons d'utiliser une version plus récente (ou de désactiver le mode de compatibilité dans Internet Explorer). Par ailleurs, afin de garantir la continuité du service, nous afficherons les sites sans styles ni JavaScript.
Le séquençage du génome a revitalisé la recherche sur les maladies infectieuses, révélant l'épidémiologie, la pathogenèse, les interactions hôte-pathogène et les processus évolutifs imposés aux pathogènes. Le complexe Mycobacterium tuberculosis (MTBC) comprend Mycobacterium bovis, un de ses membres adaptés aux animaux, responsable de la tuberculose (TB) chez les mammifères terrestres et modèle typique d'évolution bactérienne. À l'instar des autres membres du MTBC, Mycobacterium bovis est considéré comme un pathogène strictement cloné, à évolution lente, sans signe apparent de recombinaison ou de transfert horizontal de gènes. Dans cette étude, nous appliquons la génomique comparative à un jeu de données de séquençage de génome entier (WGS) composé de 70 souches bovines de M. bovis issues de différentes lignées (Europe et Afrique) afin d'étudier la diversité génétique et le potentiel évolutif de M. bovis chez les bovins. Trois méthodes différentes sont utilisées pour évaluer les signes de réorganisation. Globalement, un petit nombre d'événements de recombinaison ont été identifiés et confirmés par deux méthodes indépendantes avec un solide soutien. Néanmoins, comparée aux mutations, la recombinaison a un effet plus faible sur la diversité de *M. bovis* (r/m global = 0,037). La différence moyenne de r/m obtenue dans le complexe clonal de *Mycobacterium bovis* de notre jeu de données est cohérente avec le concept général selon lequel le degré de recombinaison peut varier considérablement entre les lignées appartenant à la même espèce taxonomique. Sur la base de ces travaux, la recombinaison chez *Mycobacterium bovis* ne peut être exclue et devrait donc faire l'objet d'études plus approfondies dans le cadre de futures recherches en génomique comparative, pour lesquelles le séquençage du génome entier (WGS) de vastes jeux de données provenant de différents contextes épidémiologiques à travers le monde est crucial. Une analyse complémentaire a ensuite été réalisée sur un jeu de données plus restreint de *Mycobacterium bovis* (n = 42) issu de l'étude de prévalence de la tuberculose multi-hôtes. Plus de 1 800 loci ont été identifiés, parmi lesquels au moins une souche présentait un polymorphisme nucléotidique unique (SNP). La plupart des SNP (87,1 %) sont situés dans la région codante, et le rapport global des mutations non synonymes (dN/dS) aux mutations synonymes dépasse 1,5, ce qui indique que la sélection positive est une force évolutive importante chez *M. bovis*. Une proportion plus élevée de SNP a été détectée dans des gènes riches en gènes impliqués dans les catégories fonctionnelles du « métabolisme des lipides », de la « paroi cellulaire et des processus cellulaires », et du « métabolisme intermédiaire et de la respiration », révélant leur importance potentielle dans la biologie et l'évolution de *Mycobacterium bovis*. Un examen plus approfondi des gènes des ancêtres du complexe *Mycobacterium tuberculosis* (MTBC) susceptibles de subir des transferts horizontaux de gènes et inclus dans le système 3R (réparation, réplication et recombinaison de l'ADN) révèle une valeur moyenne négative globale du test D neutre de Taijima, ce qui indique un balayage sélectif antérieur. Le récent goulot d'étranglement après l'expansion de la population demeure le principal moteur de l'évolution pour le pathogène obligatoire *Mycobacterium bovis* dans sa lutte contre l'hôte.
Le complexe Mycobacterium tuberculosis (MTBC) est l'un des taxons de bactéries pathogènes les plus prolifiques et un exemple typique d'évolution bactérienne. Ses membres présentent une identité nucléotidique étonnamment élevée au niveau génomique (> 99 %)1,2. Différents écotypes de MTBC peuvent provoquer la tuberculose (TB), une maladie infectieuse granulomateuse, chez un large éventail d'hôtes, des micromammifères à l'homme3,4,5. Actuellement, le complexe comprend des pathogènes humains [M. tuberculosis (Mtb), Mycobacterium africanum] et des pathogènes adaptés aux animaux (Mycobacterium bovis, Mycobacterium capitum, Mycobacterium pinnipedum, Mycobacterium microtobacter, Mycobacterium mongee, Mycobacterium miysani, Mycobacterium surika, ainsi que les espèces Bacillus chimpanzee et Bacillus dassie5,6). M. canettii (également connu sous le nom de « Nodobacter glabrata ») L'identité nucléotidique moyenne avec les mycobactéries susmentionnées est de 98 %, et des travaux de génomique comparative ont montré que M. canettii et le reste du MTBC ont récemment divergé de l'ancêtre commun.7 Compte tenu de ce concept, certains auteurs appellent M. canettii Membre du MTBC 8.
Le complexe Mycobacterium tuberculosis (MTBC) est systématiquement décrit comme un complexe clonal strict, et sa structure de population est clairement régie par une diversité réduite, des goulots d'étranglement, un balayage sélectif et une dérive génétique9,10. Dans le cadre d'une évolution clonale stricte complexe, les polymorphismes manquants ne peuvent être restaurés par recombinaison. Sur la base de cette hypothèse, les événements consécutifs de délétion génomique de la région différentielle (RD) et de TbD1 (région de délétion spécifique à Mtb 1) ont été proposés comme marqueurs moléculaires de l'évolution du MTBC2,5,11. Les études de génomique comparative et de séquençage du génome entier (WGS) confirment la division des membres adaptés à l'homme en neuf lignées (Mycobacterium tuberculosis L1 à L4, L7 et L8 ; et Mycobacterium africanum L5, L6 et L9), les lignées L2 à L4 partageant la région de délétion TbD12,11,12,13. De plus, il est proposé que les membres adaptés aux animaux partagent un ancêtre commun, qui est défini par des délétions spécifiques au clade dans RD7, RD8, RD9 et RD102, 5 et 14.
Les transferts horizontaux de gènes (THG) et les événements de recombinaison sont considérés comme rares et se produisent chez les ancêtres du complexe Mycobacterium tuberculosis (MTBC), plutôt que dans l'histoire évolutive de l'ensemble des membres de ce complexe15,16,17. Deux études préliminaires, menées par Hughes et ses collaborateurs (2002) et Gutacker et ses collaborateurs (2006), suggéraient que les événements de recombinaison pourraient contribuer à façonner les polymorphismes qui marquent des loci spécifiques chez les souches de M. tuberculosis18,19. L'absence apparente de recombinaison au sein du MTBC s'explique par : (1) le processus mécanique et la perte de la capacité de THG ; (2) la rareté des événements de THG ; (3) l'absence de conditions favorables à la recombinaison dans la niche écologique du MTBC14,17. Récemment, certaines études de séquençage du génome entier (WGS) appliquées à la souche MTBC 20 et à Mycobacterium bovis 21 ont fourni des preuves de recombinaison, les premières à montrer que les souches MTBC échangent souvent de petits fragments d'ADN, mais en raison de la variation limitée de la séquence nucléotidique, ces événements ne sont toujours pas remarqués.
Mycobacterium bovis est l'espèce de MTBC la plus fréquemment isolée chez les animaux d'élevage (principalement les bovins), bien qu'elle puisse également être isolée chez les animaux sauvages en liberté et en enclos4,22,23,24. M. bovis a évolué en cinq complexes clonaux majeurs [Européen 1 (Eu1), Européen 2 (Eu2), Européen 3 (Eu3), Africain 1 (Af1) et Africain 2 (Af2)], selon le profil de spoligotypage, des délétions spécifiques et des polymorphismes nucléotidiques (SNP)25,26,27,28,29 dans des gènes spécifiques. Ces complexes clonaux témoignent de la diversité structurale de la population de Mycobacterium bovis et de son association avec des régions géographiques. De plus, les travaux récents de séquençage du génome entier (WGS) menés par Zimpel et ses collaborateurs (2020) ont permis d'établir une phylogénie basée sur les SNP de Mycobacterium bovis, avec plus de 1 900 génomes. Cette phylogénie indique l'existence d'au moins quatre lignées différentes (nommées Lb1 à Lb4), qui ne correspondent pas entièrement au complexe clonal précédemment défini, bien qu'une spécificité géographique puisse également être confirmée30. Ces auteurs ont réalisé une analyse différentielle de la phylogénie et de la datation moléculaire, mais n'ont pas étudié la recombinaison30.
Des travaux antérieurs, utilisant différentes techniques moléculaires telles que le spoligotypage, le MIRU-VNTR (unité de répétition intercalée mycobactérienne à nombre variable de répétitions en tandem) et le typage SNP récent, ont révélé une certaine diversité génétique parmi les souches de *M. bovis* 31, 32, 33, 34 et 35. La caractérisation de cette variation génétique est devenue un outil important en épidémiologie des maladies, contribuant à une compréhension approfondie de la pathogenèse, de la virulence et de la transmission de la maladie. L'émergence du séquençage du génome entier (WGS) offre la possibilité de révéler les facteurs évolutifs à l'origine de l'adaptation et de la persistance de *Mycobacterium bovis* à différents hôtes et contextes épidémiologiques.
Dans ce travail, nous utilisons l'analyse génomique comparative sur divers ensembles de données de Mycoplasma bovis (n=70), y compris des isolats provenant de différents complexes clonaux, pour mieux comprendre le processus évolutif de Mycoplasma bovis, en particulier pour résoudre les relations phylogénétiques et les événements de recombinaison. En complément de cette analyse, un sous-ensemble de données d'isolats de *M. bovis* (n = 42) provenant d'une zone de tuberculose multi-hôtes bien caractérisée au Portugal31,36 a été étudié plus en détail afin de déterminer la non-identité des gènes. L'équilibre entre le rapport relatif des substitutions nucléotidiques sens (dN) et synonymes (dS), ainsi que la contribution évolutive de génomes spécifiques mentionnés dans la littérature,37,38 obtenus par les ancêtres du complexe *Mycobacterium tuberculosis* (MTBC) par transfert horizontal de gènes (THG), codent pour les composants du système 3R (réparation, réplication et recombinaison de l'ADN)39. Le choix des gènes obtenus par THG se justifie par leur potentiel à représenter des polymorphismes anciens, et donc à présenter une proportion plus élevée de substitutions synonymes. Les gènes inclus dans le système 3R ont été sélectionnés car des travaux antérieurs sur des souches de *M. tuberculosis* ont mis en évidence des sélections négatives/de purification générales agissant sur ces gènes, qui pourraient jouer un rôle important dans l'évolution39. Un autre objectif de ce travail est de déterminer l'existence d'événements de réorganisation. Pour cette raison, considérant que notre ensemble de données provenant du Portugal ne contient que les génomes du complexe de clones européen 2 et les souches auxquelles aucun complexe de clones n'a été attribué, nous avons décidé d'inclure des données génomiques accessibles au public afin d'obtenir enfin un représentant de tous les complexes de clones et d'améliorer la robustesse et l'étendue des résultats.
Quarante-deux génomes de Mycoplasma bovis nouvellement séquencés, provenant du contexte portugais de tuberculose multi-hôtes endémique (détails ci-dessous) et précédemment caractérisés d'un point de vue épidémiologique36, sont au cœur de cette étude. Considérant que l'ensemble de données portugais ne comprend que des représentants des complexes clonaux européens 2 et des souches sans complexe désigné, des données de séquençage de génomes entiers accessibles au public ont été ajoutées afin d'élargir cet ensemble et d'y inclure tous les représentants des complexes clonaux de M. bovis. Ainsi, trois sources de données de séquençage de génomes entiers ont été utilisées : un assemblage de génome complet ou préliminaire, jusqu'à 10 échafaudages, stockés dans la base de données NCBI (National Center for Biotechnology Information) (n = 15 isolats) ; des données stockées dans SRA (The Illumina fastq file of the sequence read archive), représentant la diversité complexe des clones de M. bovis (n = 12 isolats)30 ; et les 42 génomes nouvellement séquencés du Portugal. Mycobacterium bovis BCG (Bacille Calmette-Guérin) a été exclu de la recherche NCBI. M. bovis AF2122/97 est généralement utilisé comme génome de référence pour l'inclusion dans l'ensemble de données. En raison de l'indisponibilité publique de la séquence complète du génome représenté par le complexe de clonage African 1, et du faible nombre de génomes issus de souches représentatives d'Af2 et d'Eu1, les données de séquençage originales fournies par SRA ont été utilisées. Les travaux de Zimpel et de ses collaborateurs (2020) ont permis d'identifier le génome du complexe de clonage susmentionné et de sélectionner Mycobacterium bovis pour l'inclusion dans l'ensemble de données. Pour Eu3, un seul type de génome est décrit (Branger et al., 2020) ; le génome que nous incluons est donc un représentant distinct du complexe Eu3.
À l'échelle mondiale, cet ensemble de données comprend 70 souches de *Mycobacterium bovis* isolées chez des bovins, provenant de 8 espèces hôtes et réparties dans 12 pays entre 1985 et 2016. Parmi ces souches, 36 sont classées Eu2, 7 Eu1, 1 Eu3, 3 Af1, 4 Af2 et 19 ne sont rattachées à aucun complexe clonal (voir détails ci-dessous). Les informations détaillées (y compris le numéro d'accession) des souches de *Mycobacterium bovis* utilisées dans cette étude sont présentées dans le Tableau 1 et le Tableau supplémentaire 1.
Quarante-deux génomes complets de Mycobacterium bovis, nouvellement séquencés et provenant de zones à risque de tuberculose animale au Portugal et répartis sur plus de douze ans, sont au cœur de cette étude, les systèmes potentiels de transmission de maladies de la faune sauvage au bétail ayant fait l'objet d'une surveillance régulière31,36 (Figure supplémentaire 1). Conformément aux procédures décrites ci-après, ces souches ont été isolées de bovins (n ​​= 14), de cerfs élaphes (n = 16) et de sangliers (n = 12) entre 2003 et 2015 : la collecte et la manipulation des animaux ont été réalisées conformément aux recommandations du protocole. Les échantillons de tissus sont répertoriés dans le Manuel de l'OIE pour les animaux terrestres et ensemencés sur milieu solide et liquide au pyruvate de Stonebrink et de Löwenstein-Jensen. Les cultures sont incubées à 37 °C et la croissance est contrôlée une fois par semaine pendant au moins douze semaines. Les colonies sont conservées directement dans une solution de glycérol à -80 °C. Dans le milieu sélectif pour Mycobacterium (Middlebrook 7H9, BD Diagnostics), les échantillons archivés originaux ont subi un seul passage in vitro afin d'obtenir l'ADN nécessaire au séquençage du génome entier (WGS). Pour ce faire, la solution mère de culture congelée a été enrichie avec 5 % de pyruvate de sodium et 10 % d'ADS (50 g d'albumine, 20 g de glucose, 8,5 g de chlorure de sodium dans 1 L d'eau) sur milieu Middlebrook 7H9 à 37 °C. Après 4 semaines de culture, le milieu a été renouvelé et la croissance a été régulièrement surveillée. Les cellules ont été récoltées par centrifugation, le culot a été remis en suspension dans 500 µL de tampon phosphate salin (PBS), chauffé à 99 °C pendant 30 minutes, centrifugé, et le surnageant a été conservé à -20 °C jusqu'au WGS. Toutes les manipulations ont été effectuées en laboratoire de biosécurité de niveau 3.
La bibliothèque génomique WGS à double extrémité est préparée à partir de l'index unique de chaque échantillon d'ADN et séquencée par les technologies Illumina MiSeq (2 × 250 pb) (40 échantillons) et HiSeq (2 × 150 pb) (deux isolats) (Eurofins Genomics, Allemagne). Conformément aux instructions du fabricant, l'ADN génomique est séquencé à l'aide de l'analyseur de génome Illumina équipé d'un module double extrémité, et la bibliothèque est construite avec le kit Nextera XT DNA Library Prep d'Illumina.
Compte tenu des données extraites du SRA (n = 12), l'identification du complexe clonal peut servir de métadonnées pour les publications correspondantes 30, 41 et 43. L'analyse du génome complet, à l'exception de Mycobacterium bovis AF2122/97 et Mycobacterium bovis 3601, membres reconnus des complexes clonaux Eu1 et Eu3 respectivement 25 et 29, est identique à celle du génome complet de Mycobacterium tuberculosis H37Rv (numéro d'accès NCBI : NC_000962.3). L'alignement des génomes est réalisé à l'aide de MAFFT (programme d'alignement multiple de séquences d'acides aminés ou de nucléotides, version 7.458) avec le paramètre -addfragments48. L'absence de différents complexes clonaux et/ou la présence de SNP sont ensuite recherchées.
Les Mycobacterium bovis nouvellement séquencés (n = 42) et les lectures originales de l'ébauche de génome assemblée (n = 3) alignent le complexe avec le génome de référence Mycobacterium tuberculosis H37Rv via le pipeline vSNP et la présence de la délétion et/ou des caractéristiques SNP de différents clones. Une recherche a été menée.
Les données génomiques sont attribuées au complexe clonal correspondant grâce à l'absence de caractéristiques et/ou à la présence/absence de profils SNP et de spoligotypage. Pour les quatre assemblages préliminaires, le profil de spoligotypage ne peut être inféré ; ils sont donc inclus dans le groupe « sans complexité ».
Le flux de travail bioinformatique utilisé dans cette étude débute par un assemblage de novo et un alignement sur une stratégie de référence, dans le but d'explorer les événements de recombinaison et les polymorphismes génomiques spécifiques. La figure 1 présente un organigramme des étapes suivies. Pour l'analyse de la recombinaison, tous les génomes sont utilisés afin d'accroître la robustesse des inférences et des indicateurs associés.
Afin de réduire les erreurs lors de la génération des séquences consensus génomiques, nous avons d'abord réalisé un assemblage de novo, puis des alignements multiples de base. Le pipeline Unicycler est actuellement disponible à l'adresse https://github.com/rrwick/Unicycler49 et a été utilisé pour l'assemblage de novo de 54 génomes séquencés (42 nouvellement séquencés et 12 fichiers fastq récupérés de SRA). En résumé, avant l'assemblage initial, une analyse de la qualité des lectures a été réalisée avec FastQC version 0.11.7 (https://github.com/s-andrews/FastQC) et Trimmomatic version 0.36 (les options « Suppression des adaptateurs et autres séquences spécifiques à l'éclairage » et « Suppression des bases en fin de lecture si leur qualité est inférieure au seuil de 20 » ont été activées) (http://www.usadellab.org/cms/?page=trimmomatic)50. L'optimiseur SPAdes49 a ensuite été utilisé pour l'assemblage du génome, et Pilon version 1.1851 pour l'optimisation post-assemblage. Un mode de pontage conservateur a été sélectionné afin d'éviter les erreurs d'assemblage, et la taille des k-mers a été déterminée entre 20 % et 95 % de la longueur des lectures. Suivez les directives SPAdes et tenez compte de la taille des lectures, supprimez les contigs inférieurs à 300 pb et établissez un seuil de couverture de profondeur de lecture de 52. Dans la stratégie d'assemblage de novo, les régions génomiques telles que les paralogues hautement répétitifs proline-glutamate (PE) et proline-proline glutamate (PPE) n'ont pas été supprimées.
La qualité de l'assemblage de novo est évaluée via le pipeline QUAST (http://quast.sourceforge.net/quast.html), qui facilite le renouvellement du contig et du génome de référence M. bovis AF2122/97 (numéro d'accès NCBI LT708304.1) Mapping (voir le tableau supplémentaire 1 pour les paramètres de qualité).
Grâce au pipeline vSNP (https://github.com/USDA-VS/vSNP), le fichier FASTQ du génome de M. bovis nouvellement séquencé par Illumina est comparé au génome de référence de M. bovis AF2122/97 (LT708304.1). Conformément aux recommandations de bonnes pratiques du Genome Analysis Toolkit (GATK)53, 54, 55, des paramètres de filtrage standard ou des scores de masse des variants sont appliqués pour la recalibration. Les résultats sont filtrés avec le score de masse SAMtools le plus faible (150) et un AC de 2. Kraken (http://ccb.jhu.edu/software/kraken/) est également utilisé pour vérifier les lectures et exclure toute contamination. Le pipeline vSNP utilisé pour l'alignement sur les stratégies de séquences examine une série de SNP et de cibles définis, et exclut les scénarios d'infection mixte. La couverture du génome par les lectures est supérieure à 99 % (Tableau supplémentaire 1).
Pour éviter les erreurs de cartographie et les SNP incorrects, il convient de filtrer les variants dans les cas suivants : (1) ils sont présents dans moins de 20 lectures, (2) leur fréquence est inférieure à 0,9, (3) ils sont présents dans au moins une souche, mais présentent des lacunes dans au moins une autre souche. Le logiciel Integrated Genomics Viewer (IGV) version 2.4.19 (http://software.broadinstitute.org/software/igv/)⁵⁶ a été utilisé pour vérifier visuellement les SNP et les positions présentant des problèmes de cartographie ou d’alignement. Les gènes proline-glutamate (PE) et proline-proline-glutamate (PPE), fortement dupliqués et appartenant à une famille multigénique, sont facilement mal interprétés par le séquençage Illumina et font donc l’objet d’une analyse plus approfondie. Le flux de travail bioinformatique pour les mycobactéries a permis d’éliminer les membres du complexe de la tuberculose lors de l’utilisation de la stratégie de cartographie sur séquence pour confirmer les SNP. Par conséquent, nous avons exclu de l'analyse les gènes PE/PPE et les indels.
Selon Bovilist (http://genolist.pasteur.fr/BoviList/), tous les SNP sont classés par catégories fonctionnelles. Le pipeline SnpEff (https://pcingola.github.io/SnpEff/) est utilisé pour inférer les conséquences des SNP (modifications synonymes ou non synonymes). Une nouvelle base de données du génome de Mycobacterium bovis AF2122/97 (LT708304.1) a été créée.
L'alignement multiple du génome central a été réalisé à l'aide de Parsnp v1.2, actuellement disponible à l'adresse https://github.com/marbl/parsnp57, en utilisant 69 génomes complets/assemblages préliminaires (avec l'option -c) et *M. bovis* AF2122/97 (LT708304.1) comme référence. Quatre alignements multiples du génome central ont été effectués : un alignement portant uniquement sur les membres du complexe de clonage Eu2 (n = 37), un alignement incluant tous les membres du complexe de clonage européen (n = 44), un alignement incluant le point de jonction des complexes de clonage européen et africain (n = 51), et un alignement incluant toutes les souches de *Mycobacterium bovis* étudiées (n = 70).
L'alignement de base généré par Parsnp est utilisé pour inférer l'arbre phylogénétique de vraisemblance maximale (ML) à l'aide de CIPRES Science Gateway v3.3 (http://www.phylo.org/)58 en utilisant RAxML, et effectuer 1000 réplications guidées.
Trois algorithmes et outils bioinformatiques différents sont utilisés pour vérifier en parallèle la présence d'événements de recombinaison : le logiciel SplitsTree4, le pipeline Gubbins (lignée non biaisée par recombinaison dans les séquences nucléotidiques) et le logiciel RDP4 (programme de détection de recombinaison, version bêta 4.101).
La méthode de décomposition par division implémentée dans SplitsTree4 v4.15.1 (http://www.splitstree.org/)59 est utilisée pour calculer le réseau phylogénétique sans racine, avec une vérification statistique par le test Phi et un seuil de signification de p = 0,05. L'analyse d'alignement multiple de base de Parsnp est utilisée comme entrée, et la décomposition par division, en tant que norme de réseau, est réalisée.
Le pipeline Gubbins v2.3.1 (https://github.com/sanger-pathogens/gubbins60) est exécuté avec les paramètres par défaut afin d'évaluer l'impact de la recombinaison chez Mycobacterium bovis. L'algorithme implémenté dans le pipeline reconstruit la lignée clonale pertinente. Le génome complet/l'assemblage préliminaire de notre jeu de données et le génome de référence (bovin AF2122/97, LT708304.1) sont alignés. Le pipeline analyse ensuite la position des SNP sur chaque branche de l'arbre pour détecter le cluster de SNP représentant l'événement de recombinaison. L'origine de la branche correspond à l'absence de recombinaison, ce qui signifie que les SNP présents sur cette branche devraient être uniformément répartis. L'alignement multiple de base de Parsnp et l'arbre ML le mieux noté de RAxML sont utilisés comme fichiers d'entrée.
Enfin, afin de confirmer l'événement de réorganisation suggéré par le pipeline Gubbins, les six algorithmes implémentés dans RDP467 (RDP61, GENECONV62, Bootscan63, Maxchi64, Chimaera65 et SiScan66) sont appliqués à l'alignement multiple de base de Parsnp avec les paramètres par défaut. Nous avons déterminé qu'au moins trois des algorithmes implémentés dans RDP4 doivent systématiquement démontrer un signal significatif pour vérifier chaque événement de recombinaison.
Étant donné que les logiciels Gubbins et RDP recherchent les signaux de recombinaison en vérifiant les alignements multiples centraux dans une fenêtre allant jusqu'à 500 pb, et après avoir confirmé que l'inclusion des gènes PE/PPE lors de l'assemblage de novo n'interfère pas avec les signaux de recombinaison détectés, une analyse d'homolinéarité est effectuée au voisinage du gène identifiant l'événement de recombinaison. La carte synlinéaire, utilisant le génome complet, a été construite à l'aide de l'alignement multigénome MAUVE (http://darlinglab.org/mauve/mauve.html) afin d'exclure les translocations ou inversions génomiques locales. De plus, le génome entier a été utilisé pour réaliser une analyse d'homolinéarité sur la séquence d'acides aminés via le serveur web SyntTax (https://archaea.i2bc.paris-saclay.fr/SyntTax/).
Une analyse plus approfondie des données génomiques obtenues à partir du système portugais de tuberculose multi-hôtes vise à vérifier le polymorphisme des gènes mentionnés dans la littérature. Ces gènes sont les gènes 37 et 38, ainsi que le gène codant pour la protéine 3R, acquis par les ancêtres du complexe Mycobacterium tuberculosis (MTBC) via les composants du système de transfert horizontal de gènes (THG) (réparation, réplication et recombinaison de l'ADN) 39. Utilisez ClustalX v2.1 (http://www.clustal.org/clustal2/) et DnaSP v6.12.03 (http://www.ub.edu/dnasp/) pour calculer la diversité génique et la diversité nucléotidique (π), ainsi que les paramètres du test de neutralité D de Tajima.
Un arbre phylogénétique de vraisemblance maximale (ML) basé sur 69 isolats de Mycoplasma bovis et des génomes de référence a été obtenu (Figure 2A). Comparée aux arbres basés sur un seul gène ou sur plusieurs locus, cette stratégie permet de générer des arbres plus robustes, mais qui ne rendent pas compte de la variabilité de l'ensemble du génome et présentent donc une capacité de discrimination entre les espèces moindre68,69. La structure topologique de l'arbre ML est généralement cohérente avec la classification complexe des clones. Le génome d'Eu2 est regroupé dans une branche, et celui d'Af1 est également regroupé (Figure 2A). Ce résultat est également cohérent avec les relations évolutives connues de Mycobacterium bovis : il existe une différence significative entre le membre Eu1 et le groupe constitué de tous les autres complexes clonaux et génomes, mais le complexe clonal 30 n'est pas spécifié. La légère incohérence entre le complexe clonal et la relation observée sur l'arbre phylogénétique peut s'expliquer par le fait que le complexe clonal est décrit sur la base de régions génomiques spécifiques, tandis que l'arbre phylogénétique est basé sur de multiples alignements de génomes centraux représentant l'ensemble du génome.
L'arbre phylogénétique de vraisemblance maximale (GTR) est construit à partir de l'alignement du génome central de Mycobacterium bovis avant (A) et après (B) suppression du site de recombinaison. Les couleurs des branches représentent le complexe de clones de Mycobacterium bovis : Europe 1 est en violet, Europe 2 en rouge, Europe 3 en bleu, Afrique 1 en orange et Afrique 2 en vert. L'arbre est enraciné et dessiné à l'échelle, la longueur des branches étant mesurée en fonction du remplacement de chaque site.
Le complexe Mycobacterium tuberculosis est décrit comme ayant évolué de manière clonale, et la plupart des preuves accumulées au fil des ans soutiennent l'idée que les événements HGT et de recombinaison en cours ne se produiront pas au niveau détectable de MTBC15,17,18.
Des travaux antérieurs ont montré qu'il pourrait exister une recombinaison limitée entre les souches de MTBC20,21, tandis que d'autres n'ont pas permis d'identifier d'événements de recombinaison mesurables70,71. Il convient de réexaminer cette question en se concentrant sur Mycobacterium bovis, contrairement aux travaux précédents qui ne portaient que sur Mycobacterium tuberculosis70,71 ; ou bien de considérer le MTBC dans son ensemble, en excluant presque totalement M. bovis20 ; ou encore de se limiter aux fractions bovines. Dans cette étude, l'ensemble de données sur les mycobactéries comprend 70 souches, représentant tous les complexes clonaux, utilisées pour rechercher des recombinaisons. Cet ensemble de données est structuré selon quatre niveaux cumulatifs : (1) les membres du complexe Eu2, (2) tous les membres du complexe clonal européen (Europe), (3) le complexe clonal européen et africain (Eu + Af) et (4) l'ensemble des collections de données (y compris les génomes n'appartenant à aucun complexe clonal déjà décrit).
Pour approfondir cette hypothèse, une analyse de réseau par décomposition a été réalisée afin d'évaluer l'absence d'événements de recombinaison entre les génomes. Cette méthode permet en effet de visualiser les relations ancestrales entre les individus et de mettre en évidence les signaux phylogénétiques contradictoires. Les quatre jeux de données analysés ont confirmé l'existence de boucles dans le réseau (c'est-à-dire des zones ne convergeant pas vers un arbre unique), mais le test Phi n'a pas apporté de soutien statistique (Eu2, p = 0,0956 ; Europe, p = 0,1637 ; Eu + Af, p = 0,2774 ; l'ensemble des données, p = 0,2451), ce qui constitue une preuve insuffisante de l'existence d'événements de réorganisation (Figure 3A-D).
En Europe 2 génomes (n = 37) (A), génomes européens (n ​​= 44) (B), génomes européens et africains (n ​​= 51) (C) et l'ensemble des données (n = 70) (D).
Après cette analyse, et en tenant compte des observations cycliques dans tous les réseaux, l'algorithme de reconstruction implémenté dans le pipeline Gubbins a été appliqué pour reconstruire la lignée clonale et compléter l'estimation de l'effet de la recombinaison sur le génome de *M. bovis*. Le nombre cumulé d'événements de recombinaison, dont la plupart se sont produits dans les branches terminales (c'est-à-dire dans un seul génome), a été déduit (Tableau 2). Ces indicateurs montrent la cohérence de l'ensemble des données et indiquent que la fréquence des événements de recombinaison est 200 à 300 fois supérieure à celle des mutations. Le paramètre ρ/θ, représentant les taux relatifs de recombinaison et de mutations ponctuelles sur la branche, se situe entre 0,0037 et 0,0056 (Tableau 3). Récemment, les travaux publiés sur la souche 38 M. bovis ont démontré une valeur rho/theta plus élevée (rho/theta = 0,1) que celle obtenue dans cet ensemble de données, mais les travaux de Patané et de ses collègues ont utilisé un assemblage basé sur une référence pour déduire les paramètres de recombinaison. Un détail procédural, dû à la procédure d'assemblage, a été associé à l'abondance d'événements de recombinaison putatifs dans la branche terminale.
Ensuite, le paramètre r/m représente le rapport de diversité entre la recombinaison et l'introduction de mutations. Sa valeur moyenne se situe entre 0,025 et 0,037, indiquant que, comparativement aux mutations, la recombinaison a un impact global plus faible sur la diversité génétique de *M. bovis* (Tableau 3). À titre de comparaison, une méthode similaire a été utilisée pour estimer le paramètre r/m pour l'ensemble de données MTBC composé de 23 génomes, obtenant une valeur moyenne de 0,48620. En revanche, pour l'ensemble de données de 38 *M. bovis* de Patané et al., la valeur moyenne obtenue est de 0,98. Dans la première étude, seuls deux des 23 génomes de *M. bovis* ont été inclus (*M. bovis* BCG et la souche de référence), ce qui pourrait biaiser la valeur obtenue en raison de la surexpression du génome de *M. tuberculosis*. Dans la seconde étude, les populations de *Mycobacterium bovis* analysées provenaient principalement des États-Unis et d'animaux d'élevage. En revanche, notre jeu de données représente un plus grand nombre de localisations géographiques et d'espèces hôtes, et utilise des génomes regroupés en différents complexes clonaux présentant des caractéristiques génétiques de population distinctes, permettant ainsi une connaissance plus approfondie et plus large des populations. La différence de valeur moyenne du rapport r/m obtenue avec notre jeu de données confirme l'idée que le degré de recombinaison varie considérablement entre les lignées appartenant à la même espèce taxonomique. Ces résultats indiquent donc que le complexe clonal de *M. bovis* pourrait présenter des différences de recombinaison, comme suggéré par Didelot et Maiden72. Néanmoins, l'élargissement significatif de ce jeu de données par l'inclusion d'un plus grand nombre de génomes de *M. bovis* permettra de clarifier davantage ce point. Les paramètres r/m et ρ/θ présentent tous deux une variabilité entre les branches, un résultat cohérent avec les observations faites sur d'autres espèces bactériennes72,73.
Enfin, afin de confirmer les événements de réorganisation identifiés par le pipeline Gubbins, six algorithmes différents ont été utilisés dans le logiciel RDP4 pour tester indépendamment différentes comparaisons multiples. Globalement, moins de la moitié des événements identifiés par Gubbins ont été confirmés par RDP4 (tableaux 4 et 5). Sur l'ensemble des données, trois événements de recombinaison ont été confirmés : deux impliquant des nœuds internes et le troisième impliquant un seul génome dans une branche terminale, pour lequel aucun complexe clonal n'a pu être assigné (tableaux 4 et 5). L'identification d'événements dans les branches terminales peut indiquer que la recombinaison est toujours en cours dans les souches actuelles de *M. bovis* ou que le résultat est erroné70. Dans cette région de recombinaison hypothétique, environ 20 % des positions présentent des nucléotides indéterminés (N), affectant ainsi le signal de recombinaison (figure supplémentaire 2). De plus, cette région affecte le gène rrs, codant pour l'ARN ribosomique 16S, qui est supposé être hautement conservé. Ce signal de recombinaison putatif pourrait donc résulter d'erreurs de séquençage ou d'un mauvais alignement. Un alignement complet des génomes de Mb0003 et de Mycobacterium bovis AF2122/97 a ensuite été réalisé, confirmant la présence de nucléotides non définis et de SNP. Par conséquent, les problèmes potentiels liés à un alignement incorrect ne sont pas dus aux informations biologiques utilisées dans cette étude.
Aucune lacune ni nucléotide indéfini n'a été détectée dans les régions de recombinaison des nœuds internes (figures 4 et 5). Parmi ces événements, l'un concerne uniquement le génome d'Eu2 et affecte le gène pks12, qui code pour une possible polycétide synthase ; l'autre est enregistré dans le génome d'Eu1 et affecte le gène narX, qui code pour une possible nitrate réductase (tableau 4). De manière générale, l'analyse de recombinaison montre un nombre limité de fragments de recombinaison statistiquement significatifs, et les indicateurs obtenus suggèrent que la recombinaison a un faible impact sur la lignée de M. bovis. Le signal de recombinaison est vraisemblablement faible, mais il est important de distinguer le véritable signal d'évolution du bruit de fond, ce qui représente un défi. Afin de réduire le bruit introduit par l'assemblage basé sur un génome de référence et les problèmes d'appariement (références 70 et 71), tous les fragments, à l'exception du génome complet, ont été assemblés de zéro. La qualité de l'assemblage a ensuite été vérifiée et garantie par l'analyse du pipeline QUAST (Tableau supplémentaire 1). De plus, une série d'analyses complémentaires ont été réalisées pour confirmer la robustesse et la précision de l'étude globale. Ainsi, la qualité du séquençage des gènes narX et pks12 a été évaluée par alignement des lectures sur Mycobacterium bovis AF2122/97. La position du SNP recommandé dans la région de recombinaison a été confirmée en appliquant les critères mentionnés dans la section Méthodes (au moins 20 lectures et une fréquence de changement de 0,9). Le polymorphisme du gène narX a été pleinement confirmé dans les deux génomes (Mb1792361 et Mb7240415 ; 2,3 %) et dans les génomes du gène pks12 : gènes Mb0891, Mb1711, Mb1789, Mb1870, Mb17046, Mb1756 et Mb12. Cependant, pour le génome Mb2043, six des huit positions ne satisfont pas au critère de profondeur de lecture, car le SNP est couvert par un maximum de 17 lectures, ce qui est inférieur au seuil établi de 20. Par conséquent, la recombinaison de six génomes (8,6 %) à ce site génomique peut être confirmée (figures 4 et 5).
La visualisation détaillée de l'alignement de la région de recombinaison du jeu de données de Mycobacterium bovis affecte le gène narX, qui code pour une possible nitrate réductase. Aucune lacune ni nucléotide indéterminé n'a été détecté dans la région de recombinaison des nœuds internes. Cet événement particulier est enregistré dans le génome Eu1. La qualité du séquençage du gène narX a été évaluée en visualisant les lectures de Mycobacterium bovis AF2122/97. La localisation du SNP recommandée dans la zone de recombinaison a été confirmée en appliquant les critères mentionnés dans la section Méthodes (au moins 20 lectures et une fréquence de changement de 0,9). Le polymorphisme du gène narX a été pleinement confirmé dans les génomes de Mb1792361 et Mb7240415 (2,3 %).
Visualisation détaillée de l'alignement de la région de recombinaison du jeu de données de Mycoplasma bovis affectant le gène pks12. Aucune lacune ni nucléotide indéfini n'a été détectée dans la région de recombinaison des nœuds internes. L'événement affectant le gène pks12, codant potentiellement pour une polycétide synthase, ne contient que le génome Eu2. La qualité du séquençage de pks12 a été évaluée par alignement des lectures avec celles de Mycobacterium bovis AF2122/97. La localisation du SNP recommandée dans la zone de recombinaison a été confirmée en appliquant les critères mentionnés dans la section Méthodes (au moins 20 lectures et une fréquence de changement de 0,9). Les polymorphismes des génomes Mb0891, Mb1711, Mb1789, Mb1870, Mb1758, Mb2043 et Mb1960 ont été pleinement confirmés.
Les gènes PE et PPE possèdent des régions répétitives facilement sujettes à des erreurs de lecture et de cartographie lors du séquençage Illumina. De ce fait, ils sont généralement exclus du flux de travail bioinformatique des membres de *M. tuberculosis* uniquement lors de l'utilisation de la stratégie d'alignement sur séquence. L'inférence des événements de recombinaison appliquée dans ce travail repose sur un assemblage *de novo* sans filtrage des gènes PE/PPE. Nous pensons que la mise en œuvre de trois méthodes et algorithmes complémentaires différents, via les logiciels SplitsTree, Gubbins et RDP4, confère aux stratégies appliquées une robustesse face au traitement et au filtrage des régions réorganisées induites par les signaux d'erreur. Cependant, afin d'exclure toute interférence des gènes PE/PPE avec les logiciels Gubbins et RDP4 pour l'identification des clusters de SNP, et par conséquent l'identification des régions de recombinaison susceptibles d'affecter les gènes *narX* et *pks12*, le voisinage de ces gènes a été examiné (Figures supplémentaires 3 à 5). Chez M. bovis AF2122/97, le gène narX est séparé par narK2 et Mb1764c, tandis que pks12 est entouré par Mb2075c et Mb2073c (Figure supplémentaire 3-5). La carte générée à partir de la carte synline MAUVE du génome complet fournit des informations sur la conservation et le réarrangement des séquences géniques, révélant quatre blocs colinéaires et aucune translocation ni inversion génomique. De plus, l'analyse de complémentation avec la séquence d'acides aminés a confirmé l'homologie dans tous les génomes complets, et aucune PE/PPE n'a été trouvée dans les régions adjacentes à narX ou pks12. Concernant narX, un génome (Mb0030) présente un score de synonymie plus faible car le gène narX a été identifié comme étant composé de deux fragments (fragments 1891 et 1890). Pour pks12, en raison de similarités, Mb0030 et Mb003 ont présenté des scores de synlinéarité plus faibles, tandis que pks12 a été identifié dans deux et trois fragments, respectivement, représentant différents domaines de la protéine (Figures supplémentaires 3 à 5). Compte tenu de ces informations, et après analyse par les logiciels Gubbins et RDP4 de l'alignement multiple central sur une fenêtre de 500 pb maximum, nous avons confirmé que le gène PE/PPE n'interfère pas avec le signal de recombinaison affectant narX et pks12.
Bien que les signaux de recombinaison détectés dans cet ensemble de données puissent être considérés comme résiduels, il est vrai que la recombinaison chez M. bovis ne peut être exclue, elle devrait donc continuer à faire l'objet d'analyses plus approfondies, dans lesquelles des génomes entiers provenant de différents scénarios épidémiologiques sont séquencés.
La comparaison des arbres phylogénétiques ML obtenus avant et après correction de recombinaison (Figure 2A, B) n'a pas conduit à des changements significatifs dans la relation phylogénétique inférée, et les souches de M. bovis ont été regroupées dans le même groupe.
Après alignement de 42 séquences de M. bovis nouvellement obtenues avec le génome de référence de M. bovis AF2122/97, un alignement de SNP contenant 1816 positions polymorphes a été obtenu. La plupart des SNP (87,1 %) sont situés dans la région codante, et les gènes affectés sont caractérisés selon les catégories fonctionnelles présentées dans Bovilist (Figure 6A, B). En considérant le nombre total de gènes dans chaque catégorie fonctionnelle, les gènes de la catégorie « métabolisme des lipides » présentent le plus grand nombre de SNP, suivis par ceux de la catégorie « paroi cellulaire et processus cellulaires » et « métabolisme intermédiaire et respiration », ce qui indique leur implication dans l'évolution de M. bovis.
Analyse hiérarchique des données de M. bovis provenant du Portugal (n = 42). Nombre total de SNP et de gènes affectés enregistrés pour chaque catégorie fonctionnelle (A). Nombre total de mutations synonymes et non synonymes enregistrées par catégorie fonctionnelle (B).
À l'échelle mondiale, le rapport dN/dS moyen est supérieur à 1,5, ce qui indique que la pression évolutive globale tend à éliminer l'état ancestral et représente un scénario de choix de purification positif (diversifié ou dirigé) et/ou moins contraignant. Dans les catégories « virulence, détoxification, adaptation », « séquences d'insertion et phages » et « protéines régulatrices », plus des deux tiers des SNP sont non synonymes (figure 6B).
Dans toutes les catégories, on observe des gènes présentant de multiples SNP, ce qui se traduit par un taux de mutation moyen (c’est-à-dire le nombre moyen de SNP par gène) supérieur à 1 (Figure 6A). Les gènes Pks12 (Mb2074c), avec 15 SNP, et fas (Mb2553c), avec 8 SNP, présentent des taux de mutation plus élevés. Ces deux gènes sont impliqués dans le métabolisme des acides gras. Le gène pks code pour la polycétide synthase (PKS), une enzyme multifonctionnelle impliquée dans la biosynthèse des lipides de la paroi cellulaire mycobactérienne74,75. Ce gène code pour un polypeptide multifonctionnel impliqué dans la synthèse des mycocétides74,76. Le gène fas est impliqué dans la synthèse de l’acide mycolique. Ces deux gènes jouent un rôle important dans la biosynthèse de la paroi cellulaire au contact de l’hôte.
Afin d'étudier plus en détail l'évolution de Mycobacterium bovis, deux ensembles de gènes spécifiques ont été analysés. Des travaux antérieurs, utilisant la composition des séquences et des méthodes phylogénétiques, ont identifié des gènes acquis par les ancêtres du complexe Mycobacterium tuberculosis (MTBC) par transfert horizontal de gènes (THG) avant leur diversification37,38. Ces gènes sont listés dans le Tableau supplémentaire 2. La distribution des SNP de 77 gènes potentiellement liés au THG a été analysée, et 26 sites polymorphes ont été identifiés, présentant dans la plupart des cas (78 %) des mutations non synonymes (NS) (Tableau supplémentaire 2). Des travaux antérieurs sur le génome du MTBC ont démontré que la région putative de THG présente un ratio de SNP NS plus élevé que le reste du génome. Si l'on considère que ces régions de recombinaison ont été acquises par les ancêtres du MTBC et qu'elles surreprésentent donc les polymorphismes anciens, la proportion de mutations synonymes devrait être plus élevée. En effet, les substitutions NS sont susceptibles d'être éliminées par sélection négative, car les modifications d'acides aminés peuvent altérer la fonction de la protéine. Nos résultats indiquent donc que les conséquences fonctionnelles pourraient provenir du remplacement de gènes de type HGT, ce qui reflète leur importance pour une précieuse diversité génétique adaptative.
Parallèlement à cette analyse, les gènes codant pour les composants du système 3R (réparation, réplication et recombinaison de l'ADN) ont été examinés en détail, conformément à la liste publiée par dos Vultos et al. (2008)39. L'échange de fragments d'ADN identiques ne peut être observé directement, bien qu'il puisse s'agir d'un processus fréquent chez les bactéries étroitement apparentées, comme c'est le cas pour cet ensemble de données ; de plus, ce processus pourrait être essentiel aux mécanismes de réparation de l'ADN72 et jouer ainsi un rôle dans la recombinaison homologue. Au total, 26 positions polymorphes, réparties sur 54 gènes, ont été identifiées (Tableau supplémentaire 3). Dans cet ensemble de gènes, les mutations non spécifiques (NS) représentaient environ 65 % des variations, ce qui concorde avec les résultats antérieurs obtenus chez les souches de Mycobacterium tuberculosis.


Date de publication : 21 octobre 2021