Version classiqueVersion mobile

Initiation à la génétique des populations naturelles

 | 
Thierry De Meeûs

1. Concepts théoriques et statistiques

4. La tique Ixodes ricinus et les pathogènes (Borrelia sp.) qu’elle transmet

Texte intégral

INTRODUCTION

1Ce jeu de données, publié dans trois articles (De Meeûs et al., 2002a, 2004a, 2004b), représente un excellent exercice, car nous allons y rencontrer bon nombre de situations décrites dans le chapitre précédent. Nous allons entièrement décortiquer une nouvelle fois ce jeu de données avec les mêmes méthodes, mais aussi avec des outils plus récents que ceux qui avaient été utilisés à l’époque, ce qui sera aussi intéressant. Nous repartirons de zéro en feignant d’ignorer ce qui a déjà été fait, comme s’il s’agissait d’un jeu de données non analysé. Le jeu de données complet est téléchargeable sur mon site web.

ÉTAT DES LIEUX

2Les tiques sont des acariens hématophages qui, au cours de leur repas sanguin, peuvent transmettre des maladies à leurs hôtes vertébrés. Dans l’hémisphère nord, ce sont elles qui sont responsables de la très grande majorité des maladies à vecteur des humains et, en particulier, de la transmission de la maladie de Lyme dont l’impact économique et en santé publique est reconnu (Gubler, 1998). Encore aujourd’hui, beaucoup reste à faire pour mieux comprendre l’épidémiologie de cette maladie et la variabilité des manifestations cliniques qui la caractérise (Hubbard et al., 1998). Les tiques sont typiquement des organismes difficiles à suivre sur le terrain, et des approches par marqueur moléculaire semblent donc pertinentes dans ce cas de figure. Après une tentative peu fructueuse avec les allozymes, avec seulement deux loci peu polymorphes (Delaye et al., 1997), des microsatellites ont été développés (Delaye et al., 1998). Seuls cinq loci polymorphes avaient pu être mis au point à l’époque, ce qui était vraiment peu. Nous allons ensemble voir que, malgré cela et les problèmes rencontrés, on peut quand même recueillir beaucoup d’informations pertinentes à l’aide des méthodes décrites dans ce manuel.

3En téléchargeant le fichier “IRTotBrut.txt”, vous aurez les données brutes obtenues sur des tiques adultes échantillonnées sur la végétation (donc non gorgées), sauf pour la Tunisie où les tiques étaient fixées sur des vaches. Le fichier “IRTotBrut. txt” est un fichier texte mais que l’on peut ouvrir sous Excel si on le souhaite. Le tableau 6 donne un extrait du fichier de données brutes. Le fichier comprend neuf colonnes. La première colonne donne le nom des sites où les tiques ont été échantillonnées. Il y a huit sites en Suisse (fig. 15) et un site en Tunisie. La deuxième colonne correspond à l’année d’échantillonnage, car certains sites ont été prélevés aux printemps 1995 et 1996 et d’autres uniquement au printemps 1996. La troisième colonne correspond au sexe de la tique (F pour femelle et M pour mâle). La quatrième colonne donne le nom codé des différents individus tiques. Ce codage individuel peut être utile si on fait des analyses individus centrées telle qu’une AFC ou une construction d’arbre sur distances interindividuelles. Enfin, les cinq dernières colonnes correspondent aux génotypes (en taille d’allèles) aux cinq loci microsatellites polymorphes définis dans Delaye et al. (1998).

Image 100000000000021F0000024660664B5A.jpg

Figure 15. Localisation des sites d’échantillonnage des tiques Ixodes ricinus en Suisse et abréviations du nom des sites. Les sites marqués avec un astérisque ont été échantillonnés en 1995 et 1996

Tableau 6. Extrait du fichier de données IRTotBrut.txt

Image 10000000000005EE000007B67D2DADE5.jpg

PREMIER RECODAGE DES DONNÉES

4Certains logiciels n’aiment pas les noms longs et encore moins les accents ou autres signes cabalistiques. Par ailleurs, il est plus commode pour la lisibilité que tous les noms d’un même niveau aient le même nombre de caractères (alignement des colonnes). C’est pourquoi j’ai choisi de recoder dans IRTotBrut1.txt le nom des sites qui a été raccourci. Dans les données initiales, certains individus sont apparus avec trois ou quatre bandes à certains loci. Nous avions codé ces génotypes 333000 et 444000 pour les génotypes à trois et quatre bandes respectivement. Il convient de recoder ces données en données manquantes (000000). Nous reviendrons sur ces génotypes bizarres un peu plus tard, car ils s’avéreront utiles pour discuter des résultats des analyses de pedigrees. Nous allons procéder à une première analyse avec tous les échantillons afin de tester la panmixie locale et les déséquilibres de liaison entre loci. Nous allons pour ce faire créer un nouveau fichier où les sites et les dates seront distingués, mais aussi le sexe des tiques car on ne sait jamais à l’avance si des différences peuvent exister entre les deux sexes (Prugnolle et De Meeûs, 2002 ; Prugnolle et al., 2003), auquel cas les résultats obtenus pourraient s’en ressentir, mais surtout la discussion serait réorientée. Donc autant distinguer le sexe des individus dès le départ, quitte à ignorer ce facteur par la suite si on ne voit rien. Nous allons nommer ce fichier “IRTotTestPanmix.dat” et le mettre au format Fstat qu’il faut donc télécharger et ouvrir pour voir comment constituer un fichier à ce format. Vous pourrez aussi créer un fichier contenant le nom des sous-échantillons “IRTotTestPanmix.lab”, car un fichier de données Fstat ne contient que des chiffres. Ce fichier est constitué d’une colonne avec le nom des sous-échantillons. Vous pourrez aussi coder les données au format CREATE (qui n’existait pas au moment de réanalyser ces données) et vous servir de ce logiciel pour convertir ce fichier au format approprié.

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

5Nous allons donc tester s’il existe des déficits en hétérozygotes et des déséquilibres de liaison. Pour ce faire, il faut ouvrir Fstat. Une fois dans Fstat, il faut ouvrir le fichier “IRTotTestPanmix.dat” et cocher les cases qui vont nous être utiles ici (fig. 16). Si vous souhaitez voir apparaître les noms des sous-échantillons, il faut le spécifier par le menu “Options” de Fstat (cf. le premier recodage des données du chapitre 2 de cette deuxième partie pour une prise en main pas à pas de Create).

Image 10000000000005D8000005174C628430.jpg

Figure 16. Capture d’écran de Fstat lors de la première analyse

6Nous n’effectuons pas d’autres analyses pour le moment, car ces dernières pourraient être remises en cause par les résultats obtenus ici.

7La procédure de test de déséquilibre de liaison est assez lente, donc, si vous souhaitez que votre analyse finisse avant l’âge de la retraite, il vaut mieux dans tous les cas s’en tenir à l’option 5/100 pour le “Nominal level for multiple testing”. Mon ordinateur portable, dont l’horloge à 2.13 GHz et la mémoire vive à 2 Go témoignent d’une performance somme toute raisonnable, a mis quand même quatre heures pour effectuer cette première analyse dont le résultat est consultable dans le fichier “IRTotTestPanmix.out”. Que pouvons-nous voir dans ce fichier ?

8Les premières lignes donnent les fréquences des allèles pour chaque locus et chaque sous-échantillon, ainsi que sur l’ensemble (moyennes pondérée, W, et non pondérée, UW). Nous pouvons constater à cette occasion que chaque locus, sauf IR27, possède un très grand nombre d’allèles dont la plupart ne suivent en rien le modèle de mutation attendu de deux pas par deux pas (ce sont tous des dinucléotides). Dans ce cas, la plupart des allèles proviennent de mutations intervenues en dehors du motif microsatellite, dans les séquences flanquantes. Ce n’est pas dramatique même si non idéal. Suivent les estimateurs de Nei, en particulier ceux des diversités géniques intra-sous-échantillons (Hs) et globale (HT). Ensuite, les résultats des tests de déséquilibre de liaison sont donnés par paire de loci et par sous-échantillon et sur l’ensemble des sous-échantillons (mais toujours par paire de loci). La mention “Adjusted P-value for 5 % nominal level is : 0,000208” ne doit pas vous inquiéter. C’est le calcul du seuil de Bonferroni sur l’ensemble des tests réalisés. Comme il y a 24 sous-échantillons, cinq loci et donc 5 (5 - 1)/2 paires de loci, cela correspond à 240 tests. Le seuil corrigé par la procédure de Bonferroni à α = 0,05 est donc α’= 0,05/240 = 0,000208, seuil rarement (jamais ?) accessible, ce qui illustre une discussion que nous avons déjà eue précédemment. De toutes manières, nous ne regarderons ici que les tests multi-sous-échantillons (colonne “All”) et donc au pire, le seuil est à diviser par 10, ce qui est inutile puisque nous pouvons aussi constater qu’aucun déséquilibre de liaison n’est significatif. Les loci sont donc raisonnablement indépendants statistiquement les uns des autres. Nous pouvons donc sereinement oublier ces derniers et passer à la suite.

9Suivent les estimateurs de Weir et Cockerham dont un seul nous intéresse pour le moment, f, l’estimateur du FIS, par locus, par allèle et sur l’ensemble des allèles, sur l’ensemble des loci. Puis suivent les résultats des jackknives et bootstraps et enfin des permutations. En compilant ces résultats dans le tableau 7 et la figure 17, nous constatons de très forts et très variables déficits en hétérozygotes (tous très significatifs avec des P-values toutes inférieures à 0,0001, visibles en fin de fichier).

Tableau 7. Valeurs moyennes de f, estimateur du FIS, par locus et intervalle de confiance tels que définis par Li et Ls (limite inférieure et supérieure) obtenus pour les microsatellites d’Ixodes ricinus. Pour chaque locus, Li et Ls sont calculées à l’aide de l’erreur standard (StdErrFis) donnée par le jackknife sur les populations et la valeur du t pour 23 ddl (24 - 1) et α = 0,05 (soit 2,069, voir le tableau 3) en suivant l’équation (45). Pour la valeur globale, l’intervalle de confiance est issu du bootstrap sur les loci

Image 10000000000005AA0000019C849422F6.jpg

Image 10000000000001E300000192864C6848.jpg

Figure 17. Valeurs moyennes de f, estimateur du FIS, par locus et intervalle de confiance obtenus pour les microsatellites d’Ixodes ricinus. Pour chaque locus, les intervalles de confiance sont calculés à l’aide de l’erreur standard donnée par le jackknife sur les populations et la valeur du t pour 23 ddl (24 – 1) et a = 0,05 (soit 2,069, voir le tableau 3) en suivant l’équation (45). Pour la valeur globale, l’intervalle de confiance est issu du bootstrap sur les loci

10Ces fortes valeurs sont aberrantes étant donné qu’on sait qu’I. ricinus pratique une reproduction bi-parentale obligatoire. Des croisements systématiques entre apparentés pourraient-ils expliquer un FIS = 0,5 ? Dans la réponse 11, on décrit comment obtenir une estimation grossière du taux de croisements frère-sœur b nécessaires pour expliquer un FIS donné :

Image 10000000000002D00000004384322BAA.jpg

11Par conséquent, nous avons besoin ici de 4/5, soit 80 % de croisements frère-sœur pour expliquer nos données, ce qui est possible mais semble peu réaliste. Ixodes ricinus est en effet une tique triphasique qui change d’hôte pour chaque stade. Les adultes dont nous analysons la variabilité génétique ont donc subi deux phases de dispersion par des hôtes différents. Pour permettre un taux de 80 % de croisements frère-sœur, il faut admettre que 80 % des individus d’une même ponte restent ensemble au cours des différents stades (larvaire, nymphal et adulte) de leur vie.

12Il se pourrait, contrairement à ce qui est observé en laboratoire où aucun œuf non fécondé n’a pu éclore, que cette espèce pratique une parthénogénèse automictique d’un type qui augmente l’homozygotie (pour des descriptions des différents modes d’automixie, voir par exemple De Meeûs et al., 2007b). Seules les femelles sont en général capables de parthénogénèse. Il existe cependant une espèce de cyprès et une espèce de phasme où les mâles se reproduisent asexuellement (voir encore De Meeûs et al., 2007b) et une espèce de fourmi où mâles et femelles sont clonaux chacun de leur côté (Foucaud et al., 2010). Mais ce sont des exceptions. Si parthénogenèse il y a, les femelles devraient donc présenter de beaucoup plus gros déficits en hétérozygotes que les mâles (tous issus d’une reproduction croisée). Nous allons donc réanalyser le fichier en demandant à Fstat de nous donner les FIS par sous-échantillon, puisque nous avons fort judicieusement, il faut bien l’avouer maintenant, d’entrée de jeu distingué les deux sexes.

13Sous Fstat vous ouvrez le même fichier “IRTotTestPanmix.dat” et vous décochez toutes les cases et cochez celle qui indique “Fis” dans le cadre “Per locus and sample statistics” comme indiqué dans la figure 18. Si vous souhaitez repérer encore une fois les noms des sous-échantillons, n’oubliez pas de signaler à nouveau l’existence du fichier “IRTotTestPanmix. lab” dans le menu “Options”.

Image 10000000000005C7000003527B0A3E41.jpg

Figure 18. Capture d’écran de Fstat lors de la deuxième analyse

14Quand vous lancerez “Run”, Fstat ouvrira une boîte de dialogue avec laquelle vous pouvez décider d’écrire les résultats de cette analyse dans un nouveau fichier. Dans le cas contraire, et c’est le choix que j’ai fait, le programme écrira les résultats dans “IRTotTestPanmix.out” à la suite des analyses précédentes (fin du fichier). Qu’y découvrons-nous ? Tout d’abord que Fstat tronque les labels plus longs que six caractères. Ce n’est pas grave, car nous avons toujours le même ordre Femelles Mâles pour chaque échantillon. Et puis il suffit (sous Excel c’est facile) de faire un copier-collage spécial/transposition à partir du fichier “IRTotTestPanmix.lab”. Ensuite, comme représenté dans la figure 19, construite à partir du fichier de sortie, nous pouvons voir, qu’à part pour le locus IR08, aucune tendance claire n’apparaît. Tous ces loci présentent des déficits importants et relativement variables, mais sans lien réel avec le sexe des tiques. Ce seraient plutôt les mâles qui auraient une tendance à présenter des déficits d’hétérozygotes plus importants (nous verrons plus loin une explication possible). Pour le locus IR08 par contre, avec un FIS = 1pourlesmâles, il apparaît clairement que ce locus est situé sur le chromosome X et qu’il est donc haploïde chez les individus mâles.

Image 1000000000000316000001EFA4C13FBB.jpg

Figure 19. Estimations des FIS par locus et par sous-échantillon. Les abréviations des échantillons sont identiques à celles de la figure 15. Les échantillons de femelles sont représentés par des ronds et ceux des mâles par des carrés

15En fait, pour être précis, le locus IR08 avait été trouvé hétérozygote pour quatre individus mâles sur l’ensemble du jeu de données. Même si cela pouvait refléter des duplications toujours possibles (comme évoqué p. 118), nous avons choisi d’éliminer ces individus, car ils pouvaient correspondre à des erreurs de manipulations.

16Quoi qu’il en soit, il va donc falloir recoder les données à ce locus. Pour l’analyse des FIS, les mâles devront en effet être codés en données manquantes (000000) au locusIR08. Nous allons donc créer un nouveau fichier “IRTotTestPanmixMalManqIR08.dat” à partir du précédent et refaire l’analyse globale du FIS. Celle des déséquilibres de liaison, qui est un test génotypique, n’a aucune raison d’avoir été affectée par ce phénomène. Dans Fstat, nous cocherons donc les mêmes cases qu’en figure 16, à l’exception de celles concernant les déséquilibres de liaison.

17Dans le fichier de sortie “IRTotTestPanmixMalManqIR08.out”, nous constatons l’image suivante (voir aussi la figure 20) : rien ne change sauf pour le locus R08 qui montre les plus basses valeurs de FIS, mais qui restent très significativement (toutes les P-values sont inférieures ou égales au minimum possible 0,0001) au-dessus de la valeur nulle attendue sous panmixie. Notez au passage que je ne me sers des intervalles de confiance que pour illustration. Le FIS global reste donc très élevé (0,39), inexplicablement variable entre loci et fort variable d’un site à l’autre. Ceci suggère un rôle possible pour des allèles nuls ou de dominance d’allèles courts. L’étape suivante sera donc de mettre en évidence l’existence de tels allèles et/ou de phénomène de dominance.

Image 100000000000020B0000017946B672D1.jpg

Figure 20. Valeurs moyennes du FIS par locus et intervalle de confiance obtenus pour les microsatellites d’Ixodes ricinus, avec les mâles codés comme données manquantes pour le locus IR08. Pour chaque locus, les intervalles de confiance sont calculés à l’aide de l’erreur standard donnée par le jackknife sur les populations et la valeur du t pour 11 ddl (12 – 1) (la moitié des échantillons) et a = 0,05 (soit 2,201, voir le tableau 3) en suivant l’équation (45). Pour la valeur globale, l’intervalle de confiance est issu du bootstrap sur les loci

RECHERCHE D’ALLÈLES NULS ET DE DOMINANCE D’ALLÈLES COURTS

18Nous allons pour ce faire utiliser deux nouveaux logiciels. Micro-Checker va nous permettre d’estimer la fréquence des allèles nuls susceptibles d’expliquer, dans chaque sous-échantillon et pour chaque locus, les déficits en hétérozygotes observés. Micro-Checker permet également d’estimer si les données sont compatibles avec un bégaiement de la polymérase (stuttering) et/ou une dominance des allèles les plus courts. Pour la dominance des allèles courts, nous utiliserons également une méthode plus puissante que celle implémentée par Micro-Checker. Nous allons procéder à une régression généralisée pour la mise en œuvre de laquelle nous utiliserons le logiciel R (voir la référence complète dans la bibliographie).

