Fråga:
Transkriptkoordinatområden till genomiska koordinater
Nathan S. Watson-Haigh
2018-10-05 09:23:03 UTC
view on stackexchange narkive permalink

Jag har två GFF3-filer:

  1. Funktioner som använder transkript-ID: n som landmärken. dvs. "CDS" -funktionstyper med koordinater från transkriptionsutrymme.
  2. Funktioner som använder kromosom-ID som landmärken. dvs "exon" -funktionstyper som använder koordinater från kromosomutrymme.

Jag vill förvandla koordinaterna för funktionerna i fil 1 till koordinatutrymmet för fil 2. dvs. transkriptionsbaserade koordinater till genombaserade coodinates.

Här är ett exempel för fil1: 

  cat transcript_orfs.gff3 ## gff-version 3 ## sekvensregion Tx.1 1 4000Tx. 1 ORF_finder-gen 1 4000. +. ID = 1Tx.1 ORF_finder mRNA 1 4000. +. ID = 2; Förälder = 1Tx.1 ORF_finder exon 1501 2500. + 0 ID = 3; Förälder = 2Tx.1 ORF_finder CDS 1501 2500. + 0 ID = 4; Förälder = 2gt skiss -addintrons transcript_orfs.png transcript_orfs.gff3

coordinates in transcript space

Här är ett exempel för fil1: 

  kattgenom.gff3 ## gff-version 3 ## sekvensregion chr3A_part1 1 454103970chr3A_part1 genomet-sammansättningsgen 1001 6000. +. ID = Txchr3A_part1 genommontering mRNA 1001 6000. +. ID = Tx.1; Förälder = Txchr3A_part1 genom_montering exon 1001 3000. +. ID = Tx.1.exon1; Förälder = Tx.1chr3A_part1 genom_montering five_prime_UTR 1001 2000. +. ID = Tx.1.utr5; Förälder = Tx.1chr3A_part1 genommontering CDS 2001 3000. + 0 ID = Tx.1.cds1; Förälder = Tx.1chr3A_part1 genommontering exon 4001 6000. +. ID = Tx.1.exon2; Förälder = Tx.1chr3A_part1 genommontering CDS 4001 5000. + 2 ID = Tx.1.cds2; Förälder = Tx.1chr3A_part1 genom_montering three_prime_UTR 5001 6000. +. ID = Tx.1.utr3; Förälder = Tx.1gt skiss -addintroner genom.png genom.gff3

coordinates in genomic space

Jag vill konvertera koordinaterna för funktionerna från file1 till kromosomgenomiska koordinater. Jag förväntar mig att få ungefär följande utdata: 

  cat output.gff3
## gff-version 3 ## sekvensregion chr3A_part1 1 454103970chr3A_part1 ORF_finder gen 1001 6000. +. ID = 1chr3A_part1 ORF_finder mRNA 1001 6000. +. ID = 2; Förälder = 1chr3A_part1 ORF_finder exon 2501 3000. + 0 ID = 3.1; Förälder = 2chr3A_part1 ORF_finder CDS 2501 3000. + 0 ID = 4.1; Förälder = 2chr3A_part1 ORF_finder exon 4001 4500. + 0 ID = 3.2; Förälder = 2chr3A_part1 ORF_finder CDS 4001 4500. + 1 ID = 4.2; Förälder = 2gt skiss -addintrons output.png output.gff3

transformed coordinates

Jag har tittat på med mapFromTranscripts () från Bioconductors GenomicRanges -bibliotek men jag har gjort små framsteg när jag försöker dechiffrera manualen.

Två svar:
Daniel Standage
2018-12-13 21:11:31 UTC
view on stackexchange narkive permalink

Denna uppgift har samma "smak" som många jag har gjort tidigare, men varje fall är så subtilt olika att det är omöjligt att skriva ett generaliserat verktyg som fungerar korrekt under alla omständigheter.

R är inte mitt styrhus, men jag kunde kasta något ihop ganska snabbt med standard Python. Detta bör ge dig minst 95% av vägen dit du ska. :-) 

  #! / usr / bin / env python3från argparse import ArgumentParser, FileTypfrån samlingar importera standarddikt från reimportera sökning, subcli = ArgumentParser () cli .add_argument ('txspace', type = FileType ('r')) cli.add_argument ('chrspace', type = FileType ('r')) args = cli.parse_args () exons = defaultdict (list) för rad i args .txspace: om line.startswith ('#') eller line.strip () == '': fortsätt värden = line.split ('\ t') featuretype = värden [2] om featuretyp inte finns i ('CDS', 'exon'): fortsätt txid = värden [0] exons [txid]. lägg till (värden) för rad i args.chrspace: värden = rad.split ('\ t') om len (värden)! = 9: skriv ut ( line, end = '') fortsätt featuretype = värden [2] om featuretype inte finns i ('gen', 'mRNA'): fortsätt skriva ut (line, end = '') om featuretype == 'mRNA': txid = search ( 'ID = ([^; \ n] +)', värden [8]). Grupp (1) start = int (värden [3]) txexons = exons [txid] för exon i txexoner: exon [0] = värden [0] exon [3] = str (int (exon [3]) + start) exon [4] = str (int (exon [4]) + start) exon [8] = sub (r'Parent = [^; \ n] + ',' Parent = '+ txid, exon [8]) print (' \ t'.join (exon), end = '')
Ian Sudbery
2018-12-14 17:04:29 UTC
view on stackexchange narkive permalink

Den länkade kärnan innehåller en klass som jag skrev för att göra detta på GTF-filer:

https://gist.github.com/IanSudery/d8349c22823a475ceb489c3e8aeb448e

Den använder GTF-klassen från cgat, som finns här: https://github.com/cgat-developers/cgat-apps

Du skulle använda den för att gör den här uppgiften så här: 

  från cgat import gtffrom cgat import iotoolsimport TranscriptCoordsInterconvertertrans_coord_file = gtf.iterator (iotools.open_file (args [1]) transkriptioner = dict () för transkript i gtf.transcript_iterator (trans_coords_file): transcript_id = transcript [0] .transcript_id transcripts [transcript_id] = transcripttrans_coord_file.close () genome_coord_file = gtf.iterator (iotools.open_file) (args [2]] gtf.transcript_iterator (genome_coord_file): transcript_id = transcript [0] .transcript_id contig = transcript [0] .contig strand = transcript [0]. sträng om transcript_id inte finns i transkript: fortsätt c onverter = TranscriptCoordsInterconverter (transkript) för exon i transkriptioner [transcript_id]: genome_intervals = converter.transcript_interval2genome_intervals ((exon.start, exon.end)) för intervall i genom_intervaller: new_entry = GTF.Entry (). fromGTF (exoncontent) new_entry = contig new_entry.strand = strand new_entry.start = interval [0] new_entry.end = interval [1] print (str (new_entry) + "\ n")


ⓘ
Denna fråga och svar översattes automatiskt från det engelska språket.Det ursprungliga innehållet finns tillgängligt på stackexchange, vilket vi tackar för cc by-sa 4.0-licensen som det distribueras under.
Loading...