makeTxDbFromUCSC {GenomicFeatures} | R Documentation |
The makeTxDbFromUCSC
function allows the user to make a
TxDb object from transcript annotations available at the
UCSC Genome Browser.
Note that it uses the RMariaDB package internally so make sure that this package is installed.
makeTxDbFromUCSC( genome="hg19", tablename="knownGene", transcript_ids=NULL, circ_seqs=DEFAULT_CIRC_SEQS, url="http://genome.ucsc.edu/cgi-bin/", goldenPath_url="http://hgdownload.cse.ucsc.edu/goldenPath", taxonomyId=NA, miRBaseBuild=NA) supportedUCSCtables(genome="hg19", url="http://genome.ucsc.edu/cgi-bin/") browseUCSCtrack(genome="hg19", tablename="knownGene", url="http://genome.ucsc.edu/cgi-bin/") getChromInfoFromUCSC( genome, goldenPath_url="http://hgdownload.cse.ucsc.edu/goldenPath")
genome |
genome abbreviation used by UCSC and obtained by
|
tablename |
name of the UCSC table containing the transcript
annotations to retrieve. Use the |
transcript_ids |
optionally, only retrieve transcript annotation data for the specified set of transcript ids. If this is used, then the meta information displayed for the resulting TxDb object will say 'Full dataset: no'. Otherwise it will say 'Full dataset: yes'. |
circ_seqs |
a character vector to list out which chromosomes should be marked as circular. |
url,goldenPath_url |
use to specify the location of an alternate UCSC Genome Browser. |
taxonomyId |
By default this value is NA and the organism inferred will be used to look up the correct value for this. But you can use this argument to supply your own valid taxId here. |
miRBaseBuild |
specify the string for the appropriate build
Information from mirbase.db to use for microRNAs. This can be
learned by calling |
makeTxDbFromUCSC
is a convenience function that feeds
data from the UCSC source to the lower level makeTxDb
function.
See ?makeTxDbFromBiomart
for a similar function
that feeds data from a BioMart database.
For makeTxDbFromUCSC
: A TxDb object.
For supportedUCSCtables
: A data frame with 3 columns
(tablename
, track
, and subtrack
) and 1 row
per table known to work with makeTxDbFromUCSC
.
IMPORTANT NOTE: In the returned data frame, the set of tables associated
with a track with subtracks might contain tables that don't exist for the
specified genome.
For getChromInfoFromUCSC
: A data frame with 1 row per chromosome
(or scaffold) and with columns chrom
and length
.
M. Carlson and H. Pagès
makeTxDbFromBiomart
and
makeTxDbFromEnsembl
for making a TxDb object
from other online resources.
makeTxDbFromGRanges
and makeTxDbFromGFF
for making a TxDb object from a GRanges
object, or from a GFF or GTF file.
ucscGenomes
in the rtracklayer
package.
The supportedMiRBaseBuildValues
function for
listing all the possible values for the miRBaseBuild
argument.
The TxDb class.
makeTxDb
for the low-level function used by the
makeTxDbFrom*
functions to make the TxDb object
returned to the user.
## Not run: ## --------------------------------------------------------------------- ## A. BASIC USAGE ## --------------------------------------------------------------------- ## Use ucscGenomes() from the rtracklayer package to display the list of ## genomes available at UCSC: library(rtracklayer) ucscGenomes()[ , "db"] ## Display the list of tables known to work with makeTxDbFromUCSC(): supportedUCSCtables() ## Browse the UCSC track page for a given organism/table: browseUCSCtrack(genome="sacCer3", tablename="ensGene") ## Retrieve a full transcript dataset for Yeast from UCSC: txdb1 <- makeTxDbFromUCSC(genome="sacCer3", tablename="ensGene", circ_seqs="chrM") txdb1 ## Retrieve an incomplete transcript dataset for Mouse from UCSC (only ## transcripts linked to Entrez Gene ID 22290): transcript_ids <- c( "uc009uzf.1", "uc009uzg.1", "uc009uzh.1", "uc009uzi.1", "uc009uzj.1" ) txdb2 <- makeTxDbFromUCSC(genome="mm10", tablename="knownGene", transcript_ids=transcript_ids, circ_seqs="chrM") txdb2 ## --------------------------------------------------------------------- ## B. IMPORTANT NOTE ABOUT supportedUCSCtables() ## --------------------------------------------------------------------- ## In the data frame returned by supportedUCSCtables(), the set of ## tables associated with a track with subtracks might contain tables ## that don't exist for the specified genome: supported_tables <- supportedUCSCtables("hg38") browseUCSCtrack(genome="hg38", tablename="augustusHints") # no such table ## --------------------------------------------------------------------- ## C. RETRIEVING CHROMOSOME INFORMATION ONLY ## --------------------------------------------------------------------- chrominfo <- getChromInfoFromUCSC(genome="hg38") chrominfo ## End(Not run)