Entering edit mode
                    8.1 years ago
        longoka
        
    
        ▴
    
    40
    I'm trying to create the plot below, but with only a single transcript (uc062uto.1) displayed, not all of the variants (see below).

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
                    
                
                
Works like a charm; I think I need to verse myself a bit better on all the functions in the GenomicFeatures package. Thank you!