Convertir le fichier pour Micro-Checker et ouverture du logiciel

19Pour commencer avec Micro-Checker, nous avons besoin de transformer nos données au format Genepop qui est compatible avec ce logiciel. Ensuite, nous allons devoir créer un fichier spécial pour les données du locus IR08, lié au sexe, sans les mâles car sinon Micro-Checker risque de goûter moyennement la saveur de cette plaisanterie. Créons donc un fichier “IR08AllFem.txt” avec les données femelles pour le seul locus IR08 et un fichier “IRAutosomAll.txt” pour le reste des données. Attention, le fichier doit suivre des règles strictes sinon Micro-Checker refusera d’analyser les données. Référez-vous au fichier exemple fourni avec le logiciel et respectez les espaces et tabulations de la façon la plus scrupuleuse (ou utilisez Create).

20Lancez Micro-Checker et ouvrez “IRAutosomAll.txt” avec le menu “File”. Si tout se passe bien, vous observez l’ouverture de votre fichier avec vos données et différents menus et boutons en bas de l’écran.

Analyses des loci autosomiques du premier sous-échantillon par Micro-Checker

21Il y a un encadré en bas à gauche où il faut choisir le motif de chaque locus microsatellite. Il affiche par défaut le premier des loci (ici IR25) et un blanc pour le motif. Choisissez le motif “Mononucleotide” comme sur la figure 21.

22Nous avons déjà remarqué que nos loci microsatellites étaient peu orthodoxes. L’option mononucléotidique correspond en fait à l’option qui permet de faire face à toutes les situations. Cliquez ensuite sur le bouton “All” pour signaler que cette option est valable pour tous les loci. Cliquez ensuite sur le bouton “Analyse” (un peu plus à droite). Apparaît alors une fenêtre d’avertissement comme celle présentée en figure 22. Comme il y a des données manquantes, Micro-Checker vous demande s’il faut ou non en tenir compte. Autrement dit, les données manquantes correspondent-elles à des homozygotes nuls (blancs) et faut-il les utiliser pour le calcul des fréquences des allèles nuls par la seconde méthode de Brookfield (1996) ? La réponse étant positive, cliquez donc directement sur “Proceed” sans vous poser plus de questions.

Image 10000000000002F300000175467CA2B9.jpg

Figure 21. Capture d’écran de Micro-Checker

Image 100000000000025C0000016B3E0AE1CF.jpg

Figure 22. Cadre d’invite de commande de MicroChecker pour définir la nature des données manquantes et s’il faut en tenir compte dans le calcul des fréquences des allèles nuls

23Micro-Checker effectue plusieurs calculs et vous présente des résultats concernant le premier locus. Allez dans le menu “Tools” à “Nulls across loci” comme dans la figure 23 pour obtenir le tableau des fréquences de nuls dans le premier sous-échantillon, estimées selon différentes méthodes. Sélectionnez ce tableau avec la souris, copiez-le et sauvez-le dans un fichier (Excel, par exemple). Ensuite, regardez dans l’encadré en bas à droite (fig. 23) si le locus correspondant montre un problème de stuttering ou une dominance d’allèle court (« large allele dropout »). Si oui, notez-le dans le tableau que vous venez de créer pour sauvegarder les résultats de cette analyse puis, par le menu “Window’” (fig. 23) sélectionnez le locus suivant, etc. Vous constaterez qu’aucun locus ne présente de « stuttering » ni de dominance d’allèle court dans ce premier sous-échantillon.

Analyses des autres sous-échantillons, des autres loci autosomiques et du locus IR08

24Au centre et en bas, cliquez sur le bouton “Next Population” (voir fig. 23) pour analyser le sous-échantillon suivant en reprenant les mêmes étapes décrites en p. 125-127, jusqu’au dernier sous-échantillon. N’oubliez pas de copier le tableau des fréquences d’allèles nuls à chaque fois (dans le menu “Tools” à “Nulls across loci”, fig. 23). Ensuite, vous ferez la même chose pour le locus lié au sexe, IR08, en ouvrant le fichier correspondant “IR08AllFem.txt”.

Image 10000000000002ED000001C4BDAF9D2E.jpg

Figure 23. Sortie de MicroChecker vous indiquant, pour le locus et le sous-échantillon mentionné, la présence ou non de problèmes

Bilan des analyses avec Micro-Checker

25Nous avons constitué un fichier de résultats avec les fréquences d’allèles nuls probables, l’existence ou non de stuttering et de dominance d’allèles courts. Nous ne gardons que la méthode 2 de BROOKFIELD (1996) qui tient compte des données manquantes (blancs) comme des homozygotes nul/nul. Dans ce fichier, nous allons également insérer le nombre d’individus génotypés pour chaque locus (copiés à partir des fichiers de sortie Fstat), la fréquence attendue sous panmixie (fréquence précédente au carré) des allèles nuls pour chaque locus dans chaque sous-échantillon et sur l’ensemble des sous-échantillons, le nombre de blancs observés (compter les 000000 dans chaque sous-échantillon et sur l’ensemble), l’effectif corrigé (individus génotypés + blancs) et enfin le nombre de blancs attendus sous la double hypothèse qu’il y a panmixie et que les allèles nuls expliquent les FIS en totalité. Le tableau 8 donne un aperçu du fichier final pour le locus IR08.

Tableau 8. Synthèse des résultats de Micro-Checker pour le locus IR08 chez les femelles Ixodes ricinus. La fréquence attendue des blancs pB2² est obtenue en mettant au carré la fréquence estimée des allèles nuls selon la méthode 2 de BROOKFIELD (1996) et le nombre de blancs attendus correspondant à cette valeur multipliée par N’. N’ correspond, quant à lui, à la somme de N (individus génotypés) et des blancs observés. Pour la dernière ligne, la valeur de pB2² est obtenue en divisant le nombre total de blancs attendus par le N’ total

Image 10000000000005CF000002589DAE001C.jpg

26Pour vérifier que ces résultats expliquent correctement les FIS observés, on peut comparer la proportion de blancs observés avec celle attendue sous l’hypothèse que les allèles nuls expliquent la totalité de ces FIS. Un test binomial unilatéral avec comme fréquence attendue pB2², un nombre de réussite égal aux blancs observés pour un nombre d’essais de N’, semble ici approprié. On préfère ici un test unilatéral, car ce qui nous intéresse est de savoir si on a oui ou non moins de blancs qu’attendus. On peut facilement effectuer ce test sous R.

27Il nous faut donc lancer R et dans la fenêtre de commande taper l’instruction :

