####################################################################
## Author: Gro Nilsen, Knut Liestøl and Ole Christian Lingjærde.
## Maintainer: Gro Nilsen <gronilse@ifi.uio.no>
## License: Artistic 2.0
## Part of the copynumber package
## Reference: Nilsen and Liestøl et al. (2012), BMC Genomics
####################################################################
# Function that plots the frequency of deletions and amplifications given a threshold - by genome og chromosomes
## Input:
### segments: result from segmentation (could be a data frame with the estimates, or the segments data frame). Also possible to specify a data frame with original data
### thres.gain,thres.loss: a vector with threshold(s) to be applied for abberration calling
### pos.unit: unit used for positions (bp,kbp,mbp)
### chrom: a vector with chromosomes to be plotted. If specified the frequencies are plotted with one panel for each chromosome
### layout: number of columns and rows in plot
### ... : other optional plot parameters
## Required by: none
## Requires:
### numericChrom
### checkChrom
### genomeFreq
### chromosomeFreq
### pullOutContent
### getFreqData
plotFreq <- function(segments, thres.gain, thres.loss = -thres.gain, pos.unit = "bp", chrom = NULL, layout = c(1, 1), ...) {
# Check data input:
# Make sure data is a data frame (could be a list if it contains segmentation results or results from winsorize)
data <- segments
if (is.data.frame(data)) {
# Check if segments or data:
data <- getFreqData(data)
} else {
data <- pullOutContent(res = data, what = "estimates")
}
stopifnot(ncol(data) >= 3) # something is missing in input data
# Make sure that chromosomes in data are numeric:
data[, 1] <- numericChrom(data[, 1])
# If chrom is unspecified, the whole genome is plotted. Otherwise, selected chromosomes are plotted
type <- ifelse(is.null(chrom), "genome", "bychrom")
# Check input chrom, alt. get unique chromosomes
chrom <- checkChrom(data = data, segments = NULL, chrom = chrom)
# Check pos.unit input:
if (!pos.unit %in% c("bp", "kbp", "mbp")) {
stop("pos.unit must be one of bp, kbp and mbp", call. = FALSE)
}
# Making sure number of gain thresholds and loss thresholds are the same:
nT <- min(length(thres.gain), length(thres.loss))
thres.gain <- thres.gain[1:nT]
thres.loss <- thres.loss[1:nT]
# plot frequency, either over genome or by chromosomes:
switch(type,
genome = genomeFreq(data, thres.gain, thres.loss, pos.unit, layout, ...),
bychrom = chromosomeFreq(data, thres.gain, thres.loss, pos.unit, chrom, layout, ...)
)
}
# Function that plots frequencies by genome
## Required by: plotFreq
## Requires:
### getFreqPlotParameters
### getGlobal.xlim
### adjustPos
### framedim
### updateFreqParameters
### addToFreqPlot
### addChromlines
### addArmlines
### chromPattern
genomeFreq <- function(data, thres.gain, thres.loss, pos.unit, layout, ...) {
nr <- layout[1]
nc <- layout[2]
cr <- nr * nc
nProbe <- nrow(data)
nT <- length(thres.gain)
# Get plot parameters:
op <- getFreqPlotParameters(type = "genome", nc = nc, nr = nr, thres.gain = thres.gain, thres.loss = thres.loss, ...)
# Set global xlimits if not specified by user:
if (is.null(op$xlim)) {
op$xlim <- getGlobal.xlim(op = op, pos.unit = pos.unit, chrom = unique(data[, 1]))
}
# Adjust positions to be plotted along xaxis; i.e. get global positions, scale according to plot.unit, and get left and right pos for freq.rectangles
# to be plotted (either continuous or 1 probe long):
x <- adjustPos(position = data[, 2], chromosomes = data[, 1], pos.unit = pos.unit, type = "genome", op = op)
xleft <- x$xleft
xright <- x$xright
if (dev.cur() <= 1) { # to make Sweave work
dev.new(width = op$plot.size[1], height = op$plot.size[2], record = TRUE)
}
# Initialize:
row <- 1
clm <- 1
new <- FALSE
# Division of plotting window:
frames <- framedim(nr, nc)
# One plot for each value in thres.gain/thres.loss:
for (t in 1:nT) {
# Frame dimensions for plot t:
fig.t <- c(frames$left[clm], frames$right[clm], frames$bot[row], frames$top[row])
par(fig = fig.t, new = new, oma = c(0, 0, 1, 0), mar = op$mar)
# Calculate the percentage of samples that have estimated copy number larger than thres at given position:
freq.amp <- rowMeans(data[, -c(1:2), drop = FALSE] > thres.gain[t]) * 100
freq.del <- rowMeans(data[, -c(1:2), drop = FALSE] < thres.loss[t]) * 100
# Find default ylimits and at.y (tickmarks):
op <- updateFreqParameters(freq.del, freq.amp, op)
# Empty plot with correct limits
plot(1, 1, type = "n", ylim = op$ylim, xlim = op$xlim, xaxs = "i", main = "", frame.plot = TRUE, yaxt = "n", xaxt = "n", ylab = "", xlab = "")
# Add shifting white/grey pattern to backgroud to separate chromosomes:
chromPattern(pos.unit, op)
# main title for this plot
title(main = op$main[t], line = op$main.line, cex.main = op$cex.main)
# Add axes, labels and percentage lines:
addToFreqPlot(op, type = "genome")
# Plot frequencies:
rect(xleft = xleft, ybottom = 0, xright = xright, ytop = freq.amp, col = op$col.gain, border = op$col.gain)
rect(xleft = xleft, ybottom = 0, xright = xright, ytop = -freq.del, col = op$col.loss, border = op$col.loss)
# Add line at y=0:
abline(h = 0, lty = 1, col = "grey82", lwd = 1.5)
# Add line at y=0:
abline(h = 0, lty = 1, col = "grey82", lwd = 1.5)
# Separate chromosomes by vertical lines:
op$chrom.lty <- 1
addChromlines(data[, 1], xaxis = "pos", unit = pos.unit, cex = op$cex.chrom, op = op)
addArmlines(data[, 1], xaxis = "pos", unit = pos.unit, cex = op$cex.chrom, op = op)
# Get new page, or update column/row:
if (t %% (nr * nc) == 0) {
# Start new plot page (prompted by user)
devAskNewPage(ask = TRUE)
# Reset columns and row in layout:
clm <- 1
row <- 1
new <- FALSE
} else {
# Update column and row index:
if (clm < nc) {
clm <- clm + 1
} else {
clm <- 1
row <- row + 1
} # endif
new <- TRUE
} # endif
} # endfor
} # end function
# Function that plots frequencies by chromosomes (each chrom in separate panel)
## Required by: plotFreq
## Requires:
### getFreqPlotParameters
### adjustPos
### framedim
### updateFreqParameters
### plotIdeogram
### chromMax
### addToFreqPlot
chromosomeFreq <- function(data, thres.gain, thres.loss, pos.unit, chrom, layout, ...) {
nProbe <- nrow(data)
nT <- length(thres.gain)
nChrom <- length(chrom)
# Grid layout
nr <- layout[1]
nc <- layout[2]
# Get plot parameters:
op <- getFreqPlotParameters(type = "bychrom", nc = nc, nr = nr, thres.gain = thres.gain, thres.loss = thres.loss, chrom = chrom, ...)
# Margins for entire plot in window:
if (all(op$title == "")) {
oma <- c(0, 0, 0, 0)
} else {
oma <- c(0, 0, 1, 0)
}
mar <- c(0.2, 0.2, 0.3, 0.2)
# Adjust positions to be plotted along xaxis; i.e. scale according to plot.unit, and get left and right pos for freq.rectangles to be plotted (either continuous
# or 1 probe long):
x <- adjustPos(position = data[, 2], chromosomes = data[, 1], pos.unit = pos.unit, type = "chromosome", op = op)
xleft <- x$xleft
xright <- x$xright
# Divide the plotting window by the function "framedim":
frames <- framedim(nr, nc)
# make separate plots for each value of thres.gain/thres.loss
for (t in 1:nT) {
if (dev.cur() <= 1) { # to make Sweave work
dev.new(width = op$plot.size[1], height = op$plot.size[2], record = TRUE)
}
# Initialize row and column index:
row <- 1
clm <- 1
new <- FALSE
# For each probe, get percentage of samples amplified and deleted for this thres:
freq.amp <- rowMeans(data[, -c(1:2), drop = FALSE] > thres.gain[t]) * 100
freq.del <- rowMeans(data[, -c(1:2), drop = FALSE] < thres.loss[t]) * 100
# Find default ylimits and at.y (tickmarks):
use <- which(data[, 1] %in% chrom)
op <- updateFreqParameters(freq.del[use], freq.amp[use], op)
# Make separate plots for each chromosome:
for (c in 1:nChrom) {
# Frame dimensions for plot c:
fig.c <- c(frames$left[clm], frames$right[clm], frames$bot[row], frames$top[row])
par(fig = fig.c, new = new, oma = oma, mar = mar)
# Make list with frame dimensions:
frame.c <- list(left = frames$left[clm], right = frames$right[clm], bot = frames$bot[row], top = frames$top[row])
# Select relevant chromosome number
k <- chrom[c]
# Pick out frequencies for this chromsome
ind.c <- which(data[, 1] == k)
freqamp.c <- freq.amp[ind.c]
freqdel.c <- freq.del[ind.c]
xlim <- c(0, max(xright[ind.c]))
# Plot ideogram at below frequencies:
if (op$plot.ideo) {
# Ideogram frame:
ideo.frame <- frame.c
ideo.frame$top <- frame.c$bot + (frame.c$top - frame.c$bot) * op$ideo.frac
par(fig = unlist(ideo.frame), new = new, mar = op$mar.i)
# Plot ideogram and get maximum probe position in ideogram:
plotIdeogram(chrom = k, cyto.text = op$cyto.text, cyto.data = op$assembly, cex = op$cex.cytotext, unit = op$plot.unit)
# Get maximum position for this chromosome:
xmaxI <- chromMax(chrom = k, cyto.data = op$assembly, pos.unit = op$plot.unit)
xlim <- c(0, xmaxI)
new <- TRUE
}
# Freq.plot-dimensions:
frame.c$bot <- frame.c$bot + (frame.c$top - frame.c$bot) * op$ideo.frac
par(fig = unlist(frame.c), new = new, mar = op$mar)
# Limits:
if (!is.null(op$xlim)) {
xlim <- op$xlim
}
# Empty plot:
plot(1, 1, type = "n", ylim = op$ylim, xlim = xlim, xaxs = "i", main = op$main[c], frame.plot = FALSE, yaxt = "n", xaxt = "n", cex.main = op$cex.main, ylab = "", xlab = "")
# Add axes, labels and percentageLines:
chrom.op <- op
chrom.op$xlim <- xlim
addToFreqPlot(chrom.op, type = "bychrom")
# Plot frequencies as rectangles
rect(xleft = xleft[ind.c], ybottom = 0, xright = xright[ind.c], ytop = freqamp.c, col = op$col.gain, border = op$col.gain)
rect(xleft = xleft[ind.c], ybottom = 0, xright = xright[ind.c], ytop = -freqdel.c, col = op$col.loss, border = op$col.loss)
# Add line at y=0 and x=0
abline(h = 0, lty = 1, col = "grey90")
abline(v = 0)
# If page is full; start plotting on new page
if (c %% (nr * nc) == 0 && c != nChrom) {
# Add main title to page:
title(op$title[t], outer = TRUE)
# Start new page when prompted by user:
devAskNewPage(ask = TRUE)
# Reset columns and row in layout:
clm <- 1
row <- 1
new <- FALSE
} else {
# Update column and row index:
if (clm < nc) {
clm <- clm + 1
} else {
clm <- 1
row <- row + 1
} # endif
new <- TRUE
} # endif
} # endfor
# Add main title to page:
title(op$title[t], outer = TRUE)
if (t != nT) {
devAskNewPage(ask = TRUE)
}
} # endfor
} # endfunction
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.