De rode panda en Cserhati (6): Whole Genome K-mer Signature

Cserhati geeft in het BMC Genomics artikel aan dat hij de indeling van de de rode panda op grond van morfologie aan twijfel onderhevig vindt. Aangezien hij het meer over DNA sequenties dan over morfologie gehad heeft, zal hij bedoelen dat hij (ook) de indeling op grond van DNA aan twijfel onderhevig vindt. Nu zijn alle indelingen op grond van DNA die hij genoemd heeft op grond van vrij weinig DNA (voor 2021 begrippen). Flynn et al (2000) gebruikten bij voorbeeld de DNA sequentie van vier genen. Dat was het jaar 2000 natuurlijk, en toen waren er niet meer DNA sequenties beschikbaar.

Cserhati geeft de voorkeur aan het gebruiken van het gehele genoom voor de indeling van de rode panda: een whole genome study. Dat is ook steeds meer wat er gebeurt. Bijvoorbeeld, De Ferran et (2022) zochten niet naar orthologe genen in hun elf soorten otters, maar gebruikten genoom fragmenten als gevonden bij genoom sequencing voor het vergelijken van het DNA.

Cserhati’s voorkeur voor gebruik van het gehele genoom is dus goed te begrijpen. De methode die Cserhati gebruikt om het hele genoom te karakterisen is Whole Genome K-mer Signature, afgekort als WGKS.

Dan zijn er twee vragen: wat is WGKS? En hoe bruikbaar is WGKS?

Eerst: wat is WGKS? Dat komt in deze post aan de orde. Hoe bruikbaar WGKS is voor klassificatie van soorten komt in de volgende aflevering.

In het methoden deel van het BMC Genomics artikel schrijft Cserhati:

The WGKS algorithm that was used in the analysis is an alignment-free k-mer sequence comparison method. These methods involve the statistical comparison of k-mers between species.

A k-mer is a segment of DNA k bp long,

The k-mer signature is simply a list of all k-mers ordered in lexicographical order from AA … A to TT … T, together with their score values. For a given value k, there are 4k possible k-mers. Thus, the k-mer signature also corresponds to a vector of 4k numbers. Since octamers were analyzed, this corresponds to 65,536 possible octamers.

Cserhati zegt: k-mers van acht basen lang - octamers - turven, en dan hun scores. Ik doe dat in twee delen: octamers tellen en de scores.


1 Octamers tellen en de correlatie van de aantallen octamers vinden

DNA heeft vier basen: ACGT. Een stukje DNA van 8 basepaar lang, een octamer, kan dus alle mogelijkheden hebben van AAAAAAAA tot TTTTTTTT. Vier mogelijkheden voor de eerste plaats, vier mogelijkheden voor de tweede plaats, enzovoort. Dat is 48 = 65536 mogelijkheden. Een computer loopt langs het genoom, en leest op volgorde af welke 8 baseparen hij vindt.

In een DNA sequentie;

      gagtgggcagcactccaaataccgttaagctggagcctcggt

krijg je dus voor octamers:

vanaf base 1:                 gagtgggc

vanaf base 2:       agtgggca

vanaf base 3:       gtgggcag

en zo verder. De computer turft het aantal keren dat een volgorde van 8 basen voorkomt. In dit kleine voorbeeld van een DNA sequentie komt elke volgorde in een stukje van 8 één maal voor. Dat heet dan de k-mer signature, met 8 basen heet het octamer signature.

Bij DNA van twee verwante soorten zou je kunnen verwachten dat je om en nabij dezelfde verdeling van volgorden van 8 basen, octamers, vindt.

Alleen verschillen soorten niet alleen in belangrijk DNA maar ook in onbelangrijk DNA. Bijvoorbeeld in de lengte van een herhaling van ac - acacacac  of acacacacacacacacacaca, of een verschil in het aantal LINE1 elementen. Een huismuis bv heeft honderden LINE1 elementen (waarvan maar een paar actief als transposon (Jachowitz et al 2017)), en een andere soort muis zou duizenden LINE1 elementen kunnen hebben. Zo'n herhaling van een zelfde volgorde als in een veelvoud aan LINE1 elementen of een veelvoud aan eenvoudige repeats heeft grote invloed op de verdeling van octamers in een WGKS.


Voorbeeld bij de octamer telling:

Hier komt een lang voorbeeld, en ik heb het zo geconstrueerd dat het mogelijk moet zijn de invloed van repetitief DNA te zien.

Neem de volgende vijf sequenties. In de eerste twee sequenties is een lange of korte ac herhaling ingevoegd ten opzichte van de derde sequentie. De vierde sequentie verschilt van de eerste sequentie in het deel na de lange ac herhaling. De vijfde sequentie is na de kortere ac herhaling hetzelfde als de vierde sequentie.

>seq_1

gagtgggcagcaacacacacacacacacacacacacacacacacctccaaataccgttaagctggagcctcggt

>seq_2

gagtgggcagcaacacacacacacacacctccaaataccgttaagctggagcctcggt

>seq_3

gagtgggcagcactccaaataccgttaagctggagcctcggt

>seq_4