28binom.test (Blancs observés, N’, p = pB22, alternative = "less"

29Pour des raisons de recherche de puissance et pour limiter le nombre de tests dont la multiplication est toujours problématique (voir p. 82 en première partie), on ne fera les tests qu’avec les valeurs totales pour chaque locus. Pour le locus IR08, cela correspond aux valeurs de la dernière ligne du tableau 8. Pour ce locus, la commande devient donc :

30binom.test (4, 186, 0.0253, alternative = "less")

31Faites bien attention de respecter strictement le format (en particulier, les majuscules et minuscules sont reconnues comme des caractères différents sous R). Ici “less” signifie que le test est unilatéral dans le sens des plus petites valeurs (H1 : il y a moins de blancs observés qu’attendus) (l’instruction devient “two.sided” pour un bilatéral et “greater” pour l’autre test unilatéral). Une fois que vous avez tapé cette instruction dans R, tapez sur la touche “Entrée” et le test se fait. La P-value du test est, pour IR08, non significative (P-value = 0,4919). Les allèles nuls sont donc bien suffisants pour expliquer les déficits en hétérozygotes observés à ce locus chez les femelles, d’autant plus qu’il semble aussi exister des phénomènes de stuttering à ce locus. Pour les autres loci, on procède de la même façon. On trouve ainsi que pour les loci IR 25, IR27 et IR32, la fréquence des blancs observés est significativement inférieure à celle des blancs attendus si les allèles nuls devaient expliquer les déficits en hétérozygotes. C’est un problème car, par un phénomène de cercle vicieux, moins les allèles nuls expliquent un déficit en hétérozygotes, moins le nombre de blancs observés correspond aux attendus. Pourquoi cela ? Simplement parce que si on attend naturellement plus d’homozygotes en général, alors on devrait observer encore plus d’homozygotes nuls (blancs), en particulier (ce raisonnement ne marche cependant pas très bien s’il s’agit d’un effet Wahlund). Par ailleurs, la variance entre loci ainsi que le fait que les nuls expliquent très bien les déficits observés pour IR08 (voir plus haut), mais aussi pour IR39 (P-value = 0,312) pourraient nous inciter à exclure des causes biologiques du type régime de reproduction ou effet wahlund (voir plus loin). Notons que des phénomènes de stuttering ont été détectés pour IR25, mais seulement dans deux sous-échantillons. Pour IR32 et IR27, Micro-Checker n’a pas détecté ce phénomène pas plus qu’il n’a détecté de dominance d’allèles courts. Cependant, Micro-Checker ne travaille que dans chaque sous-échantillon de façon isolée, ce qui peut représenter une forte perte de puissance. Dans le paragraphe qui suit, nous allons utiliser une autre technique pour détecter d’éventuelles dominances d’allèles courts.

Détection de dominance d’allèles courts par la méthode de régression multiple

32Pour ce faire, nous aurons besoin de connaître, pour chaque locus et dans chaque sous-échantillon, la valeur du FIS pour chaque allèle. On peut demander à Genetix de le faire en choisissant à chaque traitement le locus et le sous-échantillon à analyser, en n’oubliant pas de zapper les mâlesaulocusIR08. On peut aussi créer autant de fichiers Fstat qu’il y a de sous-échantillons à analyser, ensuite, et parce que malheureusement Fstat ne permet pas d’analyser qu’un seul sous-échantillon, il faut créer dans chaque fichier une deuxième population fictive, de taille identique à celle à analyser et fixée à tous les loci (par exemple, tous homozygotes 170170, 150150, 123123, 235235, 129129 pour les cinq loci respectivement). Il s’agit ensuite de récupérer dans chaque sous-population les FIS de chaque allèle pour chacun des cinq loci et de créer cinq fichiers de données (un par locus) contenant pour chaque allèle son FIS, sa taille (on s’en doute), le sous-échantillon, sa fréquence allélique p dans ce sous-échantillon, le produit p (1-p), le nombre d’individus génotypés dans ce sous-échantillon N et enfin le produit p (1-p) N. Le tableau 9 donne une idée de la forme de ce fichier pour le locus IR08 que j’ai appelé “IRTotL08MalManqFisAllSizeL08.txt”. Pour fabriquer ce fichier, une feuille de calcul Excel est idéale, ensuite il suffit d’enregistrer le fichier en format texte seul.

33On peut aussi utiliser Genetix qui permet l’analyse d’un seul sous-échantillon, mais dont les sorties sont moins commodes à importer dans Excel (à vous de voir).

34La colonne Npq, qui donne en fait le résultat du produit Np (1-p), nous servira à pondérer notre régression par la taille des échantillons, mais en donnant aussi plus de poids aux allèles de fréquences proches de 0,5 (les plus polymorphes). On fait les mêmes fichiers avec les quatre autres loci. Nous allons maintenant analyser ces données avec le logiciel R.

Tableau 9. Aperçu du fichier de données pour le locus IR08 en vue de l’analyse de régression du FIS en fonction de la taille des allèles et du sous-échantillon

Image 10000000000005EA000001EAC51A3616.jpg

35Ouvrez R et dans le menu “Fichier” cliquez dans “Changer le répertoire courant...”, et allez dans le répertoire où vous avez stocké vos fichiers de données. Dans la console de travail de R, tapez la suite de commandes, chacune suivie d’un retour chariot (touche “Entrée”) :

36> data<-read.table ("IRTotL08ManqFisAllSizeL08.txt", header = TRUE)

37qui signifie que le tableau de données “data” est contenu dans le fichier nommé et que la première ligne contient le nom des colonnes. N’oubliez pas que les données manquantes se notent “NA” en majuscules et non “000000”.

38> attach(data)

  • 9 Entre temps, j’ai découvert l’existence du “Package” R-Commander ou Rcmdr qui, en quelques clics de (...)

39qui signifie que ce tableau doit être chargé en mémoire9.

40> loc8<-glm(data, formula = Fis ~ poly(Allele, 2) + Sample + Year, family = gaussian, weights = Npq)

41où loc8 est le nom d’un modèle linéaire généralisé utilisant le tableau “data” et dont la régression tente d’expliquer la valeur du FIS en fonction de la taille des allèles selon un polynôme d’ordre 2 ou quadratique (qui s’est avérée plus proche de ce qui se passe dans le cas qui nous intéresse), du sous-échantillon d’origine et de l’année. Le sexe n’a ici aucune importance puisqu’il n’y a que des femelles. Nous ne testons l’effet d’aucune interaction entre variable, car en fait je ne vois aucune raison pour qu’il en existe. Pensez à respecter les majuscules s’il y en a, car R les reconnaît comme telles. Tapez enfin :

42> anova(loc8, test = "F")

43qui renvoie à une analyse de variance utilisant la statistique F (se référer à un livre de statistique pour approfondir ces notions) et donne le résultat suivant :

Analysis of Deviance Table
Model: gaussian, link: identity
Response: Fis
Terms added sequentially (first to last)
Image 10000000000002E200000089A70170FD.jpg

44Ici, on voit que seul le sous-échantillon influence la valeur du FIS (allèles nuls, stuttering variable dans l’espace ?) qui n’explique que 14,63 % de la dispersion (100 × 3,1604/21,616), tout en étant très significatif.

45On utilise un test F, car on a supposé que la distribution des FIS suit plus ou moins une courbe de Gauss (données continues en cloche symétrique), ce qui est sûrement inexact mais ne risque guère de modifier le résultat dans un sens dramatique.

46Pour les loci suivants, nous aurons besoin de distinguer le sexe des tiques.

47Avec le locus IR25, l’analyse du fichier “IRTotSexSepFisAllSizeL25.txt” est la suivante :

> data<-read.table("IRTotFisAllSizeL25.txt",header=TRUE)
> attach(data)
> loc25<-glm(data, formula = Fis - poly(Allele, 2) + Site + Year + Sex,
family = gaussian, weights = Npq)
>anova(loc25, test="F")

48Ce qui aboutit au tableau de résultat :

Image 10000000000002F50000009F0D6D7A4E.jpg

49On aboutit à une conclusion similaire à la précédente, puisque ni le sexe ou l’année ni la taille des allèles ne comptent avec seulement environ 6,13 % de la déviance expliquée par le site qui est moins spectaculairement significatif que précédemment. Pour le locus IR27, le tableau obtenu est différent :

Image 10000000000002F3000000A363067130.jpg

50En effet, comme nous pouvons le déduire du tableau ci-dessus, le site (Sample) explique 19,96 % de la dispersion des points (5.181/25.9549) et la taille des allèles (poly (Allele, 2)) en explique 15,86 % (4.1186/25.9549) et sont tous les deux très significatifs (souligné par les trois étoiles). Ils expliquent ainsi 35,83 % de la variance. Cette valeur est conséquente eu égard à l’importante variance résiduelle attendue en général pour un estimateur de statistique F. Comme le montre la courbe décrite dans la figure 24, la relation entre FIS et taille des allèles est négative (si on exclut les trois premiers points, ce qui ne changerait rien eu égard aux intervalles de confiance), ce qui peut donc être interprété par une dominance des allèles les plus courts. Pour le locus IR32, on observe le résultat suivant :

Image 10000000000002F5000000A404667110.jpg

Image 1000000000000236000001DA1CFDB2BD.jpg

Figure 24. Relation entre taille des allèles et FIS pour le locus IR27 et sur l’ensemble des échantillons. Les intervalles de confiance à 95 % ont été obtenus avec Image 10000000000000EB0000003941B550B8.jpg Pour ce faire, les singletons (tailles d’allèles présents une seule fois comme 112 et 131) ont été réunis à la classe la plus proche

51On voit qu’en plus du site, le sexe des tiques a un effet significatif, ce qui signifie que nous avons eu raison d’en tenir compte et nous verrons ensuite pourquoi. Pour le locus IR39, le tableau obtenu est le suivant :

Image 10000000000002F3000000A26E83FB81.jpg

52Le site joue une fois encore de façon significative, mais aussi la taille des allèles, même si cette dernière n’explique même pas 2 % de la déviance et est peu significative. Par ailleurs, la figure 25 montre que la relation (augmentation globale du FIS avec la taille des allèles) n’est pas compatible avec une dominance des allèles courts. On peut donc attribuer ce résultat au hasard et au nombre de tests effectués qui augmente la probabilité d’obtenir quelque chose de significatif par hasard (revoir la première partie de ce manuel, p. 82-87). Rappelons que pour ce locus, les allèles nuls s’étaient avérés suffisants pour expliquer les déficits en hétérozygotes observés. Il est plus raisonnable ici de considérer ce résultat comme fortuit.

Image 1000000000000236000001555C2EE00F.jpg

Figure 25. Relation entre FIS et la taille des allèles au locus IR39 sur l’ensemble des échantillons. Les intervalles de confiance à 95 % ont été obtenus comme précédemment. Les sous-échantillons de moins de quatre individus ont été associés au plus proche

Bilan de l’analyse des déficits locaux en hétérozygotes

53Pour les loci IR08 et IR39, les allèles nuls semblent pouvoir expliquer les forts et variables FIS observés. Pour IR27, les allèles nuls et la dominance des allèles courts offrent conjointement une explication satisfaisante. Seul le locus IR32 offre des déficits énormes et non expliqués par les allèles nuls, le « stuttering » ou la dominance des allèles courts. Cependant, sachant que le « stuttering » n’a pu être testé que sous-échantillon par sous-échantillon (manque de puissance), que la plupart des allèles se suivent à un pas sur ce locus et compte tenu de ce que nous trouvons aux autres loci, il est possible qu’ici aussi les déficits observés proviennent d’un problème technique.

54Je peux ajouter ici qu’un module (package) de R, appelé “R-Commander”, dont je n’ai appris l’existence qu’après la rédaction de ce chapitre, permet d’accéder aux analyses effectuées dans ce paragraphe à l’aide de menus déroulants plus conviviaux que le mode commande strict.

RECHERCHE D’UNE STRUCTURE CACHÉE (EFFET WAHLUND)

Introduction

55Nous allons dans un premier temps continuer de considérer les femelles et les mâles séparément. On sait en effet qu’il y a une structure génétique spécifique pour chaque sexe dans ce jeu de données. Même si nous analyserons ceci plus tard, il n’est pas inutile de poursuivre la recherche d’explications des déficits en hétérozygotes avant d’aborder cet aspect. Nous allons donc analyser tous les sous-échantillons (mâles et femelles séparées) pour obtenir l’information sur le plus grand nombre de réplicas possibles. Ensuite, nous nous concentrerons sur 1996 en réunissant les mâles et les femelles pour faire des tests.

56Le but du jeu sera ici d’utiliser l’information multilocus de chaque individu, dans chaque sous-échantillon afin de vérifier à l’aide du logiciel BAPS (voir p. 98-101 en première partie et le tableau 1 en annexe), si certains individus peuvent être regroupés sur la base de leur ressemblance génétique. BAPS va ensuite explorer de façon itérative et répétée, en suivant plusieurs chaînes de Markhov (ou une chaîne stochastique d’optimisation suivant les versions) afin de trouver la meilleure partition (celle qui regroupe le mieux les individus) dans le sous-échantillon analysé. La partition définit un nombre donné de clusters (sous-unités) composés chacun d’un certain nombre d’individus du sous-échantillon. La qualité d’une partition se définit par un savant calcul dans le détail duquel je serai bien incapable de rentrer, mais qui dépend de la distance génétique entre les groupes définis, par rapport aux autres partitions explorées durant le processus. Il est aussi expliqué, dans les articles décrivant le logiciel, qu’une hypothèse du modèle utilisé dans l’algorithme est que les « clusters » qui composent la partition sont en équilibre de Hardy-Weinberg. Je ne suis pas certain de bien comprendre ce qui est entendu par là dans la mesure où mon expérience m’a montré que la plupart des partitions obtenues ne sont pas conformes à cet équilibre, voire même en sont très éloignées. J’ai également pu observer cela avec STRUCTURE qui fait la même hypothèse. Comme discuté dans la première partie de ce manuel, beaucoup reste à explorer concernant le fonctionnement de ces méthodes dans différentes situations. Il faudra donc vérifier si la partition obtenue (car le logiciel en donne toujours une) correspond à quelque chose de viable et pas seulement une vue de l’esprit.

57Si la partition a réellement mis en évidence des groupes cryptiques au sein des sous-échantillons susceptibles d’expliquer en partie (effet Wahlund) nos fameux déficits en hétérozygotes, il faudra ensuite trouver et explorer les hypothèses susceptibles d’expliquer le plus raisonnablement possible (mais en aveugle) ces résultats (espèces ou races d’hôtes cryptiques, sous-structures familiales, isolement par la distance entre individus sur de courtes distances).

58Il existe d’autres logiciels qui en principe font la même chose. L’avantage de BAPS réside dans sa convivialité, dans le fait qu’il accepte des fichiers de type Genepop (un peu modifiés) et qu’il m’a toujours donné de bons résultats. Le logiciel STRUCTURE est par exemple beaucoup moins commode à utiliser (et c’est un euphémisme) et, sur un même jeu de données (glossines), n’a pas offert de partitions aussi satisfaisantes que BAPS (Ravel et al., 2007). Des études comparatives de différents logiciels de clustering sont en cours, mais la longueur et la quantité des analyses font que des résultats concrets ne seront sans doute pas disponibles avant la sortie du présent ouvrage. Vous verrez aussi l’application d’un autre logiciel de même nature, Flock, plus loin dans cette partie.

Construction des fichiers BAPS

59Il faut construire un fichier pour chaque sous-échantillon. Le type est semblable à un fichier Genepop, mais avec des tabulations comme dans la figure 26 (symbolisées par des →) qui donne un exemple pour le fichier des mâles de Staadswald. On note que les mâles sont codés homozygotes pour IR08 afin que l’information multilocus soit préservée pour cinq loci. Par commodité, j’ai appelé ce fichier “IRTotBrut1Stad95M.gen”, mais vous faites comme bon vous semble.

60Ensuite, il est commode de créer un fichier texte contenant le chiffre 30 répété un grand nombre de fois (ici 50 fois), avec un espace entre chaque répétition et sur une seule ligne. Le logiciel BAPS vous demandera en effet de taper un nombre maximal probable pour les clusters. Ici, 30 m’est apparu comme largement raisonnable compte tenu des tailles de sous-échantillons. C’est à partir de ce chiffre que BAPS démarre et recherche une partition la plus probable en se limitant à ce nombre maximum de clusters. Le logiciel reprend ensuite le processus autant de fois que l’on a rentré ce chiffre (ici 50) et ne gardera que la meilleure de toutes les partitions explorées. Avoir tapé 50 fois ce chiffre dans un fichier permet de copier et coller cette séquence directement sans avoir à la retaper pour toutes les analyses. J’ai appelé ce fichier “50fois30.txt” (quelle imagination !).

Analyse des fichiers par BAPS

61Vous avez bien entendu installé BAPS sur votre machine et créé tous les fichiers nécessaires (il y en 24 normalement). Il faut maintenant lancer BAPS en cliquant sur BAPS4_RUNME.EXE. Le logiciel ouvre deux fenêtres, une fenêtre Dos dont il n’est pas vraiment nécessaire de se préoccuper maintenant et une fenêtre d’interface type Windows avec des menus que nous allons utiliser. Il est important de commencer par créer un fichier résultat. Pour ce faire, cliquez sur “File”, “Output File” et “Set” et créez un fichier en tapant son nom et en le plaçant dans le répertoire qui vous convient le mieux (là où sont vos données) (voir fig. 27).

Image 10000000000005C70000054AAF2F15C6.jpg

Figure 26. Format de fichier pour BAPS pour les tiques mâles du Staadswald en 1995. Le locus IR08 est codé homozygote chez les mâles

62Il vaut mieux garder un nom de fichier qui permette de retourner ensuite au fichier de données correspondantes. Ici, le premier fichier analysé sera “IRTotBrut1Ber95F. gen” (femelles de Berne 1995), je choisis donc ici de nommer et créer le fichier résultat “IRTotBrut1Ber95FBAPSRes.txt”. Ensuite, il faut cliquer sur le bouton “Clustering of individuals” (fig. 27). Apparaît alors une nouvelle fenêtre de dialogue qui vous propose différents formats de fichiers de données (fig. 28). Choisissez bien entendu le format Genepop en cliquant sur le bouton correspondant. Une fenêtre qui s’ouvre vous permet de naviguer vers le répertoire où se trouve IRTotBrut1Ber95F. gen que vous sélectionnez (soit en tapant son nom complet, soit en tapant *. gen et retour chariot et en double cliquant sur le fichier). Une question vous est alors posée “Do you wish to save pre-processed data ?”, cliquez sur “No”.

Image 100000000000050F000003D345038516.jpg

Figure 27. Sélection dans BAPS du fichier de résultats

Image 100000000000050F000003D3DD2F7232.jpg

Figure 28. Sélection dans BAPS du format de fichier de données à analyser

Image 100000000000049F000001725B35DC5F.jpg

Figure 29. Fenêtre de sélection de la taille maximale des clusters et du nombre de chaînes d’itérations

63C’est alors qu’apparaît une petite fenêtre permettant de sélectionner le nombre maximum de clusters, ainsi que le nombre de chaînes d’itérations à effectuer (fig. 29), comme expliqué en p. 136. Supprimez le chiffre par défaut (20) et remplacez-le par la chaîne de 30 que vous copiez à partir de “50fois30.txt”, collez cette chaîne dans la case idoine et cliquez sur “OK”.

64Les calculs démarrent et se poursuivent jusqu’à la fin où la meilleure partition est sauvée dans “IRTotBrut1Ber95FBAPSRes.txt”. Apparaissent un graphique censé représenter la partition (clusters de différentes couleurs), dont on ne va pas se servir, ainsi qu’un dialogue final vous demandant si vous souhaitez sauver ces données en vue d’une analyse ultérieure. Répondez non. Ceux qui souhaitent plus de détails sur BAPS et ses différentes possibilités et menus sont invités à consulter la documentation livrée avec le logiciel.

65Il s’agit ensuite de répéter le processus avec chacun des sous-échantillons. Ensuite, on charge le jeu de données brutes afin de le modifier. N’oubliez pas de créer un nouveau fichier de sortie à chaque fois. Dans chaque fichier de résultat BAPS sont donnés les clusters avec les individus qu’ils contiennent. Ces individus sont identifiés par leur rang d’entrée dans le jeu de données (1, 2, 3...). Par exemple, pour les femelles de Berne 1995, le fichier de résultat donne (en début de fichier) :

Image 1000000000000458000001BD31614272.jpg

Image 10000000000001EB000001BD79215715.jpg

66Le nom du fichier analysé est suivi de l’effectif de l’échantillon, du nombre de clusters dans la meilleure partition et de la valeur du Log de la valeur marginale de vraisemblance ou Log (MV) qui sert de critère à BAPS pour sélectionner la meilleure partition, c’est-à-dire celle qui présente le plus petit Log (MV). C’est bon à savoir si on souhaite relancer BAPS sur les mêmes données afin de voir s’il trouve une partition meilleure au deuxième essai. Enfin, la partition est donnée. Dans le jeu de données, il faut donc maintenant ajouter une colonne avec le numéro de cluster BAPS auquel chaque individu appartient. Il faut le faire pour tous les sous-échantillons (cf. tabl. 10). Attention, vous allez peut-être trouver des partitions légèrement différentes des miennes et avec des labels de clusters différents, c’est normal.

67Il faut ensuite créer un nouveau fichier de données où chaque sous-échantillon initial se retrouve subdivisé en autant de sous-échantillons que de clusters de BAPS qui le composent (12 pour les femelles de Berne 1995). Sous un éditeur quelconque vous fusionnez les colonnes 1, 2, 3 et 6 du tableau 10, ce qui donne pour la première ligne quelque chose du style Ber95F1. N’oubliez pas de trier les données pour que les clusters apparaissent dans l’ordre dans chaque sous-échantillon initial. Appelons le fichier contenant ces données modifiées “IRTotBAPSClustMalHomoMFSep.txt”.

68Ce n’est pas fini, car il faut maintenant coder en données manquantes le locus IR08 chez les tiques mâles. Rappelez-vous que, dans un souci de puissance, nous les avions artificiellement rendus homozygotes afin que les mâles soient pris en compte pour ce locus dans l’analyse BAPS. Maintenant, nous souhaitons calculer les nouveaux FIS de cette partition afin de voir si elle chute par rapport au jeu de données initiales. Le génotype des mâles au locus IR08 doit donc en effet être recodé 000000, car ils ne doivent pas rentrer en ligne de compte dans le calcul du FIS. Pour ce faire, il est commode soit de faire un petit programme (pour ceux qui savent), soit d’utiliser la fonction conditionnelle d’Excel. Il s’agit de créer une colonne sexe en A dans le jeu de données “IRTotBAPSClustMalHomoMFSep.txt” avec le sexe des individus (F ou M), dans une colonne libre (en H après IR39) on tape en ligne 2 (ligne du premier individu) :

69SI (A2 ="M" ;"000000" ; C2), ce qui aura pour effet d’écrire “000000” dans la case H2 si l’individu est mâle ou de recopier le génotype de la femelle au locus IR08 (contenu dans la case C2). On copie ensuite H2 et on le colle de H3 à H726 (normalement la fin du fichier). On sélectionne les cases H2 à H726, on les copie et on fait un collage spécial (on veut ne coller que la valeur et non la formule) sur C2. On supprime les colonnes H et A et on sauve en texte seul sous le nom “IRTotBAPSClustMalManqIR08MFSep.txt”. Supprimez aussi le label de la première colonne (c’est pour Genetix qui ne désire que le nom des loci).

Tableau 10. Aspect du tableau de données brutes modifiées avec l’appartenance des individus aux clusters BAPS

Image 10000000000005EE000008220A379B0D.jpg

Image 10000000000005C7000004152B9B0CF4.jpg

Figure 30. Importer les données dans Genetix

70Nous allons maintenant recalculer les FIS par locus et sur l’ensemble, avec intervalles de confiance. Nous sommes paresseux et pour ne pas avoir à supprimer les clusters d’un individu pour lequel le calcul ne se fera pas, et étant donné que Fstat ne prend pas plus que 200 sous-échantillons (avec mes partitions je me retrouve avec 368 sous-échantillons), nous allons importer notre nouveau fichier sous Genetix. Lancez Genetix et allez dans le menu “Fichier”, sélectionnez “Importer” et sélectionnez “Texte avec séparateur” et sélectionnez le fichier. Un menu apparaît et si vous avez fait comme moi, vous devez cocher les cases comme dans la figure 30. Quand cela est fait, cliquez “OK”. Si le fichier est correctement chargé, cliquez dans le menu “Fstats” et sélectionnez “Weir & Cockerham”. Cliquez OK dans la nouvelle fenêtre si vous ne changez pas le nom du fichier de sortie proposé “IRTotBAPSClustMalManqIR08MFSep.res”. Après il faut prendre les résultats dans ce fichier en prenant garde que ce qui est annoncé comme écart-type des jackknives correspond à l’erreur standard de Fstat. Il s’agit de comparer maintenant les FIS de chaque loci et leurs intervalles de confiance de jackknife sur populations (voir p. 70-74 en partie 1) avant et après clusterisation par BAPS, ainsi que les valeurs globales et leur intervalle de confiance de bootstrap sur les loci (voir p. 70-74 en partie 1). La compilation des résultats prend alors la forme de ce qui est représenté dans la figure 31. Dans cette figure, il est aisé de voir que les clusters de BAPS présentent des déficits en hétérozygotes significativement inférieurs au FIS de départ. Un test de rang de Wilcoxon pour données appariées confirme cela. Pour effectuer ce test sous R, il faut construire un fichier avec une colonne “Delta” où chaque ligne correspond à un locus.

Image 1000000000000535000002AB1159A6D9.jpg

Figure 31. Comparaison de FIS avant (données initiales : Brut) et après clusterisation par BAPS sur l’ensemble des données, par locus et sur l’ensemble (Tous)

71Chaque valeur représente la différence entre le FIS brut et le FIS BAPS au locus correspondant (ici cinq valeurs). Appelons ce fichier “DeltaFisBrutBAPS.txt”. Ensuite, sous R les commandes sont les suivantes :

>data<-read.table("DeltaFisBrutBAPS.txt", header=TRUE)
>attach(data)
>wilcox.test(Delta, alternative="greater")

72Le test est unilatéral, car ce que nous recherchons est bien un effet Wahlund. Nous attendons au départ une chute du FIS, d’où l’instruction “greater”. La P-value = 0,031 obtenue est significative. Notons aussi que la plupart des loci, mis à part IR08, gardent un fort FIS qui provient probablement des allèles nuls et autre dominance des allèles courts. Ces déficits restent très significativement au-dessus de 0 (fig. 31), ce qui rend bien compte du fait que “Hardy-Weinberg” n’est pas une nécessité pour parvenir à une partition. Par ailleurs, le FIS fait mieux qu’être faible pour IR08, il est négatif, ce qui est effectivement ce que nous attendons chez une espèce dioïque pangamique.

  • 10 Je me suis rendu compte sur le tard que Genetix contenait quelques bugs dans ce module et je consei (...)

73Il semble donc bien y avoir un effet Wahlund, contrairement à ce que la variance du FIS entre loci pouvait laisser prévoir. Reste à déterminer si cet effet provient d’une micro-structuration (en groupes familiaux, par exemple) ou de la présence d’espèces (ou races d’hôtes, ou groupes adaptatifs ou écotypes) cryptiques. Afin d’essayer d’argumenter dans un sens ou l’autre, on peut essayer de regarder l’organisation de ces différents clusters. En principe, si on a à faire à différentes espèces, ces dernières devraient apparaître clairement. Si on effectue un arbre à partir d’une matrice de distance inter-clusters, ces derniers devraient être regroupés selon l’espèce à laquelle ils appartiennent en groupes séparés par des branches relativement longues comparées aux branches séparant chaque cluster (géographique, en principe) à l’intérieur de chaque espèce. Selon Takezaki et Nei (1996), la méthode du Neighbor-Joining (NJTree) sur distances harmoniques de Cavalli-Sforza et Edwards (1967) est une bonne solution. La matrice est obtenue en important “IRTotBAPSClustMalHomoMFSep.txt” dans Genetix10, en cliquant sur le menu “Distances” puis “Sur données réelles” et en sélectionnant “Cavalli-Sforza & Edwards”. On copie la matrice obtenue afin de l’incorporer dans un fichier de type MEGA (Kumar et al., 2004) pour matrice de distances (ouvrir le fichier “IRTotBAPSClustMalHomoForNJTREENmini3CSE.meg” avec un éditeur de texte pour voir un exemple). Afin de limiter le nombre de branches et le poids des clusters ne contenant qu’un seul ou deux individus, je n’ai gardé que les clusters d’au moins 3 individus. L’arbre obtenu n’en est pas plus lisible pour autant et ce qui en ressort, c’est que les plus longues branches sont toujours celles séparant les clusters sans que se dégage une quelconque hiérarchie (on parle de râteau). Ceci plaide davantage en faveur d’une micro-structuration locale forte avec une différenciation géographique faible. On peut alors recommencer l’ensemble des opérations (BAPS-> Genetix-> MEGA) sur les échantillons de 1996 seuls et en ne séparant pas les mâles des femelles. Sur l’arbre obtenu, on ne voit pas mieux une quelconque structure sauf que les clusters tunisiens de plus de deux individus se retrouvent bien ensemble (fig. 32) avec un cluster du Tessin (Cen16 qui comporte d’ailleurs deux mâles et une femelle). Ceci ne contredit pas que l’effet Wahlund pourrait être issu de la présence dans chaque site d’individus issus des mêmes pontes. Ceci implique une forte variance du succès de survie entre pontes : de nombreux individus issus seulement de quelques pontes accèdent à l’âge adulte (voir Chevillon et al., 2007a, pour un résultat similaire sur la tique du bétail).

Image 10000000000005C7000005A53DB98AE6.jpg

Figure 32. Dendrogramme exécuté selon la méthode du NJTREE sur les distances harmoniques de Cavalli-Sforza et Edwards entre paires de clusters BAPS de taille supérieure ou égale à 3 dans chaque sous-échantillon des tiques de 1996. Les clusters tunisiens sont indiqués en rouge

Commentaires sur l’analyse des fichiers par BAPS

74Contrairement à ce qui pourrait être suggéré à la lecture du manuel d’utilisation de BAPS, les clusters obtenus ne présentent pas ici une structure panmictique, mais conservent un déficit important d’hétérozygotes sauf pour IR08. Nous verrons, avec les analyses suivantes, que ces clusters reflètent probablement en grande partie une réalité biologique de nature assez complexe (races d’hôte, structures familiales), et qui devra conduire à d’autres études. BAPS ne représente ici qu’un outil d’argumentation et d’orientation de futures investigations, pas un générateur de vérités.

CONCLUSION SUR LES DÉFICITS EN HÉTÉROZYGOTES

75À l’occasion de ces premières analyses, nous pouvons constater qu’une analyse d’un jeu de données de génétique des populations requiert de la patience, de la méthode, ainsi qu’une bonne batterie de tests. Il était cependant nécessaire d’aller jusqu’au bout avant d’aller plus loin. Nous savons maintenant que ces tiques sont structurées à une échelle locale, ce qui explique une grande partie des déficits en hétérozygotes. Cet effet Wahlund résulte probablement d’une structure en groupes familiaux. L’existence d’espèces cryptiques n’est en effet pas soutenue par nos analyses NJTREE ni par l’absence totale de déséquilibre de liaison. Nous savons également qu’une partie non négligeable de ces déficits provient de l’existence d’allèles nuls (Loci IR25, IR32, IR39). Pour ces derniers, il y a donc un risque de surestimer la différenciation entre sous-échantillons, mais seulement pour des niveaux de différenciation atteignant au moins 10 % (FST = 0,1), en dessous de quoi l’effet devient faible (Chapuis et Estoup, 2007). Nous verrons que les niveaux de différenciation entre populations d’I. ricinus se trouvent bien en dessous de cette frontière. Enfin, un locus (IR27) a montré des évidences de dominance des allèles courts. Dans la mesure où ce phénomène modifie l’hétérozygotie et l’estimation des fréquences alléliques, il faudra être constamment vigilant quant aux résultats obtenus par la suite. Si nous avions un locus de plus sans allèle nul, j’aurais même conseillé de le supprimer. Ce n’est malheureusement pas le cas. Il faudra juste vérifier que chaque résultat ultérieur n’est pas sous la dépendance de ce seul locus. L’idéal aurait été d’avoir sept loci comme IR08, mais non liés à l’X ! Mais on ne choisit pas et les problèmes de marqueurs chez les parasites et vecteurs représentent un souci récurrent.

76Une autre conclusion importante est qu’un déficit en hétérozygotes non entièrement expliqué par des allèles nuls exclut les causes endogamiques (croisements frère/sœur, autofécondation…) qui tendent à augmenter l’homozygotie et donc à dévoiler les homozygotes nuls (blancs), d’une part, et suggère, d’autre part, plutôt un effet Wahlund, qui augmente la diversité génétique Hs sans augmenter l’hétérozygotie observée (d’où augmentation du FIS, cf. équation 19 en première partie de ce manuel, p. 47). Dans le cas d’un effet Wahlund, il est donc normal que les procédures de détection d’allèles nuls ne suffisent pas à expliquer entièrement les déficits en hétérozygotes, même si ces derniers sont présents, comme l’attestent la présence fréquente d’individus blancs, ainsi que la forte variance du FIS entre loci.

STRUCTURE DES POPULATIONS ET SCHÉMAS DE DIFFÉRENCIATION

77Nous avons ici une espèce à sexes séparés. La première chose à tester est s’il n’existe pas une différence entre femelles et mâles tiques, liée par exemple à un biais de dispersion spécifique de chaque sexe (Goudet et al., 2002 ; Prugnolle et De Meeûs, 2002). En plus, on sait que c’est probablement le cas ici puisque ce signal fut détecté précédemment (De Meeûs et al., 2002a), mais aussi lors de notre recherche de dominance des allèles courts. Par ailleurs, il est intéressant de vérifier si le signal persiste en tenant compte de la microstructure en clusters, même s’il a été montré que celle-ci a peu (pas) d’effet sur la structure à plus large échelle, si la microstructure n’est pas trop forte (Fontanillas et al., 2004).

Structure génétique spécifique à chaque sexe des données brutes (sans tenir compte de BAPS)

78Comment suspecter qu’un biais de structuration existe entre mâles et femelles ? Soit en effectuant directement le test “Sex biased dispersal” de Fstat, soit, comme cela a été le cas pour les données présentes, en testant la différenciation locale entre tiques mâles et femelles. La justification de ce test est qu’un tel signal avait été suggéré chez cette espèce en Irlande pour un locus enzymatique (Healy, 1979). Nous allons donc mesurer et tester la différenciation entre mâles et femelles de chaque sous-échantillon. Pour ce faire, il faut construire un fichier Fstat (par exemple) où femelles et mâles de chaque site sont considérés comme appartenant à des échantillons différents. Appelons ce fichier “IRTotBrutSexBias.dat” et chargeons-le dans Fstat (après avoir ouvert Fstat il faut aller dans “File”, “Open”, etc.). On coche ensuite dans Fstat l’option “Fst per pair of samples” et la case “Pairwise tests of differentiation”, ainsi que la case “5/100” du “Nominal level for multiple tests”. Cette dernière case est choisie pour gagner du temps. Ici, Fstat donnera le seuil Bonferroni de significativité corrigé par le nombre de tests (276 ici). À ce seuil, une P-value sera significative si elle est inférieure ou égale à 0,05/276 = 0,00018 et Fstat ajuste le nombre de permutations nécessaires pour atteindre cette valeur, soit 5 520, ce qui est bien suffisant. Avec “1/100” on obtient 27 600, ce qui est beaucoup. En plus, à ce niveau, le Bonferroni est beaucoup trop conservateur. De toutes façons, comme nous n’allons utiliser que les résultats par paire locale de femelles et de mâles, nous n’appliquerons pas cette procédure. Après avoir cliqué sur “Run” et attendu la fin des permutations, deux fichiers sont à consulter. “IRTotBrutSexBias.fst” donne les FST par paire et “IRTotBrutSexBias-pp.pvl” donne les P-value du test de randomisation des génotypes par paire de sous-échantillons. Dans ces fichiers, il faut garder les valeurs correspondant aux paires femelle-mâle de chaque site-année. Si vous faites cela, deux probabilités sont significatives sur les 12 tests (17 %). Un test binomial peut alors être effectué sous R avec la commande suivante :

79binom.test(2, 12, p = 0.05, alternative = "greater")

80Le test est unilatéral, car on regarde si le nombre 2 n’est pas plus grand qu’attendu sous l’hypothèse nulle. Le test donne une P-value de 0,12, ce qui n’est pas vraiment significatif, mais témoigne d’un signal possible. Le test généralisé entrepris avec MultiTest et un k’= 6 donne une P-value globale seuil de 0,6015 (La notice d’utilisation de ce programme est suffisamment détaillée pour ne pas avoir à reproduire ici un tuteurage pas à pas).

81Pour effectuer le véritable test de structuration sexe-spécifique, il faut remanier quelque peu le fichier initial des données afin de le mettre au format requis par Fstat pour l’analyse du biais de dispersion sexe-spécifique (Biased dispersal menu). Référez-vous à l’aide en ligne de Fstat pour construire ce fichier. Nous allons nous focaliser sur les échantillons 1996 uniquement. Une fois ce fichier constitué, il faut lancer Fstat, cliquer sur le menu “Biased dispersal” et y charger le fichier requis. Il faut ensuite sélectionner le test “Two sided” (on n’a en principe pas de préjugé pour l’instant) et cocher toutes les options comme dans la figure 33.

82Vous remarquerez dans la figure 33 que les cases du FIS et du H0 sont cochées comme les autres, alors que cela n’a aucun sens. En effet, puisque nous avons codé les mâles homozygotes au locus IR08, il y aura nécessairement une différence mâle femelle à ce niveau. Cependant, quand cette option n’est pas cochée, on perd une partie de l’information sur Hs dans le fichier de sortie. Il conviendra donc, dans ce fichier, d’ignorer les résultats sur H0 et FIS. Le logiciel crée cinq fichiers, trois fichiers. dat au format Fstat (les données totales, les femelles, les mâles), le fichier de permutations et le fichier. res des résultats (le plus utile). Ces derniers indiquent que les femelles sont bien mieux assignées que les mâles (IAc = 0,36 et IAc = - 0,56 pour les femelles et les mâles respectivement, P-value = 0,0005) et que les femelles sont localement mois diverses génétiquement (Hs = 0,79) que les mâles (Hs = 0,81) (P-value = 0,027), ce qui va dans le sens d’un biais de dispersion femelle (les femelles disperseraient moins). Par contre, le FST et la variance d’assignement répondent en sens inverse (mais non significativement heureusement). Pourtant, ce sont ces derniers paramètres (FST et variance d’assignement) qui doivent théoriquement signaler les premiers un biais de dispersion (qui donnent les tests les plus puissants) (Goudet et al., 2002). Nous discuterons de ce paradoxe plus loin.

Image 10000000000005C70000037FFF6251F8.jpg

Figure 33. Menu et cases à cocher dans le menu “Biased dispersal”

83Afin de tester si la Tunisie n’est pas responsable seule de ce résultat, recommençons avec les données de Suisse 1996. Dans ce cas, on a des résultats comparables avec une P-value = 0,0004 pour l’assignement, mais une P-value = 0,06 marginalement significative pour Hs. Cantonnons-nous (normal pour la Suisse) au Plateau Suisse en excluant le site Monte-Ceneri du Tessin. Cette fois, les P-values tombent à 0,0002 et 0,02 pour les assignements et Hs respectivement. En restreignant l’échantillonnage aux sites du nord-ouest de la Suisse (il faut supprimer les sites Gorges-du-Trient et Dorénaz), sans oublier de le signaler en en-tête du fichier de données (il n’y a plus que cinq sites), on obtient une confirmation de ce qui était observé (tabl. 11), mais sur une échelle plus réaliste quant aux interprétations biologiques (en fin de ce chapitre). Il semble donc bien y avoir un biais de dispersion femelle (ou à tout le moins un biais de structuration génétique en faveur de ces femelles). En retirant chaque locus un à un et en recommençant l’analyse (donc cinq traitements), vous pourrez vérifier qu’aucun locus n’est responsable à lui seul du signal. On constate même, pour les données sans IR08, que le FIS est significativement supérieur chez les mâles (tabl. 11). On pourra ici se contenter de refaire ces analyse sur les échantillons du Nord-Ouest et en unilatéral pour compenser la perte de puissance. La question qui se pose ensuite est de savoir si tenir compte des résultats de BAPS (microstructuration) change cette conclusion. Pour ce faire, il faut réanalyser les données en tenant compte des clusters définis par BAPS.

Tableau 11. Résultats du test de biais de dispersion spécifique à chaque sexe sur les cinq sites du nord-ouest de la Suisse. Excepté la variance d’assignement (s² (AIc)), tous les autres paramètres plaident en faveur d’un biais de dispersion mâle (les femelles dispersent moins), avec une P-value (tests bilatéraux) très significative pour AIc et FIS et significative pour Hs. Pour le FIS, le test (unilatéral) a été réalisé en supprimant le locus IR08

Paramètres

Femelles

Mâles

P-values

AI c

0,523

-0,786

0,0002

s²(AIc)

9,970

8,611

0,3425

FST

0,001

-0,000

0,7964

Hs

0,776

0,813

0,0224

FIS

0,422

0,506

0,0081

Structure génétique spécifique à chaque sexe des données clusterisées par BAPS

84Nous prendrons ici le fichier de données 1996 de Suisse uniquement et les clusters obtenus en ne séparant pas les mâles des femelles (évidemment). Il faudra prendre garde à ne garder que les clusters contenant au moins une femelle et un mâle, car sinon Fstat va planter (comme on dit). Nous allons dans un premier temps effectuer l’analyse sur tous les clusters de tous les sites. Le label “Pop” va donc se positionner entre chaque cluster. On peut faire le test en unilatéral, mais au vu des résultats vous verrez vite qu’il convient de repartir sur une base de tests bilatéraux. Les résultats sont en effet spectaculairement divergents des précédents (tabl. 12).

Tableau 12. Résultats du test de biais de dispersion spécifique de chaque sexe d’Ixodes ricinus dans les cinq sites du nord-ouest de Suisse en tenant compte des clusters obtenus par BAPS (en ne séparant pas les mâles des femelles) et contenant au moins une femelle et un mâle. Tous les paramètres plaident fortement en faveur d’un biais de dispersion mâle (ils dispersent moins), avec des P-values (tests bilatéraux) très significatives sauf pour s² (AIc) et FIS. Pour le FIS, le test a été réalisé en supprimant le locus IR08

Paramètres

Femelles

Mâles

P-values

AIc

-0,160

0,216

0,0067

s²(AIc)

1,222

0,691

0,0706

FST

0,219

0,338

0,0012

Hs

0,584

0,524

0,0142

FIS

0,189

0,297

0,0544

85Ce résultat, très déconcertant au premier abord, est sous très forte influence du locus IR08, bien que les autres loci répondent dans le même sens (sauf peut-être IR32). Comme il s’agit peut-être d’un phénomène local, nous allons refaire les mêmes analyses, mais dans chaque site de 1996 (y compris la Tunisie). Le résultat des tests sur le FST figure dans le tableau 13. Le signal reste le même, mais semble disparaître sans le locus IR08. Il se pourrait que ce locus soit diagnostique de certains groupes de tiques. Pour vérifier cela, il faut reprendre le fichier initial de données et grouper les individus, dans chaque site, selon leur génotype au locus IR08. Ce faisant, on recalcule sur cette nouvelle partition le FIS et le FST avec Genetix, ce qui donne 0,47 et 0,02 respectivement, alors qu’on attend un faible FIS et un fort FST. IR08 n’est manifestement diagnostique de rien du tout et le fait qu’il donne les meilleurs résultats provient vraisemblablement de sa qualité (peu ou pas d’allèles nuls et très faible variance des différents estimateurs).

Tableau 13. Résultat des tests de biais de dispersion spécifique de chaque sexe sur FST, effectués dans chaque site, entre les clusters définis par BAPS et contenant au moins une femelle et un mâle. Le test global est obtenu par une procédure binomiale généralisée et les tests sans IR08 ont été effectués de façon unilatérale (les mâles dispersent moins). Utiliser le fichier d’aide de MultiTest V. 1.2. pour une description pas à pas de la procédure à suivre pour combiner les neufs tests

Sites

Cinq loci

Sans IR08

Bern

0,3250

0,2431

Monte Ceneri

0,0817

0,2827

Dorenaz

0,3199

0,3355

Eclepens

0,1306

0,2700

Gorges du Trient

0,0159

0,6392

Montmollin

0,2422

0,9079

Neuchâtel

0,0636

0,4665

Staadswald

0,0426

0,1809

Tunisie

0,1272

0,0795

Tous (Binomial)

0,0041

0,2251

86Il y a donc manifestement un effet cluster que nous essayerons d’interpréter plus loin. Afin de vérifier quand même si notre biais de dispersion spécifique femelle existe toujours même en tenant compte de l’effet Wahlund présent au sein de chaque site, la solution qui nous reste consiste à ne garder qu’un seul représentant ou une femelle et un mâle par cluster dans chaque site (nord-ouest de la Suisse 1996). On prendra le premier des individus ayant le génotype le plus complet de chaque cluster afin de conserver le plus de puissance possible. Par exemple, si dans un cluster d’un site quelconque, il n’y a que des mâles on ne prend qu’un individu, si possible génotypé aux cinq loci. Même chose pour des clusters de femelles. Pour les clusters mixtes, on prend la première femelle la plus complète et le premier mâle le plus complet. On obtient ainsi un jeu de données de cinq sites avec un nombre d’individus fortement réduit par site. C’est aussi la raison pour laquelle les tests seront unilatéraux (les femelles dispersent moins). Le résultat de cette analyse figure dans le tableau 14 où on retrouve bien le signal initial suggérant un biais de dispersion femelle, à la différence que tous les paramètres vont dans le bon sens, même si c’est toujours AIc qui donne la seule P-value significative.

Tableau 14. Résultat du biais de structuration femelle (unilatéral) sur le jeu de données réduit à un individu ou deux (une femelle et un mâle) par cluster BAPS pour les cinq sites du nord-ouest de la Suisse. Cette fois-ci, tous les paramètres vont dans le même sens (les femelles dispersent moins). Pour le FIS, le test a été réalisé sans le locus IR08

Paramètres

F

M

P-value

AIc

0,496

-0,520

0,0097

s²(AIc)

6,377

9,350

0,3341

FST

-0,008

-0,016

0,1307

Hs

0,824

0,847

0,1221

FIS

0,470

0,511

0,2220

Interpréter l’ensemble des résultats sur les biais de structuration

87Il semble bien y avoir un biais de dispersion biaisé pour les femelles (elles disperseraient moins) à l’échelle du plateau Suisse (ou même de régions plus restreintes), mais le signal est brouillé par une micro-structuration qui existe localement. Le fait que dans chaque site, les clusters trouvés par BAPS contiennent des femelles beaucoup plus hétérogènes que les mâles à l’intérieur de chaque cluster, mais beaucoup moins différentes d’un cluster à l’autre peut être interprété de deux façons. La première suggérerait que le biais de dispersion spécifique à chaque sexe s’inverse à petite échelle, mais on ne voit pas bien comment. La seconde suppose que les clusters correspondent plus ou moins à des frères et sœurs issus d’une même ponte et que les femelles ont une réussite beaucoup plus homogène que les mâles. Ne parviendraient à l’âge adulte, selon cette hypothèse, que beaucoup de mâles par ponte, mais de peu de pontes, alors que les femelles représenteraient un échantillon plus aléatoire des pontes (moins de sœurs que de frères dans chaque site). Pour confirmer cette interprétation, une approche théorique de modélisation/simulation serait nécessaire, mais dépasserait alors le cadre ambitionné par cet ouvrage. Enfin, ces clusters pourraient correspondre à des cohortes différentes (chevauchement de générations), très différenciées (dérive forte) et cela surtout pour les mâles dont beaucoup viennent d’ailleurs. Ici aussi, une approche théorique s’avérerait nécessaire. Il est cependant raisonnable d’imaginer que si les larves et les nymphes mâles sont plus souvent retrouvées sur des hôtes très dispersants, alors il y a de fortes chances que chacun de ces individus hôtes porte des mâles apparentés (surtout les larves). Une fois dispersé et gorgé, chaque groupe a une chance très inégale de trouver un habitat favorable à la mue suivante. Il en résulterait que seuls certains groupes, parfois composés d’individus très apparentés (frères), survivraient dans une zone éloignée de leur site d’éclosion, alors que beaucoup de groupes mâles seraient éliminés. Si les larves et nymphes femelles préfèrent, quant à elles, les hôtes peu dispersants (petits rongeurs), il est probable que la survie de ces femelles soit distribuée plus aléatoirement entre femelles de pontes différentes. Ceci pourrait au final expliquer notre effet Wahlund produit en majorité par les tiques mâles.

Différenciation globale et isolement par la distance

88Plusieurs éléments nous incitent ici à manquer d’optimisme. Il y a en effet de nombreux allèles nuls, un effet Wahlund local, de la dominance d’allèles courts à un locus, sans parler d’autres problèmes mis en évidence lors d’études de pedigrees (De Meeûs et al., 2004a). Si on ajoute à cela que manifestement un biais de dispersion spécifique à chaque sexe existe, supposant qu’un des deux sexes migre beaucoup (voir Goudet et al., 2002) et donc qu’une faible structuration en résulte nécessairement, la probabilité de trouver une structuration génétique devient faible, et c’est un euphémisme. Nous allons quand même tenter notre chance, et ce pour plusieurs raisons. D’abord, parce que nous ne sommes pas arrivés jusqu’ici pour se mettre à bailler aux corneilles, ensuite parce que « c’est la nuit qu’il est beau de croire à la lumière » (Rostand, 1908).

Définir différents niveaux de subdivision pour l’analyse hiérarchique

89Nous ne considérerons ici que les échantillons de 1996. Nous pouvons envisager, grâce à HierFstat (Goudet, 2005), n’importe quelle structure du moment que cette dernière reste hiérarchique. Nous allons donc dans un premier temps considérer (référez-vous au besoin à la figure 15) l’Europe-Afrique comme tout, suivi de la Tunisie versus la Suisse, puis le Tessin versus le nord des Alpes et enfin le groupe Gorges-du-Trient, Dorénaz contre le plateau Suisse (Eclepens, Montmollin, Neuchâtel, Staadswald, Bern). Référez-vous à De Meeûs et Goudet (2007) pour des détails sur la confection d’un fichier HierFstat.

Analyse hiérarchique sur données brutes (pas de cluster BAPS)

90Il faut donc créer un fichier avec quatre (hiérarchie) plus cinq (loci) colonnes. La première colonne correspond donc au continent, Cont avec 1 l’Europe (= la Suisse, et alors ?) et 2 pour l’Afrique (Tunisie). La deuxième colonne (NrdWTessin) va coder pour l’appartenance aux cantons du nord et nord-ouest de la Suisse (1), pour celle du Tessin (2) (Monte-Ceneri) ou la Tunisie (3) qui n’est pas plus subdivisée, mais doit être aussi codée dans cette colonne. La troisième colonne (NrdWNS) correspond à l’appartenance ou non au nord-ouest (1) ou au sud-ouest (Gorges-du-Trient, Dorénaz = 2) de la zone du nord des Alpes suisses. Le Tessin et la Tunisie étant codés 3 et 4 respectivement dans cette colonne. La quatrième colonne (Site) correspond aux sites eux-mêmes (1 à 9). Les cinquième à neuvième colonnes correspondent aux cinq loci, le premier, IR08, étant codé homozygote pour les tiques mâles. Appelons le fichier ainsi construit “IRTot96HierFstat. txt”. L’analyse va se faire sous HierFstat 0.04-4 (Goudet, 2006, mis à jour de Goudet, 2005) comme décrit dans De Meeûs et Goudet (2007). N’oubliez pas de remplacer les données manquantes “000000” par “NA”. Lancez le logiciel R. Chargez le package HierFstat (Menu “Package”, “Chargez le package”, “hierfstat”). Changez de répertoire pour travailler dans celui où le fichier de données “IRTot96HierFstat. txt” se trouve (Menu “Fichier”, “Changer le répertoire courant”). Dans la console R, tapez la succession de commandes (chaque ligne correspond à une commande devant être suivie d’un retour charriot), en respectant les majuscules et minuscules (distinctes en langage R) :

>data<-read.table("IRTot96HierFstat.txt", header=TRUE)
>attach(data)
>locic-data.frame(IR08,IR25,IR27,IR32,IR39)
>levels<-data.frame(Cont,NrdWTessin,NrdWNS,Site)
>varcomp.glob(levels,loci)

91Cette dernière commande produit le résultat suivant :

Image 10000000000002F6000001B311ED29D1.jpg

92Dont l’interprétation est la suivante :

93FIS = 0,459 (nous retrouvons ici un résultat ancien et sans valeur, car les mâles sont artificiellement homozygotes ici au locus IR08), FSite/NrdWNS = 0,0028, FNrdWNS/NrdWTessin = 0,0025, FNrdWTessin/Cont = 0,0057 et FCont/Total = 0,106. Toutes ces valeurs de différenciation sont très faibles sauf pour la Suisse contre la Tunisie. Il faut tester ensuite la significativité de ces différentes partitions en commençant par la plus incluse, le site :

>test.within(loci, test=Site, within=NrdWNS, nperm=1000)
$p.val
[1] 0.311

94On voit bien que le site (comme on le craignait) n’influence en rien la partition de l’information génétique. Nous allons donc supprimer ce facteur de la hiérarchie :

>levels<-data.frame(Cont, NrdWTessin, NrdWNS)
>varcomp.glob(levels, loci)
Sloc

Image 10000000000002F00000015493EB0BB2.jpg

>test.within(loci, test = NrdWNS, within = NrdWTessin, nperm = 1000)
Sp.val
[1] 0.121

95Le facteur NrdWNS, séparant les sites Dorénaz-Gorges-du-Trient de l’ensemble des sites suisses du Nord-Ouest, ne semble pas influencer davantage la structure génétique des tiques. Si nous le supprimons à son tour, nous obtenons :

>levels<-data.frame(Cont, NrdWTessin)
>varcomp.glob(levels, loci)
Sloc

Image 100000000000024000000087CDF9CEBA.jpgImage 10000000000001F5000000BCFD21B4D5.jpg

>test.within(loci, test = NrdWTessin, within = Cont, nperm = 1000)
Sp.val
[1] 0.058

96Si nous choisissons de garder le facteur NrdWTessin (marginalement significatif, P-value = 0,058) cela aboutit à :

>test.between(loci, rand.unit=NrdWTessin, test=Cont, nperm=1000)
$p.val
[1] 0.331

97Si on élimine le facteur NrdWTessin, il faut alors repasser par Fstat. Il n’y a en effet plus que trois niveaux hiérarchiques avec deux sous-populations représentées par l’ensemble des tiques suisses, d’une part et par celles de Tunisie, d’autre part. On aboutit à un FST = 0,113 très significatif (P-value < 0,0001) entre les tiques de Suisse réunies en une seule population et la Tunisie.

98Avec un Hs = 0,832, cela correspond à un FST’= FST/FSTmax = 0,113/(1 - 0,832) = 0,673, ce qui est relativement considérable et témoigne du peu de migration entre les deux pays. Par contre, à l’échelle de la Suisse, cette migration est forte et même si les Alpes apparaissent comme un facteur limitant, tout semble se passer comme si, génétiquement au moins, on avait à faire à une seule unité à cette échelle.

99Qu’en est-il si nous tenons compte des clusters trouvés par BAPS ?

Analyse hiérarchique sur données clusterisées par BAPS

100Nous allons donc utiliser le fichier de données précédent avec une colonne supplémentaire correspondant aux clusters trouvés avec BAPS. En suivant alors une procédure identique à celle décrite plus haut, nous pouvons constater que les facteurs ClusterBAPS (FClust/Site = 0,3, P-value = 0,001) qui mesurent la partition génétique entre clusters d’un même site, et Continent (FContinent/Total = 0,11, P-value = 0,001) qui mesure la différenciation entre Suisse et Tunisie, constituent les deux seuls facteurs qui structurent les sous-échantillons de façon significative.

101Si nous ne gardons qu’un mâle ou une femelle ou un individu par cluster, comme pour le tableau 14, le résultat de l’analyse par HierFstat ne montre plus aucune différenciation, à moins d’ignorer tous les facteurs sauf le continent (analyse par Fstat, FST = 0,09, P-value = 0,001).

Test d’isolement par la distance

102Nous ne travaillerons ici que sur les échantillons suisses de 1996. D’abord parce que la Tunisie est trop éloignée par rapport aux distances entre échantillons suisses. Il y aurait deux groupes de points. Procéder à un test de régression entre deux points n’a pas de sens, le plus court chemin entre eux étant nécessairement une droite, c’est dans tous les bons livres de statistiques. Or, le test d’isolement par la distance est une forme de régression où on cherche à expliquer une différence génétique croissante par un éloignement géographique. Ensuite, il n’y a pas assez d’échantillons en 1995.

103Pour le test, il faut configurer un fichier avec deux demi-matrices, l’une pour les distances géographiques entre paire de sites et l’autre pour les FST (estimés par θ) correspondants. Pour les distances géographiques, vous pouvez vous aider de la figure 15. Pour les FST, il suffit de prendre la sortie “IRTot96CH. fst” que Fstat a produit en analysant le fichier “IR96CH.dat” des données suisses 1996, si vous avez toutefois coché la case “Fst per pair of samples”. En ce qui me concerne, j’obtiens les matrices représentées dans le tableau 15. Le test va être effectué selon la méthode décrite par Rousset (1997) pour un schéma en deux dimensions. Nous allons donc effectuer un test de Mantel sur la corrélation entre le FST/(1 - FST) et le log népérien (ou naturel) de la distance géographique. Nous allons utiliser Genepop 3 pour faire ce test et donc formater les données dans ce sens et les sauvegarder dans un fichier que nous appellerons IR96CH. mig. Ce fichier doit être configuré comme présenté dans la figure 34.

Image 100000000000049F0000032C26B69CA2.jpg

Figure 34. Présentation du fichier pour tester l’isolement par la distance entre sites de prélèvement d’Ixodes ricinus en Suisse en 1996 sous Genepop

Tableau 15. Distances géographiques en km et différenciation génétique mesurée par le FST (Theta) par paire de sites d’échantillonnage d’Ixodes ricinus (abréviations comme dans la figure 15)

Image 10000000000005B20000051410F446F6.jpg

104L’étape suivante consiste à lancer Genepop 3. Éviter de double cliquer sur le fichier genepop. bat, mais préférez ouvrir une session DOS en lançant une “Invite de commandes” dans le menu “Accessoires” de Windows. Dans la fenêtre DOS, et si Genepop est dans le répertoire “Genepop” du disque D, tapez “D :”, puis “Entrée”, puis “cd Genepop”, puis “Entrée”. Vous êtes dans le répertoire Genepop. Tapez alors “isolde”, puis “Entrée” pour lancer le programme d’isolement par la distance. À l’invite, tapez le nom complet du fichier de données puis “Entrée”. Le logiciel vous demande ensuite quel type de distance (non transformée ou Log) et quel type de mesure de différenciation vous souhaitez tester (X, qui figure dans la matrice ou X/(1 - X)). À vous de choisir la méthode appropriée. Genepop vous demande ensuite la distance minimale en deçà de laquelle la mesure de corrélation ne tient plus compte des données, car en deçà d’un certain niveau la réponse a en effet tendance à ne plus suivre un modèle clair (Rousset, 1997). Réfléchissez à ce que devrait être cette distance minimale et tapez-la. Le nombre de randomisations vous est ensuite demandé. Tapez 1 000 000 pour être sûr d’obtenir une P-value suffisamment précise. Genepop vous demande, comme pour tous ses tests passant par randomisation, quatre nombres qui serviront de « graines » pour générer des nombres pseudoaléatoires qui conditionnent le départ des randomisations. Tapez ce que vous voulez entre 1 et 168, comme indiqué avec un retour charriot après chaque chiffre. Quand les randomisations sont terminées, Genepop le signale avec un bip assez désagréable, mais qui ne doit pas vous effrayer (en général, je coupe le son avant). D’après une collègue avisée (TG), il n’y a pas de son sur la version Web du programme qui, par ailleurs, est sensiblement plus conviviale. Genepop a alors créé deux fichiers, l’un porte l’extension. ISO qui vous donne le résultat du test avec les paramètres de la régression et la P-value. Le second porte l’extension. GRA et donne les coordonnées en colonne de cette régression afin de pouvoir dessiner le graphique correspondant, comme représenté dans la figure 35. On y remarque que la relation n’est pas significative. Il semble cependant qu’une tendance existe. Peut-être l’existence d’une sous-structure nuit-elle à la clarté du signal ? Pour le vérifier, il suffit de procéder au même test, mais avec les données où un ou deux (de chaque sexe) individus par cluster avaient été gardés (voir p. 151-152). Le résultat change du tout au tout puisque la relation devient très significative, comme en témoigne la figure 36. Ceci permet de calculer le voisinage Nb = 1/ b = 173 individus (Watts et al., 2007), le produit de la densité d’individus par km² par la surface de dispersion des descendants reproducteurs par rapport à leurs géniteurs, en utilisant la méthode de Rousset (1997) décrite en première partie (p. 88-90), ce qui donne² = 1/(4π0,00577) = 13,78. On peut aussi estimer le nombre d’immigrants présents dans une sous-population (Rousset, 1997), Nm = 1/2πb = 28 individus par génération. Il serait dommage de s’arrêter ici et nous allons donc essayer d’appréhender, même grossièrement, dans quelle gamme de valeurs se situe la densité de tiques afin d’en pouvoir extraire la surface de dispersion de ces tiques.

Image 10000000000001F500000193A8E2470F.jpg

Figure 35. Représentation de l’isolement par la distance entre sites suisses pour les tiques récoltées en 1996. La distance minimale de 2 km a été choisie

Image 100000000000020E0000019808397B44.jpg

Figure 36. Représentation de l’isolement par la distance entre sites suisses pour les tiques récoltées en 1996 en ne gardant qu’un ou deux individus des clusters définis par BAPS dans chaque site. La distance minimale de 2 km a été choisie

Estimation d’effectifs efficaces, extrapolation des densités et de la dispersion

Effectifs efficaces des tiques de Suisse

105Pour des raisons de commodité et de cohérence (les tiques tunisiennes n’ayant pas été échantillonnées de la même manière), nous nous focaliserons sur les échantillons de Suisse. Rappelons-nous que nous avons rencontré de gros déficits en hétérozygotes (allèles nuls et dominance d’allèles courts), ainsi que la présence d’un fort effet Wahlund. Nous ne travaillerons donc qu’à partir de méthodes indépendantes de l’hétérozygotie, telle que celle proposée par Bartley et al.(1992), basée sur les déséquilibres de liaison et implémentée par NeEstimator (Peel et al., 2004). Nous n’utiliserons que les données clusterisées par BAPS où seuls subsistent un ou deux individus par cluster dans chacun des huit sites suisses pour éviter l’effet confondant dû à l’effet Wahlund. Il faut créer un fichier par site dans un format proche de Genepop comme dans la figure 37.

Image 100000000000049F00000640D53B408B.jpg

Figure 37. Format de fichier pour NeEstimator pour les tiques de Berne 1995

106Il faut ensuite lancer le programme NeEstimator (après l’avoir installé sur votre machine, bien entendu). Une fenêtre d’avertissement sur le copyright et sur la manière idoine de citer ce logiciel apparaît. Cliquez sur OK pour accéder au programme qui apparaît dans une fenêtre comme dans la figure 38. Comme indiqué sur la figure 38, cliquez sur le menu déroulant “File” et “Open”, ce qui permet d’ouvrir la fenêtre “Analysis”.

107Dans la fenêtre “Analysis”, une série d’onglets apparaît et vous positionne sur celui du format de vos données “Data Format” où il n’y a rien à changer, car vous avez choisi le format par défaut. Allez à l’onglet “Data Files”. Là il n’y a qu’un seul bouton “Load” qui vous permet de charger votre jeux de données, ce que vous faites (fig. 39). Une fois que vous avez choisi le fichier, le logiciel vous demande à quelle génération ces données correspondent-elles. Laissez la valeur par défaut “0”, car nous n’utiliserons pas ici la méthode des moments de Waples (1989) (cf. p. 102 en première partie) et cliquez sur “OK”. Dans le menu déroulant “NeEstimator”, cliquez sur “Run” (fig. 40). Ce qui fait apparaître un message qui vous avertit qu’avec un seul échantillon, on ne peut utiliser les méthodes temporelles “Moment based” et vous demande si vous souhaitez continuer avec les méthodes à un seul échantillon. Vous répondez “Oui” bien entendu. Le résultat est affiché sous forme de tableau que je vous conseille de sauvegarder au format NeEstimator (NeA). Je conseille aussi de transcrire tous les résultats dans un tableur au fur et à mesure afin de disposer de l’ensemble dans un seul fichier. C’est ce qui est représenté dans le tableau 16.

Image 1000000000000599000003344316D062.jpg

Figure 38. Menu à l’ouverture de NeEstimator

Image 1000000000000599000004695A8A0C9B.jpg

Figure 39. Onglet “Data File” avec le bouton “Load” qu’il faut clicker

Image 10000000000002F0000001B159EF129A.jpg

Figure 40. Pour lancer l’analyse de NeEstimator

108Ici, bien que nous disposions d’échantillons espacés dans le temps (Bern, Gorges-du-Trient et Staadswald), ces échantillons ne sont séparés que d’une année, soit environ 1/3 du temps de génération d’I. ricinus. Ici, les adultes présents d’une année sur l’autre font partie de cohortes séparées et qui, même à long terme, auront du mal à échanger des gènes. La différenciation entre ces cohortes, déjà remarquée par De Meeûs et al. (2002a), va tendre à être très supérieure à celle qui existe réellement entre deux générations d’adultes reproducteurs. L’utilisation des méthodes temporelles sur nos données aboutira donc à de fortes sous-estimations des effectifs efficaces. Faites-le et vérifiez qu’effectivement, compte tenu qu’il n’y a qu’un tiers de génération séparant 1996 de 1995, les estimations obtenues par la méthode de Waples donnent des effectifs efficaces proches de 0, ce qui n’est pas très conforme à la perception que l’on peut avoir sur le terrain.

Tableau 16. Résultats synthétiques obtenus pour le calcul des effectifs efficaces (Ne) et leur intervalle de confiance à 95 % (Li et Ls) par la méthode des déséquilibres de liaison dans NeEstimator. Les valeurs infinies sont ignorées pour le calcul des moyennes. Les échantillons de 1995 sont considérés comme indépendants, car appartenant à des cohortes de tiques génétiquement isolées de celles de 1996 (le cycle d’Ixodes ricinus dure environ trois ans)

Échantillon

Ne

Li

Ls

Berne 1996

73

45

182

Berne 1995

222

79

Infini

Monte-Ceneri 1996

Infini

288

Infini

Dorénaz 1996

700

124

Infini

Eclépens 1996

Infini

81

Infini

Gorges-du-Trient 1995

177

10

601

Gorges-du-Trient 1996

75

43

219

Montmollin 1996

338

87

Infini

Neuchâtel 1996

398

93

Infini

Staadswald 1995

161

84

1 164

Staadswald 1996

Infini

374

Infini

Moyenne totale

268

119

541

109En reprenant le tableau 16, nous obtenons par conséquent un effectif efficace de 268 en moyenne sur l’ensemble des échantillons avec un intervalle de confiance à 95 % de [119, 541], avec des valeurs minimales et maximales de 73 et 700 respectivement. Ces nombres paraissent plausibles, compte tenu de l’effet Wahlund reflétant probablement un fonctionnement particulier des populations de tiques susceptible d’en réduire sensiblement l’estimation de leurs effectifs efficaces.

  • 11 Sur ces mêmes données, l’estimation avec un logiciel alternatif, LDNe (WAPLES et DO, 2008), non enc (...)

110En reprenant les données avec un ou deux individus par cluster BAPS, les valeurs obtenues sont plus grandes en moyenne (596) avec un minimum et un maximum de 75 et 1 057 respectivement11.

Extrapolation des densités et des distances de dispersion des tiques en Suisse

111Il faut dans un premier temps estimer sur quelle surface se distribuent les tiques. Ici, c’est difficile et on ne peut pas dire grand-chose de plus que les surfaces d’échantillonnage s’étendaient grossièrement sur S = 0,2km². Ceci signifie (mais vous vous en doutiez probablement) que les estimations à venir seront tout à fait approximatives. À partir de là, les densités sont faciles à calculer (Ne/S). La densité moyenne devient 1 340 tiques reproductrices/km² 95 % CI = [594, 2 706] avec un minimum et un maximum de 367 et 3 502 tiques/km² respectivement (tabl. 16). En réutilisant les résultats de la régression de l’isolement par la distance Deσ² = 13,78 (voir p. 160), on aboutit à une surface de dispersion moyenne entre adultes et leurs parents d’environ 0,01 km² [0,005, 0,023] avec un minimum et un maximum de 0,004 et 0,038 km² respectivement. Autrement dit, la distance moyenne séparant un adulte reproducteur de ses géniteurs est d’un ordre de grandeur de 100 m par génération (donc tous les trois ans environ), un intervalle de confiance à 95 % de bootstrap = [71, 152] et un maximum et un minimum de 63 à 195 m, ce qui est relativement modeste. Les données clusterisées par BAPS conduisent à une densité de 3 000 tiques par km² et une dispersion de moins de 60 m par génération. Donc, sachant que l’estimateur sans doute le moins biaisé est le produit Deσ², la dispersion par génération est, quoi qu’il en soit, extrêmement modeste à moins d’évoquer des densités (effectifs) efficaces extrêmement faibles. Il en va donc de même en ce qui concerne la propagation des maladies par les tiques.

CONCLUSIONS SUR LA BIOLOGIE ET LA GÉNÉTIQUE DES POPULATIONS D’I. RICINUS EN SUISSE

112Il existe un déficit important en hétérozygotes dans les populations d’I. ricinus (FIS = 0,39) dont une majeure partie (64 %) est expliquée par un effet Wahlund important.

113Le FIS = 0,14 résiduel correspondrait à du « stuttering », à de la dominance d’allèles courts et à des allèles nuls. Pour tester les allèles nuls dans les clusters de BAPS, on ne peut pas utiliser Micro-Checker (échantillons trop petits). Nous pouvons néanmoins tester s’il existe une relation positive entre le nombre de blancs à un locus et le FIS à ce locus. En effet, en reprenant les données clusterisées et en séparant les mâles des femelles en deux fichiers, il est facile de compter les blancs pour chaque locus avec la fonction “SI” d’Excel. Il suffit de créer autant de nouvelles colonnes qu’il y a de loci et de remplir chacune avec les instructions de type “= SI (G2 = “000000”) ; 1 ; 0)” pour inscrire “1” quand on a un blanc. À la fin de chacune de ces colonnes, on tape une instruction du type “= somme (L2 : L147)” pour obtenir la totalité des blancs à ce locus sur l’ensemble des clusters. Le FIS de chaque locus est récupérable dans les deux fichiers de sortie Fstat de l’analyse des deux jeux de données (un pour le femelles et un pour les mâles) avec les données clusterisées par BAPS, que j’ai personnellement nommés IRTotBAPSClustMalManqIR08Females.dat et IRTotBAPSClustMalManqIR08FMales.dat respectivement, et où on aura pris soin d’éliminer le locus IR08 du fichier des mâles. Quand on a fait ceci pour les femelles et les mâles, on obtient le jeu de données présenté dans le tableau 17. La corrélation entre le nombre de blancs et le FIS peut être analysée par un test de corrélation de Spearman (test non paramétrique). Ce test est facile à réaliser sous R. Si le fichier de données correspondant au tableau 17 s’appelle “AllelesNulsClustersBAPS.txt”, alors il suffit de lancer R, et de se placer dans le répertoire contenant ce fichier (menu déroulant “Fichier”, “Changer le répertoire courant”).

