makeGRangesFromDataFrame {GenomicRanges} | R Documentation |
makeGRangesFromDataFrame
finds the fields in the input that describe
genomic ranges and returns them as a GRanges object.
For convenience, coercing a data.frame or DataFrame
df
into a GRanges object is supported and does
makeGRangesFromDataFrame(df, keep.extra.columns=TRUE)
makeGRangesFromDataFrame(df, keep.extra.columns=FALSE, ignore.strand=FALSE, seqinfo=NULL, seqnames.field=c("seqnames", "chr", "chrom"), start.field=c("start", "chromStart"), end.field=c("end", "chromEnd", "stop", "chromStop"), strand.field="strand", starts.in.df.are.0based=FALSE)
df |
A data.frame or DataFrame object. |
keep.extra.columns |
|
ignore.strand |
|
seqinfo |
Either |
seqnames.field |
A character vector of recognized names for the column in |
start.field |
A character vector of recognized names for the column in |
end.field |
A character vector of recognized names for the column in |
strand.field |
A character vector of recognized names for the column in |
starts.in.df.are.0based |
|
A GRanges object with one element per row in the input.
If the seqinfo
argument was supplied, the returned object will
have exactly the seqlevels specified in seqinfo
and in the same
order.
If df
has non-automatic row names (i.e. rownames(df)
is
not NULL
or seq_len(nrow(df))
), then they will be used to
set the names of the returned GRanges object.
H. Pages, based on a proposal by Kasper Daniel Hansen
GRanges objects.
Seqinfo objects.
The makeGRangesListFromFeatureFragments
function
for making a GRangesList object from a list of fragmented
features.
The getTable
function in the
rtracklayer package for an R interface to the UCSC
Table Browser.
DataFrame objects in the IRanges package.
df <- data.frame(chr="chr1", start=11:13, end=12:14, strand=c("+","-","+"), score=1:3) makeGRangesFromDataFrame(df) gr <- makeGRangesFromDataFrame(df, keep.extra.columns=TRUE) gr2 <- as(df, "GRanges") # equivalent to the above stopifnot(identical(gr, gr2)) makeGRangesFromDataFrame(df, ignore.strand=TRUE) makeGRangesFromDataFrame(df, keep.extra.columns=TRUE, ignore.strand=TRUE) makeGRangesFromDataFrame(df, seqinfo=paste0("chr", 4:1)) makeGRangesFromDataFrame(df, seqinfo=c(chrM=NA, chr1=500, chrX=100)) makeGRangesFromDataFrame(df, seqinfo=Seqinfo(paste0("chr", 4:1))) if (require(rtracklayer)) { session <- browserSession() genome(session) <- "sacCer2" query <- ucscTableQuery(session, "Most Conserved") df <- getTable(query) ## A common pitfall is to forget that the UCSC Table Browser uses the ## "0-based start" convention: gr0 <- makeGRangesFromDataFrame(df, keep.extra.columns=TRUE) head(gr0) min(start(gr0)) ## The start positions need to be converted into 1-based positions, ## to adhere to the convention used in Bioconductor: gr1 <- makeGRangesFromDataFrame(df, keep.extra.columns=TRUE, starts.in.df.are.0based=TRUE) head(gr1) }