Hi Dr.Ivanek.I want to plot a Gene structure using GeneRegionTrack and I use the script and test gtf below.I can not plot the noncoding part(5UTR,3UTR).But I make it using the Txdb.In the Gviz Vignettes, I find the below sentence:
A nice bonus when building GeneRegionTracks from TxDb objects is that we get additional information about coding and non-coding regions of the transcripts, i.e., coordinates of the 5’ and 3’ UTRs and of the CDS regions. The class’ plotting method can use this information to distinguish between coding and non-coding regions based on the shape of the feature. All coding regions are plotted just as we have seen in the previous examples, whereas the non-coding regions are drawn as slightly thinner boxes. The distinction between coding and non-coding is made on the basis of the object’s feature values in combination with a special display parameter thinBoxFeature that enumerates all feature types that are to be treated as non-coding. Obviously this feature is available to all GeneRegionTracks, not only the ones that were build from TxDb objects. However, the coding information has to be added manually and the default value of the thinBoxFeature parameter may not be sufficient to cover all possible cases. It is up to the user to come up with a complete list of non-coding feature types depending on the source of the data.
But I did not know how to make it using the GeneRegionTracks made by GTF.
Here is the code
# loda package
library(Gviz)
library(rtracklayer)
# load gene annotation
At_gtf <- import.gff("test.gtf")
At_gtf_new <- GRanges(At_gtf,
source = At_gtf$source,
feature = At_gtf$type,
transcript = At_gtf$transcript_id,
gene = At_gtf$gene_id)
# set gene name
gene_name <- "AT5G59340"
# get gene Grange location
gene_region_select <- subset(At_gtf_new,gene == gene_name)
gene_region_select$feature <- as.character(gene_region_select$feature)
# get gene start,end
region_start <- start(gene_region_select) %>% min() - 500
region_end <- end(gene_region_select) %>% max() + 500
region_Chr <- seqnames(gene_region_select)[1] %>% as.character()
region_region_select_total_GR <- GRanges(seqnames = region_Chr,
ranges = IRanges(start = region_start,
end = region_end))
genemodel_track <- GeneRegionTrack(gene_region_select,
collapseTranscripts = "longest",
name = "gene_model",
transcriptAnnotation = "gene")
plotTracks(genemodel_track)
Chr5 Araport11 3UTR 23933338 23933407 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 exon 23933338 23933889 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 CDS 23933408 23933889 . - 2 transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 start_codon 23933887 23933889 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 stop_codon 23934327 23934329 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 CDS 23934327 23934627 . - 0 transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 exon 23934327 23935050 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";
Chr5 Araport11 5UTR 23934628 23935050 . - . transcript_id "AT5G59340.1"; gene_id "AT5G59340";