Frage:
Zusammenführen von Regionen gemäß ihrer Kennung
bapors
2018-03-20 14:03:00 UTC
view on stackexchange narkive permalink

Ich habe eine riesige Datei (20 GB ), die eine Reihe genomischer Positionen enthält, und für jede Position gibt es eine Kennung ( 4. Spalte ), was manchmal der Fall ist das gleiche.

  file1.txtchr1 10 20 ABCchr1 13 20 ABCchr1 14 21 ABCchr1 22 27 ABCchr1 29 37 ABCchr2 15 21 JJJchr2 21 31 JJJchr2 23 27 JJJchr2 35 56 JJJchr2 25 26 MMMchr3 32 42 MMMchr3 45 76 MMMchr3 88 101 MMMchr3 101 105 MMM  

Ich versuche, die Bereiche der Spalten 2 und 3 zusammenzuführen, sofern sie vorhanden sind das gleiche Chromosom und die gleiche Kennung (mit den gleichen Spalten 1st und 4th ).

Ich habe versucht, Bedtools , wie in diesem Beitrag gezeigt, der bei mir sehr ähnlich, aber nicht gleich aussieht, da sie nach unterschiedlichen Ausgaben suchen (Bereich vs. Zusammenführen)

Also die Antworten In dieser Veröffentlichung wurde groupby in bedtools2 verwendet, Quellcode von github .

Ich habe dasselbe wie folgt angewendet:

  ./groupBy -i ../file1.txt -g 1,4 -c 2,3 -o min , max | awk -v OFS = '\ t' '{print $ 1, $ 3, $ 4, $ 2}' chr1 10 37 ABCchr2 15 56 JJJchr2 25 26 MMMchr3 32 105 MMM  

Aber es verbindet nur die Bereiche von Bezeichnern. Es gruppiert nach und wird nicht zusammengeführt.

Und wenn ich die Zusammenführungsfunktion bedtools v2.26 verwende:

  >sort -k1,1 -k2,2n file1.txt > file2 .txt>cat file2.txtchr1 10 20 ABCchr1 13 20 ABCchr1 14 21 ABCchr1 22 27 ABCchr1 29 37 ABCchr2 15 21 JJJchr2 21 31 JJJchr2 23 27 JJJchr2 25 26 MMMchr2 35 56 JJJchr3 32 42 MMMchr3 45 76 MMMch3 file2.txt chr1 10 21chr1 22 27chr1 29 37
chr2 15 31chr2 35 56chr3 32 42chr3 45 76chr3 88 105  

Was zusammengeführt wird, ohne die Bezeichner zu sehen, da es mir auch diese Ausgabe gibt.

Was ich möchte, ist füge die Bereiche zusammen, wenn sie sich überlappen , wenn dieselbe 1. und 4. Spalte hat wie:

  expected_output.txtchr1 10 21 ABCchr1 22 27 ABCchr1 29 37 ABCchr2 15 31 JJJchr2 35 56 JJJchr2 25 26 MMMchr3 32 42 MMMchr3 45 76 MMMchr3 88 105 MMM  

Ich habe nur Zugriff auf 15 GB RAM, also kann ich Speichern Sie nicht die gesamte Datei im Speicher.

