Optimierung von PCR-Primern, die auf das bakterielle 16S ribosomale RNA-Gen abzielen / BMC Bioinformatics

Problembeschränkungen

Wie im vorherigen Absatz angegeben, sollte ein optimales Primer-Set-Paar gleichzeitig die Effizienz und Abdeckung maximieren und Matching-Bias minimieren. Im Folgenden beschreiben wir, wie wir diese Einschränkungen quantitativ codiert haben.

Effizienz

Die perfekten Primer-Set-Paare sollten mehrere Einschränkungen erfüllen, um die PCR-Effizienz und-Spezifität zu verbessern ., Die gleichzeitige Erfüllung aller Einschränkungen ist jedoch häufig unpraktisch, und die meisten hochmodernen Primer verletzen eine oder mehrere Einschränkungen . Wir haben uns daher entschlossen, Efficiency als Optimierungs-Score einzuführen, der viele der Einschränkungen als Fuzzy-Score-Funktionen codiert. Genauer gesagt haben wir unseren Effizienzwert als die Summe von zehn Score-Termen definiert: sieben Fuzzy-Score-Terme im Zusammenhang mit Single-Primer-Effizienzbeschränkungen, gemittelt über alle Primer in den Primer-Set-Paaren, plus drei Score-Terme im Zusammenhang mit der Effizienz der Primer-Set-Paare als Ganzes., Da alle Begriffe zwischen 0 und 1 variieren sollen, reicht der Optimierungswert von 0 (minimaler Wirkungsgrad) bis 10 (maximaler Wirkungsgrad).

Im Großen und Ganzen zählt unser Fuzzy-Score 1 für jede perfekt erfüllte Einschränkung oder alternativ einen Wert zwischen 0 und 1, je nachdem, wie nahe der Primer an der Einschränkungsgrenze liegt. Betrachten Sie als Beispiel die Primer-Schmelztemperatur Tm. Tm sollte in einer perfekten Grundierung größer oder gleich 52 Grad sein, aber 51 ist immer noch erträglich, wenn auch nicht ideal., In diesem Fall weist unsere Fuzzy-Scoring-Funktion 1 Temperaturen von 52 Grad oder mehr zu, 0 Temperaturen von 50 Grad oder weniger und betrachtet eine lineare ansteigende Funktion zwischen 50 und 52 Grad. Jeder Begriff wird im Folgenden genau beschrieben.

Die 7 Single-primer Score Terms sind:

Schmelztemperatur: Die Schmelztemperatur Tm eines Primers wird mit der Formel des nächsten Nachbarn berechnet . Der Score-Term ist 1, wenn Tm ≥ 52, 0, wenn Tm ≤ 50 und (Tm-50)/2, wenn 50 < Tm < 52.,
GC-Gehalt: Der GC-Gehalt ist der Anteil fGC von Basenpaaren in der Primersequenz, der entweder G (Guanin) oder C (Cytosin) entspricht. Der Score-Term ist 1, wenn 0,5 ≤ fGC ≤ 0,7, 0, wenn fGC > 0,7 oder fGC < 0,4 und (0,5 – fGC)/0,1, wenn 0,4 ≤ fGC < 0,5.
3′ – Endstabilität-score term 1: Bezüglich 3′-Endstabilität sind zwei Score Terms definiert. Der erste Term ist 0, wenn die letzten drei Basen des Primers vollständig aus As (Adeninen) und Ts (Thyminen) und 1 ansonsten bestehen.,
3′ – Endstabilität-Score term 2: Der zweite Score Term ist 0, wenn die letzten 5 Basen mehr als 3 Cs oder Gs enthalten, und 1 andernfalls.
Homopolymer: Ein Homopolymer ist eine Sequenz identischer Nukleotide. Der Score-Term ist 1, wenn keine Homopolymere länger als 4 nt sind, 0,5, wenn keine Homopolymere länger als 5 nt sind, und 0, wenn mindestens ein Homopolymer länger als 5 nt in der Sequenz vorhanden ist.
Selbstdimere: Das Vorhandensein von sich selbst komplementären Regionen zwischen Paaren identischer Primer kann zur Erzeugung von Selbstdimeren führen., In Anbetracht der maximalen Anzahl von Übereinstimmungen in einer spaltfreien Ausrichtung zwischen einem Primer mit seinem umgekehrten Komplement, maxM, beträgt der Score-Term 1, wenn maxM ≤ 8, 0, wenn maxM ≥ 11 und (11-maxM)/3, wenn 8 < maxM < 11.
Haarnadeln: Eine Haarnadel kann in Gegenwart von Selbstkomplementarität innerhalb der Primersequenz gebildet werden, insbesondere an ihrem 3′ – Ende., Der Score-Term ist 0, wenn für mindestens eine spaltfreie Ausrichtung zwischen dem Primer und dem umgekehrten Komplement seines 3′-Endes sowohl das letzte Nukleotid als auch 3 oder mehr der 4 vorhergehenden Nukleotide übereinstimmen, und 1 andernfalls.

