Description Usage Author(s) See Also Examples
TODO
1 2 3 4 |
H. Pagès
This man page is part of the SplicingGraphs package.
Please see ?`SplicingGraphs-package`
for an overview of the
package and for an index of its man pages.
Other topics related to this man page and documented in other packages:
The GRangesList class defined in the GenomicRanges package.
The GAlignments and GAlignmentPairs classes defined in the GenomicAlignments package.
The makeTxDbFromGFF
function and the TxDb class
defined in the GenomicFeatures package.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 | ## ---------------------------------------------------------------------
## A. LOAD THE TOY GENE MODEL AS A TxDb OBJECT AND PLOT IT
## ---------------------------------------------------------------------
toy_genes_gff()
## Note that you can display the content of the file with:
cat(readLines(toy_genes_gff()), sep="\n")
library(GenomicFeatures)
suppressWarnings(
txdb <- makeTxDbFromGFF(toy_genes_gff())
)
## Plot all the transcripts in the gene model:
plotTranscripts(txdb)
## ---------------------------------------------------------------------
## B. LOAD THE TOY READS AS A GAlignments OBJECT AND PLOT THEM
## ---------------------------------------------------------------------
## The reads are single-end reads. They are assumed to come from an
## RNA-seq experiment and to have been aligned to the exact same
## reference genome that the above toy gene model is based on.
toy_reads_sam()
toy_reads_bam()
gal <- readGAlignments(toy_reads_bam(), use.names=TRUE)
plotTranscripts(txdb, reads=gal)
plotTranscripts(txdb, reads=gal, from=1, to=320)
## ---------------------------------------------------------------------
## C. FIND THE OVERLAPS BETWEEN THE TOY READS AND THE TOY GENE MODEL
## ---------------------------------------------------------------------
grl <- grglist(gal, order.as.in.query=TRUE)
ex_by_tx <- exonsBy(txdb, by="tx", use.names=TRUE)
## Most of the times the RNA-seq protocol is unstranded so the strand
## reported in the BAM file for each alignment is meaningless. Thus we
## should call findOverlaps() with 'ignore.strand=TRUE':
ov0 <- findOverlaps(grl, ex_by_tx, ignore.strand=TRUE)
## Sort and put the overlaps in a data.frame to make them easier to
## read:
ov0 <- sort(ov0)
df0 <- data.frame(QNAME=names(grl)[queryHits(ov0)],
tx_id=names(ex_by_tx)[subjectHits(ov0)],
stringsAsFactors=FALSE)
head(df0)
## These overlaps have been manually checked and included in the
## SplicingGraphs package. They can be loaded with the toy_overlaps()
## helper:
toy_ov <- toy_overlaps()
head(toy_ov)
stopifnot(identical(df0, toy_ov[ , 1:2]))
## ---------------------------------------------------------------------
## D. DETECT THE OVERLAPS THAT ARE COMPATIBLE WITH THE GENE MODEL
## ---------------------------------------------------------------------
## First we encode the overlaps:
ovenc0 <- encodeOverlaps(grl, ex_by_tx, hits=ov0,
flip.query.if.wrong.strand=TRUE)
ovenc0
## Each encoding tells us whether the corresponding overlap is
## compatible or not with the gene model:
ov0_is_comp <- isCompatibleWithSplicing(ovenc0)
head(ov0_is_comp)
## Overlap compatibility has also been manually checked and included in
## the table returned by toy_overlaps():
stopifnot(identical(ov0_is_comp, toy_ov[ , 3]))
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.