plotData {DRIMSeq} | R Documentation |
Plot data summary
plotData(x, ...) ## S4 method for signature 'dmDSdata' plotData(x) ## S4 method for signature 'dmSQTLdata' plotData(x, plot_type = "features")
x |
|
... |
Other parameters that can be defined by methods using this generic. |
plot_type |
Character specifying which type of histogram to plot. Possible
values |
Returns a ggplot
object and can be further modified, for
example, using theme()
. Plots a histogram of the number of features
per gene. Additionally, for dmSQTLdata
object, plots a
histogram of the number of SNPs per gene and a histogram of the number of
unique SNPs (blocks) per gene.
Malgorzata Nowicka
plotPrecision
, plotProportions
,
plotPValues
# -------------------------------------------------------------------------- # Create dmDSdata object # -------------------------------------------------------------------------- ## Get kallisto transcript counts from the 'PasillaTranscriptExpr' package library(PasillaTranscriptExpr) data_dir <- system.file("extdata", package = "PasillaTranscriptExpr") ## Load metadata pasilla_metadata <- read.table(file.path(data_dir, "metadata.txt"), header = TRUE, as.is = TRUE) ## Load counts pasilla_counts <- read.table(file.path(data_dir, "counts.txt"), header = TRUE, as.is = TRUE) ## Create a pasilla_samples data frame pasilla_samples <- data.frame(sample_id = pasilla_metadata$SampleName, group = pasilla_metadata$condition) levels(pasilla_samples$group) ## Create a dmDSdata object d <- dmDSdata(counts = pasilla_counts, samples = pasilla_samples) ## Use a subset of genes, which is defined in the following file gene_id_subset <- readLines(file.path(data_dir, "gene_id_subset.txt")) d <- d[names(d) %in% gene_id_subset, ] # -------------------------------------------------------------------------- # Differential transcript usage analysis - simple two group comparison # -------------------------------------------------------------------------- ## Filtering ## Check what is the minimal number of replicates per condition table(samples(d)$group) d <- dmFilter(d, min_samps_gene_expr = 7, min_samps_feature_expr = 3, min_gene_expr = 10, min_feature_expr = 10) plotData(d) # -------------------------------------------------------------------------- # Create dmSQTLdata object # -------------------------------------------------------------------------- # Use subsets of data defined in the GeuvadisTranscriptExpr package library(GeuvadisTranscriptExpr) geuv_counts <- GeuvadisTranscriptExpr::counts geuv_genotypes <- GeuvadisTranscriptExpr::genotypes geuv_gene_ranges <- GeuvadisTranscriptExpr::gene_ranges geuv_snp_ranges <- GeuvadisTranscriptExpr::snp_ranges colnames(geuv_counts)[c(1,2)] <- c("feature_id", "gene_id") colnames(geuv_genotypes)[4] <- "snp_id" geuv_samples <- data.frame(sample_id = colnames(geuv_counts)[-c(1,2)]) d <- dmSQTLdata(counts = geuv_counts, gene_ranges = geuv_gene_ranges, genotypes = geuv_genotypes, snp_ranges = geuv_snp_ranges, samples = geuv_samples, window = 5e3) # -------------------------------------------------------------------------- # sQTL analysis - simple group comparison # -------------------------------------------------------------------------- ## Filtering d <- dmFilter(d, min_samps_gene_expr = 70, min_samps_feature_expr = 5, minor_allele_freq = 5, min_gene_expr = 10, min_feature_expr = 10) plotData(d)