Die 3 Primer-Set-Pairs score terms sind wie folgt definiert:

Schmelztemperaturbereich: Der Schmelztemperaturbereich ΔTm eines Primer-Set-Paares wird als Maximum minus dem Minimum der Schmelztemperaturen aller Primer im eingestellten Paar berechnet., Der Score-Term ist 1, wenn ΔTm ≤ 3, 0, wenn ΔTm ≥ 5 und (5-ΔTm)/2, wenn 3 < ΔTm < 5.
Dimere: Wir betrachten die maximale Anzahl von Übereinstimmungen maxM über alle möglichen Ausrichtungen zwischen allen möglichen Kombinationen von Vorwärts – und Rückwärts-Primern aus einem Primer-Set-Paar. Die Partitur Begriff ist 1, wenn maxM ≤ 8, 0, wenn maxM ≥ 11 und 11 – maxM)/3, wenn 8 < maxM < 11.,
Amplicon-Längenbereich: Aufgrund der bekannten Verringerung des PCR-Wirkungsgrades mit zunehmender Amplicon-Länge wollen wir, dass die Längen der erzeugten Amplicons in einem engen Bereich liegen. Wir möchten insbesondere vermeiden, dass Amplicons viel kürzer als die Ziellänge sind, da sie in Bezug auf die anderen überverstärkt wären. Wir wollen jedoch einen kleinen Bruchteil der Ausreißer tolerieren können, um zu vermeiden, dass potenziell wertvolle Primer-Set-Paare aufgrund einiger seltener Sequenzen bestraft werden., Bei einem repräsentativen Satz von bakteriellen 16S-Sequenzen, der von nun an als „Referenzsatz“ bezeichnet wird, betrachten wir die Differenz Δamplen zwischen dem Median und dem ersten Perzentil der Amplikonlängen über alle möglichen Amplicons hinweg, die durch Abgleich aller Kombinationen von Vorwärts-und Rückwärts-Primern aus dem eingestellten Paar mit dem Referenzsatz gebildet wird. Der Score-Term ist 1, wenn Δamplen ≤ 50 Nukleotide, 0, wenn Δamplen ≥ 100 und (100-Δamplen)/50, wenn 50 < Δamplen < 100.,

Die Auswahl der Bewertungskriterien und der Standardschwellenwert basieren auf früherer Literatur . Sowohl die Schwellenwerte als auch die Fuzzy-Toleranzintervalle können jedoch vom Benutzer unterschiedlich von der Standardeinstellung und entsprechend seinen experimentellen Anforderungen eingestellt werden, indem beim Aufruf des Befehlszeilenwerkzeugs die gewünschten Werte als Eingabeparameter angegeben werden.

Coverage

Der Coverage Score ist definiert als die Anzahl von 16S Sequenzen, die mit mindestens einem Primer übereinstimmen., In Anbetracht der Sequenzen eines Primers und eines bakteriellen 16S definieren wir die letzten 5 Nukleotide am 3′-Ende eines Primers und betrachten eine 16S-Sequenz als mit dem Primer übereinstimmend, wenn ein Bereich der 16S-Sequenz vorhanden ist, der i) genau mit dem Samen des Primers übereinstimmt; und ii) der Rest des Primers mit höchstens 2 nicht Übereinstimmungen . Eine 16S-Sequenz aus einem Referenzsatz wird als durch ein Primer-Set-Paar abgedeckt angesehen, wenn mindestens ein Vorwärts-und ein Rückwärts-Primer im Primer-Set-Paar mit der Sequenz übereinstimmen., Da die PCR-Effizienz mit der Ampliconlänge abnimmt, legen wir eine weitere Einschränkung fest: Bei einem Primer-Set-Paar und einem Referenzsatz von 16S-Sequenzen schätzen wir die Zielampliconlänge als Median der Längen aller Amplicons, die durch Abgleich aller Kombinationen von Vorwärts-und Rückwärts-Primern aus dem Primer-Set-Paar mit dem Referenzsatz erhalten werden. Wir betrachten dann nicht alle 16S-Referenzsequenzen als abgedeckt, deren Amplikonlänge mehr als 100 Nukleotide (entweder länger oder kürzer) von der Ziellänge unterscheidet.,

