een voorbeeld uitlijning met gemarkeerde MEM

indien er zowel mismatches als hiaten zijn tussen Mi en Mj, worden alle hiaten beschouwd als continu om de Gap open penalty te verminderen (Er wordt slechts één gap open penalty toegepast voor een continues gap). Dus, voor alle aangrenzende MEMs die hiaten tussen hen hebben, slechts één gat open boete wordt toegepast., De plaatsing van mismatches en de enige continue kloof is niet belangrijk, omdat het de uitlijningsscore niet zou beïnvloeden. We gaan ervan uit dat de mismatch straf constant is (dit is gebruikelijk voor DNA-sequenties).

MEM extractie

Er zijn methoden om maximale exacte overeenkomsten te extraheren tussen lange sequenties zoals een volledig genoom. Nochtans, zijn deze methodes gebaseerd op voorverwerking en indexering van één of beide opeenvolgingen die een tijdrovende verrichting is. Bijvoorbeeld, in gelezen aligner van DNA, wordt het referentie-genoom eens geïndexeerd, en de zelfde index wordt gebruikt elke keer dat een nieuwe gelezen wordt uitgelijnd., We zijn op zoek naar een snel algoritme om MEMs te identificeren tussen relatief korte sequenties die veranderen voor elke uitlijning. Een brute force methode voor dit probleem (extra bestand 1: Sectie II) is traag en inefficiënt (met de complexiteit van O(n3)). We stellen een snelle bit-level parallelle methode voor om het mem extractie proces te versnellen. Onze mem-extractiemethode is gebaseerd op de in Fig. 3b. de eerste stap is om sequenties met bit-vectoren, waar A, C, T, en G zijn gecodeerd als 00, 01, 10, en 11, respectievelijk (aanvullend bestand 1: sectie III)., Figuur 6 illustreert een voorbeeld van een sequentiepaar, samen met overeenkomstige bit-vectorrepresentaties. In een goederencomputer, is het machinewoord gewoonlijk 64 beetjes die 32 nucleotidesymbolen kunnen aanpassen. Omdat een reeks meestal groter is dan 32 symbolen, wordt de bijbehorende bit-vector opgeslagen in meerdere machine woorden. Elke bewerking op bit-vectoren van sequenties van grootte n symbolen werkt op \(\lceil \ frac {n}{32} \ rceil\) machine woorden.

representatie van sequenties met bit-vectoren. XOR-uitvoer (X) met gemarkeerde MEMs. Edges bit-vector (E) identificeert het begin en het einde van elke MEMs

Met bit-vectoren representatie van sequenties, het verschuiven van een reeks met één symbool is hetzelfde als het verschuiven van de bit-vector met twee bits, en het vergelijken van sequenties kan worden gedaan met XOR instructie (32 symbolen tegelijk). In de XOR-uitvoer (X) betekent 00 dat symbolen overeenkomen, en een reeks van 00s toont een MEM., Een set van shift en bitwise operaties zoals getoond in algoritme 1 berekent X en vervolgens de rand bit-vector (E) waarin het begin en het einde van elke MEM worden gemarkeerd met set bits (bits met een waarde van een). Figuur 6 toont de x en de e bit-vectoren met gemarkeerde MEMs. De positionering van MEMs in sequenties wordt dan berekend vanuit de rand bit-vector (extra bestand 1: sectie IV).,

Alignment algorithm

In de sectie” Approach ” laten we zien dat door verschillende combinaties van MEMs te overwegen en de uitlijningsscore voor de corresponderende uitlijningsscore te berekenen, men de combinatie van MEMs kan identificeren die resulteert in de maximale uitlijningsscore. Het onderzoeken van alle mogelijke combinatie van MEMs is echter een naïeve oplossing. Een meer systematische manier om de alignment efficiënt te vinden is het gebruik van dynamische programmering.

dynamisch programmeren is de methode om de oplossing van een probleem te benaderen door kleinere subproblemen te definiëren en op te lossen., Oplossingen voor subproblemen worden gebruikt om bij elke stap een groter probleem op te lossen. Het proces wordt herhaald totdat alle subproblemen zijn opgelost. Uiteindelijk zou de oplossing voor een van de subproblemen de oplossing zijn voor het oorspronkelijke probleem. Wanneer alle subproblemen zijn opgelost, identificeert een backtracking-proces een reeks oplossingen die bijdragen aan de uiteindelijke oplossing. In dynamische programmering, zou er een rangschikking van de inputgegevens moeten zijn waarlangs de recursieprocedure verloopt.

