1

我正在尝试使用Gviz'splotTracks函数生成轨迹图。

例如，我试图从gencode.v25.primary_assembly.annotation.gtf绘制这三个成绩单：

"ENST00000352151.9","ENST00000550488.5","ENST00000335154.9"

（全部来自基因ENSG00000161791.13：）

因此，我尝试使用仅包含以下三个成绩单记录的GeneRegionTrackdata.frame进行调用：GTF

require(rtracklayer)
#gtf.fn is the gencode.v25.primary_assembly.annotation.gtf file
gtf.df <- as.data.frame(import(gtf.fn))
my.transcripts.gtf.df <- dplyr::filter(gtf.df, transcript_id %in% c("ENST00000352151.9","ENST00000550488.5","ENST00000335154.9"))


> dim(my.transcripts.gtf.df)
[1] 171  27

并且（天真地）GeneRegionTrack尝试my.transcripts.gtf.df：

require(Gviz)
gene.region.track <-GeneRegionTrack(range=my.transcripts.gtf.df,genome="hg38",chromosome=my.transcripts.gtf.df$seqnames[1],name="Gene Model")
ideogram.track <- IdeogramTrack(genome=genome,chromosome=as.character(t.df$chromosome[1]))
plotTracks(list(ideogram.track,axis.track,gene.region.track))

我有点乱： 在此处输入图像描述

那是因为它GTF比 CDS/exon/UTR 记录有更多的记录，这将使plotTracks更有条理。

我试过了：

my.transcripts.gtf.df <- dplyr::filter(gtf.df, transcript_id %in% c("ENST00000352151.9","ENST00000550488.5","ENST00000335154.9"), type %in% c("CDS","exon","UTR"))
gene.region.track <-GeneRegionTrack(range=my.transcripts.gtf.df,genome="hg38",chromosome=my.transcripts.gtf.df$seqnames[1],name="Gene Model")
plotTracks(list(ideogram.track,axis.track,gene.region.track))

我得到： 在此处输入图像描述

这仍然是一团糟。但我注意到的是 3'UTR 似乎与所有外显子具有相同的宽度。

因此，我尝试处理 myGTF以使其Gviz与其小插图中的内容兼容：

determineGtfFeature <- function(t.gtf.df){
  if(length(which(unique(t.gtf.df$type) == "CDS")) > 0){
    res.gtf.df <- dplyr::filter(t.gtf.df,type =="CDS")
    res.gtf.df$type <- NA
    res.gtf.df$type <- "protein_coding"
    if(length(which(unique(t.gtf.df$type) == "UTR")) > 0){
      utr.gtf.df <- dplyr::filter(t.gtf.df,type == "UTR")
      utr.gtf.df <- utr.gtf.df[order(utr.gtf.df$start),]
      if(utr.gtf.df$strand[1] == "+"){
        utr5.idx <- which(utr.gtf.df$end < min(res.gtf.df$start))
        utr3.idx <- which(utr.gtf.df$start > max(res.gtf.df$end))
      } else{
        utr3.idx <- which(utr.gtf.df$end < min(res.gtf.df$start))
        utr5.idx <- which(utr.gtf.df$start > max(res.gtf.df$end))
      }
      utr.gtf.df$type <- NA
      utr.gtf.df$type[utr5.idx] <- "utr5"
      utr.gtf.df$type[utr3.idx] <- "utr3"
      res.gtf.df <- rbind(res.gtf.df,utr.gtf.df)
    }
  } else{
     if(grepl("pseudogene",t.gtf.df$transcript_type)==T){
       res.gtf.df <- dplyr::filter(t.gtf.df,type =="exon")
       res.gtf.df$type <- "pseudogene"
     } else{
       res.gtf.df <- dplyr::filter(t.gtf.df,type =="exon")
       res.gtf.df$type <- "lincRNA"
     }
  }
  res.gtf.df$type <- factor(res.gtf.df$type,levels=unique(res.gtf.df$type))
  res.gtf.df <- res.gtf.df[,-which(colnames(res.gtf.df) == "transcript_type")]
  return(res.gtf.df)
}

接着：

my.transcripts.gtf.df <- do.call(rbind,lapply(c("ENST00000352151.9","ENST00000550488.5","ENST00000335154.9"),function(t)
  determineGtfFeature(dplyr::select(dplyr::filter(gtf.df, transcript_id == t),seqnames,start,end,width,strand,type,gene_id,exon_id,transcript_id,gene_name,transcript_type))))
colnames(my.transcripts.gtf.df)[c(1,6,7,8,9,10)] <- c("chromosome","feature","gene","exon","transcript","symbol")
tid.starts <- sort(sapply(unique(my.transcripts.gtf.df$transcript),function(t) min(dplyr::filter(my.transcripts.gtf.df,transcript == t)$start)))
my.transcripts.gtf.df$transcript.index <- sapply(my.transcripts.gtf.df$transcript,function(t) which(names(my.transcripts.gtf.df) == t))
my.transcripts.gtf.df <- my.transcripts.gtf.df[order(my.transcripts.gtf.df$transcript.index,my.transcripts.gtf.df$start),]
my.transcripts.gtf.df <- my.transcripts.gtf.df[,-which(colnames(my.transcripts.gtf.df) == "transcript.index")]
my.transcripts.gtf.df$gene <- factor(my.transcripts.gtf.df$gene,levels=unique(my.transcripts.gtf.df$gene))
my.transcripts.gtf.df$exon <- factor(my.transcripts.gtf.df$exon,levels=unique(my.transcripts.gtf.df$exon))
my.transcripts.gtf.df$transcript <- factor(my.transcripts.gtf.df$transcript,levels=unique(my.transcripts.gtf.df$transcript))
my.transcripts.gtf.df$symbol <- factor(my.transcripts.gtf.df$symbol,levels=unique(my.transcripts.gtf.df$symbol))

所以：

gene.region.track <- GeneRegionTrack(range=my.transcripts.gtf.df,genome=genome,chromosome=as.character(t.df$chromosome[1]),name="Gene Model")
plotTracks(list(ideogram.track,axis.track,gene.region.track))

现在给出： 在此处输入图像描述

只有utr3中间的成绩单具有适当的较小宽度，而其他两个则没有。

据我检查my.transcripts.gtf.df是兼容的：

 data(geneModels)

用在Gviz's小插图中。

所以我的问题是：

  1. 有没有一种简单的方法可以从标准文件Gviz中绘制请求，而不必经历编辑它的麻烦。可能不是，因此需要与. 那么知道为什么utr3的行为不端吗？transcriptGTFGTFgeneModels data.frame

  2. 哪个参数displayPars(gene.region.track)控制轨道之间的间距？如果我将displayPars(gene.region.track)$stackHeight <- 0.1轨道设置得更窄，从而在它们之间留下很大的空间，我想缩小那个空间。

4

0 回答 0