gagtgggcagcaacacacacacacacacacacacacacacacacactcttctggtccccacagactcagagaga

>seq_5

gagtgggcagcaacacacacacacacacactcttctggtccccacagactcagagaga

 Als ik wat punten toevoeg om de ac herhaling beter uit te laten komen:

>seq_1

gagtgggcagca...acacacacacacacacacacacacacacacac....ctccaaataccgttaagctggagcctcggt

>seq_2

gagtgggcagca...acacacacacacacac....................ctccaaataccgttaagctggagcctcggt

>seq_3

gagtgggcagca.......................................ctccaaataccgttaagctggagcctcggt

>seq_4

gagtgggcagca...acacacacacacacacacacacacacacacacac...tcttctggtccccacagactcagagaga

>seq_5

gagtgggcagca...acacacacacacacacac...................tcttctggtccccacagactcagagaga

Een heel eenvoudig programma, CLUSTAL, geeft een alignment van de sequenties, en twee fylogenetische bomen:

CLUSTAL O(1.2.4) multiple sequence alignment

 

 

seq_4      gagtgggcagcaacacacacacacacacacacacacacacacacactcttctggtccc--        58

seq_5      Gagtgggcagcaacacacaca----------------cacacacactcttctggtccc--        42

seq_1      gagtgggcagcaacacacacacacacacacacacacacacacacctcc---aaataccgt        57

seq_2      gagtgggcagcaacacacac----------------acacacacctcc---aaataccgt        41

seq_3      gagtgggcagca--------------------------------ctcc---aaataccgt        25

           ************                                   *      * ** 

 

seq_4      -----cacagactcagagaga   74

seq_5      -----cacagactcagagaga   58

seq_1      taagctggagcctcggt----   74

seq_2      taagctggagcctcggt----   58

seq_3      taagctggagcctcggt----   42

                   ** *** *    

Figuur 1. Grafische weergave als Phylogram van de overeenkomst tussen de vijf sequenties


Figuur 2.  Grafische weergave als Phylogenetic tree van de overeenkomst tussen de vijf sequenties

In het fylogram worden de eerste drie sequenties bij elkaar geplaatst, in de fylogenetische boom op grond van Neighbour-joining krijgen de sequenties 4 en 5 gezelschap van de eerste sequentie. In dat geval stuurt de herhaling sequentie 1 naar sequenties 4 en 5. 

Hoe reageert WGKS op herhalingen?

Van dezelfde vijf sequenties is een octamer signature opgesteld. Er zijn 83 verschillende octamers gevonden over de vijf sequenties. Daarvan komen 81 nul of één keer voor in een sequentie. De overige twee zijn acacacac en cacacaca. De octamer signatures beginnen:

string

seq 1

seq 2

seq 3

seq 4

seq 5

aaataccg

1

1

1

0

0

aacacaca

1

1

0

1

1

aagctgga

1

1

1

0

0

aataccgt

1

1

1

0

0

acacacac

13

5

0

14

6

acacacct

1

1

0

0

0

acacactc

0

0

0

1

1

acacctcc

1

1

0

0

0

acactctt

0

0

0

1

1

acagactc

0

0

0

1

1

 

De correlatiematrix laat zien dat sequentie 3 behoorlijk afwijkt van de andere vier sequenties:

        

 

seq 1

seq 2

seq 3

seq 4

seq 5

seq 1

1

0.905

0.020338

0.915432

0.758678

seq 2

0.905

1

0.262668

0.693422

0.511118

seq 3

0.020338

0.262668

1

-0.27067

-0.43957

seq 4

0.915432

0.693422

-0.27067

1

0.935843

seq 5

0.758678

0.511118

-0.43957

0.935843

1

Als het aantal repeats in sequenties 1 en 4, en in sequenties 2 en 5, omhoog gegooid wordt krijg je een dergelijke correlatiematix:

 

seq 1

seq 2

seq 3

seq 4

seq 5

seq 1

1

0.992977

-0.07366

0.985514

0.964082

seq 2

0.992977

1

-0.00813

0.963535

0.936946

seq 3

-0.07366

-0.00813

1

-0.1939

-0.25839

seq 4

0.985514

0.963535

-0.1939

1

0.992437

seq 5

0.964082

0.936946

-0.25839

0.992437

1

Bij octamer tellingen en hun correlaties overweegt de invloed van de herhalingen op de overeenkomst in de sequenties. De hoeveelheid informatief DNA is te laag ten opzichte van de sequenties om nog enig verschil te maken.


2 Octamer scores berekenen

Cserhati gebruikt niet de octamer tellingen, maar een 'octamer score'.

Nu komt eerst wat de score is, dan hoe je de score berekent, dan wat de gevoeligheid van de score is.

De octamer score is volgens het Python Script motief programma :

Score Sc = (O-E)/(O+E)

waarin het waargenomen aantal octamers gelijk is aan O en verwachte aantal octamers gelijk is aan E. Voor waargenomen aantal gelijk aan verwacht aantal is de score Sc gelijk aan nul. Als het waargenomen aantal nul is, O=0, is de score gelijk aan -1.