we Sorteren alle MEMs op basis van de positie van hun einde in query sequence (EQ)., MEMs die op dezelfde positie eindigen, worden willekeurig geordend. Het j-de subprobleem is het vinden van de uitlijning van de onderdelen T en Q die eindigen op j-de MEM Mj (respectievelijk T en Q). We zullen laten zien dat deze volgorde van MEM voldoende is om de juiste recursie te ondersteunen.

in de gesorteerde lijst van MEMs geeft EQi=EQj aan dat een van Mi of Mj de andere MEM volledig overlapt in de query-reeks. Aangezien in Φ2B van algoritme 2 het overlappingsgebied is uitgesloten, kunnen Mi en Mj niet in dezelfde uitlijning zijn., Zo worden ith en jth subproblemen onafhankelijk van elkaar opgelost en kan de volgorde van i en j in de gesorteerde lijst willekeurig zijn. Als EQk>EQj (k>j in de gesorteerde lijst), kan Mk geen deel uitmaken van de uitlijning die eindigt op Mj. Zo kunnen de JTH subproblemen onafhankelijk van de oplossing voor het kth subprobleem worden opgelost. Merk op dat het ook mogelijk is om MEMs te sorteren op basis van hun eindpositie in de doelreeks (ET) met behulp van een soortgelijke rechtvaardiging.

ons voorgestelde dynamische programmeeralgoritme (DP-MEM) is uitgewerkt in algoritme 2., Voor het voorbeeld MEMs geëxtraheerd in Fig. 3b, de dynamische programmeertabel en de tussenwaarde berekend in het algoritme worden weergegeven in Fig ‘ s. Respectievelijk 7 en 8. De input voor DP-MEM is de lijst van MEMs waar elke MEM (Mj) een triplet van gehele getallen is . De tweede invoer n is het aantal MEMs in de lijst. De output S is de uitlijningsscore voor de sequenties. Het algoritme drukt de indices af van alle MEMs die de uitlijning vormen, waarbij de eerste en de laatste afgedrukte getallen respectievelijk de indices zijn van de meest rechtse en de meest linkse MEMs in de uitlijning., Alle stappen van algoritme 2 worden als volgt becommentarieerd:

dynamische programmeertabel gebruikt in algoritme 2 om geëxtraheerde MEMs in Fig. 3b. cel i en j vertegenwoordigen de waarde van \(S_{i}^{j}\). Lege cellen worden niet geëvalueerd in Φ2. Evaluatie van cellen met kruisteken worden overgeslagen in Φ2A. initiële waarde van Sj wordt berekend in Φ1. De uiteindelijke waarde van Sj en de bron ervan (wat Sj maximaliseert) worden voor elke rij gemarkeerd. De hoogste Sj (S13) is de uitlijningsscore., M13 is de laatste MEM in de uitlijning en de Mem daarvoor is MW = M3. Aangezien W = -1,M3 is de eerste MEM in de uitlijning. Het scoresysteem voor deze uitlijning is Rm=2,Px=3,Po=4 en Pe=1

tussenwaarden om \(s_{i}^{j}\) in Fig. 7. Merk op dat Sij in deze figuur verwijst naar \(S_{i}^{j}\)

• Φ1: elke MEM scoort voor alle bijbehorende symbolen., Merk op dat er LJ overeenkomende symbolen in Mj. Sj vertegenwoordigt de hoogste uitlijningsscore voor de uitlijning die eindigt op Mj. Het initialiseren van Sj in deze stap is vergelijkbaar met het berekenen van de partiële uitlijningsscore wanneer alleen Mj is opgenomen in de uitlijning. W wordt gebruikt voor backtracking. De waarde -1 geeft aan dat de huidige Sj wordt verkregen door alleen MJ in de uitlijning te beschouwen.

• Φ2: berekening van Sj voor elke MEM (Mj)., Om SJ te berekenen, voor elke MEM Mi waar Mi verschijnt voor Mj in de lijst, voegt het algoritme MJ toe aan de uitlijning eindigend op Mi (uitbreiding van eerder gevonden uitlijningen) en zoekt naar de extensie die SJ maximaliseert (het oplossen van een groter subprobleem met eerder opgeloste subproblemen).

• Φ2A: extensie overslaan als dit niet mogelijk is. Als ETi>ETj dan bevat Mi een deel van de doelreeks die voorbij de uitlijning eindigt op Mj en de uitbreiding is niet mogelijk. Als EQi = EQj of ETi = ETj of SQi≥SQj of STi≥STj dan overlapt een van de MEMs de andere MEMS volledig., In dit geval kunnen Mi en Mj niet in een uitlijning samen zijn.