Tableau 17. Données pour la régression entre le nombre de données manquantes (génotypes « blancs ») et la valeur des FIS pour les différents loci (chez les mâles et les femelles pris séparément)

Sexe

Locus

Blancs

FIS

Femelles

IR08

10

-0,030

IR25

50

0,256

IR27

22

0,201

IR32

47

0,253

IR39

45

0,076

Mâles

IR25

51

0,368

IR27

21

0,010

IR32

74

0,473

IR39

30

0,115

114Ensuite, il faut taper les instructions suivantes :

>data<-read.table("AllelesNulsClustersBAPS.txt",header=TRUE)
>attach(data)
>cor.test(data$NBlancs, data
$FIS, alternative="two.sided",
method="spearman")

Image 1000000000000448000002CE8C6F8DFA.jpg

Figure 41. Relation entre le nombre de génotypes blancs trouvés par locus et le FIS mesuré à ce locus sur l’ensemble des clusters de BAPS de l’ensemble des données microsatellites d’Ixodes ricinus (1995-1996, Tunisie et Suisse)

115ce qui renvoie au résultat :

Spearman's rank correlation rho
data: data$NBlancs and data$FIS
S = 8, p-value = 0.0007496
alternative hypothesis: true rho is not equal to 0
sample estimates:
rho
0.9333333