Matching-bias

Bei einem Referenzsatz von 16S-Sequenzen und einem Primer-Set-Paar misst der dritte Optimierungswert die Variabilität der Anzahl der Kombinationen von Vorwärts-und Rückwärts-Primern, die jeder 16S-Referenzsequenz entsprechen. Die Abdeckungsvariabilität aufgrund einer übereinstimmenden Verzerrung sollte minimiert oder zumindest berücksichtigt werden, wenn die Studie die relative Häufigkeit der verschiedenen Bakterienarten quantifizieren soll, da die Amplifikations-Verzerrung gegenüber den Arten, die durch mehr Kombinationen von Vorwärts-und Rückwärts-Primern abgedeckt werden, größer ist., Als Maß für Matching-Bias nutzen wir den Variationskoeffizienten der Abdeckung über die Zielsequenzen hinweg, berechnet als Standardabweichung über den Mittelwert der Anzahl der Kombinationen, die jeder Sequenz entsprechen.

Referenzsatz von 16S-Sequenzen, Vorbereitung und Annotation

Zur Optimierung der drei obigen Scores verlassen wir uns auf einen repräsentativen Satz bakterieller 16S-Sequenzen, die aus einer öffentlichen 16S-Sequenzdatenbank, GreenGenes, extrahiert wurden ., Die GreenGenes 16S-Sequenzdatenbank ist in Operational Taxonomic Units (OTUs) organisiert, bei denen es sich um verschachtelte Sequenzcluster in der Datenbank handelt, die auf verschiedenen Ebenen der Interclusterähnlichkeit organisiert sind. Für jede Ebene der Ähnlichkeit wird jedem Cluster eine Referenzsequenz zugeordnet, die allen anderen Sequenzen im selben Cluster maximal ähnlich ist . Der Satz von Referenzsequenzen kann somit als repräsentative Teilmenge der gesamten Sequenzdatenbank betrachtet werden und wird immer genauer, um die Ähnlichkeit zwischen den Clustern (und damit die Anzahl der Referenzsequenzen) zu erhöhen., Wir wählten eine 85% ige Inter-Cluster-Ähnlichkeitsstufe als einen guten Kompromiss zwischen Repräsentativität und Komplexität, entsprechend einer Reihe von 5088 repräsentativen Sequenzen, die zur Bewertung der Optimierungskriterien verwendet werden sollen.

Obwohl die GreenGenes-Taxonomie bei der Annotation der Bakterien-und Archaea-Domänen sehr empfindlich ist, dient sie nicht zur Unterscheidung von Sequenzen, die zu Eukaryoten oder Viren gehören., Aus diesem Grund haben wir beschlossen, 16S-Bakteriensequenzen mithilfe der ursprünglichen NCBI-Taxonomie neu zu kommentieren, um unter den repräsentativen Sequenzen nur diejenigen zu identifizieren, die zur Bakteriendomäne gehören. Da Domäneninformationen in der NCBI-Annotation für etwa 20% der Sequenzen fehlen, haben wir ein Ad-hoc-Verfahren entwickelt, um Bakteriensequenzen unter diesen zu identifizieren. Das Verfahren wird in den ergänzenden Materialien ausführlich beschrieben (siehe zusätzliche Datei 1)., Wir haben uns konservativ dafür entschieden, nur die als Bakterien kommentierten Sequenzen sowohl in unserer kuratierten, NCBI-basierten Annotation als auch in der ursprünglichen GreenGenes-Annotation zu berücksichtigen. Dies führte zu einem Satz von 4573 repräsentativen 16S-Sequenzen, die zur Bakteriendomäne gehören.,

Optimierungsalgorithmus

