Question:
Inférer les fonctionnalités UTR manquantes dans le fichier GFF3
l0110
2017-12-29 23:49:41 UTC
view on stackexchange narkive permalink

Je travaille sur un fichier GFF qui ne contient pas les informations 5'UTR et 3'UTR. Par exemple:

  ctg123. gène 1050 9000. +. ID = gene00001; Nom = EDEN ctg123. ARNm 1050 9000. +. ID = mRNA00001; Parent = gene00001; Nom = EDEN.1 ctg123. exon 1050 1500. +. ID = exon00002; Parent = ARNm00001 ctg123. exon 3000 3902. +. ID = exon00003; Parent = ARNm00001 ctg123. exon 5000 5500. +. ID = exon00004; Parent = ARNm00001 ctg123. exon 7000 9000. +. ID = exon00005; Parent = ARNm00001 ctg123. CDS 1201 1500. + 0 ID = cds00001; Parent = mRNA00001; Nom = edenprotein.1 ctg123. CDS 3000 3902. + 0 ID = cds00001; Parent = mRNA00001; Nom = edenprotein.1 ctg123. CDS 5000 5300. + 0 ID = cds00001; Parent = mRNA00001; Name = edenprotein.1  

Existe-t-il un moyen d'ajouter des lignes dans le GFF avec les plages 5'UTR et 3'UTR?

Devrions-nous également nous attaquer au volet négatif? Je veux dire, pouvez-vous avoir des lignes gff faisant référence au brin négatif et, si oui, pouvons-nous supposer qu'elles seront correctes? Que les positions seront données par rapport au brin `+` et que la position finale sera "avant" celle de départ?
oui, nous devons également faire face au volet négatif. Voici juste un exemple simplifié.
Ensuite, veuillez [modifier] votre question et inclure également un exemple du volet négatif. Nous aurions également besoin de savoir si chaque gff aura une seule entrée ou si vous pouvez avoir plusieurs gènes dans la même gff. Plus vous nous fournissez d'informations précises, plus il est probable que nous pourrons vous apporter une réponse utile. Cependant, gardez à l'esprit que les UTR peuvent être épissés, ce n'est donc pas un bon moyen d'obtenir les informations dont vous avez besoin. Obtenir une annotation complète de l'endroit où vous avez téléchargé ce serait mieux.
Trois réponses:
Daniel Standage
2019-01-09 20:36:20 UTC
view on stackexchange narkive permalink

Entre autres choses, le programme canon-gff3 du AEGeAn Toolkit déduit et imprime les fonctionnalités manquantes comme les UTR.

  $ canon -gff3 < genes.gff3 > genes-with-utrs.gff3  
Devon Ryan
2017-12-30 19:10:19 UTC
view on stackexchange narkive permalink

Malheureusement, vous pouvez charger un tel fichier GFF3 dans R (en utilisant le package GenomicFeatures) et accéder aux UTR, mais ils ne sont pas ensuite enregistrés si vous utilisez la fonction export de rtracklayer (cela aurait été la solution la plus simple).

Comme je ne sais rien pour faire ça de façon désinvolte, j'ai passé quelques minutes à en écrire une. L'utilisation serait addGFF3UTRs.py input.gff output.gff . Vous pouvez modifier un certain nombre d'étiquettes pour qu'elles apparaissent comme vous le souhaitez. Cela n'a pas été testé de manière approfondie, mais il semble produire les bons résultats dans un petit test que je viens d'exécuter.


Le code est hébergé sur GitHub au lien que j'ai donné ci-dessus et j'inclus également ici:

  #! / usr / bin / env pythonimport argparsedef parseAttributes (kvps): d = dict () pour kvp dans kvps: k, v = kvp.split ("=") d [k] = v return ddef parseGFF3 (fname, topLevel = "gene", secondLevel = "mRNA", CDS = "CDS", exon = "exon"): f = open (fname) genes = dict () transcripts = dict () g2t = dict () # Un dictionnaire avec les ID de gène associés à leurs ID de transcription pour la ligne dans f: if line.startswith ("#"): continue cols = line.strip () .split ("\ t") attribs = parseAttributes (cols [8] .split (";")) ID = attribs ["ID"] if cols [2] == topLevel: genes [ID] = ligne g2t [ID] = [] elif cols [2] == secondLevel: parent = attribs ["Parent"] g2t [parent ] .append (ID) # Génère une liste de [ligne, [entrées exon], [entrées CDS]] si l'ID ne figure pas dans les transcriptions: transcriptions [ID] = [ligne, [], []] sinon: transcriptions [ID] [0] = line else: parent = attribs ["Parent"] si le parent n'est pas dans les transcriptions: transcriptions [parent] = [None, [], []] if cols [2] == exon: transcripts [parent] [1 ] .append ((int (cols [3]), int (cols [4]), ligne))