• Φ2B: berekening van de lengte van de overlapping tussen Mi en Mj. Als \({MO}_{i}^{j}\) kleiner of gelijk is aan nul, dan bestaat er geen overlap.

• Φ2C: een kopie van Mj bewaren voordat overlapping wordt uitgesloten.

• Φ2D: als er overlap bestaat, exclusief overlap van Mj

• Φ2E: berekening van de eindpositie van Mj in T en Q.

• Φ2F: berekening van de afstand (aantal symbolen) tussen Mi en Mj in T en Q.

• Φ2G: berekening van het aantal mismatches en hiaten tussen Mi en Mj.,

• Φ2H: berekening van de boete voor de mismatches en hiaten tussen Mi en Mj (\(P_{i}^{j}\)). Als er een kloof bestaat, wordt slechts één open kloof boete afgetrokken.

• Φ2I: Computeruitlijningsscore \(\left (s_{i}^{j} \ right)\) wanneer Mj wordt toegevoegd aan de uitlijning die eindigt op Mi. De score voor alle overeenkomende symbolen in Mj (Lj×Rm) wordt toegevoegd aan de uitlijningsscore voor de uitlijning die eindigt op Mi (Si). Dan wordt de boete voor de gaten en mismatches tussen Mi en Mj\(\left (P_{i}^{j}\right)\) afgetrokken.,

• Φ2J: als de uitbreiding van Mj tot de uitlijning eindigend op Mi resulteert in een score \(\left (s_{i}^{j}\right)\) hoger dan de huidige score voor Mj (Sj) dan wordt de nieuwe score opgeslagen in Sj. Ook W is ingesteld op i om de Mi bij te houden die de score voor Mj maximaliseert.

• Φ2K: de waarde van Mj vóór uitsluiting herstellen zodat Mj in andere uitlijningsuitbreidingen kan worden gebruikt.

• Φ3: zoekt naar de MEM met de hoogste Sj. Dit MEM is de laatste MEM in de uitlijning (Me)., De hoogste score (Se) wordt geretourneerd als S, wat de hoogste uitlijningsscore is voor de gegeven sequenties. De index van de MEM die SJ maximaliseert wordt opgeslagen in e om terug te gaan van mij.

• Φ4: in de uitlijning is het directe voorgaande MEM voor mij degene die de uitlijningsscore voor mij maximaliseert. De index van dergelijke MEM wordt opgeslagen in W. als gevolg daarvan, de iteratie van f W W bezoekt de index van alle MEMs in de uitlijning. Als W gelijk is aan -1, is Mf de eerste MEM in de uitlijning en wordt de iteratie gestopt.,

in ons algoritme straffen we geen mismatches en hiaten voor de eerste MEM en na de laatste MEM in de uitlijning. Dit resulteert in een lokaal uitlijningsalgoritme. Door deze sancties te overwegen genereert het algoritme een globale uitlijning (aanvullend bestand 1: sectie V).

de vergelijking om \(p_{i}^{j}\) te berekenen in Φ2H van algoritme 2 gaat ervan uit dat er geen passend symbool is tussen T en Q in het gebied tussen Mi en Mj (alle symbolen worden geteld als mismatches of hiaten)., Hoewel deze aanname niet waar is voor alle Mi, is het altijd waar voor de Mi die leidt tot maximum \(S_{i}^{j}\) wat het effect van de aanname dat onjuist is voor andere Mi overstijgt. Als bewijs, neem aan dat er een passend symbool in het gebied tussen Mi en Mj. Het bijpassende symbool zou een MEM (Mk) zijn. Mk is al uitgebreid tot de uitlijning die eindigt op Mi. Bij uitbreiding van Mj naar MK wordt dus een hogere score behaald in vergelijking met uitbreiding van MJ naar Mi.

Chaining colinear seeds zoals besproken in zijn op grote schaal gebruikt in de uitlijning van grote sequenties, dat wil zeggen genome-to-genome uitlijning., Het is ook gebruikt om kandidaat-regio ‘ s te identificeren voor een read gegeven een reeks MEMs in BWA. Chaining algoritmes met scoring zijn vergelijkbaar met het dynamic programming algoritme dat we hebben voorgesteld (DP-MEM). Nochtans, zijn er verschillen die DP-MEM geschikt maken voor paarsgewijze uitlijning van korte opeenvolgingen. DP-MEM houdt er rekening mee dat alle MEMs binnen een bepaalde gapgrootte in de invoer worden verstrekt en optimaliseert het aantal iteratie in het algoritme. DP-MEM implementeert ook een heuristische benadering om het effect te compenseren van korte MEMs verwijderd uit de invoerlijst resulterende hiaten tussen MEMs.,