Da das Problem der optimalen Primerwahl die gleichzeitige Optimierung verschiedener konkurrierender Scores erfordert, kann es als ein multiobjektives Optimierungsproblem gegossen werden, bei dem der Suchraum die Menge aller möglichen Primer-Set-Paare ist und eine Bewertungsfunktion oder ein Optimierungskriterium definiert werden kann, um Effizienz und Abdeckung zu maximieren und Matching-Bias zu minimieren., Wenn mehr als ein Kriterium gleichzeitig optimiert werden muss, die zu optimierenden Ziele jedoch widersprüchlich sind, interessiert man sich in der Regel nicht für eine einzige Lösung, sondern für den Satz Pareto optimaler Lösungen, d. H. Für den Satz von Lösungen, für die keines der Ziele verbessert werden kann, ohne mindestens ein anderes Ziel zu opfern ., Das Ergebnis der Multi-Objective-Optimierung ist nicht mehr ein eindeutiges optimales Primer-Set-Pair, wie bei der Single-Objective-Optimierung, sondern eine Sammlung von Primer-Set-Paaren, die nicht schlechter als jedes andere Primer-Set-Pair und streng besser sind nach mindestens einem der Kriterien. Genauer gesagt, für das tri-objektive Optimierungsproblem der Maximierung der Effizienz (E)-und Abdeckung (C)-Optimierungswerte und der Minimierung der Matching-Bias (M) – Punktzahl, wie im vorherigen Abschnitt definiert, werden Kandidatenprimer-Set-Paare gemäß einem objektiven Funktionsvektor f = (f E ; f C ; fM) ausgewertet., Bei zwei Primer-Set-Paaren p und p‘ sagen wir, dass p p‘ (p ≺ p‘) dominiert, wenn und nur wenn f (p) ≠ f (p‘), fE (p) ≥ fE (p‘), fC (p) ≥ fC (p‘) und fM (p) ≤ fM (p‘). Wenn kein p‘ existiert, so dass p ‚ ≺ p, wird das Primer-Set-Paar p als Pareto-optimal bezeichnet. In diesem Zusammenhang besteht das Ziel der optimalen Primerwahl darin, die Menge aller pareto-optimalen Primer-Set-Paare zu bestimmen (oder zu approximieren), deren Bild im Tri-Objective-Raum als Pareto-Front bezeichnet wird .

Um nach der optimalen Pareto-Front zu suchen, verlassen wir uns auf den von Dubois-Lacoste et al. vorgeschlagenen zweiphasigen iterierten Best improvement Local Search-Ansatz., und effektiv genutzt in Sambo et al. und Borrotti et al. für die optimale multi-objektive Gestaltung von Experimenten.

Die lokale Suche beginnt mit einer ersten Lösung und verfeinert sie iterativ, indem sie kleine lokale Änderungen anwendet und jedes Mal ihre Auswirkungen auf die Lösungsqualität bewertet.es hört auf, wenn keine weiteren lokalen Änderungen die Lösung verbessern können. Der Prozess wird von mehreren verschiedenen Startpunkten aus iteriert und die beste Lösung, die jemals gefunden wurde, wird als Annäherung an das unbekannte Optimum zurückgegeben ., Eine gängige Erweiterung der lokalen Suche auf den multiobjektiven Fall besteht darin, von einer Reihe anfänglicher Pareto-Lösungen auszugehen, eine Lösung von vorne abzutasten, mit lokaler Suche eine zufällige Skalarisierung des Problems zu optimieren, dh eine lineare Kombination der Optimierungswerte mit Gewichten, die gleichmäßig zufällig von der Einheit Simplex abgetastet wurden, die Pareto-Front zu aktualisieren und zu iterieren, bis eine Beendigungsbedingung erfüllt ist .,

Die Prozedur MULTI-OBJECTIVE-SEARCH, deren Pseudocode im Folgenden angegeben wird, erhält als Eingaben den gewünschten Bereich von Ampliconlängen (rangeamplen), einen repräsentativen Satz von 16S-Sequenzen (repset), einen anfänglichen Satz von (möglicherweise degenerierten) Primerpaaren (init) und die Anzahl der Neustarts (nres). Das Verfahren beginnt mit der Auswahl aller möglichen Primerpaare mit der gewünschten Ampliconlänge, Primerlänge (zwischen 17 und 21 Nukleotiden) und Zieldomäne (Bakterien oder Universal) aus init.,

Degenerierte Primerpaare werden in nicht degenerierte Primer-Set-Paare konvertiert und einem Archiv hinzugefügt. Das Verfahren iteriert dann nrest-Zeiten, wobei jedes Mal ein zufälliges Primer-Set-Paar pstart von der Pareto-Front und ein zufälliger Vektor α von relativen Gewichten für die Optimierungswerte abtastet werden, wobei die Gewichte gleichmäßig von der Einheit simplex abgetastet werden; Das Verfahren löst dann eine Skalarisierung des Multi-Objective-Problems, dh ein Single-Objective-Problem, bei dem eine lineare Kombination der drei Ziele mit relativen Gewichten α maximiert wird, und fügt das Ergebnis dem Archiv hinzu., Zu diesem Zweck werden Effizienz -, Abdeckungs-und Matching-Bias – Scores auf ihr Maximum normalisiert, so dass jeder normalisierte Score zwischen 0 und 1 liegt und Matching-Bias als 1-Matching-Bias neu definiert wird, so dass es maximiert werden kann als die anderen Scores., amplikon-Länge in rangeamplen