elif cols [2] == CDS: transcriptions [parent] [2] .append ((int (cols [3]), int (cols [4]), ligne)) f.close () retourne les gènes, transcriptions, g2tdef findUTRs (transcripts, fivePrime = True): pour k, t dans transcripts.items (): # S'il n'y a pas d'entrées CDS alors sautez si len (t [2]) == 0: continue # Récupère le brin ("." sera traité comme "+") strand = t [0] .split ("\ t") [6] # Pour 3 'UTR, permutez simplement le brin sinon cinqPrime: strand = "+" if strand == "- "else" - "exons = [(s, e) pour s, e, _ dans t [1]] CDSs = [(s, e) pour s, e, _ dans t [2]] exons.sort () CDSs.sort () UTRs = [] if strand! = "-": final = CDSs [0] [0] for s, e in exons: if e < final: UTRs.append ((s, e)) elif s < final: UTRs.append ((s, final - 1)) else: break else: final = CDSs [-1] [- 1] pour s, e dans les exons: if e < final: continue elif s > final: UTRs.append ((s, e)) else: UTRs.append ((final + 1, e) ) t.append (UTRs) def saveGFF (oname, genes, transcripts, g2t, fivePrime, threePrime): o = open (oname, "w") o.write ("## gff-version 3 \ n") pour geneID , geneLine dans genes.items (): o.write (geneLine) # Parcourez chaque transcrit pour transID dans g2t ​​[geneID]: t = transcripts [transID] si t [0] est None: continue o.write (t [0 ]) cols = t [0] .strip (). split ("\ t") # Ecrire les exons, puis CDS, puis 5'UTR, puis 3'UTR pour l'exon dans t [1]: o.write (exon [2]) pour CDS dans t [2]: o.write (CDS [2]) pour idx, UTR dans enumerate (t [3]):
o.write ("{} \ t {} \ t {} \ t {} \ t {} \ t. \ t {} \ t. \ t" .format (cols [0], cols [1], fivePrime , UTR [0], UTR [1], cols [6])) o.write ("ID = {}. FivePrimeUTR {}; Parent = {} \ n" .format (transID, idx, transID)) pour idx , UTR dans enumerate (t [4]): o.write ("{} \ t {} \ t {} \ t {} \ t {} \ t. \ T {} \ t. \ T" .format ( cols [0], cols [1], threePrime, UTR [0], UTR [1], cols [6])) o.write ("ID = {}. threePrimeUTR {}; Parent = {} \ n". format (transID, idx, transID)) o.close () parser = argparse.ArgumentParser (description = "Analyser un fichier GFF3 sans entrées UTR et les ajouter. Notez que ce programme suppose un fichier GFF3 raisonnablement formaté") parser.add_argument ("--fiveUTRname", default = "five_prime_UTR", help = "Le libellé des 5 entrées UTR (par défaut:% (default) s)") parser.add_argument ("- threeUTRname", default = "three_prime_UTR", help = "Le libellé pour 3 'entrées UTR (par défaut:% (default) s)") parser.add_argument ("- topLevelID", default = "gene", help = "Le' type 'désignant l'entrée de niveau supérieur (par défaut:% (par défaut) s) ") parser .add_argument ("- secondLevelID", default = "mRNA", help = "Le 'type' désignant l'entrée de deuxième niveau, généralement quelque chose comme 'mRNA' ou 'transcript' (par défaut:% (par défaut) s)") parser.add_argument ("- exonID", default = "exon", help = "Le 'type' désignant les exons (par défaut:% (default) s)") parser.add_argument ("- CDSID", default = "CDS ", help =" Le 'type' désignant CDS (par défaut:% (valeur par défaut) s) ") parser.add_argument (" input ", metavar =" input.gff3 ", help =" Input file name ") parser.add_argument ( "output", metavar = "output.gff3", help = "Output file name") args = parser.parse_args () genes, transcripts, g2t = parseGFF3 (args.input, topLevel = args.topLevelID, secondLevel = args.secondLevelID , exon = args.exonID, CDS = args.CDSID) # Ajouter les UTRsfindUTRs (transcriptions) findUTRs (transcripts, fivePrime = False) # Ecrire les sortiesaveGFF (args.output, genes, transcripts, g2t, args.fiveUTRname, args.threeUTRname )  
vkkodali
2019-03-28 20:56:39 UTC
view on stackexchange narkive permalink

Vous pouvez utiliser le script add_utrs_to_gff.py fourni par NCBI RefSeq à cette fin. L'exécution du script sur votre exemple produit ce qui suit:

  ./add_utrs_to_gff.py example.gff3ctg123. gène 1050 9000. +. ID = gene00001; Nom = EDENctg123. ARNm 1050 9000. +. ID = mRNA00001; Parent = gene00001; Nom = EDEN.1ctg123. exon 1050 1500. +. ID = exon00002; Parent = ARNm00001ctg123. exon 3000 3902. +. ID = exon00003; Parent = ARNm00001ctg123. exon 5000 5500. +. ID = exon00004; Parent = ARNm00001ctg123. exon 7000 9000. +. ID = exon00005; Parent = ARNm00001ctg123. CDS 1201 1500. + 0 ID = cds00001; Parent = mRNA00001; Nom = edenprotein.1ctg123. CDS 3000 3902. + 0 ID = cds00001; Parent = mRNA00001; Nom = edenprotein.1ctg123. CDS 5000 5300. + 0 ID = cds00001; Parent = mRNA00001; Nom = edenprotein.1ctg123. cinq_prime_UTR 1050 1200. +. ID = utr00002; Parent = ARNm00001ctg123. three_prime_UTR 5301 5500. +. ID = utr00004; Parent = ARNm00001ctg123. three_prime_UTR 7000 9000. +. ID = utr00005; Parent = mRNA00001  


Ce Q&R a été automatiquement traduit de la langue anglaise.Le contenu original est disponible sur stackexchange, que nous remercions pour la licence cc by-sa 3.0 sous laquelle il est distribué.
Loading...