optimalisatie

gegeven sequenties van lengte n, vereist het algoritme om MEMs te extraheren (gegeven in de sectie “mem extraction”) 2(n−1) shift−en 2n-1-vergelijkingsbewerkingen op bit-vectoren (elke handeling op \(\lceil \frac {n}{32} \rceil \) machinewoorden) die resulteren in een algoritme met complexiteit van O(n2) om randbit-vectoren te produceren voor het gegeven paar sequenties. De complexiteit van de functie die positionering van MEMs uit de rand bit-vector berekent en ze sorteert op basis van EQ moet nog worden toegevoegd. Verder, als M MEMs worden geëxtraheerd, Φ2 van algoritme 2 (DP-MEM) heeft de complexiteit van O(m2)., Gezien de complexiteit van mem extractie en DP-MEM, MEM-Align is veel langzamer dan de beschikbare uitlijningsalgoritmen. Om het proces te versnellen, presenteren we verschillende optimalisaties voor MEM-Align die snellere runtime bereikt door nauwkeurigheid op te offeren. Aan de andere kant introduceren we methoden om de nauwkeurigheid te verbeteren met een minimaal prestatieverlies.

Banded alignment

Banded alignment is een bekende heuristische methode om het uitlijningsproces te versnellen. Deze techniek beperkt het patroon van de gaten in de uitlijning (aanvullend bestand 1: sectie VI)., Als de uitlijning tussen twee sequenties dit patroon niet volgt, zal het algoritme de uitlijning niet identificeren. In traditionele dynamische programmering, wordt de uitlijning verkregen na het berekenen van de waarde van alle cellen in de tabel. Echter, met de banded alignment optimalisatie, alleen cellen op de diameter en dicht bij diagonaal worden geëvalueerd. gl is de breedte van de band in gestreepte uitlijning waar cellen verder dan gl aan de diameter niet worden geëvalueerd. Donkere cellen in Fig. 1 toon het geval waarin gl=1.,

In tegenstelling tot de traditionele dynamische programmeerbenadering heeft MEM-Align geen vergelijkbare tabel om gestreepte uitlijning toe te passen. We ontdekten echter dat we hetzelfde effect kunnen simuleren door het aantal shift operaties in het mem extractie proces te beperken. Als we bijvoorbeeld de query-reeks naar gl naar rechts en naar links verschuiven, bereiken we banded alignment met de band van gl. Banded-alignment verminder de complexiteit van mem-extractie van O (n2) naar O (n. (2gl+1)) waarbij gl een kleine en vaste waarde is. De complexiteit van mem-extractie is dus O (n) wanneer banded alignment wordt toegepast., Ook, met de genoemde gestreepte uitlijning, is het waarschijnlijk dat minder MEMs worden geëxtraheerd die de daaropvolgende algoritmische stappen versnelt.

Short mem removal

We merkten dat de meerderheid van de geëxtraheerde MEMs kort zijn en het resultaat zijn van willekeurig overeenkomende symbolen. Om MEM-Align te versnellen, worden MEMs korter dan sl uitgefilterd tijdens MEM-extractieproces. Dit vermindert het aantal MEMs worden geëxtraheerd en verwerkt (vervolgens versnellen van het algoritme). Het filteren van korte MEM wordt gedaan door algoritme 1 uit te breiden met een set shift-en bitwise-bewerkingen (aanvullend bestand 1: sectie VII).,

aan de andere kant zijn er zeldzame gevallen waarin korte MEMs deel uitmaken van de uitlijning; bijvoorbeeld een overeenkomend symbool omgeven door mismatches. Zonder dat alle MEMs in de invoerlijst staan, is DP-MEM niet in staat om dezelfde uitlijning te vinden als wanneer alle MEMs in de invoerlijst staan. Om te compenseren voor verloren korte MEMs in de invoer, wij wijzigen Φ2H van DP-MEM om de mogelijkheid van het hebben van korte wedstrijden tussen MEMs overwegen (aanvullend bestand 1: sectie VIII).