116La figure 41 illustre la relation positive forte entre les deux variables.

117Nous pouvons également tester de nouveau la dominance d’allèles courts au locus IR27 en prenant les FIS par allèle donnés par la sortie Fstat sur les mâles et les femelles séparément. Nous ne distinguerons en effet pas les clusters, car il y en a trop. Sous R, la procédure est comparable à celle utilisée en p. 129-134, sauf que nous n’utiliserons que le polynôme d’ordre deux de la taille des allèles et le sexe comme variables explicatives du FIS. Le résultat est de nouveau très significatif et on explique donc toujours une forte proportion du FIS par ce phénomène de dominance des allèles courts, comme illustré par la figure 42. Suivent les instructions R :

>datac-read.table("FISAlleleSizeIR27ClustersBAPS.txt",header=TRUE)
>attach(data)
>loc27<-glm(data, formula = Fis
~ poly(Allele, 2) + Sexe, family = gaussian)
>anova(loc27, test="F")

Image 1000000000000440000003069DCAECE8.jpg

Figure 42. Régression entre taille des allèles et FIS observés au locus IR27 dans les clusters d’I. ricinus définis par BAPS. Il n’y avait pas assez de données pour calculer des intervalles de confiance

118ce qui renvoie au résultat suivant :

Analysis of Deviance Table
Model: gaussian, link: identity
Response: Fis
Terms added sequentially (first to last)

