Cuffdiff merged.gtf to UCSC bed format with gene_id

Question

[ashuseq2015]

Dear all,

I have the following format from Cufflinks/cuffmerge, I would like to convert to a bed format, by merging the transcripts (transcript_id). in other words, I would like to have the files converted to bed file with gene_id only.

the input is:

chr1 Cufflinks exon 295 1580 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000001"; exon_number "1"; oId "CUFF.4.1"; tss_id "TSS1";
chr1 Cufflinks exon 3851 4424 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000001"; exon_number "2"; oId "CUFF.4.1"; tss_id "TSS1";
chr1 Cufflinks exon 7276 7377 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000001"; exon_number "3"; oId "CUFF.4.1"; tss_id "TSS1";
chr1 Cufflinks exon 8527 8720 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000001"; exon_number "4"; oId "CUFF.4.1"; tss_id "TSS1";
chr1 Cufflinks exon 11556 13757 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000001"; exon_number "5"; oId "CUFF.4.1"; tss_id "TSS1";
chr1 Cufflinks exon 1518 1557 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000002"; exon_number "1"; oId "CUFF.4.2"; tss_id "TSS2";
chr1 Cufflinks exon 3851 4424 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000002"; exon_number "2"; oId "CUFF.4.2"; tss_id "TSS2";
chr1 Cufflinks exon 7276 7377 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000002"; exon_number "3"; oId "CUFF.4.2"; tss_id "TSS2";
chr1 Cufflinks exon 8527 8720 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000002"; exon_number "4"; oId "CUFF.4.2"; tss_id "TSS2";
chr1 Cufflinks exon 11556 13757 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000002"; exon_number "5"; oId "CUFF.4.2"; tss_id "TSS2";
chr1 Cufflinks exon 1746 2079 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000003"; exon_number "1"; oId "CUFF.4.3"; tss_id "TSS3";
chr1 Cufflinks exon 3851 4424 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000003"; exon_number "2"; oId "CUFF.4.3"; tss_id "TSS3";
chr1 Cufflinks exon 7276 7377 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000003"; exon_number "3"; oId "CUFF.4.3"; tss_id "TSS3";
chr1 Cufflinks exon 8527 8720 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000003"; exon_number "4"; oId "CUFF.4.3"; tss_id "TSS3";
chr1 Cufflinks exon 11556 13757 . + . gene_id "XLOC_000001"; transcript_id "TCONS_00000003"; exon_number "5"; oId "CUFF.4.3"; tss_id "TSS3";
chr1 Cufflinks exon 31100 33382 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000004"; exon_number "1"; oId "CUFF.5.1"; tss_id "TSS4";
chr1 Cufflinks exon 36218 36411 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000004"; exon_number "2"; oId "CUFF.5.1"; tss_id "TSS4";
chr1 Cufflinks exon 37561 37662 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000004"; exon_number "3"; oId "CUFF.5.1"; tss_id "TSS4";
chr1 Cufflinks exon 40514 41087 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000004"; exon_number "4"; oId "CUFF.5.1"; tss_id "TSS4";
chr1 Cufflinks exon 42859 43146 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000004"; exon_number "5"; oId "CUFF.5.1"; tss_id "TSS4";
chr1 Cufflinks exon 31100 33382 . - . gene_id "XLOC_000002"; transcript_id "TCONS_00000005"; exon_number "1"; oId "CUFF.5.3"; tss_id "TSS5";

the output I would like to have is:

chr1 Cufflinks exon 295 1580 . + . gene_id "XLOC_000001"
chr1 Cufflinks exon 31100 33382 . - . gene_id "XLOC_000002"

with the chr start (column 4) and chr end (column 5), representing start and end of all other isofoms or longest.

Thank you very much for your help!

Regards,
Ashu

[GouthamAtla]

the output I would like to have is:

chr1 Cufflinks exon 295 1580 . + . gene_id "XLOC_000001"
chr1 Cufflinks exon 31100 33382 . - . gene_id "XLOC_000002"

This just shows the first exon of each XLOC. Is that what you want?

Answer 1

[GouthamAtla]

There is a python script dexseq_prepare_annotation.py in DEXSeq package that can be used get the start and end coordinates of a gene (gene_id).

python dexseq_prepare_annotation.py cufflinks.gtf out.gtf
grep aggregate_gene out.gtf > desired_out.gtf

or to a bed file:

grep aggregate_gene out.gtf | sed 's/"//g' | awk -v OFS="\t" '{ print $NF,$4,$5}'​

Output:

XLOC_000001    295    13757
XLOC_000002    31100    43146

[bioinfo17]

This is very informative! thanks