Er kunnen moeilijkere gevallen zijn waarin in de uitlijning meerdere korte MEMs bestaan tussen twee MEMs (zie Fig. 9)., De enige manier om de score voor het gebied tussen Mi en Mj in Φ2H correct te identificeren is door een globale afstemming op dit gebied toe te passen. Φ2H is echter een frequente operatie en moet snel blijven. Daarom hebben we besloten om het probleem gedeeltelijk op te lossen door mogelijke gevallen te beperken (een heuristische methode).

een voorbeeld dat meerdere korte MEM toont in het kleine gebied tussen Mi en Mj in een uitlijning

als er gaten zijn in het gebied tussen Mi en Mj, gaan we ervan uit dat er alleen een continue opening aan de linkerkant of aan de rechterkant van het gebied. Dan zijn slechts twee uitlijningen mogelijk voor het gebied., Het aantal overeenkomende symbolen wordt voor beide gevallen op een sequentiële manier geteld en degene die resulteert in maximale overeenkomsten wordt genomen als het aantal overeenkomsten tussen Mi en Mj (aanvullend bestand 1: sectie IX). De sequentiële vergelijking is een dure operatie en we bedenken een methode om de sequentiële vergelijking waar mogelijk te vermijden (aanvullend bestand 1: sectie X).

elk ander geval dat niet aan bovenstaande aanname voldoet, resulteert in een uitlijning met een lagere score., Gezien het lage percentage hiaten en mismatches is de mogelijkheid van meerdere hiaten en mismatches in een klein gebied echter gering.

efficiënte uitlijning extensie

In Φ2 van DP-MEM breidt Mj alle uitlijningen uit die eindigen op {M1 … Mj-1} (indien mogelijk). Echter, voor elke Mj is er een kleinere deelverzameling Ωj⊆{M1 … Mj-1} zodanig dat door MJ uit te breiden tot alle uitlijningen eindigend in een Mi ω Ωj de uitlijning die eindigt in Mj wordt gevonden (Eq. 2). Met andere woorden, er zouden minder \(s_{i}^{j}\) moeten worden geëvalueerd. Definitie van de set Ωj en het bewijs van Eq. 2 worden verstrekt in aanvullend bestand 1: sectie XI., De definitie van Ωj wordt beïnvloed wanneer korte mem verwijdering optimalisatie wordt toegepast (aanvullend bestand 1: sectie XII).

$$ \max\limits_{M_{i} \in \omega_{j}}{S_{i}^{j}} = \max\limits_{1 \leq i \leq j-1}{S_{i}^{j}} $$

(2)

hybride uitlijning

om de nauwkeurigheid van het algoritme te behouden, hebben we besloten een hybride methode gebruiken die een combinatie is van mem-align en Smith-Waterman algoritme. We definiëren drie gevallen waarin MEM-Align onnauwkeurig kan zijn., Als de uitlijning van een paar sequenties in één van deze gevallen valt, gebruiken we het Smith-Waterman algoritme om sequenties uit te lijnen. Deze gevallen zijn:

• wanneer de sequenties repetitief zijn en het aantal geëxtraheerde MEMs de drempelwaarde TM overschrijdt. We ontdekten dat MEM-Align waarschijnlijk onnauwkeurige uitlijning veroorzaakt bij het uitlijnen van repetitieve sequenties. Een passende TM-waarde vermindert de kans op het melden van onnauwkeurige afstemming met een verwaarloosbare toename van de gemiddelde verwerkingstijd.,

• wanneer de berekende uitlijningsscore voor de uitlijning gegenereerd door MEM-Align lager is dan een drempelwaarde TS. Dit geval doet zich meestal voor wanneer er een hiaat in de alignment is die niet kan worden geïdentificeerd als gevolg van banded alignment.

• wanneer er geen MEM langer dan sl bestaat om te worden geëxtraheerd (een zeldzaam geval)., Als sl is ingesteld op een hoge waarde en de gelijkenis tussen sequenties laag is,

hoewel het verzenden van sequentieparen naar een extern algoritme resulteert in extra berekening, blijft het aantal sequenties dat naar het externe algoritme wordt verzonden klein indien geschikte waarden worden gekozen voor TM en TS.

MEMs op afstand overslaan

wanneer de afstand tussen Mi en Mj groot is, is het niet waarschijnlijk dat Mi en Mj als aangrenzende MEMs in de uitlijning staan., Daarom slaat het algoritme de uitbreiding over als de afstand tussen Mi en Mj langer is dan een drempelwaarde TD (waardoor het aantal \(S_{i}^{j}\) dat moet worden geëvalueerd verder wordt verminderd). Deze optimalisatie verbetert de prestaties lichtjes met een verwaarloosbaar neveneffect op de nauwkeurigheid.