INTRODUCTION

1Ce jeu de données a fait l’objet d’un article (Bouyer et al., 2009). Il permettra de réviser plusieurs notions et de mettre en pratique de nouvelles méthodes d’analyse telles que l’isolement par la distance entre individus. Comme pour les tiques, ce jeu de données est téléchargeable, mais dans un format différent. Ce fichier s’appelle “TsetseJerCoordGeo&Trap&SexTotData.xls”. Dans ce chapitre, nous ne détaillerons que les analyses nouvelles et irons la plupart du temps très vite sur les notions déjà illustrées dans le précédent chapitre.

ÉTAT DES LIEUX

2Les trypanosomoses africaines figurent parmi les plus sérieuses des maladies tropicales négligées (Schofield et Kabayo, 2008). L’OMS estime que le nombre total de cas de maladie du sommeil avoisine les 300 000 personnes (WHO, 2006a). Par ailleurs, d’après la FAO, le coût économique des trypanosomoses animales (nagana) atteint 4,74 milliards de dollards US par année (FAO, 2000). En 2001, plusieurs pays africains ont lancé le Pan African Tsetse and Trypanosomosis Eradication Campaign (PATTEC) afin d’établir une lutte concertée contre cette plaie à l’origine de nombreux problèmes de faim, d’appauvrissement et de frein au développement d’une agriculture durable dans les zones rurales d’Afrique subsaharienne (http://www.africa-union.org/Structure_of_the_Commission/depPattec.htm). Glossina palpalis s. l. (une des espèces de mouches tsé-tsé) est un des plus importants vecteurs de trypanosomoses humaine et animales en Afrique de l’Ouest. En Guinée, la sous-espèce Glossina palpalis gambiensis (Gpg) transmet la maladie du sommeil avec une prévalence relativement élevée (Camara et al., 2005). Au Burkina Faso, c’est un vecteur majeur de nagana, en particulier dans le bassin de la rivière Mouhoun où se situent les échantillons que nous allons analyser (Bouyer et al., 2006). La connaissance des schémas de dispersion et de tailles de populations est un pré-requis nécessaire au développement d’une lutte raisonnée pour le contrôle des populations de vecteurs (Tabachnick et Black, 1995). Pour les mouches tsé-tsé, comme pour les autres espèces de vecteurs, les estimations directes par marquage-recapture sont fastidieuses et coûteuses et pas nécessairement très fiables dans le cas des mouches tsé-tsé (Terblanche et Chown, 2007). Comme nous allons le voir, les marqueurs génétiques et les outils de la génétique des populations peuvent apporter une solution très efficace.

3Le jeu de données concerne des échantillons de Gpg prélevées dans quatre zones le long de la rivière Mouhoun (fig. 52).

PREMIER RECODAGE DES DONNÉES

4Les données brutes se présentent comme dans le tableau 23. La première colonne indique le site de prélèvement (A, H, C ou D, comme dans la figure 52). Les deuxième et troisième colonnes correspondent aux coordonnées GPS des pièges suivies du nom du piège en quatrième colonne, du sexe et du nom des individus glossines génotypés en colonnes cinq et six respectivement. Suivent les génotypes des allèles aux sept loci étudiés avec une colonne par allèle et donc 14 colonnes (colonnes 7 à 21). Vous remarquerez que les loci liés à l’X possèdent cette lettre dans leur nom (comme pour PgpX11, par exemple) et que les mâles ont été codés homozygotes pour ces loci, ce dont il faudra se souvenir au moment de tester la panmixie. Les données manquantes sont, quant à elles, codées par des “0”. Pour tous les tests liés à l’hétérozygotie locale, il faut créer un second fichier “TsetseJerCoordGeo&Trap&SexTotDataMalManq.xls” où les mâles sont manquants aux loci liés à l’X.

5Ensuite, nous allons utiliser un nouveau logiciel très pratique qui peut convertir facilement nos deux fichiers dans des formats variés, y compris pour les programmes dont nous avons besoin. Ce programme s’appelle Create v 1.1 (Coombs et al., 2008). Vous lancez Create et remplissez la fiche comme dans la figure 53.

6Quand vous sélectionnez le fichier Excel, le programme vous demande dans quelle fiche Excel12 se trouvent les données. Cliquez sur celle qui convient (la 1 en principe). Cliquez ensuite sur “Proceed”. Le programme vous demande de vérifier qu’il a bien pris en compte ce qu’il fallait en vous montrant l’exemple du premier individu. Répondez oui si ça colle. Un second menu apparaît qu’il vous faut remplir comme en figure 54. Vous obtenez ainsi quatre fichiers, deux pour les données en format Genepop et Fstat, et deux pour le nom des populations. Faites la même chose pour “TsetseJerCoordGeo&Trap&SexTotDataMalManq.xls”. Nous allons dans un premier temps tester les déséquilibres de liaison avec “TsetseMouhouMalHomo-FSTAT.dat” et les F_IS avec “DataTsetseMouhoun MalManq-FSTAT.dat”, fichiers Fstat que vient de créer Create. Vous pouvez renommer ces fichiers avec des noms moins longs. Vous pouvez aussi éditer les fichiers *.lab et supprimer les colonnes supplémentaires qu’a créé Create (je ne sais pas pourquoi il fait ça) et qui risquent de générer des problèmes ensuite. Ne gardez que la première colonne de ces fichiers, qui correspond à l’identifiant des sous-populations.

Tableau 23
Extrait du jeu de données brutes des génotypes des individus Gpg capturées le long de la Mouhoun. Le tableau est tronqué pour les derniers loci. Notez qu’un locus occupe deux colonnes

Tableau 23 (suite)

Figure 53
Fiche Menu pour Create pour convertir le fichier de données brutes de mouches tsé-tsé de la Mouhoun au format désiré

PREMIÈRES ANALYSES : INDÉPENDANCE ENTRE ALLÈLES DANS ET ENTRE LOCI

Déséquilibres de liaison au sein des quatre zones

7Lancez Fstat et chargez le fichier “TsetseMouhouMalHomo”. Testez les déséquilibres de liaison en demandant le test “for each pair of loci in each population” et au “nominal level” 1/100 afin d’avoir assez de précision. Dans le fichier de sortie correspondant, nous constatons que seul un test est significatif entre les loci 1 et 2 (c’est-à-dire entre PgpX11 et PgpX13) avec une P-value = 0,0044. Cette P-value ne reste pas significative après correction de Bonferroni (0,0044×21 = 0,09) et un test significatif sur 21 représente environ 5 % des tests, ce qui est la proportion attendue sous l’hypothèse nulle. Avec la procédure "binom. test" sous R, nous pouvons calculer la probabilité avec laquelle nous pouvons observer une fois un test significatif au seuil α = 0,0044 sous l’hypothèse nulle H₀. Cette probabilité est P-value = 0,0889. On peut donc considérer qu’à l’échelle de chaque zone, il y a indépendance entre loci.

Test de la panmixie dans les quatre zones d’échantillonnage

8Chargez “DataTsetseMouhounMalManq. dat” dans Fstat et demandez le F_IS par locus et population, les estimations de Weir et Cockerham et testez Hardy-Weinberg dans les sous-échantillons avec 10 000 permutations d’allèles entre individus. Le résultat peut être résumé dans le tableau 24. On y constate un fort déficit en hétérozygotes très significatif, mais aussi une forte variance du F_IS entre loci. Une recherche d’allèles nuls, de « stuttering » ou de dominance d’allèles courts s’avère nécessaire.

ANALYSE PAR MICRO-CHECKER

9Il faut ici traiter les femelles seules pour les loci hétérosomaux du site A (un fichier de plus), pour les autres sites il n’y pas de problème (pas de mâle). En passant par Create, vous transformez vos deux fichiers Excel en format Micro-Checker, le premier pour les loci liés à l’X en A, le second pour toutes les données (et on ne regardera pas le résultat des loci liés au sexe en A). Lancez micro-Checker. N’oubliez pas de préciser le pas de mutation correct. Bon je vous aide, mis à part BX104 et GpCag (mononucléotides) et C102 (trinucléotides), tous les loci sont dinucléotidiques. Les analyses montrent que les allèlles nuls expliquent très bien tous les résultats, y compris pour les mâles. En effet, pour les femelles et loci autosomaux de A, et pour tous les loci en H, C et D, il y a plus de blancs observés qu’attendus par la méthode de Brookfield. Pour les loci liés à l’X chez les mâles du site A, les différences ne sont pas significatives. Il semble même y avoir du « stuttering » pour le locus C102 en D. Cependant, l’effet Wahlund ne peut non plus être totalement écarté, ainsi que nous allons pouvoir le vérifier. Vous pourrez aussi vérifier qu’il ne semble pas exister de dominance d’allèles courts ici.

Tableau 24
Résultat du test de Hardy-Weinberg sur le F_IS dans les différentes zones de capture des tsé-tsé, par locus et sur l’ensemble et résultat des tests de permutation

MISE EN ÉVIDENCE D’UNE SOUS-STRUCTURATION À L’INTÉRIEUR DES ZONES A, H, C ET D

10Vous pouvez tester en zone A s’il existe un biais de structuration spécifique au sexe entre pièges, à titre d’exercice, et constater qu’il n’y a aucune signature d’un tel phénomène dans ces données. Nous allons rechercher un possible effet Wahlund comme une cause possible d’excès d’homozygotie chez les tsét-tsé d’une même zone : d’abord en analysant le F_IS à une échelle plus réduite (piège), ensuite par analyse bayésienne de clusterisation comme pour les tiques et enfin en recherchant un isolement par la distance entre individus le long du cours d’eau.

Analyse par piège

11En prenant chaque piège comme une sous-population potentielle et en recalculant le F_IS, on obtient une valeur plus faible de 0,144, significativement inférieure à la précédente (test de Wilcoxon pour données appariées comme pour les tiques, P-value = 0,0391), mais toujours significativement supérieure à 0 (P-value = 0,0001). Il semble donc bien que chaque piège recèle, au moins en partie, des mouches plus apparentées que des mouches prises au hasard dans chaque zone. À cause de la faiblesse des échantillons, Micro-Checker ne peut être utilisé ici. Nous devons donc trouver une méthode alternative afin de rechercher si les allèles nuls peuvent contribuer à expliquer les déficits en hétérozygotes rencontrés. Une méthode pratique consiste à regarder s’il existe une relation entre le nombre de blancs par locus et le F_IS effectivement mesuré à ce locus, dans chaque sous-échantillon. Nous obtenons ainsi les données du tableau 25.

12On lance ensuite une analyse de corrélation. Pour plus de sécurité, on utilisera une analyse dite non paramétrique à l’aide du coefficient de corrélation de Spearman. Dans R, les commandes seront (en respectant les majuscules et minuscules, test unilatéral car on a un préjugé de la direction du signal) :

13Le résultat est un coefficient de corrélation de Spearman ρ = 0,46 très significatif (P-value = 0,0073) (fig. 55).

14Vous remarquerez que le coefficient de corrélation est légèrement différent de celui publié dans l’article de Molecular Ecology (ρ = 0,499 et P-value = 0,0048). La différence provient de trois pièges de l’échantillon en zone A (le seul où il y avait des mâles) où le F_IS est différent. Cela provient certainement du recodage des mâles ou plus probablement du fait que je n’ai éliminé aucun sous-échantillon ici, même ceux de taille 1. De toutes manières, cela ne change pratiquement rien. Ce genre de petits problèmes est fréquent quand le nombre d’analyses différentes à effectuer est très grand, comme cela a été le cas ici. C’est pour cela que j’ai choisi d’en parler, car cela arrive et il ne faut pas le cacher. Ce genre d’erreurs (assimilables à celles éventuelles associées au génotypage/sexage, etc.), inévitables à la longue, n’est cependant pas en mesure de générer un signal quelconque, et va plutôt contribuer à masquer les signaux de faibles amplitudes. Ici, les allèles nuls expliquent donc bien en partie les F_IS. En mettant au carré le coefficient de corrélation trouvé, on réalise qu’environ 21 % seulement de la variance de ces derniers est expliquée par ce phénomène (16 % si on utilise le coefficient de détermination normal). Il est donc raisonnable de rechercher d’autres facteurs responsables de ces déficits en hétérozygotes.

Tableau 25
F_IS et nombre de blancs (homozygotes nuls supposés) par piège-site et par locus

Clusters BAPS

15Ici, en ce qui concerne mon analyse, les clusters trouvés par BAPS semblent expliquer une très grande partie du déficit en hétérozygotes, voire la totalité, puisque nous passons d’un F_IS = 0,175 à un F_IS = 0,031 non significativement différent de 0 cette fois (~ panmixie locale). Il semble donc bien que l’effet Wahlund soit responsable de la plus grande part du déficit en hétérozygotes. Il semble aussi que les pièges eux-mêmes capturent des mouches issues de voisinages différents puisque le F_IS intra-piège, même s’il baisse, reste fortement positif. La correspondance entre les clusters BAPS et pièges est à cet égard mauvaise (vérifiez-le), même s’il arrive fréquemment que des mouches du même piège se retrouvent dans le même cluster BAPS. Les allèles nuls jouent peut-être, quant à eux, un petit rôle également, comme semblent le montrer les analyses de MicroChecker et de corrélation avec le nombre de blancs. Cependant, une régression du nombre de blancs trouvés dans les différentes zones n’explique que peu la dispersion des F_IS par loci et piège-zone (R² = 0,16). Cependant, en cas de pangamie, c’est un F_IS légèrement négatif qui est attendu. Le F_IS des clusters de BAPS étant légèrement positif, il est possible que la contribution des allèles nuls, même modeste, soit réelle. Mais c’est bien l’effet Wahlund qui explique le mieux les données.

Isolement par la distance entre individus

16La plupart de ces pièges contiennent trop peu d’individus génotypés pour mettre en œuvre la même procédure que pour les tiques. Par ailleurs, nous savons que les pièges ne représentent qu’approximativement des voisinages (si voisinage il y a) puisque pièges et clusters BAPS ne sont pas en très bon accord. Nous pouvons cependant utiliser ici la procédure d’isolement par la distance entre individus (Rousset, 2000 ; Watts et al., 2007) implémentée par le logiciel Genepop 4 (Rousset, 2008) téléchargeable à partir du site http://kimura.univ-montp2.fr/~ rousset/Genepop.htm. Il faut recoder les données pour chaque zone (un fichier par zone A, H, C et D) de telle sorte que chaque individu est considéré comme une sous-population comme dans la figure 56.

17Il faut ensuite copier le logiciel Genepop. exe dans le répertoire où se trouvent les quatre fichiers que nous venons de créer. En ce qui me concerne, et n’écoutant que mon imagination débordante, j’ai nommé les quatre fichiers A. txt, H. txt, C. txt et D. txt. Ici, si vous êtes cardiaque, je vous conseille de couper le son de votre ordinateur, car quand Genepop est content il le manifeste par un bruit effroyable à réveiller un mort. On clique deux fois sur Genepop. exe et une fenêtre apparaît où le nom du fichier vous est demandé. À l’invite, tapez “A. txt” puis “Entrée”. Lisez les informations et si vous êtes d’accord retapez “Entrée”. Un menu apparaît. C’est l’option 6 qui nous intéresse. Tapez donc “6”. Un sous-menu apparaît dont l’option 5 est celle qu’il faut implémenter. Tapez “5”. On vous demande si vous souhaitez effectuer le test avec la statistique â (un équivalent du F_ST/(1 - F_ST) pour la différenciation entre individus) ou ê. D’après Watts et al. (2007), pour une structure en une dimension, comme c’est le cas le long de la rivière Mouhoun, la statistique ê est meilleure quand la taille de voisinage (Nb pour neighbourhood) Nb = 4Dσ² > 10 000 individus et â est plus performant quand Nb < 10 000 individus. Commençons par â et nous prendrons ê ensuite. Tapez “a”. On vous demande si vous souhaitez faire le test avec le logarithme népérien des distances géographiques ou non. Tapez “d” car nous sommes dans un contexte unidimensionnel (cf. p. 89 en première partie). On vous demande la distance minimale à considérer pour la régression. Comme le test n’en tiendra pas compte, que le biais ne risque pas d’être important (Watts et al., 2007, voir aussi le commentaire de Rousset dans la documentation de Genepop 4) et qu’il n’y a pas de log, tapez “0”. Le nombre de randomisations à effectuer pour le test de Mantel vous est demandé. Tapez “1 000 000”. En fonction de l’ordinateur le processus Markovien prend plus ou moins de temps. Le programme vous demande te taper “Return” (soit “Entrée”). Le résultat est disponible dans A.txt.ISO. On fait de même avec H, C et D. Pour changer de fichier de données, il faut taper “C” dans le menu général de Genepop. Ne soyez pas étonnés si, à partir de H, le test de Mantel démarre sans vous demander votre avis. C’est comme ça. C’est Genepop. Vous vous apercevez que le calcul ne se fait pas pour D. En fait, cela ne se termine jamais, car il y a un problème dans le fichier et un bug dans Genepop. Ouvrez D. txt. Il faut supprimer le 12^e individu (000000 partout), sauvez puis recommencez, ça marche ! Ensuite, on enregistre les quatre fichiers de données sous un autre nom pour les analyses avec ê, par exemple A_e.txt, H_e.txt, C_e.txt et D_e.txt (quelle imagination ! mais où va-t-il les chercher ?). Nous nous retrouvons donc avec huit fichiers *.ISO que nous pouvons ouvrir avec n’importe quel éditeur de texte. Vous pouvez aussi remettre le son. Nous souhaitons savoir si 4Dσ² > 10 000. D’après ce que nous avons vu en p. 89 de la première partie de ce manuel, le voisinage est égal à Nb = 1/b = 4Dσ². Nous souhaitons vérifier si Nb > 10 000 afin de décider si c’est le paramètre â ou ê qu’il vaut mieux utiliser. C’est le cas uniquement pour le site A avec la statistique â. En outre, vous remarquerez que la statistique â donne de bien meilleurs résultats avec ces données de glossines. On sait par ailleurs que ce type de tests est très conservateur et que â n’est pas biaisé alors que ê l’est (Watts et al., 2007). Nous ne considérerons donc que les résultats obtenus avec â.

18Les résultats pour les quatre zones et la moyenne sur l’ensemble figurent dans le tableau 26. Il y a donc bien un isolement par la distance, mais les pentes sont très faibles. Cela signifie que les voisinages sont très lâches (beaucoup d’échange entre voisins, σ grand) et/ou de grande taille (D grand). Pour visualiser cette relation, nous pouvons utiliser les sorties *.GRA de genepop qui contiennent deux colonnes, la première avec les distances géographiques et la seconde avec la distance génétique a. Nous pouvons les charger sous Excel en précisant que les colonnes sont délimitées par des espaces et tracer le graphique de la figure 57.

19Nous avons maintenant besoin d’estimer des effectifs efficaces pour essayer d’obtenir une idée (mais ce sera à la louche) des densités.

Effectifs efficaces

20Les seules méthodes disponibles ici sont celles basées sur l’hétérozygotie et les déséquilibres de liaison, où nous serons obligés de considérer l’absence d’allèles nuls et d’effet Wahlund. Ces phénomènes (que nous savons probables pour ces échantillons de mouches tsé-tsé) auront tendance à produire des surestimations de N_e pour les méthodes basées sur l’hétérozygotie, et des sous-estimations pour les méthodes basées sur les déséquilibres de liaison. Nous utiliserons trois méthodes. La méthode d’identité intra et inter locus de Vitalis et Couvet (Vitalis et Couvet, 2001a, b, c) est implémentée par le logiciel Estim qui accepte le format Genepop pour le fichier de données. La méthode des déséquilibres de liaison de Bartley et al. (1992) est modifiée comme décrit dans l’aide du logiciel NeEstimator. La méthode des excès en hétérozygotes se fait très simplement en utilisant l’estimateur de Weir et Cockerham du F_IS dans la formule N_e = 1/(- 2F_IS) - F_IS/(1 + F_IS) (Balloux, 2004) qui ne donne bien entendu un résultat valide que si F_IS < 0.

21Pour ce faire, les effectifs par piège étant bien insuffisants, nous allons devoir regrouper les mouches de différents pièges en fonction de leur proximité (voir fig. 52). Ceci ne va pas arranger l’effet Wahlund, mais nous n’avons pas le choix. Ces regroupements sont synthétisés dans le tableau 27 (trois premières colonnes). Les pièges isolés ne contenant qu’une seule mouche ne sont pas représentés dans ce tableau. Vous pouvez essayer avec une autre stratégie de regroupement pour vérifier si on retrouve des valeurs équivalentes. Pour Estim, il faut recoder les données de départ au format Genepop, avec données manquantes pour les mâles aux loci hétérosomaux, et les regroupements du tableau 27. Pour NeEstimator (déséquilibres de liaison), il faut autant de fichiers qu’il y a de groupes de pièges définis dans le tableau 27 avec les mâles codés homozygotes pour les loci liés à l’X. Lesmême fichiers que pour Estim pourront être utilisés pour l’estimation des F_IS par groupe de pièges, soit en utilisant Genepop directement, soit en traduisant les fichiers pour un autre logiciel (Fstat, Genetix). Notez que NeEstimator donne aussi une estimation basée sur les excès d’hétérozygotes (Luikart et Cornuet, 1999), mais contenant des inexactitudes corrigées par la méthode de Balloux. On peut aussi utiliser le fichier recodé Genepop pour une analyse par LDNe (Waples et Do, 2008), qui implémente une méthode basée sur les déséquilibres de liaison non biaisée (ou beaucoup moins) pour les petits échantillons, alors qu’on sait que la méthode de Bartley est biaisée quand la taille des échantillons est inférieure à la taille efficace des populations étudiées (England et al., 2006 ; Waples, 2006). Ceux qui s’en rappellent constateront que les méthodes implémentées par Estim et celles basées sur les excès d’hétérozygotes n’avaient pas été utilisées pour les tiques (p. 160-164 de la seconde partie). Chez les tiques, la forte présence d’allèles nuls en plus de la dominance des allèles courts au locus IR27 rendaient caduque toute approche basée sur les corrélations d’allèles intra-individuelles. Ici, il n’y a pas de dominance d’allèles courts et les allèles nuls sont peu influents, même si on ne peut totalement exclure leur impact (voir plus haut).

Tableau 27
Stratégie de regroupements par piège de Gpg le long du Mouhoun, en se basant sur la figure 52 et estimation des effectifs efficaces. Les résultats sont donnés pour les trois méthodes utilisées pour des résultats autres que l’infini, 0 ou NA (not available) (cases vides). Le nombre d’individus génotypés par piège est donné (N _Traps)

Tableau 27 (suite)

22Le logiciel Estim (http://www.ecoanthropologie.cnrs.fr/spip.php?article296) utilise un fichier au format Genepop. Dans la mesure où Estim utilise les identités intra-individuelles, interindividuelles, inter-échantillons et leur corrélation entre loci, et que par ailleurs l’hypothèse d’un modèle en îles est faite, il est clair que nous ne sommes pas tout à fait dans les critères orthodoxes de cette méthode. Il vaut mieux considérer chaque zone (A, H, C, D) séparément, car cela influence les résultats (comme vous pourrez le vérifier). Notez que la stratégie de regroupement diffère quelque peu de celle de l’article de Bouyer et al. (2009) avec des résultats légèrement différents. On peut donc charger le fichier contenant tous les groupes de pièges de la zone A dans Estim. Mon fichier s’appelle “TsetseMouhounAllMalManqNearestPooledA.gen”. N’oubliez pas de supprimer les pièges isolés ne contenant qu’une mouche, car Estim ne va pas apprécier. On lance donc Estim et on charge son fichier. Pour qu’il apparaisse, on tape *.gen dans la case appropriée comme dans la figure 58 ou alors on change l’extension du fichier de. genà. txt et on clique dessus deux fois.

23Ensuite, on appelle la commande “Identity measures” du menu “Analysis” (fig. 59). On obtient alors une fenêtre résultat dont on fait descendre le curseur pour pouvoir enregistrer (cliquer sur “Save”) (fig. 60). Je l’ai enregistrée sous le nom “NeEstimA.txt” Cliquez ensuite sur l’option “Ne inferences” du menu “Analysis” et sauver en gardant le même nom, car ces nouveaux résultats sont écrits à la fin de la fenêtre précédente. Vous pouvez ouvrir le fichier résultat avec un éditeur de texte. Recommencez la même opération pour chacune des zones restantes. Les résultats sont que seules deux collections de pièges donnent des valeurs exploitables : le piège 8 de la zone A (N_e = 3,19, m = 0,27) et le groupe de pièges (16, 17) de la zone H (N_e = 2,08, m = 0,55) (tabl. 27). La méthode de Waples et Do, pour laquelle il suffit de charger le fichier Genepop avec toutes les données en appuyant sur le bouton “Search”, ne donne aucun résultat ici (pas d’estimation possible) comme c’est très souvent le cas, mais la plupart des limites inférieures paramétriques disponibles indiquent de très faibles valeurs de N_e. Néanmoins, et sans garde-fou solide pour la méthode implémentée, il faudra bien garder à l’esprit qu’on aura peut-être ici des valeurs très sous-estimées par la méthode des déséquilibres de liaison de Bartley. La méthode des déséquilibres de liaison de Bartley est implémentée par NeEstimator comme pour les tiques (un fichier par groupe de piège ici) (voir p. 160-164 dans la seconde partie de ce manuel). Enfin, la méthode de Balloux est très simple puisqu’il suffit de calculer les F_IS par groupe de pièges (avec Fstat, par exemple) et d’utiliser la formule N_e = 1/(- 2F_IS) - F_IS/(1 + F_IS) et de ne garder que les valeurs de N_e positives. Tous les résultats sont compilés dans le tableau 27. Nous avons maintenant besoin de transformer ces effectifs en densités.

Densités efficaces

24Nous allons utiliser une autre méthode que celle de Bouyer et al. (2009), pour changer. Nous allons simplement pour chaque méthode prendre l’effectif efficace moyen par piège (pondéré pour le nombre d’individus et de pièges) et diviser cette valeur par la distance minimale (en m) entre deux pièges (tabl. 28). Cette distance minimale se trouve en zone A (facile à trouver dans le fichier que vous avez créé pour la figure 56) et est de D_mini = 71 m. Le calcul de pondération est assez particulier. Pour ce faire, j’ai multiplié le N_e par le nombre d’individus N_ind capturés dans les pièges correspondants : N_eP = N_e × N_ind. Pour chaque N_e, j’ai calculé le produit du nombre de mouches par le nombre de pièges correspondant n_pièges : N_P = N_ind × n_pièges. J’ai ensuite fait la somme des N_eP = ∑N_eP et des N_P = ∑N_P. L’effectif efficace moyen est ensuite calculé par le rapport de ces deux valeurs.

25La densité D_c est ensuite obtenue en divisant cette valeur par D_mini = 71, ce qui permet ensuite de déduire les dispersions σ à partir des valeurs de Dσ² du tableau 26 :

26On comprend bien que les valeurs obtenues (tabl. 28) ne pourront être que très approximatives.

27Par conséquent, à partir des données génétiques et de leur analyse (isolement par la distance) et de calculs de densités efficaces, nous pouvons inférer que le long du Mouhoun les tsé-tsé ont des densités comprises entre 12 et 176 mouches par km et une dispersion (distance entre adultes reproducteurs et leurs parents) comprise entre 131 m et 1 620 m. Ces inférences sont remarquablement convergentes avec celles du papier de Bouyer et al. (2009) (tabl. 29) et donc avec les données issues de marquage-relâchage et recapture (MRR) de tsé-tsé marquées en zone A. Dans la mesure où la stratégie de regroupement fut ici légèrement différente, de même que le choix de calcul des densités, ce résultat confirme la robustesse des résultats du papier. La convergence avec les données de marquage-recapture souligne également l’efficacité des outils de génétique des populations, en particulier la méthode de Rousset (1997) pour estimer Dσ².

Tableau 28
Calcul détaillé des densités (en mouches par m) et de la dispersion des glossines (distance entre adultes reproducteurs et leurs parents en m) le long du Mouhoun

Conclusions : isolement par la distance intra-zone (rolling on the river)

28Nous avons bien mis en évidence une sous-structure au sein des zones A, H, C et D. Le fait que les déficits en hétérozygotes persistent au sein de chaque piège, auquel s’ajoute la non-correspondance parfaite entre clusters BAPS et pièges alors qu’un isolement par la distance existe bel et bien, plaide pour deux interprétations complémentaires. Il semble bien y avoir quelques allèles nuls, mais ces derniers n’expliquent qu’une faible partie des déficits en hétérozygotes observés. L’effet Wahlund explique probablement la majeure partie des déficits. Il provient de deux causes. La première est inhérente aux systèmes d’isolement par la distance, d’une nature plus ou moins continue, et de la nature nécessairement discrète du piégeage des tsé-tsé. La seconde raison, qui dépend de la première, provient de la mauvaise correspondance entre dispersion trophique, plus large, et dispersion reproductrice (accouplements et larvipositions) plus restreinte (homing). Cette information est capitale si nous parvenons un jour à déterminer avec précision les micro-conditions écologiques qui poussent les tsé-tsé à revenir se reproduire et larviposer à l’endroit où elles ont émergé. Il reste aussi à déterminer quelle influence la densité (compétition) a sur la dispersion de reproduction afin d’évaluer si nos estimations restent valables dans le cadre de campagnes de contrôle et/ou d’élimination.

DIFFÉRENTIATION ENTRE LES QUATRE ZONES

Analyse HierFstat du jeu de données total partitionné par BAPS

29Nous savons qu’une différenciation existe bien à une mini (voire micro) échelle à l’intérieur de chacune des zones A, H, C et D. Nous devons donc tenir compte de ce niveau de structuration à micro-échelle avant d’estimer et tester l’existence d’une différenciation entre zones. Cette information est utile, car elle pourrait permettre d’estimer le temps nécessaire à une recolonisation d’une zone éliminée par la zone la plus proche. Nous allons devoir utiliser HierFstat une nouvelle fois. Considérant que les regroupements définis par BAPS pourraient mieux regrouper les individus de la même unité populationnelle par rapport aux pièges et pour faire autre chose que dans le papier initial, nous allons prendre comme niveau le plus imbriqué les clusters BAPS de chaque zone définis en p. 199. Le niveau suivant sera la zone (A, H, C, D) et enfin la totalité. Il y aura ainsi quatre niveaux définis avec leur F, l’individu (F_IS), le sous-groupe défini par BAPS dans la zone (F_SZ), la zone dans le tout (F_ZT), auxquels s’ajoutent bien sûr les F_IZ, F_IT et F_ST moins intéressants pour nous.

30En procédant comme pour les tiques (voir p. 156 dans la seconde partie de ce manuel), et en prenant soin de recoder les mâles homozygotes pour les locus hétérosomiques, on obtient une forte valeur pour F_SZ ≈ 0,22 et une valeur négative pour F_ZT ≈ - 0,03. Il semble que toute l’information soit contenue à l’intérieur des zones et qu’il ne reste plus assez de variation pour distinguer les zones entre elles. Le « supplementary information » de Rougeron et al. (2009), présenté ci-dessous, permet de mieux comprendre ce problème inhérent aux statistiques F hiérarchiques et renforcé par l’homoplasie des microsatellites.

Comprendre le manque de structure inter-zones avec un peu de théorie

31Quand nous avons, comme c’est le cas ici, quatre niveaux hiérarchiques (individus, sous-populations, archipels et totalité), quatre paramètres d’identité peuvent être définis : Q_i, la probabilité que deux allèles d’un locus d’un individu pris au hasard soient identiques ; Q_s, la probabilité que deux allèles à un locus, de deux individus pris au hasard dans la même sous-population soient identiques ; Q_a, la probabilité que deux allèles à un locus, de deux individus pris au hasard dans deux sous-populations différentes dans un même archipel soient identiques ; et Q_T, la probabilité que deux allèles à un locus, de deux individus pris au hasard dans deux sous-populations différentes et deux archipels différents soient identiques (cf. fig. 61).

32Nous pouvons définir six indices de fixation : F_IS (consanguinité individuelle relative à celle des sous-populations), F_SA (consanguinité des sous-populations relative à celle des archipels), F_IA = 1-(1 - F_IS)(1 - F_SA) (consanguinité individuelle relative à celle des archipels), F_AT (consanguinité des archipels relative à celle de la population totale), F_ST = 1-(1 - F_SA)(1 - F_AT) (consanguinité des sous-populations relative au total) et F_IT = 1-(1 - F_IS)(1 - F_ST) (consanguinité des individus relative à la population totale). Ces indices peuvent être exprimés, en suivant la méthode proposée par Cockerham (1969, 1973), en fonction des probabilités d’identité définies plus haut dans ce paragraphe (on peut aussi consulter les p. 40-48 de la première partie de ce manuel) :

33Si nous nous concentrons maintenant sur les indices de fixation qui reflètent la différenciation génétique entre sous-populations du même archipel et entre archipels, il n’y a alors plus que F_SA et F_AT qui nous intéressent. Si pour une raison quelconque, Q_a est très petit (migration très faible entre sous-populations), il est alors facile de voir par l’équation (66) que F_SA sera très grand (~ Q_s si Q_a ~ 0). Silamigration est très faible entre sous-populations d’un même archipel, il est alors probable que celle entre archipels soit au moins aussi faible et donc que Q_T ~ 0 et F_AT ~ Q_a. À partir de là, il est facile de voir que, quand la différenciation est extrême entre les sous-populations celle-ci sera nécessairement faible (en apparence) entre archipels. En fait, cela veut juste dire que la différenciation entre sous-populations est très forte, que ce soit entre sous-populations du même archipel ou de deux archipels différents, et la distinction entre archipels n’apporte pas suffisamment d’information avec ces outils. Ce dernier point peut être illustré mathématiquement par le fait que dans ce cas F_ST ~ F_SA.

34C’est donc probablement ce phénomène qui empêche partiellement de détecter un quelconque signal entre zones chez les glossines du Mouhoun. Le fait que la zone soit perturbée et que donc l’isolement puisse être récent entre les différentes zones peut également contribuer à brouiller l’image. En effet, alors que l’isolement par la distance est un phénomène qui se met très vite en place et devient détectable en quelques générations, comme le montrent nos simulations (Bouyer et al., 2009), le F_ST met un certain nombre de générations à atteindre l’équilibre migration, mutation, dérive. Reprenons ces simulations.

Comprendre le manque de structure inter-zones avec un peu de simulations

35Pour effectuer ces simulations, nous allons utiliser Easypop v 2.0.1 (Balloux, 2006, mise à jour de Balloux, 2001). Le problème avec les simulations, c’est de choisir un jeu de paramètres pertinent, car il y a une infinité de combinaisons possibles. Pour limiter notre travail, nous allons utiliser encore une fois les résultats de l’article de Rousset (1997) où l’on peut lire que dans un dispositif en une dimension on peut démontrer que, si N est le nombre d’individus d’un dème, m la proportion de migrants, D_e la densité efficace d’individus par km², σ la dispersion efficace (distance entre individus reproducteurs et leurs parents) et ε la distance entre deux dèmes adjacents :

36Nmε = D_eσ² (67)

37La distance entre deux dèmes sera :

38(68)

39En explorant les possibles (qui collent le mieux aux données), on peut obtenir D_eσ² = 700, N = 30 et m = 0,5 et donc ε = 50. Ce qui voudrait dire que la distance entre deux dèmes serait de l’ordre de 50 m. On va donc supposer que nos estimations étaient les meilleures dans la zone A et que nous y avons sous-estimé la taille des sous-populations. La distance entre les zones A et H est d’environ 70 km. Ces zones font environ 3 km chacune et nous allons simuler deux zones de même nature de 3 000 m, soit 3 000/50 (ε = 50, distance entre deux dèmes), 60 sous-populations chacune environ, séparées de 70 km, soit 1 400 sous-populations environ. Nous avons donc besoin de simuler 1 520 populations de taille 30 et échangeant 0,5 proportion de migrants dans un « stepping-stone » en une dimension.

40Il faut maintenant lancer la simulation. Copiez Easypop dans le répertoire où vous souhaitez travailler, et double-cliquez dessus (sur le fichier programme pas sur le répertoire). Il faut ensuite répondre à toutes les questions. Vous souhaitez simuler des diploïdes à sexes séparés (dioïques) qui se croisent au hasard (on ne va pas se compliquer la vie). Nous voulons 1 520 populations de 30 individus avec un sexe-ratio équilibré (je dis bien UN sexe-ratio, car sexe est masculin en français et ratio de même en latin, quoiqu’en disent de nombreux écologues mal instruits), soit 15 femelles et 15 mâles. Vous souhaitez simuler un « steppingstone » en une dimension tout au long de la simulation avec un taux de migration de 0,5 pour les deux sexes. Vous allez simuler 10 (plus rond que 7) loci indépendants avec un modèle de mutation KAM, 99 allèles possibles et un taux de mutation (le même pour tous) de 0,0001 qui correspond à un taux raisonnable, mais vous pourrez essayer avec 10^-3 (consultez Ellegren, 2000 ; Balloux et Lugon-Moulin, 2002 ; Ellegren, 2004). Nous allons commencer avec une variabilité maximale (99 allèles équitablement répartis dans les 100 sous-populations), car cela fait gagner du temps (démarrer avec un seul allèle requiert un nombre important de générations avant d’obtenir quelque chose d’utilisable). Nous allons simuler 1 000 générations et récupérer tous les individus des 1 520 sous-populations. Nous ne voudrons pas connaître l’ascendance de nos individus. Nommez les fichiers résultats comme bon vous semble et ne demandez qu’une réplication (cela suffira ici). N’oubliez pas de valider chacun de vos choix par un retour chariot, sinon vous risquez d’attendre longtemps. En fonction de la puissance de votre ordinateur, la simulation durera plus ou moins longtemps (11 minutes avec ma double CPU 2.2 GHz avec 3.5 Go de RAM). Quand la simulation est terminée, Easypop vous demande le nom du fichier de sauvegarde des paramètres de la simulation (très utile !) et crée trois fichiers résultats : le détail de l’évolution de la simulation au cours des générations au niveau de divers paramètres (nombre d’allèles, diversité génétique, F-statistiques de Wright, etc.) (*.equ), un fichier de données Fstat (*.dat) et un fichier au format Genepop (*.gen). Il va falloir créer des fichiers pour tester l’effet « sous-structuration », l’effet Wahlund et l’interaction entre les deux. Pour l’effet sous-structuration, on ouvre le fichier. dat avec un bon éditeur de texte et on ne garde que six sous-populations dans les deux zones extrêmes « A-like », c’est-à-dire qu’on ne garde que les sous-populations 5-15-25-35-45-55 et 1 465-1 475-1 485-1 495-1 505-1 515 que l’on recode de 1 à 12 en gardant bien à l’esprit que 1-6 = A1 et 7-12 = A2 (les deux zones extrêmes). Les populations marginales sont à éviter13, c’est pourquoi on exclut les sous-populations 1 et 1 520. On crée un fichier de type HierFstat avec en première colonne la zone (1 et 2), en deuxième la sous-population (1 à 12) et en colonnes 3 à 12 les génotypes aux 10 loci. On lance R et on charge le package hierfstat et on se met dans le répertoire où on a créé ce fichier que j’ai personnellement appelé MouhounA-001HierFstat. txt. Puis après les commandes habituelles :

41on obtient le résultat suivant :

42Nous pouvons remarquer que le F_IS est négatif, ce qui est normal puisque nous avons des sexes séparés. Ensuite, on a un fort F_SZ et un F_ZT beaucoup moins fort, comme prévu, mais cependant beaucoup plus important que celui observé entre zones pour nos données réelles de tsé-tsé. Nous allons refaire un fichier HierFstat en prenant deux zones adjacentes : sous-populations 5-15-25-35-45-55 (Zone 1) et 65-75-85-95-105-115 (Zone 2).

43On diminue certes de beaucoup le F_ZT, mais pas autant que celui avec les sous-groupes définis par BAPS. Les clusters BAPS ne sont donc sans doute pas très au point. Il reste d’ailleurs un F_IS > 0 dans les clusters au lieu d’un excès. Il se peut donc que d’une part un effet Wahlund, combiné à des allèles nuls ainsi qu’à un isolement encore trop récent entre zones, expliquent le manque de différenciation observé entre zones, alors que l’on sait que les passages de mouches d’une zone à l’autre sont quasi impossibles (en principe). En fait, une analyse HierFstat en prenant les pièges comme sous-populations donne :

44La différenciation entre sous-échantillons (pièges) est fortement diminuée par rapport à l’attendu selon le modèle d’isolement par la distance, la différenciation entre zones semble plus substantielle alors que la corrélation intra-individuelle (F_IS) est forte en raison d’un effet Wahlund et des allèles nuls et dans une moindre mesure du codage homozygote des mâles (peu nombreux en Zone A) aux loci hétérosomaux (dans les pièges le vrai F_IS = 0,14, voir plus haut).

CONCLUSIONS

45Il y a une forte micro-structuration que le maillage des pièges ne permet pas de rendre avec précision, en particulier en ce qui concerne la taille des dèmes et leur distance entre eux, à moins que le voisinage soit strictement continu. Il n’en reste pas moins qu’un isolement par la distance a pu être mis en évidence et que l’on sait que la distance entre deux pièges doit être réduite, si l’on souhaite affiner notre vision de la structure des populations de G. palpalis gambiensis le long du Mouhoun (soit moins que 71 m entre pièges). Les estimations de densités et de dispersions convergent avec celles des données MRR, ce qui incite à une certaine confiance malgré les effets Wahlund, les petits échantillons et les allèles nuls qui ont très certainement considérablement brouillé les signaux. À l’avenir, pour les tsé-tsé de forêt galerie (le Mouhoun fait actuellement l’objet d’une campagne d’éradication par le PATTEC et ne sera donc pas étudiable avant longtemps), des pièges distants de 20 m, et le génotypage de davantage d’individus par piège sur des loci de meilleure qualité devraient permettre des estimations beaucoup plus précises encore. En particulier, la différenciation entre zones qui est apparue très incertaine pourra ainsi davantage être précisée et, par conséquent, une probabilité de dispersion entre ces zones estimée plus clairement.