2
$\begingroup$

I making a bed file for RSeQC, so it can do things like compute the number of reads from exons, introns, 5"UTRs, etc.

I want to use a bed file that corresponds to my GTF file, so I use gtf2bed to make a bed file, like this:

awk '{ if ($0 ~ "transcript_id") print $0; else print $0" transcript_id \"\";"; }' Homo_sapiens.GRCh38.gtf | gtf2bed - > $Homo_sapiens.GRCh38.bed

This file throws some errors, among other when using it with junction_annotation.py:

Reading reference bed file:  /rst1/2017-0205_illuminaseq/data/references/Reference_Genomes/GRCh38.95/Homo_sapiens.GRCh38.87.bed  ... Traceback (most recent call last):
  File "/rst1/2017-0205_illuminaseq/tools/snakemake_conda_env/f4c16777/bin/junction_annotation.py", line 125, in <module>
    main()
  File "/rst1/2017-0205_illuminaseq/tools/snakemake_conda_env/f4c16777/bin/junction_annotation.py", line 109, in main
    obj.annotate_junction(outfile=options.output_prefix,refgene=options.ref_gene_model,min_intron=options.min_intron, q_cut = options.map_qual)
  File "/rst1/2017-0205_illuminaseq/tools/snakemake_conda_env/f4c16777/lib/python2.7/site-packages/qcmodule/SAM.py", line 3762, in annotate_junction
    exon_starts = map( int, fields[11].rstrip( ',\n' ).split( ',' ) )
ValueError: invalid literal for int() with base 10: 'gene_version'

and when I compare this bed file to the RSeQC provide (dated) one these are the differences:

My bed file I created:

1       11868   12227   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000456328"; transcript_version "2"; exon_number "1"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-202"; transcript_source "havana"; transcript_biotype "processed_transcript"; exon_id "ENSE00002234944"; exon_version "1"; tag "basic"; transcript_support_level "1";
1       11868   14409   ENSG00000223972 .       +       havana  gene    .       gene_id "ENSG00000223972"; gene_version "5"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_id "";
1       11868   14409   ENSG00000223972 .       +       havana  transcript      .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000456328"; transcript_version "2"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-202"; transcript_source "havana"; transcript_biotype "processed_transcript"; tag "basic"; transcript_support_level "1";
1       12009   12057   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; exon_number "1"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; exon_id "ENSE00001948541"; exon_version "1"; tag "basic"; transcript_support_level "NA";
1       12009   13670   ENSG00000223972 .       +       havana  transcript      .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; tag "basic"; transcript_support_level "NA";
1       12178   12227   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; exon_number "2"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; exon_id "ENSE00001671638"; exon_version "2"; tag "basic"; transcript_support_level "NA";
1       12612   12697   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; exon_number "3"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; exon_id "ENSE00001758273"; exon_version "2"; tag "basic"; transcript_support_level "NA";
1       12612   12721   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000456328"; transcript_version "2"; exon_number "2"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-202"; transcript_source "havana"; transcript_biotype "processed_transcript"; exon_id "ENSE00003582793"; exon_version "1"; tag "basic"; transcript_support_level "1";
1       12974   13052   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; exon_number "4"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; exon_id "ENSE00001799933"; exon_version "2"; tag "basic"; transcript_support_level "NA";
1       13220   13374   ENSG00000223972 .       +       havana  exon    .       gene_id "ENSG00000223972"; gene_version "5"; transcript_id "ENST00000450305"; transcript_version "2"; exon_number "5"; gene_name "DDX11L1"; gene_source "havana"; gene_biotype "transcribed_unprocessed_pseudogene"; transcript_name "DDX11L1-201"; transcript_source "havana"; transcript_biotype "transcribed_unprocessed_pseudogene"; exon_id "ENSE00001746346"; exon_version "2"; tag "basic"; transcript_support_level "NA";

The RSeQC provided bed file:

1       11868   14409   ENST00000456328 0       +       11868   14409   0       3       359,109,1189,   0,744,1352,     DDX11L1
1       12009   13670   ENST00000450305 0       +       12009   13670   0       6       48,49,85,78,154,218,    0,169,603,965,1211,1443,        DDX11L1
1       17368   17436   ENST00000619216 0       -       17368   17436   0       1       68,     0,      MIR6859-1
1       14403   29570   ENST00000488147 0       -       14403   29570   0       11      98,34,152,159,198,136,137,147,99,154,37,        0,601,1392,2203,2454,2829,3202,3511,3864,10334,15130,   WASH7P
1       29553   31097   ENST00000473358 0       +       29553   31097   0       3       486,104,122,    0,1010,1422,    MIR1302-2
1       30266   31109   ENST00000469289 0       +       30266   31109   0       2       401,134,        0,709,  MIR1302-2
1       30365   30503   ENST00000607096 0       +       30365   30503   0       1       138,    0,      MIR1302-2
1       34553   36081   ENST00000417324 0       -       34553   36081   0       3       621,205,361,    0,723,1167,     FAM138A
1       35244   36073   ENST00000461467 0       -       35244   36073   0       2       237,353,        0,476,  FAM138A
1       52472   53312   ENST00000606857 0       +       52472   53312   0       1       840,    0,      OR4G4P

Why are the files so different? How would I make a simpler bed file form my GTF file?

$\endgroup$
1
  • $\begingroup$ Please show us an excerpt of your GTF file (a few lines) and i) the output you get from it and ii) the output you want to get from it. As it stands, I don't understand what you mean by "simpler bed file" nor why the files are different (can't know that if we don't know what data you are starting from). $\endgroup$
    – terdon
    Feb 21 '19 at 16:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.