plotGeneSummary: Plot JCC score summary for a given gene

In csoneson/jcc: Calculate JCC Scores

Plot JCC score summary for a given gene

Description

This function generates a multi-panel plot, including (i) the observed coverage and gene model, (ii) the estimated transcript abundances and (iii) the observed and predicted junction coverages within the gene.

Usage

plotGeneSummary(gtf, junctionCovs, useGene, bwFile, txQuants,
  txCol = "transcript_id", geneCol = "gene_id", exonCol = "exon_id")

Arguments

gtf	Path to gtf file with genomic features. Preferably in Ensembl format.
junctionCovs	A data.frame with junction coverage information, output from calculateJCCScores.
useGene	A character string giving the name of the gene to plot.
bwFile	Path to a bigwig file generated from a genome alignment BAM file.
txQuants	A data.frame with estimated transcript abundances, output from combineCoverages.
txCol	The column in the gtf file that contains the transcript ID.
geneCol	The column in the gtf file that contains the gene ID.
exonCol	The column in the gtf file that contains the exon ID.

Value

A ggplot object with a multi-panel plot. The size is optimized for display in a device of size approximately 12in wide and 10in high.

Author(s)

Charlotte Soneson

References

Soneson C, Love MI, Patro R, Hussain S, Malhotra D, Robinson MD: A junction coverage compatibility score to quantify the reliability of transcript abundance estimates and annotation catalogs. bioRxiv doi:10.1101/378539 (2018)

Examples

## Not run: 
gtf <- system.file("extdata/Homo_sapiens.GRCh38.90.chr22.gtf.gz",
                   package = "jcc")
bam <- system.file("extdata/reads.chr22.bam", package = "jcc")
biasMod <- fitAlpineBiasModel(gtf = gtf, bam = bam,
                              organism = "Homo_sapiens",
                              genome = Hsapiens, genomeVersion = "GRCh38",
                              version = 90, minLength = 230,
                              maxLength = 7000, minCount = 10,
                              maxCount = 10000, subsample = TRUE,
                              nbrSubsample = 30, seed = 1, minSize = NULL,
                              maxSize = 220, verbose = TRUE)
tx2gene <- readRDS(system.file("extdata/tx2gene.sub.rds", package = "jcc"))
predCovProfiles <- predictTxCoverage(biasModel = biasMod$biasModel,
                                     exonsByTx = biasMod$exonsByTx,
                                     bam = bam, tx2gene = tx2gene,
                                     genome = Hsapiens,
                                     genes = c("ENSG00000070371",
                                               "ENSG00000093010"),
                                     nCores = 1, verbose = TRUE)
txQuants <- readRDS(system.file("extdata/quant.sub.rds", package = "jcc"))
txsc <- scaleTxCoverages(txCoverageProfiles = predCovProfiles,
                         txQuants = txQuants, tx2gene = tx2gene,
                         strandSpecific = TRUE, methodName = "Salmon",
                         verbose = TRUE)
jcov <- read.delim(system.file("extdata/sub.SJ.out.tab", package = "jcc"),
                   header = FALSE, as.is = TRUE) %>%
 setNames(c("seqnames", "start", "end", "strand", "motif", "annot",
           "uniqreads", "mmreads", "maxoverhang")) %>%
 dplyr::mutate(strand = replace(strand, strand == 1, "+")) %>%
 dplyr::mutate(strand = replace(strand, strand == 2, "-")) %>%
 dplyr::select(seqnames, start, end, strand, uniqreads, mmreads) %>%
 dplyr::mutate(seqnames = as.character(seqnames))
combCov <- combineCoverages(junctionCounts = jcov,
                            junctionPredCovs = txsc$junctionPredCovs,
                            txQuants = txsc$txQuants)
jcc <- calculateJCCScores(junctionCovs = combCov$junctionCovs,
                          geneQuants = combCov$geneQuants)
bwFile <- system.file("extdata/reads.chr22.bw", package = "jcc")
plotGeneSummary(gtf = gtf, junctionCovs = jcc$junctionCovs,
                useGene = "ENSG00000070371", bwFile = bwFile,
                txQuants = combCov$txQuants)

## End(Not run)

csoneson/jcc documentation built on March 26, 2024, 1:24 a.m.