Ontwikkeling van een monitoringsprotocol voor heikikker Bepaling van de aanwezigheid en verhouding bruine en heikikker via (e)DNA-barcoding en opstart genetische monitoring Rein Brys, Karen Cox & Sabrina Neyrinck
Auteurs: Rein Brys, Karen Cox & Sabrina Neyrinck Instituut voor Natuur- en Bosonderzoek Reviewers: Jeroen Speybroeck Instituut voor Natuur- en Bosonderzoek Het INBO is het onafhankelijk onderzoeksinstituut van de Vlaamse overheid dat via toegepast wetenschappelijk onderzoek, data- en kennisontsluiting het biodiversiteitsbeleid en -beheer onderbouwt en evalueert. Vestiging: INBO Geraardsbergen Gaverstraat 4 9050 Geraardsbergen www.inbo.be e-mail: rein.brys@inbo.be Wijze van citeren: Brys, R., Cox, K. & Neyrinck S. (2020). Ontwikkeling van een monitoringsprotocol voor heikikker: bepaling van de aanwezigheid en verhouding bruine en heikikker via (e)DNAbarcoding en opstart genetische monitoring. Rapporten van het Instituut voor Natuur- en Bosonderzoek 2020 (8). Instituut voor Natuur- en Bosonderzoek, Brussel. DOI: doi.org/10.21436/inbor.17686189 D/2020/3241/060 Rapporten van het Instituut voor Natuur- en Bosonderzoek 2020 (8) ISSN: 1782-9054 Verantwoordelijke uitgever: Maurice Hoffmann Foto cover: Yves Adams / Vildaphoto
Š 2020, Instituut voor Natuur- en Bosonderzoek
ONTWIKKELING VAN EEN MONITORINGSPROTOCOL VOOR HEIKIKKER Bepaling van de aanwezigheid en verhouding bruine en heikikker via (e)DNA-barcoding en opstart genetische monitoring
Rein Brys, Karen Cox & Sabrina Neyrinck
doi.org/10.21436/inbor.17686189
www.inbo.be
Dankwoord We danken Jeroen Speybroeck, Koen Van Muylem en Maurice Hoffmann voor hun hulp bij de staalname en/of het nalezen van het rapport, en Quentin Mauvisseau voor de ondersteuning bij de ontwikkeling van de ddPCR assays. Graag willen we ook An Van Breusegem hartelijk danken voor de hulp bij de DNA-extractie van eitjes en larven. Tot slot zijn we Afdeling Natuur en Bos erkentelijk voor het verschaffen van de toegang en toelating (nr. ANB/BL-FF/V17000035) om op hun terreinen het voorliggende onderzoek te kunnen uitvoeren.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 2 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Samenvatting Heikikker (Rana arvalis) is een Habitatrichtlijnsoort die binnen Vlaanderen op gewestelijk niveau gemonitord moet worden. Dit brengt met zich mee dat de aanwezigheid van de soort in het landschap accuraat in kaart moet worden gebracht en dat populatiegroottes van de soort in de tijd dienen opgevolgd te worden. Hierbij stellen zich echter praktische problemen omdat de soort erg moeilijk is waar te nemen en bepaling van eilegsels een verkeerd beeld kan geven doordat ook bruine kikker (Rana temporaria) op deze locaties vaak eitjes afzet. Binnen deze studie onderzochten we in eerste instantie in hoeverre (e)DNA-barcoding technieken zouden kunnen helpen bij de bepaling van de aanwezigheid van heikikker in een waterpartij. Hiervoor werden enkele vijvers en vennen onderzocht in het Nationaal Park Hoge Kempen en het Midden Limburgs vijvergebied. Ook gingen we na of deze methodiek kan helpen bij het inschatten van de verhouding heikikker versus bruine kikker eilegsels op plaatsen waar de soort aan reproductie doet en populatiegroottes via eitellingen worden opgevolgd. Dit gebeurde op verschillende locaties in het Klein Schietveld. Uit de eDNAresultaten afkomstig van waterstalen genomen doorheen het jaar in verschillende vijvers en vennen waar heikikker al is waargenomen, blijkt dat we met eDNA ddPCR-analysen de aanwezigheid van heikikker vroeg in het voorjaar perfect kunnen aantonen (eind maart – begin april). Op locaties waar aan reproductie wordt gedaan is op basis van dezelfde eDNA methodiek het echter veel minder vanzelfsprekend om de verhouding aan heikikker versus bruine kikker eilegsels te bepalen. Deze resultaten matchen ook niet altijd even goed met de werkelijke verhouding aan beide soorten. Deze informatie kan veel efficiënter en accurater worden bekomen door van elk eilegsel één eitje te bemonsteren, dit te mixen en van het resulterende DNA extract via ddPCR-analyse de verhouding aan heikikker en bruine kikker DNA te bepalen. In tweede instantie gingen we binnen deze studie na of een geselecteerd panel aan genetische merkers bruikbaar zou kunnen zijn voor de genetische monitoring van de Vlaamse heikikker populaties. Dit werd getest op dezelfde heikikkerpopulatie van het Klein Schietveld. Hierbij bleken de geselecteerde merkers voldoende variabel om heikikkerpopulaties genetisch te karakteriseren. Zo bleek de omvangrijke populatie van het Klein Schietveld een grote genetische variatie te herbergen met een groot aantal effectieve ouders. We kunnen besluiten dat de gebruikte (e)DNA barcoding methodieken en de geselecteerde genetische merkers beide bruikbare instrumenten zijn om op termijn heikikkerpopulaties in het landschap te detecteren, en hun grootte, genetische diversiteit en onderlinge connectiviteit adequaat te kunnen bepalen.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 3 van 40
Aanbevelingen voor monitoring en beheer Dit onderzoek heeft aangetoond dat eDNA-gebaseerde technieken erg geschikt kunnen zijn om de aanwezigheid van heikikker in een waterpartij te bepalen. Veel moeilijker is het een exacte maat voor populatiegrootte te bepalen. Hiervoor blijkt het tellen van eilegsels nog steeds een goede optie, tenminste als het aandeel legsels van bruine kikker mee in rekening wordt genomen. Uit de geteste methodes komt naar voor dat dit het meest efficiĂŤnt bepaald kan worden door van elk legsel (of een representatief aandeel van het totaal aantal legsels) op een welbepaalde locatie ĂŠĂŠn eitje te nemen. Op basis van het resulterende DNA-extract van het mengsel van deze eitjes, kunnen we in een volgende stap via ddPCR analysen de verhouding tussen beide soorten accuraat achterhalen. Het aantal heikikkermoeders op een welbepaalde locatie kan op die manier erg nauwkeurig bepaald worden door de telling van de eilegsels met de bekomen DNA-verhouding te vermenigvuldigen. Hierbij dient opgemerkt te worden dat dit uiteraard geen volledig non-invasieve methode is, waardoor steeds de afweging moet worden gemaakt of het voor een welbepaalde situatie wenselijk is om poelen en eilegsels op deze manier (beperkt) te verstoren. Verder heeft deze studie aangetoond dat met behulp van de geteste genetische merkers een genetische karakterisatie en/of monitoring van heikikkerpopulaties verwezenlijkt kan worden. Zeker voor populaties waar men bijvoorbeeld twijfelt aan de conservatiestatus, waar recent ingrepen in de omgeving gebeurden, of waar kennis van de onderlinge connectiviteit van populaties ontbreekt, kan dit zeker een meerwaarde hebben.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 4 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
English abstract Moor frog (Rana arvalis) is listed in the Habitats Directive, which implies that its populations must be monitored over time. Therefore, the species distribution range needs to be regularly mapped, in combination with an accurate temporal monitoring of the size of its present populations. In case of Moor frog this is, however, very difficult given the hidden behaviour of this species. Although counting the number of egg clutches can offer a good proxy of the population size of Moor frog, it is frequently documented that also Common frog (Rana temporaria) lays its egg clutches at the same locations. Within this study, we initially tested to what extent eDNA barcoding based on water samples can help to determine the presence of Moor frog within a waterbody and which period of the year is most optimal to do so. For this purpose we regularly sampled (from winter until summer) several ponds and fens located at the Hoge Kempen National Park and the Midden Limburg pond area. Based on these analyses, it appeared that we can detect the presence of Moor frog most accurately early in the spring (late March - early April). Secondly, we determined whether eDNA or DNA barcoding can additionally help to obtain a ratio of the density of Moor frogs versus Common frogs, when sexual reproduction takes place and egg clutches are deposited at a certain spot. This part of the study was performed at Klein Schietveld. The eDNA results obtained from water sampled around the egg clutches revealed that this method is quite inaccurate in predicting the actual ratio of Moor versus Common frog eggs. In contrast, ddPCR analyses on tissue material of a pooled sample of eggs, through the collection of one egg per clutch, showed that the DNA ratio of both species perfectly predicts the actual number of clutches that are present per species. In other words, combining the count number of egg clutches with the ratio of both species by means of ddPCR analyses on a pooled sample of eggs thus offers an adequate method to monitor population sizes of both species. Finally, we investigated whether a selected panel of genetic markers could prove useful for the genetic monitoring of the Flemish Moor frog populations. The selected markers (microsatellites) were found to be sufficiently variable to genetically characterize Moor frog populations located at the area of Klein Schietveld. Based on this result we can conclude that this population is characterized by high genetic variation together with a large effective number of breeders. In conclusion, both the DNA barcoding methodologies and the selected genetic markers are useful instruments for appropriately detecting Moor frog populations in the landscape, determining the size of their (sub)populations and finally examining genetic diversity and the level of connectivity.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 5 van 40
Inhoudstafel 1
Inleiding ................................................................................................................................. 9
2
Materiaal en methode ........................................................................................................ 11 2.1
2.1.1
Temporele detectie van bruine en heikikker via eDNA ...................................... 11
2.1.2
Soortidentficatie en populatiegenetica............................................................... 12
2.2
3
Studiegebieden ........................................................................................................... 11
Staalname.................................................................................................................... 13
2.2.1
Temporele detectie ............................................................................................. 13
2.2.2
Soortidentificatie en populatiegenetica.............................................................. 14
2.2.3
eDNA-afgifte onder gecontroleerde condities .................................................... 15
2.3
eDNA-extractie uit water ............................................................................................ 15
2.4
DNA-extractie uit eitjes en larven ............................................................................... 15
2.5
ddPCR-analysen ........................................................................................................... 16
2.6
Soortidentificatie van de eilegsels via RFLP ................................................................ 16
2.7
Populatiegenetische analyses ..................................................................................... 17
2.7.1
Genotyping .......................................................................................................... 17
2.7.2
Genetische structuur ........................................................................................... 18
2.7.3
Genetische variatie en effectief aantal ouders ................................................... 19
Resultaten en discussie ....................................................................................................... 21 3.1
eDNA afgifte van hei- en bruine Kikker ....................................................................... 21
3.1.1
Temporele detectie onder veldcondities ............................................................ 21
3.1.2
eDNA-afgifte onder gecontroleerde condities .................................................... 22
3.2
Bepaling van de verhouding bruine kikker- en heikikkereilegsels .............................. 24
3.2.1
Soortidentificatie van de eilegsels via RFLP ........................................................ 24
3.2.2
Bepaling van de soortenverhouding op basis van mengstalen van eitjes .......... 25
3.2.3
Bepaling van de soortverhouding op basis van eDNA uit water ......................... 27
3.3
Populatiegenetisch onderzoek .................................................................................... 28
4
Conclusies ............................................................................................................................ 35
5
Referenties .......................................................................................................................... 36
6
Bijlagen ................................................................................................................................ 38
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 6 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Lijst van figuren Figuur 1: Parende heikikkers (Rana arvalis) vroeg in het voorjaar (foto: Remo Schulze – met licentie onder CC BY-SA 3.0)........................................................................................ 9 Figuur 2: Voorbeeld van (a) een koor van roepende heikikker mannetjes (foto: Christian Fischer, met licentie onder CC BY-SA 3.0) en (b) eilegsels die door heikikker worden afgezet in een ondiepe plas (b, foto: Christian Fischer, met licentie onder CC BY-SA 3.0). ........................................................................................................................... 10 Figuur 3. Locaties van de vijvers in het Vijvergebied Midden-Limburg, waar doorheen het jaar staalname werd uitgevoerd voor eDNA-detectie van heikikker en bruine kikker. ... 11 Figuur 4. Locaties van de vennen in het Nationaal Park Hoge Kempen, waar doorheen het jaar staalname werd uitgevoerd voor eDNA detectie van heikikker en bruine kikker. ... 11 Figuur 5 Locaties van vennen en poelen met eilegsels van hei- of bruine kikker. * is een poel waar referentiestalen van heikker werden genomen voor de soortidentificatie en die ook in het mesocosmexperiment werden ingezet (zie verder). ......................... 12 Figuur 6. Waterpartijen die doorheen het seizoen werden bemonsterd in functie van de temporele detectie van heikikker en bruine kikker via eDNA in respectievelijk het Nationaal Park de Hoge Kempen (a) en het Vijvergebied Midden-Limburg (b). ...... 13 Figuur 7. Illustraties van de staalnamecampagne op 23/03/2017 op het Klein Schietveld. (a) Typische cluster van eilegsels van heikikker en/of bruine kikker, (b) staalname van een drietal eitjes uit elk van de verschillende legsels per cluster, en c) staalname van water rondom de cluster van legsels voor het eDNA-onderzoek. ............................ 14 Figuur 8 Resultaat van een RFLP-analyse met duidelijke amplificatie en opsplitsing van de verschillende DNA-fragmentjes in gescheiden bandjes na gelelectroforese. Elke verticale kolom is een staal. ...................................................................................... 17 Figuur 9. Temporeel eDNA-profiel (gemeten via ddPCR) (gemiddelde van twee analysen per staal ± SE) van bruine kikker en heikikker in het water van twee vennen en vijvers in respectievelijk het Nationaal Park de Hoge Kempen (a-b) en het Vijvergebied Midden-Limburg (c-d). Waterstalen werden driewekelijks genomen vanaf eind februari t.e.m. eind juli 2018..................................................................................... 21 Figuur 10. Temporeel profiel van de watertemperatuur (C°) gemeten op elk van de eDNA bemonsteringen in twee vennen 3a en 3b en vijvers 4A en 4b in het Nationaal Park Hoge Kempen en het Vijvergebied Midden-Limburg. Stalen werden driewekelijks genomen vanaf eind februari t.e.m. eind juli 2018. ................................................. 22 Figuur 11. Relatie tussen het aantal bruine kikker- (a) en heikikkerlarven (b) in elk van de 11 mesocosmen en het eDNA signaal gemeten via ddPCR. .......................................... 23 Figuur 12. Relatie tussen het aandeel (%) heikikker legsels en de pH gemeten in het water rondom de clusters van eilegsels. ............................................................................. 25 Figuur 13. Verband tussen a) het totaal aantal eDNA-kopieën van bruine kikker en heikikker (bepaald via ddPCR) (gemiddelde van twee analyses ± SE) en b) de proportionele verhouding bepaald op basis van het gemiddelde van het eDNA-signaal tussen beide soorten binnen elk staal, en dit in functie van een voor bruine kikker aflopend en voor heikikker oplopend aandeel eitjes op een totaal van 10 eitjes per staal. ... 26 Figuur 14. Verband tussen het percentage heikikker eDNA (t.o.v. de totale eDNA-concentratie van bruine kikker en heikikker, gemeten via ddPCR) en het aandeel heikikkereilegsels op de onderzochte locaties in het Klein Schietveld. .................... 27 Figuur 15. Verband tussen het totaal aantal heikikker eilegsels die op elk van de locaties in het Klein Schietveld aanwezig waren het totaal aantal heikikker eDNA moleculen gemeten in het omliggende water (via ddPCR). ....................................................... 28 Figuur 16. De eerste twee assen van de principale coördinaten analyse, die respectievelijk 30 en 23% van de variatie verklaren. De analyse is gebaseerd op de Fst -waarden tussen ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 7 van 40
poelen met meer dan 10 stalen van heikikker, exclusief de genotypes van broers/zussen. De nummers in de codes komen overeen met deze in Figuur 5. ..... 29 Figuur 17. Grafische weergave van de sPCA resultaten met (a) de globale structuur en (b) de lokale structuur t.o.v. de x- en y-coĂśrdinaten van de locaties. De blokjes geven de waarde van lagged scores volgens grootte en kleur: zwart is positief, wit is negatief. ................................................................................................................................... 30 Figuur 18. Spatiale autocorrelatie van genetische verwantschap en geografische afstand voor (a) de totale dataset en (b) de dataset zonder broers/zussen.................................. 31 Figuur 19 Verwachte heterozygositeit (gendiversiteit) Âą standaarddeviatie per poel. .............. 32 Figuur 20 Visualisatie van het GLM model met het briggse logaritme van de LD schattingen van het effectief aantal ouders als respons. De voorspelde waarden worden weergegeven door de rechte. Het logaritme van de zuurstofconcentratie wordt weergegeven op de x-as en de oppervlakte van de poelen in kleurschakering van donker- naar lichtblauw (van klein naar groot). ....................................................... 34
Lijst van tabellen Tabel 1. Lijst van aantal getelde en bemonsterde eilegsels per poel op het Klein Schietveld. De code van de poel komt overeen met de codes gegeven in Figuur 5. Het percentage legsels van heikikker werd berekend rekening houdend met de soortbepaling van larven (poelen 1 en 24). Ook het aantal legsels waarbij de RFLP-analyse geen resultaat gaf zijn gegeven. ........................................................................................... 24 Tabel 2 Schattingen van effectief aantal ouders (Nb) volgens twee methodes: de linkage disequilibriummethode (LD) en de verwantschapsmethode met bijhorende betrouwbaarheidsintervallen tussen haakjes. De aangepaste waarde van Nb LD is Nb adjust. Het aantal ouders bepaald met COLONY wordt gegeven. Door het aantal bemonsterde legsels (dus het aantal moerders) te delen door het aantal ouders, kan de verhouding moeders/vaders berekend worden..................................................... 33 Tabel 3 Modelsamenvatting van het glm met log10(Nb) als responsvariabele. ......................... 34
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 8 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
1 INLEIDING Sinds 2016 voert Natuurpunt soortmonitoring uit in opdracht van INBO en ANB om te voorzien in de noden voortvloeiend uit de Habitatrichtlijn. Eén van de soorten die in deze richtlijn in de Bijlagen II en IV staat vermeld, is de heikikker (Rana arvalis) (Figuur 1). Dit houdt in dat deze soort op gewestelijk niveau dient te worden gemonitord. Voor de kikkers en padden die binnen deze context onder de aandacht vallen, stelden De Bruyn et al. (2015) een monitoringsprotocol op. Doordat heikikker (net zoals poelkikker) een erg moeilijke soort is om op te volgen, werd het betreffende meetnet voor deze soort tot dusver nog niet opgestart, omdat men van oordeel was dat deze soort te moeilijk zou zijn om te laten monitoren door het vrijwilligersnetwerk van Natuurpunt.
Figuur 1: Parende heikikkers (Rana arvalis) vroeg in het voorjaar (foto: Remo Schulze – met licentie onder CC BY-SA 3.0)
Daarom werd voorgenomen om binnen het INBO een pilootstudie uit te voeren, waarin we onderzoeken in hoeverre, met het inzetten van de technieken van (e)DNA-barcoding, een oplossing kan worden gevonden om voor deze soort een accuratere inschatting te kunnen maken van 1) zijn aanwezigheid en 2) de populatiegrootte. Traditioneel wordt de aanwezigheid en inschatting van de populatiegroottes van kikkers en padden gemaakt op basis van tellingen van het aantal roepende mannetjes (Figuur 2a), legsels en/of larven (Figuur 2b). Bij heikikker stelt zich hierbij het probleem dat de kooractiviteit van deze soort zeer kort is (vaak maar één of enkele dagen op uitzonderlijk mooi weer in de periode tussen 15 maart en 15 april), vaak op locaties gebeurt die erg moeilijk bereikbaar zijn en dat de kikkers dan erg gevoelig zijn aan verstoring. Heikikker komt vaak samen voor met bruine kikker (Rana temporaria), een soort die morfologisch sterke gelijkenissen vertoont met ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 9 van 40
heikikker, waardoor ook zichtwaarnemingen moeilijk bruikbaar zijn om correcte kwantitatieve inschattingen van de aanwezigheid van heikikker te bekomen. De bruine kikker deponeert eilegsels bovendien vaak in dezelfde poelen en vennen als heikikker. Omdat het onderscheiden van de legsels van beide soorten ook al moeilijk is, is de inschatting van populatiegroottes op basis van de eilegsels evenmin een betrouwbare methode.
Figuur 2: Voorbeeld van (a) een koor van roepende heikikker mannetjes (foto: Christian Fischer, met licentie onder CC BY-SA 3.0) en (b) eilegsels die door heikikker worden afgezet in een ondiepe plas (b, foto: Christian Fischer, met licentie onder CC BY-SA 3.0).
Om toch een gestandaardiseerd meetnet te kunnen ontwikkelen wilden we daarom nagaan of het mogelijk zou zijn de monitoring uit te voeren aan de hand van genetische methodes. In eerste instantie werd besloten om na te gaan in hoeverre we via eDNA de aanwezigheid van heikikker en bruine kikker in waterpartijen accuraat kunnen aantonen. We willen hierbij eveneens onderzoeken wanneer tijdens het jaar we dit via deze methodiek het best onderzoeken. Daarnaast willen we met dit onderzoek nagaan of er een afgeleide methodiek/protocol kan worden bekomen door eilegseltellingen te combineren met DNA barcoding-methodieken op basis van weefsel (eitjes) en/of eDNA. Hierbij willen we onderzoeken of: 1) we op basis van het nemen van een reeks eitjes (1 eitje per legsel, waarbij een legsel rond de 500 Ă 3000 eitjes kan bevatten), via moleculaire technieken kunnen komen tot een soortbepaling van deze eitjes, en dus de verhouding aan bruine- en heikikkereitjes in een cluster van legsels, en/of 2) deze informatie even goed uit de eDNA (environmental DNA) patronen van het omringende water kan worden gehaald, met het oog op een zo eenvoudig, minst invasieve/verstorende en betaalbaar mogelijke monitoring. Tot slot willen we binnen deze pilootstudie de heikikkereitjes die worden bemonsterd eveneens benutten voor een genetische karakterisering van de onderzochte populatie. Omdat nog geen panel van bruikbare genetische merkers werd samengesteld voor de genetische monitoring van de Vlaamse heikikkerpopulaties, werd deze set van stalen daarvoor ingezet. Aan de hand van de populatie-genetische analyses probeerden we specifiek voor de onderzochte populatie een idee te krijgen over de spatiale structuur en de genetische variatie. Deze waarden kunnen dan als verdere referentie benut worden bij toekomstige genetische monitoring van Vlaamse heikikkerpopulaties. ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 10 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
2 MATERIAAL EN METHODE 2.1 STUDIEGEBIEDEN 2.1.1 Temporele detectie van bruine en heikikker via eDNA Om te onderzoeken of we via eDNA de aanwezigheid van heikikker en/of bruine kikker in het water kunnen aantonen en wanneer dit doorheen het seizoen (lente – zomer) het best wordt uitgevoerd, werd in 2018 een selectie gemaakt van vier vijvers/vennen die verspreid lagen over twee deelgebieden waarvan geweten is dat beide soorten er zich ophouden.
Figuur 3. Locaties van de vijvers in het Vijvergebied Midden-Limburg, waar doorheen het jaar staalname werd uitgevoerd voor eDNA-detectie van heikikker en bruine kikker.
Figuur 4. Locaties van de vennen in het Nationaal Park Hoge Kempen, waar doorheen het jaar staalname werd uitgevoerd voor eDNA detectie van heikikker en bruine kikker.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 11 van 40
Twee vijvers waren gelegen in Vijvergebied Midden-Limburg (Figuur 3) en de twee andere (vennen) bevonden zich in het Nationaal Park Hoge Kempen, nabij de bron van de Kikbeek (Figuur 4).
2.1.2 Soortidentficatie en populatiegenetica Voor de inschatting van de populatiegrootte en populatiegenetica, werd het Klein Schietveld (gelegen op grondgebied van de gemeentes Brasschaat, Kalmthout en Kapellen) als studiegebied geselecteerd, omdat van dit gebied geweten is dat ze een grote heikikker (meta)populatie herbergt en het een overzichtelijk en redelijk toegankelijk heidegebied is. Hier werd half maart 2017 een gebiedsdekkende screening ondernomen om zo veel mogelijk locaties met eilegsels van heikikker of bruine kikker in kaart te brengen. In totaal werden op deze manier 26 locaties langsheen poelen en vennen geselecteerd, waar op het moment van de opname ĂŠĂŠn of meerdere eilegsels in clusters waren afgezet (Figuur 5). Locaties 22a, 22b en 22c zijn delen van eenzelfde poel waar aparte clusters van eilegsels te vinden waren, die duidelijk in ruimte van elkaar gescheiden lagen.
Figuur 5 Locaties van vennen en poelen met eilegsels van hei- of bruine kikker. * is een poel waar referentiestalen van heikker werden genomen voor de soortidentificatie en die ook in het mesocosmexperiment werden ingezet (zie verder). ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 12 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
2.2 STAALNAME 2.2.1 Temporele detectie Voor de temporele screening van heikikker en bruine kikker via eDNA werden de daarvoor geselecteerde vijvers en vennen vanaf eind februari (20/02/2018) tot eind juli (27/07/2018) om de drie weken bezocht en bemonsterd. Van elk waterlichaam werd telkens een “gebiedsdekkende” bemonstering uitgevoerd. Dit bestond erin een groot aantal sub-stalen (0,5 l) langsheen de oever (elke 2 meter) te nemen d.m.v. een lange stok waaraan op het einde een steriel zakje werd bevestigd. Het bemonsterde water van elke ven of vijver werd voor elke bemonsteringsronde samengevoegd tot één groot mengstaal. Op deze manier werd elk van de geselecteerde vijvers en vennen van winter tot zomer 10 maal bemonsterd (Figuur 6). Van elk mengstaal werd dan nadien een deel van het water gefilterd. Hiervoor werd een filter (type Sterlitech) gebruikt, waarbij het water met behulp van een peristaltische pomp onder grote druk door de filter wordt geperst zodat alle mogelijke eDNA partikels op de filter achterblijven. (zie verder). Na saturatie van de filters werden ze droog geblazen en onmiddellijk bij -20°C ingevroren tot verdere analyse.
Figuur 6. Waterpartijen die doorheen het seizoen werden bemonsterd in functie van de temporele detectie van heikikker en bruine kikker via eDNA in respectievelijk het Nationaal Park de Hoge Kempen (a) en het Vijvergebied Midden-Limburg (b). ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 13 van 40
2.2.2 Soortidentificatie en populatiegenetica Voor de bepaling van de populatiegroottes via soortidentifcatie en genetische screening van de populaties in het Klein Schietveld werden enkele dagen na de eerste screening van het gebied eitjes bemonsterd op al de locaties waar eilegsels van hei- en of bruine kikker werden waargenomen. Per locatie werden hiervoor telkens een drietal eitjes per legsel en minimum 4 tot maximaal 51 legsels genomen, afhankelijk van het aantal aanwezige legsels (Figuur 7a-b). Deze werden gekoeld in het water bewaard waarin ze zich bij staalname bevonden. Op enkele plaatsen waar de larven al uitgekomen waren, werden willekeurig minimaal 30 individuen bemonsterd. Voor de bemonstering van eitjes en/of larven werd voordien een afwijking op het Soortenbesluit (ANB/BL-FF/17-000035) verkregen van het Agentschap voor Natuur & Bos. Tijdens deze tweede bemonsteringsronde werd een telling van het totaal aantal aanwezige eilegsels uitgevoerd, en indien verschillend van de vorige telling (in sommige gevallen waren er in tussentijd nog enkele legsels bijkomend afgezet) werd dit aantal bijgesteld. Tot slot werd voor een selectie van 20 locaties een liter water bemonsterd rondom de legsels (Figuur 7c), om dit in een later stadium te kunnen onderwerpen aan eDNA analysen. Alle waterstalen werden zo snel mogelijk na staalname ingevroren (-20°C), om degradatie van het eDNA zoveel mogelijk te beperken. Drie weken later, wanneer alle eitjes in het veld ontwikkeld waren tot kikkervisjes, werd deze procedure van waterstaalname op exact dezelfde locaties herhaald. Tegelijkertijd werden enkele parameters van het water gemeten zoals: pH, temperatuur, zuurstofgehalte en conductiviteit. De ingevroren waterstalen zijn later gefilterd over een Sterivex filter.
Figuur 7. Illustraties van de staalnamecampagne op 23/03/2017 op het Klein Schietveld. (a) Typische cluster van eilegsels van heikikker en/of bruine kikker, (b) staalname van een drietal eitjes uit elk van de verschillende legsels per cluster, en c) staalname van water rondom de cluster van legsels voor het eDNA-onderzoek.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 14 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
2.2.3 eDNA-afgifte onder gecontroleerde condities Om onder gecontroleerde condities in het water de eDNA afgifte van heikikker en bruine kikker op het stadium van kikkervisjes te kunnen bepalen, werden in het voorjaar van 2017 in het veld 120 eitjes per soort verzameld. Voor heikikker werden hiervoor eitjes genomen op het Klein Schietveld, op een locatie waar we 100% zeker waren dat het om heikikker ging. Voor bruine kikker werden eitjes bemonsterd in een polderlandschap (Bazel), zodat er eveneens 100% zekerheid was dat we met bruine kikker te maken hadden. Deze kikkervisjes werden vervolgens in de volgende ratio’s per soort verdeeld over 11 bakken: 0:10, 1:9, 2:8, 3:7, 4:6, 5:5, 6:4, 7:3, 8:2, 9:1 en 10:0. Elk van deze bakken (70 x 70 x 50 cm) werd twee dagen voor aanvang van het experiment gevuld met 20 l leidingwater en buiten opgesteld zodat de temperatuur min of meer overeenkwam met de veldcondities, en mogelijke chloorresten uit het water konden verdampen, waarna de kikkervisjes onder de verschillende ratio’s in de bakken werden uitgezet. Nadat de kikkervisjes 48 uur in de bakken vertoefden werd per bak 0,5 l water genomen en onmiddellijk gefilterd over een filter (type Sterivex). Deze procedure werd per bak twee maal herhaald zodat in totaal 22 filters werden bekomen. Na filtering werden de filters onmiddellijk ingevroren bij -20°C tot verdere analyse (eDNA-extractie en uiteindelijke detectie).
2.3 eDNA-EXTRACTIE UIT WATER Alle eDNA-filters en afgeleide eDNA-stalen werden bewaard en behandeld in het labo in onze vestiging in Geraardsbergen dat speciaal voor eDNA-onderzoek, en dus lage DNAconcentraties, is ingericht. De eDNA-filters werden allen volgens het licht gewijzigde SXCAPSULE-extractieprotocol geëxtraheerd. Dit protocol is ontwikkeld voor filters zonder preservatiebuffer, zoals beschreven in Spens et al. 2017. Doorheen de eDNA workflow in het labo bouwen we een reeks testen of controlestappen in die moeten toelaten de betrouwbaarheid van het volledig proces en de uiteindelijke kwaliteit van de resultaten te toetsen. Dit geeft ons onder andere inzichten over de efficiëntie van de eDNA-extractie en de kwaliteit van het geëxtraheerde eDNA, voordat tot effectieve interpretatie van de resultaten wordt overgegaan. Om vals-negatieve detectie uit te sluiten, zowel bij eDNA-extractie als op het niveau van de PCR amplificatie, werd aan de Sterlitech filters, genomen voor de temporele screening, systeemvreemd DNA toegevoegd.
2.4 DNA-EXTRACTIE UIT EITJES EN LARVEN Omdat extractie op enkel het embryo een merkelijk betere DNA-opbrengst opleverde bij voorafgaande testen, werd eerst het uitwendige geleimembraan van elk eitje met een scalpel van het embryo gescheiden. Per eilegsel werden twee individuele eitjes geëxtraheerd. De overgebleven eitjes werden op ethanol en koel bewaard. De extractie gebeurde met de Blood & Tissue Kit (Qiagen) met een lysisstap van 1 uur. Het verkregen extract werd geëlueerd in 70 µl buffer (elutie twee maal uitgevoerd). De extractie van de larven gebeurde ook met de Blood & Tissue Kit en elutie in 140 µl. DNA-extracten die bestemd waren voor ddPCR-analyse werden, na DNA concentratiebepaling, verdund naar 1 ng /µl vooraleer binnen te brengen in het labo voor eDNA-onderzoek. Naast het extraheren van individuele eitjes, werd een bijkomende test uitgevoerd waarbij mengstalen van eitjes werden genomen. Dit waren elf mengstalen met telkens 10 eitjes, met ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 15 van 40
variabele samenstelling van heikikker en bruine kikker. Hiervoor werd hetzelfde protocol voor DNA-extractie gebruikt als voor de afzonderlijke eitjes. Voorafgaand zijn de mengstalen van 10 eitjes verbrijzeld met een TissueLyzer (3 min bij 30 Hz).
2.5 DDPCR-ANALYSEN Voor de detectie en kwantificatie van zowel het eDNA uit de waterstalen als deze afkomstig van de gemengde weefselstalen (d.z. mengstalen van eitjes van beide soorten) gebruikten we de digital-droplet-PCR (of kortweg ddPCR) methode. Deze methodiek maakt gebruik van een verstuiving van de PCR-reactie, waarbij een groot aantal druppels of droplets wordt aangemaakt. Binnen deze droplets zitten dan zowel de (e)DNA-template van het extract als de reagentia vervat, en hierop kan de PCR-reactie plaatsvinden. Doordat de PCR-reactie in deze mini-droplets gebeurt (< 1 nl), is deze methode veel minder onderhevig aan eventuele inhibitoren die de amplificatie van het target-DNA (d.i. een specifieke regio in het DNA van de doelsoort) kunnen belemmeren of verhinderen. Bijkomend geeft deze methode, door het grote aantal droplets dat per analyse wordt gegenereerd en dus onderzocht (> 20000), een exacte en erg betrouwbare inschatting van het totaal aantal eDNA partikels dat per μl aanwezig is, dit in tegenstelling tot de indirecte relatieve waarde die men bekomt bij traditionele qPCR-analyses (kwalitatieve PCR). Voor de detectie van eDNA afkomstig van heikikker en bruine kikker werden twee soortspecifieke primer/probe-sets ontwikkeld. Hiervoor werd het Geneious Software pakket gebruikt (versie 11.1.5; https://www.geneious.com). Voor heikikker betrof dit een primer/probe-set dat een 113 bp groot amplicon amplificeert op het Cytochrome C Oxidase subunit 1 gen (COI), in het mitochondriaal genoom. Voor bruine kikker was dit een primer/probe-set dat een 130 bp amplicon overspant dat eveneens gelegen is op het COI-gen. Als referentiemateriaal werd DNA gebruikt afkomstig van individuen van beide soorten uit de betreffende studiegebieden, alsook sequenties uit de NCBI databank (https://www.ncbi.nlm.nih.gov/nucleotide/) van individuen afkomstig uit de ons omringende landen. Specificiteit van beide assays werd in silico getoetst t.o.v. een groot aantal in Vlaanderen voorkomende amfibie- en vissoorten, en dit d.m.v. de ‘multiple alignment functie’ in Geneious (zoals rugstreeppad, boomkikker, knoflookpad, etc., zie Bijlage 1 en 2). Deze resultaten werden bijkomend visueel gescreend op hun nauwkeurigheid. Tot slot werd voor het bekomen primer/probe-assay van beide soorten een in vitro validatie uitgevoerd met DNA uit weefselstalen van andere amfibiesoorten die we eveneens in de betreffende vijvers en vennen kunnen aantreffen, om mogelijke interferentie voor 100% te kunnen uitsluiten.
2.6 SOORTIDENTIFICATIE VAN DE EILEGSELS VIA RFLP Per legsel werd één eitje geanalyseerd. Via de RFLP-methode (d.i. restriction fragment length polymorphism) werd de soortbepaling van elk van de eitjes uitgevoerd (Palo & Merilä, 2003). Deze methode benut de variatie die bestaat tussen soorten in een gedeelte (605 bp) van het cytochroom b (cytb) gen in het mitochondriaal genoom. Door gebruik te maken van het restrictie-enzyme HaeIII wordt de DNA-sequentie in kleinere stukjes geknipt die door verschillen in hun lengte kunnen worden gescheiden en zichtbaar gemaakt d.m.v. horizontale agarose gelelectroforese (AGE). Door deze methode toe te passen op het DNA van de ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 16 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
bemonsterde eitjes kan een duidelijk onderscheid worden gemaakt tussen heikikker en bruine kikker (Figuur 8). De PCR werd uitgevoerd op 2 µl onverdund DNA in een totaal volume van 20 µl (1x PCR buffer, 250 µM dNTPs, 0,1 µM van elke primer, 1 U Taq polymerase en 1,5 mM MgCl2). De PCRamplificatie werd gecontroleerd door middel van horizontale agarose gelelectroforese (2% agarose). De PCR-condities beschreven door Palo & Merilä (2003) werden aangepast als volgt: initiële stap van 3 min bij 94 °C, gevolgd door 35 cycli van 30 s bij 94°C, 30 s bij 49 °C waarbij de primers aan het DNA zullen hechten en 30 s bij 72 °C (annealing temperatuur).
Figuur 8 Resultaat van een RFLP-analyse met duidelijke amplificatie en opsplitsing van de verschillende DNA-fragmentjes in gescheiden bandjes na gelelectroforese. Elke verticale kolom is een staal.
Tenslotte volgde een finale elongatiestap van 5 min bij 72 °C. Aan 10 µl PCR product werd 1 U HaeIII, 2 µl 10x Buffer R (Thermo Scientific) en 7 µl ultra-puur water (Milli-Q van Millipore Corporation) toegevoegd. Deze oplossing werd gedurende 2 uur geïncubeerd bij 37 °C. De restrictie werd gecontroleerd door middel van horizontale agarose gelelectroforese (2% agarose). In totaal werd 11% van de stalen een tweede keer geanalyseerd ter controle.
2.7 POPULATIEGENETISCHE ANALYSES 2.7.1 Genotyping Verschillende gepubliceerde microsatellietmerkers werden uitgetest op een beperkt aantal stalen met als doel een voldoende aantal polymorfe merkers te verkrijgen voor deze studie, die bovendien in een optimaal aantal multiplexsets gelopen kunnen worden. Van de uiteindelijke set van 19 merkers werden er twee in simplex en de overige in twee multiplexsets gelopen (Bijlage 3). Naast het specifieke volume aan primers met een concentratie van 10 µM, zoals vermeld in Bijlage 3, werd in elke multiplex- en simplexmix 1 µl DNA met een concentratie van c. 5 ng/µl en 5 µl Multiplex Master Mix (Qiagen) toegevoegd en verder aangevuld met ultra-puur water tot een totaal volume van 10 µl. ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 17 van 40
De PCR-condities waren voor de eerste multiplex als volgt: 15 min bij 95°C gevolgd door 35 cycli met telkens 30 s bij 95°C, 45 s bij 50°C (annealing temperatuur) en 45 s bij 72°C; na een elongatiestap van 10 min bij 72°C wordt de temperatuur op 4°C gebracht gedurende 15 min waarna een temperatuur van 15°C aangehouden wordt. Voor multiplex 2 werden dezelfde condities gebruikt maar met een annealing temperatuur van 55°C. Voor merker RA13 werd een touch-down procedure gevolgd met het volgende programma: 15 min bij 94°C, 22 cycli van 30 s bij 94°C, 30 s beginnend bij 60°C met een afname van een halve °C per cyclus, 30 s bij 72°C gevolgd door 15 cycli van 30 s bij 94°C, 30 s bij 50°C en 30 s bij 72°C, eindigend met 10 min bij 72°C en 15 min bij 4°C. Het PCR product wordt dan op 15°C bewaard. Voor merker RtU4 gelden dezelfde condities met uitzondering van de touch-down cycli. Hier zijn 35 cycli gebruikt met telkens 30 s bij 94°C, 30 s bij 46°C en 30 s bij 72°C. De genotyping gebeurde op een ABI 3500 Genetic Analyzer (Applied Biosystems) waarbij een verdunning van de PCRproducten van 1/10 werd gehanteerd voor beide multiplexsets. De simplexsets werd in poolplex geanalyseerd met een verdunning van 1/50. Het identificeren van de allelen gebeurde met het programma GeneMapper 4.1 (Applied Biosystems). Doorgaans werden twee eitjes per legsel geanalyseerd. Circa 3% van de stalen werd twee tot vijf maal herhaald binnen en/of tussen platen. Wanneer eitjes een slechte amplificatie gaven, werden ze vervangen door genotypes van nieuwe eitjes van hetzelfde legsel. De merkers werden gecontroleerd op afwijkingen van het Hardy-Weinberg evenwicht, op de aanwezigheid van nulallelen en op linkage disequilibrium (d.i. de niet-toevallige associatie tussen allelen van verschillende loci binnen een populatie) met genepop 4.6 (Rousset, 2008). Voor de beoordeling van het percentage nulallelen kozen we de methode van grootste aannemelijkheid (maximum likelihood) volgens het Expectation Maximization algoritme van Dempster et al. (1977). Ter vergelijking werd ook het aandeel nulallelen bepaald met MLNullFreq (Kalinowski & Taper, 2006). Omdat er een sterke mogelijkheid bestaat op het overbemonsteren van sterk verwante individuen (broers of zussen) en deze een invloed kunnen hebben op de voornoemde analyses, werden de tests herhaald met enkel één per koppel broers/zussen over te houden in de dataset. De broers/zussen werden geïdentificeerd met behulp van het programma COLONY 2.0.6.5 (Jones & Wang, 2010). Enkel de combinaties met een probabiliteit hoger dan 0,9 werden als zeker beschouwd. De analyse gebeurde met hoge precisie, drie runs en met een polygaam paringssysteem. De informatie dat eitjes afkomstig van hetzelfde legsel afkomstig waren van dezelfde moeder en omgekeerd kon in de analyse mee opgenomen worden.
2.7.2 Genetische structuur We hebben de spatiale structuur binnen de populatie onderzocht met behulp van paarsgewijze differentiatie volgens Weir & Cockerham (1984) (Fst; met 1000 bootstraps) tussen poelen en met een spatiale principale componenten analyse (sPCA) met respectievelijk de R pakketten diveRsity 1.9.90 (Keenan et al., 2013) en adegenet 2.1.0 (Jombart, 2008) in de R omgeving versie 3.5.2 (R Core Team, 2018). De Fst resultaten werden gevisualiseerd in een principale coördinaten analyse (PCoA). Een sPCA maakt gebruik van allelfrequenties en brengt hun genetische variabiliteit en spatiale autocorrelatie in rekening met Moran’s I (Moran, 1948; Moran, 1950). Als netwerk gebruikten we een neighbourhood-by-distance netwerk met 0 m als minimale afstand (binnen een poel) en 1 km als maximale afstand, een afstand die als maximale afstand tijdens de lentemigratie werd gevonden voor de soort (Kovar et al., 2009). We onderzochten of er een relatie bestaat tussen genetische en geografische afstand (isolation by distance; IBD) via een Manteltest, geïmplementeerd in adegenet, met 9999 permutaties. Daarnaast werd de spatiale autocorrelatie geanalyseerd van de genotypes met ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 18 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
SPAGeDi 1.5a (Hardy & Vekemans, 2002). Vooreerst werd nagegaan welke maat van verwantschap een goede schatter zou zijn via simulaties in COANCESTRY 1.0.1.9 (Wang, 2011). De gesimuleerde populatie bestond uit 100 paren voor elk van vier specifieke relaties: ouderkind, broer/zus, halfbroer/halfzus en neef/nicht. Tegelijkertijd krijg je paren van onverwanten gesimuleerd. De simulaties gebeuren op basis van de beschikbare genotypes, rekening houdend met de gevonden fout per merker en mogelijke allel drop-out. Dat laatste werd bepaald met het programma Pedant op basis van de dubbele stalen (Johnson & Haydon, 2007a; Johnson & Haydon, 2007b). Als schatter voor verwantschap werd deze van Lynch & Ritland (1999) gekozen, die tegelijkertijd relatief goed scoorde in de simulaties (zie resultaten) en ook door SPAGeDi berekend kon worden. We lieten het programma een optimale verdeling van afstanden over acht klassen bepalen (10 klassen voor de set zonder broers/zussen), waarbij rekening wordt gehouden met het aantal paren, de proportie van alle individuen en het aantal keer dat elk individu in elke klasse opgenomen is. We stelden 1000 permutaties in.
2.7.3 Genetische variatie en effectief aantal ouders Genetische variatie onder de vorm van geobserveerde en verwachte heterozygositeit (Ho en He), allelische rijkdom (Ar; met rarefactie) en de fixatie-index (Fis) werden geschat met het R pakket hierfstat 0.04-22 (Goudet, 2005). Daarnaast werd de genetische effectieve populatiegrootte bepaald via een paar momentschatters. Omdat we enkel eitjes bemonsterd hebben binnen eenzelfde voortplantingsseizoen, bepalen we hier eerder het effectief aantal ouders (Nb) dan de effectieve populatiegrootte. De ene methode baseert zich op de verwantschap tussen de individuen via het programma COLONY en de andere op de mate van linkage disequilibrium (LD) genaamd LDNe (Waples & Do, 2008) met behulp van de software NeEstimator 2.01 (Do et al., 2014). Zeldzame allelen worden best buiten beschouwing gelaten in deze laatste schatting (hier met grenswaarde 0,02). We hebben Nb bepaald per poel en over de hele populatie heen. De schattingen bekomen met de LD-methode werden aangepast zoals voorgesteld door Waples et al. (2014). Met behulp van twee levensgeschiedeniskenmerken, volwassen leeftijd en volwassen levensduur, kan een betere schatting van Nb van een enkele cohorte gemaakt worden wanneer we met overlappende generaties te maken hebben. We gebruikten daarvoor een leeftijd van twee jaar waarop heikikker seksueel rijp is (Lardner & Loman, 2009; Söderman, 2006) en een maximale leeftijd van 11 jaar (Ishchenko, 2005). We noemen deze aangepaste schatting Nb adjust. De analyse met COLONY bepaalde naast effectieve grootte ook het aantal ouders van de onderzochte cohorte. Het totaal aantal legsels van heikikker is daarnaast een indicatie van het aantal voortplantende vrouwtjes. Dit aantal werd per poel berekend aan de hand van de resultaten van de soortbepaling, waarbij de proportie bemonsterde legsels van heikikker werd vermenigvuldigd met het totaal aan getelde legsels. Door het aantal bemonsterde legsels te vergelijken met het aantal ouders kan de verhouding tussen mannetjes en vrouwtjes benaderd worden en bijgevolg het totaal aantal reproducerende adulten. De potentiële invloed van de in het veld gemeten abiotische kenmerken van de poelen (pH, temperatuur, zuurstofgehalte en conductiviteit) op de schattingen van Nb en op het aantal genkopieën/µl gemeten via ddPCR op waterstalen werd bepaald met generalised linear models (GLM). Als extra covariate werd de geschatte oppervlakte van de poel toegevoegd. Deze werd berekend met behulp van luchtfoto’s van 2017 (Informatie Vlaanderen) en ArcGIS ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 19 van 40
10.4.1 (Esri). In het model met ddPCR-resultaten als respons werd tevens het aantal reproducerende adulten van heikikker als verklarende variabele toegevoegd. Via achterwaartse eliminatie werd op basis van het Akaike Informatie Criterium (AIC) een finaal model geselecteerd. Vooraleer de volledige modellen werden gevormd, controleerden we of er geen collineariteit bestond tussen variabelen met behulp van Variantie Inflatie Factor (VIF) met R pakket usdm 1.1-18 (Naimi et al., 2014), waarbij we een grenswaarde van drie hanteerden. De analyses werden uitgevoerd met R pakket nlme 3.1-137 (Pinheiro et al., 2016).
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 20 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
3 RESULTATEN EN DISCUSSIE 3.1 EDNA AFGIFTE VAN HEI- EN BRUINE KIKKER 3.1.1 Temporele detectie onder veldcondities De temporele eDNA-patronen in de vier verschillende vijvers en vennen met gekende heikikker- en bruine kikker populaties, tonen aan dat we via het nemen van een mengstaal in staat zijn om de aanwezigheid van beide soorten duidelijk aan te tonen (Figuur 9). Het temporele eDNA profiel van beide soorten geeft een duidelijke eerste piek vroeg in het voorjaar (vanaf half maart tot begin april). Deze valt samen met het actief worden van de adulte individuen na de winter en de periode van reproductie wanneer de temperaturen van het water in het voorjaar stijgen (zie Figuur 10).
Figuur 9. Temporeel eDNA-profiel (gemeten via ddPCR) (gemiddelde van twee analysen per staal Âą SE) van bruine kikker en heikikker in het water van twee vennen en vijvers in respectievelijk het Nationaal Park de Hoge Kempen (a-b) en het Vijvergebied Midden-Limburg (c-d). Waterstalen werden driewekelijks genomen vanaf eind februari t.e.m. eind juli 2018.
Een vergelijking van het temperatuurverloop van het water (gemeten op elk van de bemonsteringen) met de temporele eDNA-profielen geeft aan dat wanneer de ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 21 van 40
watertemperatuur na de winter boven de 10°C stijgt, dit beide soorten lijkt aan te zetten om tot reproductie over te gaan (Figuur 9 en 10).
Figuur 10. Temporeel profiel van de watertemperatuur (C°) gemeten op elk van de eDNA bemonsteringen in twee vennen 3a en 3b en vijvers 4A en 4b in het Nationaal Park Hoge Kempen en het Vijvergebied Midden-Limburg. Stalen werden driewekelijks genomen vanaf eind februari t.e.m. eind juli 2018.
Drie weken na het bereiken van de eerste eDNA-piek, valt het eDNA-signaal even snel weer terug als dat het verscheen, wat allicht te verklaren is door het feit dat de adulte individuen zich na de paring en het afzetten van de eitjes niet langer in het water ophouden. Verder valt bij heikikker ook op dat de sterke terugval van de eerste piek gevolgd wordt door een graduele eDNA-toename in het water met een duidelijke tweede piek rond eind mei – begin juni. Deze wordt mogelijk veroorzaakt door de opgroeiende juvenielen, die daarna ook in metamorfose gaan en het water verlaten. Opmerkelijk hierbij is echter wel dat een tweede piek zich niet of nauwelijks voordoet bij bruine kikker. Dit opmerkelijke verschil tussen beide soorten is mogelijk te verklaren doordat: bruine kikker een snellere metamorfose van zijn larven vertoont en juvenielen sneller het water verlaten (Loman et al., 2002), grotere sterfte van bruine kikker eitjes/larven omdat de habitat voor de soort mogelijk suboptimaal is, bruine kikkerlarven een lagere eDNA-afgifte vertonen dan heikikker, wat eveneens naar voor kwam in het mesocosmexperiment (zie verder, Figuur 11).
3.1.2 eDNA-afgifte onder gecontroleerde condities De eDNA-afgifte van heikikker- en bruine kikkerlarven in het mesocosmexperiment toont voor beide soorten een significant positief verband aan tussen de graduele toename van elke soort en concentratie eDNA van beide soorten (r = 0,969; P = 0,0001 voor bruine kikker en r = 0,915; P = 0,003 voor heikikker; Figuur 11a-b).
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 22 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Figuur 11. Relatie tussen het aantal bruine kikker- (a) en heikikkerlarven (b) in elk van de 11 mesocosmen en het eDNA signaal gemeten via ddPCR.
Verrassend is echter wel het opmerkelijk significante verschil in eDNA-afgifte tussen beide soorten (t = -4,001; P = 0,0007), waarbij heikikker gemiddeld genomen 25 maal meer eDNA in het water lijkt af te geven dan bruine kikker. Dit is opmerkelijk aangezien de larven die voor dit experiment werden gebruikt zich allen in hetzelfde levensstadium bevonden en ook in grootte erg gelijkaardig waren. Of dit sterke verschil een natuurlijk verschil in eDNA-afgifte is dat zich ook onder natuurlijke condities voordoet, is niet duidelijk. Het zou bijvoorbeeld kunnen dat dit verschil deels in de hand zou zijn gewerkt door een hogere stress bij heikikkerlarven door in dit experiment gebruik te maken van leidingwater. Het is namelijk niet ondenkbaar dat het sterke verschil in pH (in vergelijking met de condities in hun natuurlijke leefgebied), een stressfactor heeft veroorzaakt die een verhoogde eDNA afgifte in de hand kan werken.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 23 van 40
3.2 BEPALING VAN DE VERHOUDING BRUINE KIKKER- EN HEIKIKKEREILEGSELS 3.2.1 Soortidentificatie van de eilegsels via RFLP In Tabel 1 worden de resultaten van de RFLP-analyse gegeven. De proportie legsels van heikikker per bemonsterde poel varieert van 0 tot 100% met een gemiddelde van 76%. In bepaalde gevallen was er geen amplificatie of amplificatie van slechte kwaliteit. Omdat deze in aantal beperkt waren, werden ze niet herhaald.
Tabel 1. Lijst van aantal getelde en bemonsterde eilegsels per poel op het Klein Schietveld. De code van de poel komt overeen met de codes gegeven in Figuur 5. Het percentage legsels van heikikker werd berekend rekening houdend met de soortbepaling van larven (poelen 1 en 24). Ook het aantal legsels waarbij de RFLP-analyse geen resultaat gaf zijn gegeven.
Code poel
Aantal legsels
Aantal legsels bemonsterd (aantal larven)
Aantal legsels zonder soortbepaling
01 13 13 0 02 70 31 1 1 75-85 9 (30a) 0 10 37 17 1 11 16 12 1 13 47 18 0 14 35 14 0 16 60-65 37 0 17 35 21 1 18 50-55 51 0 19 20-25 21 1 20 4 4 0 21 20-25 18 1 22a 10013 0 120 22b 70-80 17 0 22c 38 15 0 24 500 0 (30a) 25 50-55 25 0 3 14 14 0 30 42 36 0 31 50 19 2 32 3 3 1 4 40 21 1 5 9 9 0 6 85 41 1 9 17 13 0 a Alle bemonsterde larven waren bruine kikkers
Aantal legsels bruine kikker 0 27 6 0 0 2 0 27 13 1 0 0 10 1 1 0
Aantal legsels heikikker
% heikikker
13 3 3 16 11 16 14 10 7 50 20 4 7 12 16 15
0 0 36 1 0 0 0 0 0
25 14 0 16 2 20 9 40 13
100 10 4 100 100 89 100 27 33 98 100 100 41 92 94 100 0 100 100 0 89 100 100 100 100 100
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 24 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
De RFLP bleek een zeer eenvoudige methode te zijn om de soort van de eitjes te identificeren. In het Klein Schietveld waren de meeste bemonsterde poelen bezet door eitjes van heikikker, met een gemiddelde van 76% legsels van heikikker (standaard fout = 7,4). Slechts in zes poelen van de 24 (waarbij 22a, b en c als ĂŠĂŠn poel beschouwd worden) zijn er meer legsels van bruine kikker.
Figuur 12. Relatie tussen het aandeel (%) heikikker legsels en de pH gemeten in het water rondom de clusters van eilegsels.
Wanneer deze resultaten gecombineerd worden met de pH-metingen, uitgevoerd tijdens het nemen van de waterstalen voor het eDNA-onderzoek, blijkt er zich een sterk verband voor te doen tussen het percentage heikikkerlegsels en de gemeten pH (Figuur 12). Over alle gemeten vennen was de pH gemiddeld 5,24, met een mediaan van 4,98. In de vennen en plassen waar de pH < 5 was, bleek het aandeel heikikkerlegsels 100%, maar eens de pH boven de 5 stijgt, blijkt bruine kikker eveneens op deze locaties aanwezig te zijn. De twee locaties waar een pH > 6 werd gemeten, bleken ook de locaties waar bruine kikker in aandeel eilegsels het sterkst domineerde (> 90% van de eilegsels). Dit patroon geeft aan dat bruine kikker zich ofwel niet in en om de typische zure vennen ophoudt (die zich min of meer centraal in het studiegebied bevonden) ofwel dat de soort overal wel aanwezig is maar geen eilegsels deponeert in water boven een bepaalde zuurtegraad.
3.2.2 Bepaling van de soortenverhouding op basis van mengstalen van eitjes Binnen de ddPCR-resultaten afkomstig van het experiment waarin een DNA-extract werd bekomen van 10 eitjes met variabele samenstelling van heikikker en bruine kikker, blijkt er zich een significant verband voor te doen tussen de DNA-concentratie van elke soort en hun respectievelijk aandeel in elk van de 11 mengstalen (zie Figuur 13a). Het valt hierbij wel op dat er voor beide soorten zich een aanzienlijke variatie voordoet tussen de verschillende mengstalen. Opvallend is echter dat wanneer de concentraties van beide soorten binnen elk staal proportioneel t.o.v. elkaar worden uitgezet (% heikikker-DNA t.o.v. % bruine kikker-DNA), ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 25 van 40
deze variatie tussen de verschillende stalen bijna volledig verdwijnt en we erg accurate inschattingen krijgen die bijna perfect overeenstemmen met de verhouding aan heikikker en bruine kikker eitjes (Figuur 13b). Dit geeft dus aan dat deze methode op termijn erg efficiĂŤnt zou kunnen zijn voor de bepaling van het aandeel heikikker versus bruine kikker eilegsels in functie van de monitoring van de heikikkerpopulaties.
Figuur 13. Verband tussen a) het totaal aantal eDNA-kopieĂŤn van bruine kikker en heikikker (bepaald via ddPCR) (gemiddelde van twee analyses Âą SE) en b) de proportionele verhouding bepaald op basis van het gemiddelde van het eDNA-signaal tussen beide soorten binnen elk staal, en dit in functie van een voor bruine kikker aflopend en voor heikikker oplopend aandeel eitjes op een totaal van 10 eitjes per staal.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 26 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
3.2.3 Bepaling van de soortverhouding op basis van eDNA uit water De eDNA-concentraties van heikikker en bruine kikker die werden aangetroffen in het water dat werd bemonsterd rondom de eilegsels in het Klein Schietveld, bleken voor beide soorten hoog te zijn. Over alle onderzochte locaties heen werden er van heikikker gemiddeld 6541,2 kopieĂŤn per Âľl eDNA-extract gemeten, terwijl dit er voor bruine kikker 607,8 waren (op de locaties waar deze soort aanwezig was). Wanneer we kijken naar het aandeel (%) heikikker eDNA (t.o.v. van de totale eDNA-concentratie van beide soorten), dan blijkt dit significant gerelateerd te zijn (r = 0,820; P = 0,0001) met het werkelijk aandeel heikikkereilegsels op deze locaties (bepaald via de RFLP-analyses zie paragraaf 3.2.1.) (zie Figuur 14).
Figuur 14. Verband tussen het percentage heikikker eDNA (t.o.v. de totale eDNA-concentratie van bruine kikker en heikikker, gemeten via ddPCR) en het aandeel heikikkereilegsels op de onderzochte locaties in het Klein Schietveld.
Wanneer we kijken naar de absolute heikikker eDNA-concentratie die wordt opgepikt op de locaties waar heikikker eilegsels heeft afgezet, lijkt er zich een zwakke trend voor te doen met het totaal aantal heikikkereilegsels dat op deze plaatsen werd afgezet, hoewel niet significant (r = 0,440; P = 0,0874; Figuur 15). In beide gevallen, zowel bij de proportionele als de absolute data, valt echter op dat er een aanzienlijke variatie bestaat tussen de opgemeten eDNAconcentraties en het werkelijke aandeel heikikkereilegsels bepaald via de RFLP-analyses. Het is met andere woorden wel mogelijk om een ruwe inschatting te bekomen van het aandeel heikikkerlegsels op basis van eDNA uit het omringende water, maar een accurate inschatting kan via deze methode niet worden bekomen. Hiervoor kunnen verschillende factoren worden aangehaald die hier mogelijk verantwoordelijk voor zijn. Zo is het niet evident, zeker op plaatsen waar veel eilegsels zijn afgezet, om op een gehomogeniseerde manier water rondom de eilegsels te bemonsteren en dragen ingesloten legsels allicht weinig bij aan de totale eDNA concentratie. Belangrijker is allicht de impact die adulte individuen op de eDNA patronen veroorzaken en die niet kan worden uitgesloten noch ingeschat worden. Dit alles doet dan ook de vraag rijzen of eDNA-bepalingen wel de meest nauwkeurige en efficiĂŤnte manier is om het aandeel heikikker versus bruine kikkereilegsels te kunnen bepalen (zie verder).
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 27 van 40
Figuur 15. Verband tussen het totaal aantal heikikker eilegsels die op elk van de locaties in het Klein Schietveld aanwezig waren het totaal aantal heikikker eDNA moleculen gemeten in het omliggende water (via ddPCR).
3.3 POPULATIEGENETISCH ONDERZOEK Enkele loci bleken op de totale dataset geen variatie te vertonen binnen het Klein Schietveld: Rtempµ4, Rtempµ5 en Rt2Ca2-22. De loci RECALQ en RlatCa41 zijn niet in Hardy-Weinbergevenwicht (P < 0,0001), ook na exclusie van broers/zussen. Bovendien bevatten ze vele nulallelen, respectievelijk gemiddeld 22% en 17% volgens de genepop-resultaten. Hoge proporties nulallelen komen ook voor in RCIDII (gemiddeld 57%) en RNTYR2 (gemiddeld 24%), maar niet wanneer berekend met ML-NullFreq. Tegelijkertijd vertonen RECALQ en RlatCa41 per poel problematische hoeveelheden nulallelen. De proporties zakken slechts licht wanneer we één individu per koppel broers/zussen overhouden in de dataset. Er waren enkele paren loci in LD (vijf combinaties) in maximaal twee poelen. Door exclusie van broers/zussen zakte het aantal significante combinaties naar twee: RA11 versus R13 in poel 3, RlatCa18 versus RA13 in poel 6. Vanwege het Hardy-Weinberg-onevenwicht en de hoge proportie nulallelen, lieten we loci RECALQ en RlatCa41 achterwege voor de verdere analyses. De fout op de resterende merkers ligt op 0,6%. Door stalen met slechte amplificatie (> 50% ontbrekende waarden) te herhalen of te vervangen door stalen afkomstig van dezelfde legsels, zijn alle bemonsterde legsels vertegenwoordigd in de uiteindelijke dataset met minimaal een staal. De range van fragmentgroottes voor de merkers RRD590 en RtµP bleek te overlappen. De multiplexset zal dus in de toekomst aangepast moeten worden wanneer een bredere staalname zou gebeuren, bijvoorbeeld in verschillende populaties, met een verhoogde kans op meer allelen voor beide loci. Paarsgewijze Fst-waarden varieerden tussen -0,0037 en 0,1113 (exclusief broers/zussen: 0,0023 en 0,1057). De poelen 32, 1 en 02 waren het sterkst gedifferentieerd (resultaten niet weergegeven). Hoewel de waarden significant verschillen van nul; hoewel dit logisch lijkt vanwege hun ligging aan de uiterste rand van de populatie, blijkt dit minder sterk op te gaan voor poelen 01 en 31. Vermoedelijk heeft het lage aantal stalen van poelen 32, 1 en 02 hiermee te maken. De Fst-waarden voor deze poelen hadden dan ook zeer ruime ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 28 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
betrouwbaarheidsintervallen. Als we bijgevolg de poelen met staalgroottes kleiner dan 10 weglaten, zijn poelen 01, 17, 16, 21 en 31 het sterkste gedifferentieerd (Figuur 16), met een maximale Fst-waarde van 0,057. Poelen 01 en 31 liggen opnieuw aan de respectievelijke zuidelijke en noordelijke rand van de populatie. Bijgevolg lijken toch eitjes in de zuidelijke en noordelijke poelen te verschillen van elkaar en van de centraal gelegen eilegsels. De individuen in poelen 16, 17 en 21 wijken van dit patroon af.
Figuur 16. De eerste twee assen van de principale coördinaten analyse, die respectievelijk 30 en 23% van de variatie verklaren. De analyse is gebaseerd op de Fst -waarden tussen poelen met meer dan 10 stalen van heikikker, exclusief de genotypes van broers/zussen. De nummers in de codes komen overeen met deze in Figuur 5.
De Mantel r waarde van 0,11 was significant (P = 0,0001). De aanwezigheid van IBD werd ook vastgesteld uit de sPCA resultaten (Figuur 17a), waarbij de globale structuur (positieve spatiale autocorrelatie) een noord-zuid trend vertoonde. Deze werd bepaald door een enkele principale component. Enkel poel 18 leek hiervan af te wijken. Deze poel is niet zichtbaar in de Figuur 17a omdat hij een zeer kleine positieve lagged score heeft. Ook de lokale structuur (negatieve spatiale autocorrelatie) werd voornamelijk bepaald door één principale component (Figuur 17b). In plaats van een trend vertoonde de structuur eerder vier groepen elk bestaande uit nabijgelegen poelen, al lijken de centraal gelegen poelen ook hier meer op elkaar. Opnieuw wijkt poel 18 af en in mindere mate ook 19 en 21. Samenvattend kunnen we stellen dat er globaal een IBD-structuur is, met enkele poelen die lokaal lijken af te wijken. Mogelijk heerst daar een iets andere sociale structuur die spatiaal beïnvloed wordt door het landschap of andere omgevingsfactoren.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 29 van 40
Figuur 17. Grafische weergave van de sPCA resultaten met (a) de globale structuur en (b) de lokale structuur t.o.v. de x- en y-coördinaten van de locaties. De blokjes geven de waarde van lagged scores volgens grootte en kleur: zwart is positief, wit is negatief.
De beste schatter van verschillende types verwantschappen uit de simulaties van COANCESTRY bleek TrioML te zijn. Deze leverde de meest nauwe betrouwbaarheidsintervallen voor de gesimuleerde waarden en de hoogste correlatie tussen gesimuleerde en echte waarden (r = 0,11, P < 0,0001). Deze schatter komt niet voor in SPAGeDi en werd dus vervangen door de coëfficiënt van Lynch & Ritland (1999) die de hoogste correlatie met TrioML had (r = 0,756) en een gelijkaardige correlatie tussen de gesimuleerde en geschatte waarden (r = 0,086, P < 0,0001). De resultaten van de spatiale autocorrelatie gaven aan dat tot en met een maximale afstand van 305m de verwantschappen significant verschilden van random paring (Figuur 18a). Bij grotere afstanden werden ze niet-significant. In de laatste afstandsklassen van 900-1950m was de relatie significant negatief, maar deze moet met voorzichtigheid behandeld worden vanwege een hoge coëfficiënt van variatie voor het aantal keer dat elk individu is vertegenwoordigd in die klasse (> 1). Wanneer broers/zussen verwijderd werden uit de dataset is de relatie tussen verwantschap en afstand nog significant tot een afstand van 219m, maar is opnieuw net significant tussen 343m en 431m (P = 0,035) (Figuur 18b). ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 30 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Figuur 18. Spatiale autocorrelatie van genetische verwantschap en geografische afstand voor (a) de totale dataset en (b) de dataset zonder broers/zussen.
Het IBD-patroon wordt dus hoogst waarschijnlijk bepaald door de gemiddelde dispersieafstand die heikikker aflegt binnen het Klein Schietveld. Hier lijkt dat rond 300m te liggen, al is er waarschijnlijk spatiale variatie door de omgeving, lokale dichtheden van individuen en door de ligging van de poelen (centraal of aan de periferie). De genetische variatie verandert weinig na verwijdering van de broers/zussen. Ze is bovendien gelijkaardig in alle poelen (Ho = 0,51 – 0,66; He = 0,51 – 0,62 (Figuur 19); Ar = 3,1 – 3,7). Hoewel Fis varieert van -0,09 tot 0,13 is ze nergens significant verschillend van nul. Hoewel we hier grotendeels ander microsatellieten gebruiken, lijken de waarden van genetische variatie gelijkaardig aan deze in andere populaties (bijv. Arens et al., 2007).
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 31 van 40
Figuur 19 Verwachte heterozygositeit (gendiversiteit) Âą standaarddeviatie per poel.
Zoals eerder vermeld kunnen een abnormaal aantal broers/zussen in de steekproef een effect hebben op de schattingen, ook op Nb. Soms wordt aangeraden om die stalen van broers/zussen te verwijderen, maar uit recent onderzoek lijkt dit niet altijd tot de juiste resultaten te leiden (Waples & Anderson, 2017). Omdat we hier telkens twee eitjes per legsel hebben bemonsterd en er dus een dominantie aan ĂŠĂŠn of enkele families niet mogelijk is, hebben we voor de schattingen van Nb alle stalen behouden. Het effectief aantal ouders werd per poel maar ook over de hele dataset heen bepaald (Tabel 2). In enkele gevallen kon geen eindig betrouwbaarheidsinterval bekomen worden voor de schattingen. Dit is vooral wanneer de staalgrootte laag is. Als we een som nemen van de schattingen per poel, dan liggen de totalen dicht bij elkaar en relatief dicht bij de schatting van Nb met de verwantschapsmethode voor de totale populatie. De LD-schatting van Nb voor de totale populatie ligt beduidend lager (409 na aanpassing). Omdat er IBD aanwezig is in de populatie, kan de geschatte Nb gebaseerde op de gehele dataset een onderschatting zijn van de echte waarde door de aanwezigheid van mixture-LD wanneer je staalname meerdere nabuurschappen bevat (Neel et al., 2013; Waples & England, 2011). Met een genetische nabuurschap bedoelen we hier een groep van individuen die sterk interageren. Hier bleek sterke interactie of verwantschap aanwezig te zijn binnen een afstand van ongeveer 300 m. De schattingen van Nb via de verwantschapsmethode zullen daarentegen steeds blijven stijgen naarmate meer legsels bemonsterd worden. Zo is gebleken dat de schattingen van Nb met de verwantschapsmethode vooral dicht bij de echte waarde ligt als het aantal stalen de werkelijke Nb benaderd of overstijgt (Ackerman et al., 2017). We hebben niet alle legsels bemonsterd, laat staan alle poelen. Op basis van de getelde legsels zouden we een totaal aantal van c. 725 legsels van heikikker geteld hebben. De verhouding tussen aantal moeders en vaders bleek gemiddeld 0,41 te zijn per poel. Wanneer we over de totale populatie heen kijken, ligt deze ratio op 0,59. Als we bijgevolg een gemiddelde ratio van 0,5 nemen, zouden er minstens 1450 reproducerende adulten aanwezig zijn geweest in Klein Schietveld in 2017.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 32 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Tabel 2 Schattingen van effectief aantal ouders (Nb) volgens twee methodes: de linkage disequilibriummethode (LD) en de verwantschapsmethode met bijhorende betrouwbaarheidsintervallen tussen haakjes. De aangepaste waarde van Nb LD is Nb adjust. Het aantal ouders bepaald met COLONY wordt gegeven. Door het aantal bemonsterde legsels (dus het aantal moerders) te delen door het aantal ouders, kan de verhouding moeders/vaders berekend worden.
Poelcode
Nb LD
Nb adjust
Nb verwantschap
Aantal ouders
Aantal eilegsels
Moeders/ vaders
01 02 1 10 11 13 14 16 17 18 19 20 21 22A 22B 22C 25 3 31 32 4 5 6 7 Som KS
16,1 (10,3-28) 6,4 (2,4-32,7) 2,9 (2,1-6,9) 31,8 (21,7-51,5) 24,6 (16,8-40,2) 29,8 (20,9-46,1) 22,3 (14,7-37,7) 21 (12,9-41,6) 8,2 (5,2-13) 54 (40,7-74,1) 48,7 (31,5-88,6) 2,8 (2-7) 15 (9,5-26,5) 55,7 (31,1-172,2) 31 (22-47,1) 18,5 (13,2-27,2) 41,6 (32,3-55,8) 19,6 (12,7-33,6) 34,8 (23,3-58,1) 16,9 (3,8- ∞) 40,3 (24,5-80,5) 38,3 (18,3-278,8) 75,4 (50,2-128,2) 21,2 (15,6-30,3) 676,9 381,3 (290,1-519,1)
17,3 6,9 3,1 34,1 26,4 32 23,9 22,5 8,8 58 52,3 3 16,1 59,8 33,3 19,9 44,6 21 37,3 18,1 43,3 41,1 80,9 22,8 726,5 409,2
26 (14-52) 7 (2-41) 10 (4- ∞) 29 (17-54) 15 (7-34) 31 (18-55) 24 (13-47) 18 (10-39) 17 (9-39) 95 (69-129) 33 (20-56) 9 (4-28) 17 (8-41) 19 (10-39) 29 (17-51) 36 (20-71) 45 (29-69) 28 (14-60) 38 (24-64) 5 (2-32) 41 (26-68) 12 (6-28) 64 (46-95) 22 (12-43) 670 615 (541-691)
28 7 7 42 30 43 38 26 17 112 57 10 20 33 37 38 58 34 42 8 55 26 91 27
13 3 3 17 12 16 14 10 7 50 21 4 8 12 16 15 25 14 18 3 21 9 41 13 Gemidd. 365
0,46 0,43 0,43 0,4 0,4 0,37 0,37 0,38 0,41 0,45 0,37 0,4 0,4 0,36 0,43 0,39 0,43 0,41 0,43 0,38 0,38 0,35 0,45 0,48 0,41 0,59
619
Om Nb te benutten als monitoringstool, zouden we om praktische redenen het aantal stalen moeten kunnen reduceren. Met de verwantschapsmethode zou dit echter leiden tot onderschatte waarden enerzijds, met de LD methode tot schattingen nog steeds beïnvloed door mixture-LD anderzijds. De optie om volwassen dieren van verschillende leeftijden te bemonsteren voor een schatting van generationele effectieve populatiegrootte moet nog onderzocht worden voor deze soort. Niettemin zal je ook daar mogelijk een effect ondervinden van mixture-LD op de schatting. Eventueel kunnen temporele methodes additionele schattingen van Ne geven. Er werden geen modellen bekomen met significante covariaten die het aantal gedetecteerde eDNA-genkopieën/µl konden verklaren. Er is met andere woorden geen duidelijke link tussen het aantal genkopieën en het aantal reproducerende adulten. Dit aantal werd immers afgeleid van het aantal eilegsels van heikikker en daarmee kon ook geen significant verband gevonden ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 33 van 40
worden (zie 3.2.3 en Figuur 15), zelfs niet met inbegrip van de poeloppervlakte in het model. Deze laatste is mogelijk niet nauwkeurig genoeg, want vele poelen waren moeilijk te identificeren op luchtfoto en bleken op terrein soms in verbinding te staan met andere wateren bij een stijgend waterpeil. Bovendien leggen heikikkers vaak bij voorkeur eitjes in (vanuit de lucht slecht zichtbare) ondiepe, grasrijke zones af. Met Nb (LD schatting) als responsvariabele leek zuurstofgehalte van belang (AIC= 25,14; Tabel 3, Figuur 20). Ook de oppervlakte van de poel zat nog in het model, maar was net niet significant. De outlier in Figuur 20 wordt veroorzaakt door het lage aantal eilegsels van heikikker (3). Voor de Nb schattingen met de verwantschapsmethode kwam het model met het zuurstofgehalte als verklarende variabele als beste uit de bus, maar dan met een net nietsignificante coĂŤfficiĂŤnt (P = 0,0833). Het aantal eilegsels van heikikker per poel is natuurlijk sterk gecorreleerd met de schattingen van Nb. Het is dus mogelijk dat zuurstofrijkere poelen worden verkozen voor ovipositie door heikikkervrouwtjes (Pintar & Jr, 2017).
Tabel 3 Modelsamenvatting van het glm met log10(Nb) als responsvariabele.
Waarde Intercept
S.E.
t-waarde
P
-2,438684 1,5393259 -1,584254 0,1372
log10(O2) 3,595641
1,4940956 2,406567
0,0317
Opp. poel 0,000984
0,0004682 2,100860
0,0557
Figuur 20 Visualisatie van het GLM model met het briggse logaritme van de LD schattingen van het effectief aantal ouders als respons. De voorspelde waarden worden weergegeven door de rechte. Het logaritme van de zuurstofconcentratie wordt weergegeven op de x-as en de oppervlakte van de poelen in kleurschakering van donker- naar lichtblauw (van klein naar groot). ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 34 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
4 CONCLUSIES Uit het uitgevoerde onderzoek blijkt dat we via eDNA-analysen de aanwezigheid van zowel heikikker als bruine kikker in een waterpartij eenduidig kunnen aantonen. Het temporele eDNA-patroon vertoont hierbij een sterke piek vroeg in het voorjaar (eind maart – begin april), die voor beide soorten sterk overlappend is. Dit is dan ook de beste periode om in ongekende systemen via eDNA na te gaan of de soorten er al dan niet aanwezig zouden zijn. Om een idee te bekomen van de verhouding heikikker versus bruine kikker eilegsels op locaties waar er zich natuurlijke afzet van eitjes voordoet, noodzakelijk voor een adequate inschatting van de populatiegroottes van heikikker, blijkt dat een indirecte benadering via eDNA uit het omliggende water voor deze doelstelling een omslachtige en weinig accurate methode is. Een veel efficiëntere manier bestaat erin om van elk eilegsel 1 eitje te bemonsteren deze te mengen tot één staal per waterpartij, en op het DNA extract (afkomstig van dit mengstaal) een (ddPCR) meting te doen van de concentratie DNA-moleculen van heikikker versus bruine kikker. Dit onderzoek heeft namelijk aangetoond dat de proportionele verhouding van deze concentraties een nagenoeg exacte weerspiegeling is van het aandeel heikikker- versus bruine kikkereitjes. Ook de individuele soortidentificatie van elk eitje afzonderlijk via de RFLPmethode is erg nauwkeurig, maar uiteraard veel arbeidsintensiever en ook qua analysekost duidelijk duurder aangezien alle eitjes in het veld afzonderlijk dienen opgeslagen te worden en verder afzonderlijk dienen te worden geanalyseerd in het labo. Voor zowel de RFLP-methode per eitje als de ddPCR-methode op een mengstaal van eitjes moet een representatieve steekproef van de legsels genomen worden om een degelijke uitspraak te kunnen maken over de verhouding tussen de soorten binnen een poel. De genetische merkers zijn hier samengesteld om verder op een routineuze wijze te gebruiken voor de genetische monitoring van Vlaamse heikikkerpopulaties. Ze blijken variabel genoeg om ook een enkele populatie genetisch te karakteriseren. Klein Schietveld bevat zoals verwacht een gezonde heikikkerpopulatie met een relatief hoge genetische variatie en een groot aantal effectieve ouders voor het jaar 2017. Ook al kan dat aantal in de tijd variëren, het lijkt groot genoeg om jaren van mindere reproductie op te vangen. Hoe het zit met de overleving van de larven werd hier niet onderzocht. De populatie in het Klein Schietveld lijkt op basis van waarnemingen een eerder geïsoleerde populatie, maar ook dat kan enkel besloten worden na onderzoek van omringende populaties. Dezelfde merkers kunnen benut worden om connectiviteit te onderzoeken. Voor genetische monitoring van heikikkerpopulaties zijn de maten voor genetische variatie hier gebruikt, zoals verwachte heterozygositeit, zeker nuttig en bruikbaar ook wanneer de staalname minder dekkend zou gebeuren. Uit onze analyses blijkt immers dat die parameters niet sterk variëren per poel binnen een populatie waar voornamelijk IBD de spatiale structuur bepaalt (ofwel geografische afstand) en voor het overige de poelen gelegen zijn in goed verbonden habitat. Dus ook bij bemonstering van een selectie van poelen zouden we dezelfde waarden verkregen hebben. Voor het effectief aantal ouders lijkt dit anders te liggen en zal (samen met het gebruik van generationele effectieve populatiegrootte) verder onderzocht moeten worden.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 35 van 40
5 REFERENTIES Ackerman M.W., Hand B.K., Waples R.K., Luikart G., Waples R.S., Steele C.A., Garner B.A., McCane J., Campbell M.R. (2017). Effective number of breeders from sibship reconstruction: empirical evaluations using hatchery steelhead. Evol Appl 10(2):146160. Arens P., van der Sluis T., van’t Westende W.P.C., Vosman B., Vos C.C., Smulders M.J.M. (2007). Genetic population differentiation and connectivity among fragmented Moor frog (Rana arvalis) populations in The Netherlands. Landscape Ecol 22(10):1489-1500. Berlin S., Merilä J., Ellegren H. (2000). Isolation and characterization of polymorphic microsatellite loci in the common frog, Rana temporaria. Mol Ecol 9(11):1938-1939. De Bruyn L., Speybroeck J., Maes D., De Knijf G., Onkelinx T., Piesschaert F., Pollet M., Van Calster H., Westra T., Quataert P. (2015). Monitoringsprotocol Amphibia: kikkers en padden. INBO.R.2015.11336466. Dempster A.P., Laird N.M., Rubin D.B. (1977). Maximum likelihood from incomplete data via the EM algorithm. J Roy Stat Soc B Met:1-38. Do C., Waples R.S., Peel D., Macbeth G.M., Tillett B.J., Ovenden J.R. (2014). NeEstimator v2: reimplementation of software for the estimation of contemporary effective population size (Ne) from genetic data. Mol Ecol Resour 14(1):209-214. Garner T.W.J., Tomio G. (2001). Microsatellites for use in studies of the Italian Agile Frog, Rana latastei (Boulenger). Conserv Genet 2(1):77-80. Goudet J. (2005). hierfstat, a package for r to compute and test hierarchical F-statistics. Mol Ecol Notes 5(1):184-186. Hardy O.J., Vekemans X. (2002). SPAGeDi: a versatile computer program to analyse spatial genetic structure at the individual or population levels. Mol Ecol Notes 2(4):618-620. Ishchenko V.G. Growth of brown frogs of fauna of Russia: some problems of study of growth in amphibians. In: N. A., Tsinenko O., (editors); (2005) August 12-16 2003; St. Petersburg, Russia. Russ J Herpetol. p 153-157. Johnson P.C.D., Haydon D.T. (2007a). Maximum-likelihood estimation of allelic dropout and false allele error rates from microsatellite genotypes in the absence of reference data. Genetics 175(2):827-842. Johnson P.C.D., Haydon D.T. (2007b). Software for quantifying and simulating microsatellite genotyping error. Bioinform Biol Insights 1:71-75. Jombart T. (2008). adegenet: a R package for the multivariate analysis of genetic markers. Bioinformatics 24(11):1403-1405. Jones O.R., Wang J. (2010). COLONY: a program for parentage and sibship inference from multilocus genotype data. Mol Ecol Resour 10(3):551-555. Kalinowski S., Taper M. (2006). Maximum likelihood estimation of the frequency of null alleles at microsatellite loci. Conserv Genet 7(6):991-995. Keenan K., McGinnity P., Cross T.F., Crozier W.W., Prodöhl P.A. (2013). diveRsity: an R package for the estimation and exploration of population genetics parameters and their associated errors. Methods Ecol Evol 4(8):782-788. Kovar R., Brabec M., Bocek R., Vita R. (2009). Spring migration distances of some Central European amphibian species. Amphibia-Reptilia 30(3):367. Lardner B., Loman J. (2009). Does landscape and habitat limit the frogs Rana arvalis and Rana temporaria in agricultural landscapes? A field experiment. Appl Herpetol 6(3):227. Loman, J. (2002). Temperature, genetic and hydroperiod effects on metamorphosis of brown frogs Rana arvalis and R. temporaria in the field. Journal of Zoology 258: 115-129. ///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 36 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Lynch M., Ritland K. (1999). Estimation of pairwise relatedness with molecular markers. Genetics 152(4):1753-1766. Moran P.A.P. (1948). The interpretation of statistical maps. J Roy Stat Soc B Met 10(2):243-251. Moran P.A.P. (1950). Notes on continuous stochastic phenomena. Biometrika 37(1/2):17-23. Naimi B., Hamm N.A.S., Groen T.A., Skidmore A.K., Toxopeus A.G. (2014). Where is positional uncertainty a problem for species distribution modelling? Ecography 37(2):191-203. Neel M.C., McKelvey K., Ryman N., Lloyd M.W., Short Bull R., Allendorf F.W., Schwartz M.K., Waples R.S. (2013). Estimation of effective population size in continuously distributed populations: there goes the neighborhood. Heredity 111(3):189-199. Palo J.U., Merilä J. (2003). A simple RFLP method for identification of two ranid frogs. Conserv Genet 4(6):801-803. Pidancier N., Gauthier P., Miquel C., Pompanon F. (2002). Polymorphic microsatellite DNA loci identified in the common frog (Rana temporaria, Amphibia, Ranidae). Mol Ecol Notes 2(3):304-305. Pinheiro J., Bates D., DebRoy S., Sarkar D., Team R.C. (2016). Linear and nonlinear mixedeffects models. R package version 3.1-127 ed. Pintar M.R., Jr W.J.R. (2017). Out with the old, in with the new: oviposition preference matches larval success in Cope's gray treefrog, Hyla chrysoscelis. J Herpetol 51(2):186-189. R Core Team. (2018). R: A language and environment for statistical computing. Vienna, Austria: R Foundation for Statistical Computing. Rousset F. (2008). genepop’007: a complete re-implementation of the genepop software for Windows and Linux. Mol Ecol Resour 8(1):103-106. Rowe G., Beebee T.J.C. (2001). Polymerase chain reaction primers for microsatellite loci in the common frog Rana temporaria. Mol Ecol Notes 1(1‐2):6-7. Söderman F. (2006). Comparative population ecology in moor frogs with particular reference to acidity. Uppsala: Acta Universitatis Upsaliensis. Spens J., Evans A.R., Halfmaerten D., Knudsen S.W., Sengupta M.E., Mak S.S.T., Sigsgaard E.E., Hellström M. (2017) Comparison of capture and storage methods for aqueous macrobial eDNA using an optimized extraction protocol: advantage of enclosed filter. Methods in Ecology and Evolution, 8: 635-645. Vos C.C., Jong A.D., Goedhart P., Smulders M. (2001). Genetic similarity as a measure for connectivity between fragmented populations of the moor frog (Rana arvalis). Heredity 86(5):598-608. Wang J. (2011). COANCESTRY: a program for simulating, estimating and analysing relatedness and inbreeding coefficients. Mol Ecol Resour 11(1):141-145. Waples R.S., Anderson E.C. (2017). Purging putative siblings from population genetic data sets: a cautionary view. Mol Ecol 26(5):1211-1224. Waples R.S., Antao T., Luikart G. (2014). Effects of overlapping generations on linkage disequilibrium estimates of effective population size. Genetics 197(2):769-780. Waples R.S., Do C. (2008). LDNE: a program for estimating effective population size from data on linkage disequilibrium. Mol Ecol Resour 8(4):753-756. Waples R.S., England P.R. (2011). Estimating contemporary effective population size on the basis of linkage disequilibrium in the face of migration. Genetics 189(2):633-644. Weir B.S., Cockerham C.C. (1984). Estimating F-statistics for the analysis of population structure. Evolution 38(6):1358-1370.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 37 van 40
6 BIJLAGEN Bijlage 1 Schematisch overzicht van de positie van de primers en probe voor amplificatie van heikikker (Rana arvalis) DNA op het COI fragment, met de respectievelijke forward en reverse primer (donkergroen), en de locatie van de regio waarop de probe bindt (rood), t.o.v. andere amfibiesoorten die in Vlaanderen kunnen worden aangetroffen.
Forward primer
Probe
Reverse primer
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 38 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be
Bijlage 2 Schematisch overzicht van de positie van de primers en probe voor amplificaite van bruine kikker (Rana temporaria) DNA op het COI-fragment, met de respectievelijke forward en reverse primer (donkergroen), en de locatie van de regio waarop de probe bindt (rood), t.o.v. een heel aantal andere amfibiesoorten die in Vlaanderen kunnen worden aangetroffen.
Forward primer
Reverse primer
Probe
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// www.inbo.be
doi.org/10.21436/inbor.17686189
Pagina 39 van 40
Bijlage 3 Informatie over de gebruikte microsatellietmerkers. Naam merker
Multiplexnr.
Fluorescent label
Volume (µl)
Originele bron/ GenBank acc. nr.
RCIDII 1 PET 0,2 (Vos et al., 2001) RNTYR2 1 PET 0,6 (Vos et al., 2001) RA14 1 NED 0,2 EU871714 RtsB14 1 VIC 0,2 (Berlin et al., 2000) RlatCa41 1 FAM 0,2 (Garner & Tomio, 2001) RECALQ 1 PET 0,2 (Vos et al., 2001) RC08604 1 FAM 0,4 (Vos et al., 2001) RRD590 2 NED 0,2 (Vos et al., 2001) RlatCa18 2 VIC 0,2 (Garner & Tomio, 2001) Rtempµ9 2 FAM 0,05 (Rowe & Beebee, 2001) RA11 2 PET 0,2 EU871712 Rt2Ca2-22 2 FAM 0,2 (T. Garner, niet gepubliceerd) RtµP 2 NED 0,2 (Pidancier et al., 2002) Rtempµ4 2 NED 0,05 (Rowe & Beebee, 2001) Rtempµ5 2 PET 0,4 (Rowe & Beebee, 2001) RA03 2 FAM 0,4 EU871710 RA04 2 VIC 0,2 EU871711 RA13 3 (simplex) VIC 0,4 EU871713 RtU4 4 (simplex) FAM 0,4 (Berlin et al., 2000) N.B.: de primers werden verdund volgens volgend recept: 10 µl forward primer + 2,5 µl reverse primer + 37,5 µl ultra-puur water.
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////// Pagina 40 van 40
doi.org/10.21436/inbor.17686189
www.inbo.be