Image 10000000000002E90000005CD1F8F018.jpg

119Comme cela a été vu au début de ce paragraphe, la majeure partie (64 %) du FIS initial est expliquée par un effet Wahlund. Cet effet Wahlund est plus prononcé chez les mâles dont les clusters se trouvent plus différents entre eux que les femelles. Dans l’hypothèse de spécificités d’hôtes différentes des larves et/ou nymphes femelles et mâles, des groupes de larves ou nymphes mâles fortement apparentés seraient transportés ensemble sur le même hôte très dispersant (oiseau) avec de forts risques de tomber sur des sites défavorables lors du détachement, à la fin du repas sanguin. Les mâles retrouvés adultes dans nos échantillons correspondraient alors aux quelques groupes d’apparentés ayant eu la chance de tomber ensemble dans un site favorable. Les larves ou nymphes femelles seraient, quant à elles, plus souvent retrouvées sur des hôtes très peu dispersants, comme des petits rongeurs très territoriaux. Il en résulterait un apparentement réparti beaucoup plus aléatoirement pour les femelles dans chaque site. Il y a un fort biais de dispersion spécifique à chaque sexe (les femelles dispersent très peu). Ce biais est partiellement masqué par l’effet Wahlund, et il est plus facilement visible quand cet effet est corrigé (données réduites), et l’indice d’assignement corrigé AIc semble à cet égard beaucoup plus robuste que sa variance vAIc et le FST.

120Cet effet Wahlund nuit considérablement à l’image perçue au niveau de la structuration à l’échelle de la Suisse. Quand cet effet est contrôlé (au moins en grande partie), on observe un isolement par la distance très significatif, et les adultes non gorgés d’I. ricinus paraissent distribués en populations locales de tailles importantes (plus de 1 000 tiques par km²) et se dispersant difficilement à plus de 200 m par génération.

121Il reste cependant bien d’autres questions et toutes ces hypothèses doivent être testées sur le terrain. Cette étude ouvre de nombreuses et prometteuses perspectives de recherche que je vous laisse le soin de discuter.

INTERACTIONS AVEC LES MICROPATHOGÈNES TRANSMIS

Introduction

122La tique I. ricinus transmet un très grand nombre de pathogènes à ses multiples hôtes, dont la borréliose de Lyme qui, dans les régions boréales, représente un poids économique et en santé publique important (Gubler, 1998). Les agents de la borréliose de Lyme appartiennent au complexe d’espèces Borrelia burgdorferi sl. Il existe actuellement 12-13 espèces (ou génoespèces) reconnues dans le complexe dont 4, B. burgdorferi, B. garinii, B. afzelii et B. spielmanii, sont pathogènes pour l’homme (borréliose de Lyme) et deux autres, B. Valaisiana et B. lusitaniae, sont de pathogénicité incertaine (Richter et al., 2006 ; Postic et al., 2007). Ces différentes espèces ne sont d’ailleurs pas responsables de symptômes identiques et présentent des spécificités d’hôtes réservoirs différentes (De Meeûs et al., 2004b). En Europe de l’Ouest, B. burgdorferi est préférentiellement retrouvée chez l’écureuil roux, B. afzelii chez des campagnoles, des mulots et aussi l’écureuil roux, B. garinii plutôt chez des oiseaux et B. spielmanii uniquement chez le loir (Richter et al., 2006). Quant à B. valaisiana, elle se montre spécifique des oiseaux, alors que B. lusitanaiae n’a été trouvée, en dehors de la tique I. ricinus, que chez un lézard (Postic et al., 2007). L’épidémiologie de ces pathogènes reste largement mal connue et les résultats obtenus précédemment par nos analyses suscitent un certain nombre de questions. S’il y a spécificité différente des tiques immatures, sachant que les borrélies sont spécifiques des hôtes, les tiques des deux sexes devraient présenter des prévalences différentes pour les différentes espèces de borrélies. En particulier, les femelles devraient porter davantage de borrélies d’hôtes peu mobiles (B. burgdorferi, B. afzelii) et les mâles celles d’hôtes plus mobiles (B. garinii, B. valaisiana). Ensuite, il est possible que l’infection par les borrélies puisse modifier le schéma de migration. Enfin, dans la mesure où un conflit/coopération pourrait exister au sein des tiques, existe-t-il une corrélation entre la présence des différentes espèces de borrélies au sein de tiques ?

Présentation des données

123Toutes les tiques échantillonnées en Suisse pour cette étude avaient été coupées en deux, et une moitié envoyée à l’Institut de zoologie de Neuchâtel pour détermination de présence de borrélies et détermination de l’espèce (sondes moléculaires). L’autre moitié a été gardée dans l’alcool et un grand nombre utilisé pour génotypage microsatellite. Les données sont contenues dans le fichier TotBrutBorIR.txt où toutes les informations nécessaires sont disponibles. La présence ou l’absence de chaque espèce de borrélie trouvée est notée par un 1 ou un 0 dans la colonne correspondante. Un grand nombre de borrélies n’ont pu être déterminées au niveau de l’espèce (colonne “Bbundet”) et seules trois espèces ont été trouvées : B. burgdorferi (Bbss), B. afzelii (Bba) et B. garinii (Bbg, trouvée trois fois).

Distribution des différentes borrélies dans les femelles et mâles d’I. ricinus

124Pour cette analyse, nous allons devoir effectuer une régression logistique pour chaque espèce de borrélie (Bbundet, Bbss, et Bba). Bbg, trop rare sera laissée de côté. On va chercher à expliquer la présence de telle ou telle autre espèce de borrélie par le site, l’année et le sexe de la tique, ainsi que les interactions. Nous allons donc avoir besoin de R une fois de plus. Comme c’est le sexe que l’on souhaite tester ici, nous allons mettre ce facteur en premier (l’ordre compte dans les modèles de R). Après avoir lancé R et s’être positionné dans le répertoire approprié, on tape les commandes suivantes :