Ik geef nu twee benaderingen om naar de gevoeligheid van de score voor afwijkingen van verwacht te kijken: de relatieve afwijking van verwacht en de absolute afwijking van verwacht.


i) Ten eerste kijken we naar de relatieve afwijking van verwacht, waarbij de genoomgrootte geen rol speelt.

Als O=xE, met x ≥ 0, wordt Sc = (xE-E)/(xE+E) = (x-1)/(x+1) .

Voor x=0 is Sc=-1; voor x=1 is Sc =0; voor x nadert tot oneindig is de limiet van Sc = 1. De score Sc is sterk asymmetrisch.

De score Sc als functie van x ziet er als volgt uit:

Lineaire x-as:

Figuur 3 Score Sc = (x-1)/(x+1) als functie van x

Logarithmische x-as

Figuur 4 Score Sc = (x-1)/(x+1) als functie van x

De score lijkt nogal gevoelig voor afwijkingen tussen waargenomen aantal O en verwacht aantal E voor vrij kleine afwijkingen van verwacht, en ongevoelig voor grote afwijkingen van verwacht.



ii) Ten tweede kijken we naar de absolute afwijking van verwacht, waarbij de genoomgrootte een rol speelt.

Er zijn 48 = 65536 verschillende octamers. Stel de genoomgrootte op N basepaar. Bij eerste ruwe benadering is het verwachte aantal van een octamer E=Nx4-8. Noem het waargenomen aantal van een octamer O=E+n

De score wordt nu Sc= (E+n-E)/(E+n+E) = n/(2E+n). De totale genoomgrootte en de absolute afwijking spelen nu een rol in de score.

In de figuur staan scores voor genoomgroottes van N=106, N=107, N=108, N=109, als verschillende lijnen; en verschil tussen waargenomen en verwacht van 101 tot 105 op de x-as.

Bij een groot genoom, in de orde van 109 bp, moet de absolute afwijking tussen waargenomen en verwacht groot zijn wil de score noemenswaard veranderen. Bij een kleiner genoom leiden vrijwel alle veranderingen tot hoge scores.

Figuur 5. De score als functie van verschil tussen waargenomen en verwacht voor verschillende genoomgroottes

Kleine duplicaties hebben daarmee nauwelijks effect op de score, terwijl transposonfrequenties grote effecten hebben. Verschillen in frequentie van transposons als LINE1 met 15-20% van het genoom gaan een grote invloed op het genoom hebben, omdat ze in zulke grote aantallen voorkomen.

Cserhati zegt:

Even if the genome is partially or completely duplicated, then the score value will not change. This is because both the Observed and Expected values will increase by the proportion that the duplicated genome is compared to the pre-duplication genome

Dit klopt niet. Als een deel van het genoom gedupliceerd is, zal het verwachte aantal E van een octamer minder afhangen van de aanwezigheid van duplicaties dan het waargenomen aantal O. De invloed van het niet-gedupliceerde deel van het genoom zal in E overwegen.



3 Moraal van het voorbeeld

De octamer score Sc = (O-E)/(O+E) is dermate niet-lineair dat het te betwijfelen valt of dit een werkbare maat voor enige genoom eigenschap is.

In een correlatiematrix op grond van WGKS gaan verschillen in grote hoeveelheden repetitief DNA een grote invloed hebben.

Als we werken met octamer signatures van verwante soorten, hebben we een vrij overeenkomstig octamer patroon afkomstig van informatief DNA, en heeft repetitief DNA tegen die achtergrond een grote invloed op de correlaties tussen de signatures van de soorten.

Als we werken met octamer signatures van soorten die ver van elkaar staan in hun fylogenie, verwachten we dat er veel verschil in octamer patroon afkomstig is van informatief DNA. Tegen die achtergrond van veel verschil door informatief DNA kan de invloed van repetitief DNA verminderen, verdwijnen, of bij gebruik van een groot aantal soorten van een groep, uitmiddelen.

De verwachting is daarom dat WGKS octamer patronen redelijk de hoofdlijnen van de fylogenie zullen volgen, maar bij vergelijkingen op kleinere schaal - tussen soorten van een familie bijvoorbeeld, of bij superfamilies - een verkeerd beeld zullen geven.


***

Cserhati, M., 2021, A tail of two pandas – whole genome k-mer signature analysis of the red panda (Ailurus fulgens) and the Giant panda (Ailuropoda melanoleuca), BMC Genomics 22: 228

https://bmcgenomics.biomedcentral.com/articles/10.1186/s12864-021-07531-3

Jachowicz, J.W., Bing, X., Pontabry, J., Bošković,A., Rando, O.J., & Torres-Padilla, M-J (2017) LINE-1 activation after fertilization regulates global chromatin accessibility in the early mouse embryo. Nature Genetics 49: 1502

de Ferran, V., Figueiró, H., de Jesus Trindade, F, en 17 anderen , & Eizirik, E. (2022) Phylogenomics of the world’s otters. Current Biology 32; 3650–3658,

The Python script motif_analysis_k-1.py at github.com/csmatyi/motif_analysis was used to generate WGKS profiles

https://github.com/csmatyi/motif_analysis


Reacties