2 zu Archiv Hinzufügen der entsprechenden nicht-degenerierte primer-set-Paare

3 für r = 1 bis nrest

4 pf = PARETO-FRONT(Archiv)

5 Probe pstart von pf

6 Sample α von 3, mit Σi ai = 1

7 p = LOKALE-SUCHE(pstart , α , repset)

8-Add p to archive

9 return Archiv

Einzel-Ziel-Optimierung erhält man mit die Beste Verbesserung Lokale-Suche-Algorithmus : beginnend von einer ersten primer-set-Paares, der LOKALE-SUCHE-Algorithmus durchläuft die Primer der set-paar und für jeden primer, scannt die Nachbarschaft, ich.,e. der Satz aller möglichen lokalen Störungen der Grundierung. Lokale Störungen bestehen in allen möglichen Flips eines Nukleotids (Beurteilung der drei anderen möglichen Basen) und allen möglichen Zusätzen und Entfernungen eines Nukleotids an den Extremitäten., Die Suche im Lösungsraum wird mit dem bestmöglichen lokalen Suchansatz durchgeführt: Nachdem die gesamte Nachbarschaft wie oben erläutert generiert wurde, wählt der Algorithmus die beste Nachbarstörung aus, beginnt damit, die nächste Nachbarschaft zu generieren, und iteriert, bis sie eine Lösung erreicht, für die keine bessere Nachbarstörung gefunden werden kann. Das Verfahren wird beendet, wenn keine weiteren lokalen Verbesserungen auf eine Grundierung im Primer-Set-Paar angewendet werden können., Die Funktion WEIGHTED-SCORE berechnet die drei Optimierungswerte aus einem Primer-Set-Paar und der repräsentativen Menge, multipliziert die Werte mit den relativen Gewichten α und gibt die Summe der Ergebnisse zurück.

Wir haben ein Softwaretool entwickelt, das unseren Ansatz umsetzt und unter der GNU General Public License als mopo16S Software tool (Multi-Objective Primer Optimization for 16S experiments) unter veröffentlicht., mopo16S ist als Multithreading-C++ – Befehlszeilentool implementiert; das Software-Tool stützt sich auf die effizienten Algorithmen und Datenstrukturen aus der SeqAn-Bibliothek und verwendet die OpenMP-Bibliothek für Multithreading.,>

4 for i = 1 to |pcurr|

5 pri = i-th-primer von pcurr

6 pnew = pcurr mit allen möglichen Ergänzungen und Löschungen von einer Basis an den Extremitäten und Ersatz von eine Basis von pri –

7 scorenew = WEIGHTED-SCORE(pnew , α , repset)

8, wenn scorenew > scorebest

9 pbest = pnew

10 scorebest = scorenew

11 pcurr = pbest

12 zurück pcurr

State-of-the-art-primer-Paare, die als erste Antworten

Wir haben die online-Datenbank probeBase als eine Quelle von Kandidaten-primer-set-Paare verwendet werden, die als erste Lösungen, die von mopo16S., Die Datenbank enthält mehr als 500 Paare von (möglicherweise degenerierten) Primern und berichtet für jeden Primer seine Sequenz, den Strang und die Position, an der er mit dem Referenz-16S-Escherichia-coli-Gen übereinstimmt, und die Zieldomäne, für die er ausgelegt ist (entweder Bakterien, Archaeen oder universell).,

Angesichts eines gewünschten Bereichs für die Zielampliconlänge als Eingabe von mopo16S haben wir alle Primerpaare aus der probeBase-Datenbank ausgewählt, die alle folgenden Eigenschaften erfüllen:

Ampliconlänge im gewünschten Bereich;
Länge beider Primer größer oder gleich 17 nt und kleiner oder gleich 21 nt;
oder Universelle Zieldomäne beider Primer.,

Da unser Ansatz darin besteht, mit Sätzen nicht degenerierter Primer zu arbeiten, ersetzen wir im Falle von Degenerationen entweder im vorderen oder im umgekehrten Primer den degenerierten Primer durch einen entsprechenden Satz nicht degenerierter Primer, erhalten durch Zuweisen aller möglichen Kombinationen von Werten zu den degenerierten Nukleotiden im Primer. Ein Beispiel für dieses Verfahren ist in Tabelle 1 angegeben.

Wir berechneten die drei Werte für jedes Primer-Set-Paar und identifizierten unter diesen die Primer-Set-Paare, die die anfängliche Pareto-Front bildeten.