> data<-read.table("TotBrutBorIR.txt", header=TRUE)
> attach(data)

125afin de faire lire l’ensemble du jeu de données à R (NB le > est automatiquement inséré par R). On spécifie ensuite le modèle en tapant la commande (sur une ligne) :

> Bba<-glm(data, formula =Bba ~ Sex + Site + Year + Sex:Site + Sex:Year +
Sex:Site:Year, family = binomial(link = logit))

126On remarque que l’interaction entre facteurs est codée avec un “ :” et que la régression est logistique, car on spécifie bien qu’elle appartient à la famille binomiale avec un lien “logit” de la moyenne. Le lien logit signifie juste que la fonction qui relie la probabilité moyenne de la variable à expliquer (PBba probabilité de trouver une Bba) est du type log (PBba/(1 - PBba)) et la variance égale à PBba/(1 - PBba). Dans notre cas, la variance est en fait inférieure à cette valeur et il y a sous-dispersion, ce dont nous discuterons plus loin.

127Ensuite, il s’agit de tester le modèle par la commande :

> anova (Bba, test ="Chi")

128Le test est en effet un Chi2, car nous comparons des fréquences. Cette commande renvoie au résultat suivant :

Analysis of Deviance Table
Model: binomial, link: logit
Response: Bba
Terms added sequentially (first to last)

Image 10000000000002D6000000B0242E26BC.jpg

