This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Plotting a single transcript with ggbio

I'm trying to create the plot below, but with only a single transcript (uc062uto.1) displayed, not all of the variants (see below).

enter image description here

What am I missing?

The code I'm using to make this image:

# simple genomic/transcript map

rm(list = ls()) # clear workspace

Gene <- "HTT"

require(biomaRt)
require(ggbio)
require(Biostrings)
require(ShortRead)
require(BSgenome)
require(GenomicFeatures)

## get gene ranges, strand, based on gene name (biomaRt functions)
mart <- useMart("ENSEMBL_MART_ENSEMBL", host = "useast.ensembl.org")
datasets <- listDatasets(mart)
mart <- useDataset("hsapiens_gene_ensembl",mart)

r1 <- getBM(filters="external_gene_name",
            values = Gene,
            attributes=c("external_gene_name",
                         "chromosome_name",
                         "start_position",
                         "end_position",
                         "strand",
                         "refseq_mrna",
                         "entrezgene",
                         "ucsc",
                         "ensembl_gene_id",
                         "hamap"),
            mart = mart)
r1 <- r1[r1$refseq_mrna != "",]

CHR <- paste("chr",r1$chromosome_name[1], sep = "")
START <- r1$start_position[1]
END <- r1$end_position[1]
if (r1$strand[1] == 1) {
  STRAND = "+"
  } else {
    STRAND = "-"
  }

library(TxDb.Hsapiens.UCSC.hg38.knownGene) # knowngene db for hg38
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene
gr <- GRanges(CHR, IRanges(START, END), strand = STRAND)
p <- autoplot(txdb,
              which = gr,
              ylab = paste(Gene," [",CHR,":",START,"-",END,"]", sep = ""),
              geom = "full",
              coord = "genome",
              layout = "linear",
              strand = STRAND)
p
ggbio genomicfeatures
• 4,597 views • link
written by longoka • longoka

Works like a charm; I think I need to verse myself a bit better on all the functions in the GenomicFeatures package. Thank you!

• 0 views • link
written by longoka • longoka

1 answer

There may be better ways, but one way would be creating your own TxDB object:

library(ggbio)
library(GenomicFeatures)
library(TxDb.Hsapiens.UCSC.hg38.knownGene)
txdb <- TxDb.Hsapiens.UCSC.hg38.knownGene

gr <- subset(transcripts(txdb), tx_name == "uc062uto.1")

tx <- makeTxDbFromUCSC(
  genome="hg38",
  tablename="knownGene",
  transcript_ids="uc062uto.1",
  url="http://genome.ucsc.edu/cgi-bin/",
  goldenPath_url="http://hgdownload.cse.ucsc.edu/goldenPath" )

p <- autoplot(tx,
              which = gr,
              geom = "full",
              coord = "genome",
              layout = "linear" )
p
• 0 views • link
written by h.mon h.mon

Log in to answer this question.

More posts like this