#' lightPlot
#'
#' Function to plot a particular region from a single chromosomes for a sample
#' @param cnv_file File to plot. [Required]
#' @param from Start coordinate of plotting window [Default from=2950000 ]
#' @param to End coordinate of plotting window [Default to=3400000]
#' @param chrom Specify the chromosome to plot [Default 'X']
#' @param ylim Adjust y axis limits in plot [Default ylim=c(-5,5)]
#' @param tick Adjust spacing of ticks [Defualt 1Mb]
#' @param bp1 Draw a vertical line at given coordinate to mark bp1
#' @param bp2 Draw a vertical line at given coordinate to mark bp2
#' @param title Plot title
#' @import scales ggplot2 dplyr RColorBrewer plyr RcppRoll
#' @keywords plot region
#' @export
#' @examples regionPlot(cnv_file="data/w500/test.window-500.cnv", from=3050000, to=3450000, chrom="X", ylim=c(-7,7), bp1=3129368,bp2=3352041, tick=100000, title="222Kb DEL on X")
lightPlot <- function(cnv_file, from=NULL, to=NULL, chrom=NULL, ylim=c(-5,5),
tick=100000, bp1=NULL, bp2=NULL, title=NULL, position=NULL) {
if (missing(cnv_file) ) {
if (!missing(position) || missing(from) || missing(to) || missing(chrom)){
stop("Required arguments: `cnv_file` `from`, `to`, `chrom`", call. = FALSE)
}
}
if (!missing(position)){
vars <- unlist(strsplit(position, ":|-"))
chrom <- vars[1]
bp1 <- as.integer(vars[2])
bp2 <- as.integer(vars[3])
length <- bp2 - bp1
if(length>2e5){
from <- bp1 - 50000
to <- bp2 + 50000
tick = 5e5
rounder = 1e5
} else if(length>1e4) {
from <- bp1 - 10000
to <- bp2 + 10000
tick = 1e4
rounder = 1e4
}
else{
from <- bp1 - 5000
to <- bp2 + 5000
tick = 1e3
rounder = 1e4
}
}
# cat("Processing", cnv_file, "\n")
# cat("Specified region", from, "-", to, "on", chrom, "\n")
if(to > 2500000 && from < 3750000 && chrom == "X"){
# cat("Specified Notch locus\n")
notch <- TRUE
} else {
notch <- FALSE
}
if ( !missing(bp1) && !missing(bp2) || !missing(position) ){
# cat("Drawing lines for breakpoints: bp1:", bp1, " bp2:", bp2, "\n")
draw_bps <- TRUE
} else {
draw_bps <- FALSE
}
base = basename(cnv_file)
parts <- strsplit(base, "[.]")[[1]]
sample <- parts[1]
# cat("Light-plotting region ", chrom, ":", from, "-", to, " for ", sample, "\n", sep="")
read_file_in <- data.table::fread(cnv_file, header = T, select=c('chromosome', 'position', 'log2'))
clean_file <- cleanR(read_file_in, region=T)
if (missing(title)) title <- paste0(sample, chrom)
from <- plyr::round_any(from, rounder, floor)
to <- plyr::round_any(to, rounder, ceiling)
region <- clean_file %>%
dplyr::filter(chromosome == chrom) %>%
dplyr::filter(position >= from,
position <= to) %>%
dplyr::mutate(mavlog2 = RcppRoll::roll_mean(log2, 10, fill=0))
ylim <- plyr::round_any(max(abs(region$mavlog2)), 1, ceiling)
if (ylim < 3)
ylim = 3
ylim <- c(-ylim, ylim)
p <- ggplot(region, aes(x = position, y = mavlog2))
p <- p + geom_line(aes(colour = chromosome), size=2, show.legend = FALSE)
p <- p + scale_x_continuous("Mb", expand = c(0.001, 0.001))
p <- p + scale_y_continuous("Log2 FC ratio", limits=c(min(ylim),max(ylim)), expand = c(0.02, 0.02), breaks = seq(min(ylim), max(ylim),by=1))
p <- p + cleanTheme() +
theme(axis.text.x = element_text(angle = 90, hjust=1,vjust = 0.5))
# if(notch){
# kirreStart = 2740384
# dunceEnd = 3343000
# if (kirreStart < from){
# kirreStart = from
# }
# if (dunceEnd > to){
# dunceEnd = to
# }
# p <- p + annotate("rect", xmin=kirreStart, xmax=3134000, ymin=(min(ylim)+0.5), ymax=min(ylim), alpha=.75, fill="#CFAEAEFE")
# p <- p + annotate("rect", xmin=3134000, xmax=3172000, ymin=(min(ylim)+0.5), ymax=min(ylim), alpha=.75, fill="#8FBD80FE")
# p <- p + annotate("rect", xmin=3176000, xmax=dunceEnd, ymin=(min(ylim)+0.5), ymax=min(ylim), alpha=.75, fill="#A9D0DEFE")
#
# p <- p + annotate("text", x = kirreStart+100000, y = (min(ylim)+0.25), label="Kirre", size=6)
# p <- p + annotate("text", x = 3153000, y = (min(ylim)+0.25), label="Notch", size=6)
# p <- p + annotate("text", x = dunceEnd-100000, y = (min(ylim)+0.25), label="Dunce", size=6)
# }
if(draw_bps){
p <- p + geom_vline(xintercept = bp1, colour="slateblue", alpha=.7, linetype="dotted")
p <- p + geom_vline(xintercept = bp2, colour="slateblue", alpha=.7, linetype="dotted")
}
# scale_bar_start <- (from+tick)
# scale_bar_end <- (scale_bar_start + tick)
#
# scale_text <- paste(tick/1000000, "Mb")
# p <- p + annotate("rect", xmin=scale_bar_start, xmax=scale_bar_end, ymin=(min(ylim)+0.6), ymax=(min(ylim)+0.7))
# p <- p + annotate("text", x=(scale_bar_start + (tick/2)), y = (min(ylim)+0.9), label=scale_text, size=8)
p <- p + ggtitle(title)
p
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.