Mögliches Duplikat von [Bettdatensätze basierend auf Name zusammenführen] (https://bioinformatics.stackexchange.com/questions/2265/merging-bed-records-based-on-name)
Dieser Beitrag ist bereits in meinem Eintrag verlinkt und sie suchen nach unterschiedlichen Ausgaben, also nicht nach einem Duplikat
Vier antworten:
#1
+6
dariober
2018-03-20 15:38:19 UTC
view on stackexchange narkive permalink

Wenn ich das richtig verstehe, können Sie Dummy-Chromosomen erstellen, die durch Zusammenführen von Chromosom und Identifikator hergestellt, mit Bettwerkzeugen zusammengeführt, Chromosomen und Identifikatoren zurückgespalten werden. Zum Beispiel

  awk -v OFS = '\ t' '{print $ 1 "_" $ 4, $ 2, $ 3}' file1.txt \ | mergeBed \ | sed 's / _ / \ t /' \ | awk -v OFS = '\ t' '{print $ 1, $ 3, $ 4, $ 2}'  

(Angenommen, "_" ist nicht Teil des Chromosomen- oder Identifikatornamens, passen Sie dies entsprechend an ).

Vielen Dank für die Antwort, es funktioniert gut. Würden sich "Bedtools" jedoch beschweren, wenn die Spalte "$ 4th" Sonderzeichen enthält?
@bapors warum sollte die 4. Spalte Sonderzeichen haben?
@bapors Wenn Sie meinen, dass sich Bedtools beschweren würden, wenn das geänderte 1. Feld Sonderzeichen enthalten würde, bearbeiten Sie bitte Ihre Frage, um ein repräsentatives Beispiel für Ihre Eingabedaten zu geben. Das heißt, ich habe gerade versucht, die Zeichen "& ^% * () @ # ^" zu verwenden, und Bettwerkzeuge hatten kein Problem.
Nein, meine Eingabe enthält keine Sonderzeichen, meine Frage war ausschließlich theoretisch. Perfekt!
#2
+4
user1141
2018-03-20 19:06:24 UTC
view on stackexchange narkive permalink

Nur um dies hinzuzufügen, können Sie mit bedtools jeden Bezeichner separat verarbeiten, was angesichts Ihrer großen Datei möglicherweise praktisch ist oder nicht.

  für ID in $ (perl -ane 'print $ F [3]. "\ n";' file1.txt | sort | uniq); grep $ ID file1.txt | \ perl -ane 'chomp; print join ("\ t", @ F [0 .. $ # F]). "\ n"; ' | \ sort -V | \ bedtools merge -i stdin -c 4 -o differentdone | sortiere -V  
#3
+3
Devon Ryan
2018-03-20 14:42:12 UTC
view on stackexchange narkive permalink

Das folgende Python-Skript macht das, was Sie wollen, und sollte relativ speichereffizient sein. Es verarbeitet jeweils ein einzelnes Chromosom. Sortieren Sie also entweder die Eingabe oder stellen Sie sicher, dass mindestens Einträge in einer Gruppe nebeneinander liegen.

  # ! / usr / bin / env pythonimport sys # Verarbeitet eine sortierte Liste von Tupeln (z. B. [(0, 10), (10, 20), ...]) # Überlappende Elemente werden zusammengeführt und eine Liste der zusammengeführten Regionen wird # zurückgegeben (z. B. [(0, 20), ...]) def merge (Werte): Start, Ende = Werte [0] o = [] für v in Werten [1:]: wenn v [0] < = Ende: wenn v [1] > Ende: Ende = v [1] sonst: o.append ((Start, Ende)) Start, Ende = v o.append ((Start, Ende)) return o # Prozess a einzelnes Chromosom, eine Gruppe bei einem zeitgesteuerten ProzessChrom (d, chrom): für k, v in d.items (): v.sort () tuples = merge (v) für t in tuples: print ("{} \ t { } \ t {} \ t {} ". Format (chrom, t [0], t [1], k)) def main (): currentChrom = Keine d = dict () für Zeile in sys.stdin: cols = Linien plit () wenn cols [0] == currentChrom: wenn cols [3] nicht in d: d [cols [3]] = [] d [cols [3]]. append ((int (cols [1]), int (cols [2]))) else: wenn currentChrom: processChrom (d, currentChrom) d = dict () d [cols [3]] = [] d [cols [3]]. append ((int (cols [) 1]), int (cols [2]))) currentChrom = cols [0] wenn currentChrom: processChrom (d, currentChrom) wenn __name__ == "__main__": main ()  

Die Verwendung ist python foo.py < foo.txt und die Ausgabe lautet:

  chr1 10 21 ABCchr1 22 27 ABCchr1 29 37 ABCchr2 25 26 MMMchr2 15 31 JJJchr2 35 56 JJJchr3 32 42 MMMchr3 45 76 MMMchr3 88 105 MMM  

Sie werden feststellen, dass dies nicht unbedingt sortiert ist. Daher möchten Sie die Ausgabe möglicherweise über bedtools sort oder sort weiterleiten.

Vielen Dank für Ihre Antwort. Meinen Sie damit, dass die Eingabedatei nach 1. oder 4. Spalte sortiert werden soll?
Sie können es auch nach der ersten Spalte sortieren. Sie können das Programm alternativ neu schreiben, um Dateien zu verarbeiten, die nach der vierten Spalte sortiert sind, aber das aktuelle Skript ist nicht dafür geschrieben.
Sicher. Der Code ist sauber und funktioniert gut. Ich bin froh, dass es anpassungsfähig ist
#4
+3
Chris_Rands
2018-03-20 20:04:13 UTC
view on stackexchange narkive permalink

Hier ist eine ziemlich einfache (und hoffentlich lesbare) native Python-Lösung. Es wird davon ausgegangen, dass die Bettdatei vor dem Parsen sortiert wird:

  aus itertools import groupbyclass BedRecord (Objekt): def __init __ (self, chr, start, end, feature): self.chr = chr self.start = start self.end = end self.feature = featuredef parse_records (in_file): mit open (in_file) als f: für Zeile in f: Ausbeute BedRecord (* line.strip (). split ()) def merge_groups ( group): output = [next (group)] für Datensatz in Gruppe: if record.start < = output [-1] .end: output [-1] .end = max (output [-1] .end, record. end) else: output.append (record) return outputdef main (in_file): für _, group in groupby (parse_records (in_file), lambda x: (x.chr, x.feature)): für r in merge_groups (group) : print ('\ t'.join ([r.chr, r.start, r.end, r.feature])) if __name__ ==' __main__ ': main (' test.bed ')  
Vielen Dank für Ihre Antwort. Es wird davon ausgegangen, dass die Eingabe nach welcher Spalte sortiert wird.
@bapors effektiv nach Chromosom, dann Feature dann Startwert, können Sie leicht Code hinzufügen, um nach diesen Kriterien zu sortieren, wenn diese Annahme verletzt wird


Diese Fragen und Antworten wurden automatisch aus der englischen Sprache übersetzt.Der ursprüngliche Inhalt ist auf stackexchange verfügbar. Wir danken ihm für die cc by-sa 3.0-Lizenz, unter der er vertrieben wird.
Loading...