Warning message:
In method(x = x[, varseq <= i, drop = FALSE], y = object$y, weights =
object$prior.weights, :
des probabilités ont été ajustées numériquement à 0 ou 1

129Nous constatons que seuls les termes “Site” et “Year” semblent importer et que le logiciel n’est apparemment pas très satisfait de la qualité des données. Pour simplifier ce modèle, une commande pratique est la commande “step” qui permet d’analyser la qualité de différents modèles plus simples en retirant et ajoutant des termes l’un après l’autre en commençant par les interactions d’ordre supérieur (celles faisant appel au plus grand nombre de facteurs). Ceci est évalué à l’aide d’un critère appelé AIC (Akaike Information Criterion) (Akaike, 1974) dont la valeur, qui doit être minimisée, est une mesure de la qualité d’ajustement du modèle statistique estimé par rapport aux données. Il ne s’agit pas d’un test, mais d’un outil d’aide à la sélection du modèle le plus simple permettant d’expliquer au mieux les données, le modèle doté du plus petit AIC étant le meilleur (cf. réponse 12 pour plus de précisions). En tapant donc la commande :

> step (Bba)

130nous obtenons les résultats pour une série de différents modèles de plus en plus simples où les différents termes sont retirés un à un en commençant par l’interaction la plus complexe (Sex : Site : Year), qui est éliminée, l’AIC obtenu (338,69) s’avérant inférieur à celui du modèle complet (343,81), puis les interactions plus simples (Sex : Site et Sex : Year), jusqu’à ce que le retrait des facteurs conduisent à une augmentation de l’AIC par rapport au précédent. Ci-dessous sont présentés le début et la fin du processus :

Image 10000000000002F30000011E89ECDA89.jpg

131etc.

Image 10000000000002F0000000DC537A1CDE.jpg

132La dernière ligne présentée ci-dessus donne le meilleur modèle. Suivent des informations sur les coefficients associés aux différents facteurs que nous n’allons pas utiliser, ainsi que des messages d’alertes sur la mauvaise qualité des données (on ne fait pas de miracles). Il s’agit maintenant d’analyser en détail ce meilleur modèle avec la série d’instructions (pour gagner du temps on peut copier le modèle ci-dessus et le coller après avoir tapé "Bba2 < -") :

> Bba2<-glm(formula = Bba - Site + Year, family = binomial(link = logit),
data = data)
> anova(Bba2, test="Chi")

133qui renvoie au résultat :

Analysis of Deviance Table
Model: binomial, link: logit
Response: Bba
Terms added sequentially (first to last)

Image 10000000000002370000005AC67F7CB8.jpg

134La conclusion est donc qu’en ce qui concerne Bba, seuls le site et l’année importent. Ils expliquent respectivement 100 × 35,72/358,68 = 10 % et 100 × 8,86/358,68 = 2 % de la déviance totale. En procédant d’une manière identique pour Bbg, nous observons qu’aucune des variables n’explique les données alors que pour Bbss, en plus du site qui explique 28 % de la déviance totale (P-value < 0,001), le sexe des tique explique 3 % de la déviance (P-value = 0,007). Enfin, pour Bbundet le site seul explique 15 % de la déviance totale (P-value < 0,001).

135Comme je l’ai déjà signalé plus haut, la dispersion des résidus ne suit probablement pas une loi binomiale et la variance est probablement différente de P/(1 - P). Pour vérifier cela, il faut calculer le paramètre φ = Var (µ)×(1-µ)/µ qui est ici inférieur à 1 (sous-dispersion) en particulier pour Bbss. On peut le calculer facilement avec la fonction "quasibinomial" (voir réponse 13). Comme seul Bbss a donné quelque chose de significatif pour le sexe des tiques, nous allons vérifier cela sur cette bactérie. Sous R, après avoir chargé le fichier de données si ce n’est déjà fait, nous allons taper les instructions suivantes :

> Bbss<-glm(data, formula =Bbss ~ Sex + Site, family =quasibinomial(link

= "logit"))
> summary(Bbss)

136ce qui renvoie au résultat suivant (je ne garde que ce qui est le plus utile) :

Image 10000000000005D300000285643CB4CE.jpg

137Nous voyons donc que le le coefficient de dispersion est petit (0,52), il y a donc bien sous-dispersion (pour Bba ϕ = 0,93, il n’y a pratiquement pas de sous-dispersion pour cette espèce-ci). Nous voyons également que le sexe des tiques est important (significatif) avec une estimation négative pour les mâles (les données partielles, corrigées des autres effets, sont centrées sur 0). Ceci est vérifiable en tapant la commande anova (Bbss, test = "F") (les modèles quasi se testent avec un F), ce qui donne :

Image 10000000000002F200000077DFAAB433.jpg

138Sachant que le comportement des modèles quasi en régression logistique peut s’avérer étrange quand l’événement étudié (présence de Bbss) est rare, ce qui est notre cas, on est en droit de chercher à renforcer ce résultat. En fin de compte, nous cherchons juste à vérifier si nous n’avons pas plus de Bbss chez les tiques femelles que chez les mâles, puisque ces borrélies sont spécifiques de petits rongeurs peu dispersants, supposés être davantage parasités par les larves et nymphes femelles que mâles, quel que soit le site ou l’année. On peut donc calculer parmi les tiques infectées par Bbss, la proportion de tiques femelles et mâles et comparer cette proportion à ½ par un test binomial. Sur 34 tiques infestées par Bbss, 26 étaient femelles, ce qui conduit à la P-value du test binomial (sous R, binom. test (26, 34, p = 0.5, alternative ="greater")) Pbino_26/34,0.5 = 0,0015, ce qui est équivalent aux résultats précédents. Vous vous demandez alors pourquoi vous ai-je cassé les pieds avec toutes ces régressions, alors qu’il était si simple de commencer par le test binomial ? La réponse est simple. D’abord, il n’est pas inutile d’apprendre à taquiner les régressions linéaires généralisées qui servent très souvent et, ensuite, dans une publication, une régression logistique en « quasi-likelihood » va avoir beaucoup plus de classe (en apparence) qu’un petit test binomial et impressionner beaucoup plus facilement ces referees désobligeants qui empoisonnent si souvent nos soumissions d’articles.

139Donc Bbss, borrélie d’écureuil en Suisse, est plus fréquente chez les tiques adultes femelles que mâles, suggérant ainsi une prédisposition de ces femelles à se nourrir sur cet hôte quand elles sont aux stades larvaire et/ou nymphal.

Co-occurrence des différentes espèces de borrélies

140Les différentes espèces de borrélies peuvent se retrouver en compétition, car elles partagent la même espèce de vecteur. Elles sont donc potentiellement en conflit et on pourrait s’attendre à un évitement. Au contraire, il pourrait y avoir association positive si les intérêts convergent ou si l’une des deux espèces immunodéprime ses hôtes et favorise ainsi l’entrée d’autres pathogènes. Il est donc intéressant de tester si ces borrélies se rencontrent au hasard ou non. La problématique est identique à une recherche d’association statistique entre deux états (infecté/non infecté) de deux caractères (espèce x, espèce y). On peut donc simplement appliquer la même procédure que pour un test de déséquilibre de liaison. Il suffit donc de coder la présence de chaque borrélie comme un locus et l’absence par 11 et la présence par 22. Il y a donc quatre loci (Bba ; Bbg, Bbss, Bbundet) avec chacun deux allèles (1 ou 2), toujours homozygotes (ou haploïdes). Pour ce faire, il suffit d’ouvrir le fichier “TotBrutBorIR.txt” et d’y remplacer, dans l’ordre, tous les 1 en 22 et tous les 0 en 11 et de fusionner les colonnes Site year sex pour obtenir quelque chose de la forme (fig. 43).

Image 100000000000033D0000021165AD2EC9.jpg

Figure 43. Début du fichier « TotBrutBorIRCoOccur. txt »

141Enregistrons ce fichier en le nommant “TotBrutBorIRCoOccur.txt” et importons-le dans Genetix afin de le convertir au format Fstat. Cliquez sur Fichier, Importer. Choisissez l’option fichier texte et double-cliquez sur “TotBrutBorIRCoOccur.txt”. Choisissez les options séparateur tabulation, un chiffre par allèle, décochez la case de l’identifiant des individus et cliquez sur OK.

142Il faut ensuite cliquer sur le menu Link. Dis et choisir Black & Kafsur comme sur la figure 44, ce qui aura pour effet de lancer une fenêtre de choix que vous devrez rendre comme dans la figure 45.

143Cliquez ensuite sur OK et les résultats s’affichent dans TotBrutBorIRCoOccur.lkd.

144Cliquez ensuite sur Outils, Conversion et FSTAT et nommez le fichier “TotBrutBorIRCoOccur.dat”. Genetix construira donc un fichier où seront considérées comme appartenant à des populations différentes les tiques de sites, d’années et de sexes différents. Ouvrons ce fichier sous Fstat et sélectionnons les mêmes options qu’en figure 46.

145Constatez que nous ne gardons que les fréquences alléliques (cela pourrait servir) et ce qui nous intéresse, le test de déséquilibre de liaison. On choisit dans un premier temps le niveau 5/100 pour aller plus vite. Cliquez sur “Run” et ensuite ouvrez le fichier “TotBrutBorIRCoOccur.out”. Vous constatez que seulement 2 640 permutations ont été effectuées. Recommencez donc avec le niveau 1/100 pour le menu “Nominal level for multiple tests”. Le résultat peut être synthétisé dans le tableau 18.

Image 100000000000057C000002E5406F3451.jpg

Figure 44. Menu Link. Dis

Image 100000000000057C000003944F8F8EC3.jpg

Figure 45. Choix à faire dans le menu de Black et Krafsur

Image 10000000000005C7000003F479D57C2A.jpg

Figure 46. Menu Fstat pour tester l’association entre borrélies

146On y voit clairement une association positive entre Bbss, Bba et Bbg, même si les associations avec Bbg sont marginalement significatives, on peut considérer que le signal existe eu égard à la grande rareté de Bbg (puissance très faible du test). Il est intéressant de noter pour information que Bbundet, vraisemblablement composée d’une mixture de Bbg (très largement sous-représentée ici) et Bbv (B. valaisiana curieusement absente de l’échantillon) donnent des valeurs essentiellement négatives pour R (IJ), ce qui fait regretter plus encore que les déterminations de l’époque aient connu autant de problèmes. Il n’en reste pas moins qu’une forte corrélation positive lie Bbss, Bbg et Bba, qui est confirmée si on teste la co-occurrence des trois espèces dans la même tique rencontrée une fois dans l’échantillon des 73 tiques femelles de Neuchâtel en 1996, et pas à Bern comme annoncé dans la partie résultat de l’article de De Meeûs et al. (2004b) (on ne relit jamais assez ses épreuves). Il y a N = 73 observations, une fréquence observée de 4/73, 6/73 et 1/73 pour Bbss, Bbg et Bba respectivement, donc une fréquence attendue de p = (4 × 6 × 1)/(73)3 pour l’événement de co-occurrence des trois borrélies dans la même tique, événement observé avec la fréquence k = 1. Cette fréquence observée peut être comparée à l’attendue par un test binomial. Sous R, tapez “binom. test (1, 73, p = 0.00006169, alternative =" two. sided ")”, ce qui donne une P-value = 0,0045 très significative. Cette P-value est en fait égale à la probabilité de l’événement lui-même puisqu’il n’y en a pas de plus rare possible. Elle est donc égale à la probabilité (dans une loi binomiale) de tirer une seule fois Bbss-Bbg-Bba dans 73 tirages et où la probabilité de tirer Bbss-Bbg-Bba une fois (un tirage aléatoire) est de 0,00006169, soit (cf. n’importe quel manuel de statistiques) :

Image 100000000000051800000106CA82D599.jpg

Tableau 18. Valeurs (moyennes pondérées par les tailles de sous-échantillons) de corrélation entre la présence de chacune des deux bactéries considérées (ligne 2-2 dans la sortie de Genetix, colonne R (IJ)) et P-value (test G multi-échantillon de Fstat) correspondante. Une valeur de R (IJ) positive indique une association positive, alors qu’une valeur négative indique une répulsion

Borrélies (I x J)

R (IJ)

P-value

Bbss × Bba

0,292

0,00008

Bbss × Bbg

0,496

0,05311

Bbss × Bbundet

-0,069

1

Bba × Bbg

0,109

0,09348

Bba × Bbundet

-0,017

0,91598

Bbg × Bbundet

-0,030

1

147Cette corrélation est donc très forte. Elle peut être due au fait que les tiques infectées correspondent à des individus sensibles et que les autres individus sont résistants. Cette corrélation peut également provenir du fait qu’être infecté par une des trois borrélies tend à favoriser l’infection par les deux autres (par immunosuppression, par exemple). Ceci peut être testé en ne regardant que les tiques infectées. La corrélation existe-t-elle toujours ?

148Nous allons donc recommencer la même opération, mais avec un fichier ne contenant que les tiques infectées par au moins une borrélie (donc sensibles). Le résultat de cette analyse est présenté dans le tableau 19. L’image qui en ressort est totalement différente puisque les associations sont presque toutes négatives, témoignant d’une exclusion entre espèces de borrélies au sein des tiques avec deux exclusions très significatives entre Bbundet (vraisemblablement composée de Bbg et Bbv) et les deux borrélies les plus fréquentes Bbss et Bba (R (IJ) ≈ - 0,6, P-value ≤ 0,003). Il y a donc dans les populations d’I. ricinus des tiques plus sensibles ou plus exposées que d’autres aux borrélies de toutes espèces où ces dernières s’accumulent. Par contre, au sein de ces tiques plus sensibles, ou davantage exposées, les différentes espèces de borrélies tendent à « s’éviter », en particulier Bbundet, que l’on peut raisonnablement soupçonner d’être essentiellement composée de borrélies d’oiseaux (Bbg et Bbv) et les borrélies de micromammifères (Bbss et Bba). Ceci relance le débat sur l’existence de spécificités différentielles des tiques et donc de races d’hôtes, en particulier vis-à-vis des oiseaux et des micromammifères, au moins aux stades larve et nymphe. Il pourrait s’agir aussi d’une manipulation des tiques par les borrélies qui « inciteraient » une tique infectée par telle ou telle autre borrélie à faire son prochain repas sur la même espèce hôte que celle l’ayant infectée. Nous aurons l’occasion de discuter davantage cela dans le paragraphe suivant.

Tableau 19. Valeurs (moyennes pondérées par les tailles de sous-échantillons) de corrélation entre la présence de chacune des deux bactéries considérées (ligne 2-2 dans la sortie de Genetix, colonne R (IJ)) et P-value (test G multi-échantillon de Fstat) correspondante, mais avec les tiques infectées par au moins une borrélie uniquement (tiques saines exclues du test). Une valeur de R (IJ) positive indique une association positive, alors qu’une valeur négative indique une répulsion

Borrélies

R (IJ)

P-value

Bbss × Bba

-0,471

0,11049

Bbss × Bbg

0,412

1

Bbss × Bbundet

-0,622

0,00304

Bba × Bbg

-0,212

0,67559

Bba × Bbundet

-0,632

0,0001

Bbg × Bbundet

-0,290

1

Occurrence des différentes espèces de borrélies et génétique des tiques

149Dans cette partie, nous rechercherons s’il existe une relation entre la génétique des tiques et leur probabilité d’infection par chacun des quatre types de borrélies. On peut répondre à cette question de trois manières. Soit en testant la différenciation génétique entre tiques infectées et non infectées dans chaque sous-échantillon, soit en testant la différenciation, dans chaque sous-échantillon, entre tiques infectées par des borrélies différentes, enfin en procédant à un test de biais de structuration, comme nous l’avons fait pour le sexe des tiques, mais avec le statut infecté/non infecté à la place.

Différenciation entre tiques infectées et non infectées

150Il faut construire un fichier par espèce de bactérie Bbss, Bba et Bbundet (il n’y a pas assez de Bbg). On doit changer de nom de population pour chaque site, année et sexe. La figure 47 donne un exemple de fichier pour Bbss.

Image 1000000000000264000001A68A9A2D77.jpg

Figure 47. Exemple des premiers individus du fichier de données pour tester la différentiation entre tiques infectées (I) et non infectées (U) par Bbss. Nous voyons ici les tiques de Dorénaz 1996 femelles

Image 10000000000005C7000003A95CBBE4C2.jpg

Figure 48. Exemple des cases à cocher pour une analyse de différenciation par paire d’échantillons, exemple des tiques infectées ou non par Bbss

151Il suffit ensuite de convertir ce fichier au format Fstat (en passant par Genetix, par exemple) et de procéder sous Fstat au calcul des FST par paire de sous-échantillons et au test de différenciation par paire, comme indiqué dans la figure 48. Vous constatez que j’ai coché la case 1/1000 pour le nominal level afin d’obtenir au moins 10 000 permutations et donc d’obtenir des P-values assez précises. Le fichier de données s’appelle “ForPairedBbss.dat” et les fichiers de sortie qui nous intéressent sont “ForPairedBbss.fst” pour récupérer les valeurs de FST par paire qui nous intéressent et “ForPairedBbss-pp.pvl” où nous allons récupérer les P-values correspondantes. Attention, dans ces fichiers, seules les comparaisons entre tiques infectées et non infectées du même sexe, de la même année et du même site nous intéressent. Le résultat pour Bbss est présenté dans le tableau 20.

Tableau 20. Compilation des résultats obtenus lors de l’analyse de la différenciation entre paires de sous-échantillons infectés et non infecté par Bbss. La combinatoire est obtenue par la moyenne non pondérée des FST et un test binomial généralisé pour les P-values

Sous-échantillon

FST

P-value

Dorénaz 1996 femelles

-0,008

0,6477

Dorénaz 1996 mâles

-0,030

0,3226

Eclepens 1996 femelles

0,008

0,1206

Eclepens 1996 mâles

0,027

NA

Gorges-du-Trient 1996 femelles

-0,034

0,9171

Montmollin 1996 mâles

-0,027

NA

Neuchâtel 1996 femelles

-0,001

0,7250

Combinatoire

-0,009

0,5179

152Vous remarquerez que la combinaison des cinq tests disponibles a été effectuée à l’aide de la procédure binomiale généralisée de Teriokhin et al. (2007) effectuée à l’aide du logiciel MultiTest (De Meeûs et al., 2009). En effet, à partir de quatre tests, je préfère utiliser cette procédure plutôt que le test Z de Stouffer (Whitlock, 2005). Pour effectuer ce dernier, chaque P-value individuelle est transformée en son équivalent de la distribution Z centrée sur 0 et d’écart-type 1. Sous Excel, on tape =si (b2="na"; " ";si (b2>0 . 99 99 ; loi. normale. inverse (0.999 9 ; 0 ;1) ; loi .normale. inverse(b2 ;0;1)) ). B2 correspond aux coordonnées de la case du tableau Excel où la P-value à transformer se trouve. Cette commande renvoie une absence de résultat quand “NA” est rencontré et tient compte du fait qu’une P-value de 1 n’est pas transformable et la P-value = 0,9999 est choisie comme limite supérieure. Enfin, l’équivalent de la P-value en Z centrée réduite de moyenne 0 et d’écart-type 1 est calculé. Les valeurs Zi obtenues sont ensuite combinées dans la formule (Whitlock, 2005) :

Image 100000000000006E00000048299F810B.jpg, où k est le nombre de tests (=0.3266 ici).

153La P-value globale s’obtient ensuite par un retour à la loi normale, soit sous Excel : = loi. normale. standard (Zs)(= 0,628 ici). Vous trouverez un argumentaire plus détaillé dans De Meeûs et al. (2009) pour les situations où la procédure binomiale généralisée ou le test Z doivent ou peuvent être utilisés.

154Si on procède de la même façon pour Bba et Bbundet, le même type de résultat émerge, même quand on ne distingue pas le sexe des tiques (échantillons plus grands) puisque pour ces deux catégories de borrélies, nous avons vu que le sexe des tiques n’importait pas. Ce résultat est rassurant car, étant donné que les marqueurs sont non codants (donc neutres) et indépendants, il eut été difficile d’interpréter une différenciation entre tiques infectées et non infectées, à moins d’évoquer l’existence d’espèces cryptiques de tiques et une spécificité des borrélies.

Différenciation entre tiques infectées par différentes borrélies

155Ici, il faut ne garder que les tiques infectées et définir comme sous-population les tiques du même sexe, échantillonnées la même année, dans le même site et ayant le même statut infectieux. Notons qu’une tique infectée par Bba et Bbss ne fera pas partie de la même sous-population qu’une tique infectée par Bbss seule. On met ensuite le fichier au format Fstat et on lance la procédure de FST par paire. Ce faisant, vous constaterez que la plupart des tests sont infaisables, c’est normal. Les résultats sont compilés dans le tableau 21. En toute rigueur les tests, qui ne sont pas tous indépendants, devraient subir la correction de Bonferroni. Cependant, étant donné la faiblesse des échantillons (manque total de puissance), nous nous abstiendrons de le faire. Le seul FST positif est obtenu entre Bba et Bbundet, mais il n’est pas significativement plus grand que 0. Eu égard à la faiblesse des tailles de sous-populations ici, nous décidons que rien ne permet d’affirmer l’existence d’une différence génétique entre tiques infectées par différentes bactéries et rien ne permet de l’exclure formellement au moins pour ce qui concerne le couple Bba/Bbundet. S’il existe des races d’hôtes chez I. ricinus, ce n’est pas avec ces données qu’on peut le montrer.

Biais de structuration spécifique associé au pathogène

156Ici, il faut reprendre les données pour chaque espèce de bactérie et créer un fichier de type Genepop comme ce qui a été fait en p. 147-150, sauf qu’ici les tiques sont distinguées en fonction de leur statut infectieux et non par leur sexe, tel que dans la figure 49. Notons que nous ne traitons que les sites prélevés en Suisse et où au moins une tique infectée est trouvée. Parce qu’il y a un biais de structuration sexe-spécifique, ainsi que des différences d’infection, les femelles et les mâles sont analysés séparément. Cependant, parce que la taille des échantillons est très faible (peu de borrélies trouvées et identifiées), nous combinerons le tout dans un seul fichier (gain de puissance). On prendra soin de distinguer les tiques d’années et de sexe différents comme appartenant à des populations différentes (séparées par un “pop” dans le fichier).

Tableau 21. Compilation des résultats des tests de différenciation, parmi les tiques infectées, par paire en fonction de l’espèce de bactérie présente et pour les paires effectivement trouvées. Quand plusieurs tests indépendants sont disponibles ils sont combinés : les FST sont des moyennes non pondérées, alors que les P-values ont été obtenues par la procédure Z (il y a en effet systématiquement moins de quatre tests ici)

Borrélies

Sous-échantillon

FST

P-value

Bbss/Bba

Dor96F

-0,0095

0,8577

Gor96F

0,0000

0,6628

Combinés

-0,0047

0,8540

Bba/Bbundet

Mon96F

-0,0357

1

Sta96F

0,1025

0,0662

Sta96M

0,0454

0,1687

Combinés

0,0374

0,7657

Bba/Bbss + Bba

Dor96F

-0,0501

0,8560

Bbss/Bbss + Bba

Dor96F

0,0004

0,5998

Bba/Bba + Bbundet

Sta96M

0,0269

0,0676

Bbundet/Bba + Bbundet

Sta96M

-0,0394

0,8043

+ signifie la co-occurrence de deux espèces de borrélies

157Quand le fichier est constitué, il faut ensuite lancer Fstat et cliquer sur le menu “Biased dispersal”. La fiche correspondante apparaît alors. Il faut ensuite charger le fichier à analyser en cliquant le menu “File” et “Open” et cocher les cases comme en figure 50 puis sur le bouton “Go !”. Pour une raison que j’ignore, il faut cocher tous les paramètres si on souhaite obtenir le résultat du test sur Hs, en particulier FIS et H0 qui ne sont guères utiles ici, car nous avons codé les mâles homozygotes pour IR08.

158Le résultat est contenu dans un fichier de type nom du fichier.res (un fichier par espèce de borrélie). Le résultat principal concerne le test du FST (et aussi la relatedness, ce qui est normal si on regarde sa définition dans la documentation de Fstat) et est présenté dans le tableau 22.

Image 10000000000001E2000002537084DEA5.jpg

Figure 49. Type de données pour le test de biais de structuration pathogène spécifique. Les tiques infectées sont notées avec un I et les saines avec un U. Se référer à l’aide de Fstat pour plus de détails sur la constitution d’un tel fichier. Il est important de ne pas oublier que les allèles doivent être à deux chiffres et que les colonnes sont séparées par des espaces et non des tabulations

Tableau 22. Résultat du test basé sur le FST de biais de structuration génétique pathogène spécifique des tiques pour les différentes espèces de borrélies pour lesquelles assez de données étaient disponibles (Bbg exclue). On remarque une structuration significativement plus forte pour les tiques infectées (I) par Bba par rapport aux tiques non infectées par cette borrélie (U)

Bbss

Bba

Bbundet

U

0,001

0,002

0,000

I

-0,015

0,076

-0,045

P-value

0,4998

0,0033

0,1764

Image 10000000000005C70000037B9F3243D6.jpg

Figure 50. Cases à cocher pour l’analyse du biais de dispersion pathogène spécifique. Le test demandé est bilatéral, car nous n’avons en principe pas d’à priori. Toutes les cases sont cochées, même les cases “Fis” et “Ho” (inutiles à cause du codage de IR08), car on souhaite obtenir le résultat pour Hs

159Il y a donc bien un biais de structuration dû à l’infection par Bba. Plusieurs hypothèses peuvent expliquer ce résultat. La première hypothèse implique que certaines tiques, plus sensibles à l’infection par Bba sont aussi pléiotropiquement moins mobiles. Les marqueurs utilisés étant des microsatellites non codants, cela impliquerait un déterminisme génomique peu vraisemblable. Par ailleurs, l’absence totale de différenciation entre tiques infectées et non infectées (montrée en p. 179-182) discrédite cette interprétation.

160La deuxième hypothèse implique l’existence d’au moins deux espèces cryptiques dont l’une, moins mobile que la seconde, serait plus sensible à l’infection par Bba. Notons que nous n’avons noté aucun déséquilibre de liaison (attendu en pareil cas). Par ailleurs, si on calcule avec Fstat le FIS des tiques en séparant celles infectées par Bba de celles qui ne le sont pas, on ne retrouve aucune diminution du FIS (~ 0,45 pour les infectées et ~ 0,44 pour les non infectées par Bba). Il n’existe pas de différenciation significative entre tiques infectées et non infectées. Cette interprétation n’est donc aucunement soutenue.

161La troisième interprétation possible impliquerait l’existence d’une adaptation locale des borrélies qui infecteraient plus facilement les tiques locales (résidentes) que les immigrantes. Deux arguments vont à l’encontre de cette hypothèse.

162La première est que les tiques mâles et femelles qui en principe n’ont pas la même dispersion (les femelles dispersent en principe peu ou pas, cf. p. 147-153) ne sont pas infectées différemment par Bba (p. 173). Par ailleurs, c’est le partenaire le plus mobile des deux qui doit en théorie être le mieux adapté localement (Gandon et al., 1996 ; Gandon, 2002). Or ici, les tiques sont modestement structurées alors que l’on pense que les borrélies le sont beaucoup plus (Qiu et al., 1997). C’est donc l’hôte (la tique) qui devrait être adapté localement et non l’inverse.

163La quatrième hypothèse implique une survie plus faible des tiques migrantes quand ces dernières sont infectées par Bba. Comme les tiques femelles sont moins mobiles que les mâles, ce sont ces derniers qui devraient être les plus affectés par ce phénomène. Ceci est testable en refaisant l’analyse sur les tiques femelles et mâles séparément. Cela suppose une survie au stress moins bonne des larves et/ou nymphes infectées par Bba.

164La cinquième hypothèse est la plus séduisante. Elle implique une manipulation des larves et nymphes par la borrélie. Cette borrélie est spécifique de petits rongeurs. Il est donc plus intéressant pour elle d’être injectée dans un petit rongeur, peu dispersant, que dans un oiseau ou un grand mammifère, hôtes beaucoup plus mobiles. Les Bba capables de manipuler les tiques qu’elles infectent de sorte que ces dernières préfèrent se fixer sur un petit rongeur plutôt que sur d’autres hôtes seraient donc avantagées. Cette hypothèse est testable en laboratoire, mais cela n’a malheureusement jamais été fait. Cela implique aussi, comme pour l’hypothèse précédente, que les femelles, déjà très peu mobiles, seront moins affectées par le biais de structuration Bba-spécifique que les mâles.

Biais de structuration spécifique au pathogène et au sexe

165Nous allons utiliser la même procédure que précédemment, mais en divisant le fichier en deux : un fichier pour les tiques femelles et un autre pour les tiques mâles. Cette fois, les tests seront faits de manière unilatérale avec I (infectés) comme catégorie la plus philopatrique. Il y a deux raisons à cela. La première est que l’on connaît d’avance le sens du signal. La seconde raison est que les échantillons étant encore plus petits, nous aurons besoin d’encore plus de puissance dans le test. Nous ne nous occuperons que du test sur le FST. Les tests sont tous les deux significatifs avec P-value = 0,0497 pour les tiques femelles et P-value = 0,0123 pour les tiques mâles et une apparente très forte différence de signal entre les deux, comme indiqué dans la figure 51.

166Nous pouvons également constater la formidable variance chez les mâles infectés (très peu nombreux). Nous pouvons effectuer un test unilatéral de Wilcoxon pour données appariées comme en p. 143 quand nous avions comparé les FIS des données brutes avec ceux des données clusterisées par BAPS. Ici, l’unité d’appariement reste le locus (donc cinq données), mais la statistique est la différence de FST entre tiques infectées et non infectées chez les femelles et les mâles. Le fichier à tester contiendra donc les différences des différences appariées : (FSTMI-FSTMU) - (FSTFI-FSTFU). Le test unilatéral (FSTMI -FSTMU > FSTFI -FSTFU) montre que la différence n’est pas significative, même si la P-value reste relativement faible (0,17). Ceci illustre les limites de notre jeu de données (beaucoup trop de données manquantes).

Image 100000000000042F0000036073A50AD1.jpg

Figure 51. Différence comparée entre tiques mâles et femelles du FST mesuré entre tiques infectées par Bba et celles qui ne le sont pas

CONCLUSIONS SUR LES BORRÉLIES ET I. RICINUS EN SUISSE

167Au cours de nos analyses, nous avons constaté que Bbss, borrélie d’écureuil, était plus souvent retrouvée chez les tiques mâles que femelles, ce qui est attendu si, comme le suggérait le biais de dispersion sexe-spécifique détecté chez ces tiques, les larves et nymphes femelles préfèrent se nourrir sur des rongeurs (peu dispersants). Rien de tel n’a pu être trouvé pour Bba pour laquelle ceci était attendu également, peut-être parce qu’une certaine quantité de tiques infectées par cette borrélie fait partie du stock Bbundet. Quant à Bbg, trop rarement détectée, d’autres études seront requises afin de déterminer si, comme attendu, elle est plus souvent retrouvée chez les tiques mâles.

168Certaines tiques sont plus sensibles ou plus exposées à l’infection par les borrélies en général, comme l’attestent les fortes corrélations positives observées sur les co-occurrences des trois espèces Bbss, Bba et Bbg. En se concentrant sur ces tiques sensibles (infectées par au moins une borrélie), il y a un évitement manifeste. Les corrélations deviennent toutes négatives, exception faite de l’association Bbss × Bbg, pour qui le faible nombre de Bbg détectées rend les choses difficiles à interpréter, et très significatives pour les couples Bba × Bbundet et Bbss × Bbundet. Cette dernière observation peut laisser spéculer que ces borrélies indéterminées soient majoritairement des borrélies d’oiseaux (Bbg et Bbv) très déficitaires dans notre jeu de données. Dans ce cas, nous pourrions proposer que les larves et nymphes sensibles se subdivisent en tiques ne se nourrissant que sur une gamme limitée d’hôtes réservoirs de borrélies spécifiques. Tout dépend de l’identité spécifique de ces Bbundet. Les données ne permettent pas d’exclure l’existence d’une telle spécificité en races d’hôtes. La manipulation de la spécificité des tiques par les borrélies ne peut pas non plus être exclue. C’est aussi cette manipulation qui expliquerait le biais de structuration des tiques infectées par Bba. D’une manière générale, on ne peut que regretter le nombre de données manquantes qui limite nos conclusions mais aussi remarquer que, malgré cela, de nombreuses perspectives nouvelles de recherche ont émergé qui illustrent la puissance des outils offerts par la génétique des populations.

Notes

9 Entre temps, j’ai découvert l’existence du “Package” R-Commander ou Rcmdr qui, en quelques clics de souris, permet d’effectuer ces commandes automatiquement.

10 Je me suis rendu compte sur le tard que Genetix contenait quelques bugs dans ce module et je conseillerai d’utiliser plutôt MSA pour le calcul de distances, bien qu’ici cela n’ait pas changé grand-chose, raison pour laquelle j’ai laissé l’analyse telle qu’elle. Pour l’utilisation de MSA, se référer à la seconde partie de ce manuel, p. 254.

11 Sur ces mêmes données, l’estimation avec un logiciel alternatif, LDNe (WAPLES et DO, 2008), non encore connu au moment de la rédaction de ce chapitre et dont l’utilisation est détaillée plus loin, donne une moyenne de Ne = 223.

© IRD Éditions, 2012

Conditions d’utilisation : http://www.openedition.org/6540

Rechercher dans OpenEdition Search

Vous allez être redirigé vers OpenEdition Search