% \VignetteIndexEntry{Make genomic plots}
% \VignetteDepends{circlize}
% \VignetteKeywords{circos R}
% \VignetteKeywords{circlize}
% \VignettePackage{circlize}
%\VignetteCompiler{knitr}
%\VignetteEngine{knitr::knitr}

<<echo = FALSE>>=
library(knitr)
opts_chunk$set(fig.pos = "", fig.align = "center")

library(circlize)
circos.genomicInitialize = function(...) {
    circos.par(unit.circle.segments = 200)
    circlize::circos.genomicInitialize(...)
}

circos.initializeWithIdeogram = function(...) {
    circos.par(unit.circle.segments = 200)
    circlize::circos.initializeWithIdeogram(...)
}
@

\documentclass{article}

\usepackage[T1]{fontenc}
\usepackage[sc]{mathpazo}
\renewcommand{\sfdefault}{lmss}
\renewcommand{\ttdefault}{lmtt}

\setcounter{tocdepth}{2}

\title{An introduction of making genomic plots with \textbf{circlize}}

\author{Zuguang Gu {\tt <z.gu@dkfz.de>}}
\usepackage[margin=1in, a4paper]{geometry}

\usepackage{hyperref}

\begin{document}

\maketitle 


Since circos plots are mostly used in genomic research, the \textbf{circlize} package 
particularly provides functions which focus on genomic plots. These functions 
are synonymous to the basic circos graphic functions but expect special format of input data:

\begin{itemize}
  \item {\tt circos.genomicTrackPlotRegion}: create a new track and add graphics.
  \item {\tt circos.genomicPoints}: low-level function, add points
  \item {\tt circos.genomicLines}: low-level function, add lines
  \item {\tt circos.genomicRect}: low-level function, add rectangles
  \item {\tt circos.genomicText}: low-level function, add text
  \item {\tt circos.genomicLink}: add links
\end{itemize}

The genomic functions are implemented by basic circos functions (e.g. {\tt circos.trackPlotRegion}, {\tt circos.points}),
thus, you can customize your own plots by both genomic functions and basic circos functions.

\section{Input data}

Genomic circos functions expect input data as a data frame or a list of data frames 
in which there are at least three columns. The first column is genomic category 
(e.g. chromosome), the second column is the start positions in the genomic category and 
the third column is the end positions. Following columns are optional where numeric values 
or other related values are stored. Such data structure is known as {\it BED} format 
and is broadly used in genomic research.

\textbf{circlize} provides a simple function {\tt generateRandomBed} which can generate 
random genomic data. Positions are uniformly generated from human genome and the number 
of regions on chromosomes approximately proportional to the length of chromosomes. 
In the function, {\tt nr} and {\tt nc} are number of rows and numeric columns that 
users want. Please note {\tt nr} are not exactly the same as the number of rows which 
are returned by the function. {\tt fun} is a self-defined function to generate random values.

<<>>=
library(circlize)
set.seed(999)

bed = generateRandomBed()
head(bed)
bed = generateRandomBed(nr = 200, nc = 4)
nrow(bed)
bed = generateRandomBed(nc = 2, fun = function(k) runif(k))
head(bed)
@

\section{Initialize the layout}

\subsection{Initialize with cytoband data}

\subsubsection{Basic usage}

Similar as general circos plots, the first step is to initialize the plot with genomic categories.
In most situations, genomic categories are measured by chromosomes. The easiest way is to used
{\tt circos.initializeWithIdeogram} (figure \ref{fig:genomic_initialize_ideogram} A):

<<genomic_initialize_ideogram_1, eval = FALSE, echo = 1:2>>=
par(mar = c(1, 1, 1, 1))
circos.initializeWithIdeogram()
text(0, 0, "default", cex = 0.7)
text(-0.9, 0.9, "A", cex = 1.5)
@

By default, the function will initialize the plot with cytoband data of hg19. You can also 
use your own cytoband data by specifying the path of your cytoband file (no matter 
compressed or not) or providing your cytoband data as a data frame.  An example for 
cytoband file is \url{http://hgdownload.soe.ucsc.edu/goldenpath/hg19/database/cytoBand.txt.gz}.

<<echo = 2:6>>=
pdf(NULL)
cytoband.file = paste0(system.file(package = "circlize"), "/extdata/cytoBand.txt")
circos.initializeWithIdeogram(cytoband.file)
cytoband.df = read.table(cytoband.file, colClasses = c("character", "numeric",
    "numeric", "character", "character"), sep = "\t")
circos.initializeWithIdeogram(cytoband.df)
circos.clear()
invisible(dev.off())
@

If you want to read cytoband data from file, please explicitly specify {\tt colClasses} 
arguments and set the class of position columns as {\tt numeric}. The reason is since 
positions are represented as integers, {\tt read.table} would treat those numbers as 
{\tt integer} by default. In initialization of circos plot, \textbf{circlize} needs to 
calculate the summation of all chromosome lengths. The sumation of such large integers 
would throw error of data overflow.

For simple use, users can also specify abbreviation of the species and the function will 
download cytoband file from UCSC server automatically (If it exists in UCSC. As you can guess,
the URL template we use is {\tt http://hgdownload.soe.ucsc.edu/goldenpath/\$species/database/cytoBand.txt.gz}).

<<eval = FALSE>>=
circos.initializeWithIdeogram(species = "hg18")
circos.initializeWithIdeogram(species = "mm10")
@

By default, the function will use all chromosomes which are available in cytoband data to 
initialize the circos plot. Users can also choose a subset of chromosomes by specifying 
{\tt chromosome.index}. This argument is also for 
ordering chromosomes (figure \ref{fig:genomic_initialize_ideogram} B).

<<genomic_initialize_ideogram_2, eval = FALSE, echo = 1>>=
circos.initializeWithIdeogram(chromosome.index = paste0("chr", 10:1))
text(0, 0, "subset of chromosomes", cex = 0.7)
text(-0.9, 0.9, "B", cex = 1.5)
@

\subsubsection{Order chromosomes}

Initialization step is important for circos plot. It controls orders of chromosomes
which are going to be put on the circle. There are several ways to control the order:

\begin{itemize}
  \item If {\tt chromosome.index} is set, the order of {\tt chromosome.index} is taken as 
        order of chromosomes.
  \item If {\tt chromosome.index} is not set and {\tt sort.chr} is set to {\tt TRUE}, chromosomes will be sorted 
      first by numbers then by letters (figure \ref{fig:genomic_initialize_ideogram} C).
  \item If {\tt chromosome.index} and {\tt sort.chr} are not set neither:
  \begin{itemize}
    \item If {\tt cytoband} is provided as a data frame, and if the first column is not a factor, the order of chromosomes would be 
        {\tt unique(cytoband[[1]])} (figure \ref{fig:genomic_initialize_ideogram} D). 
    \item If {\tt cytoband} is provided as a data frame, and if the first column is a factor, 
        the order of chromosome would be {\tt levels(cytoband[[1]])}  (figure \ref{fig:genomic_initialize_ideogram} E). 
    \item If {\tt cytoband} is specified as a file path, or {\tt species} is specified, the order 
        of chromosomes depends on the original order in the source file.
  \end{itemize}
\end{itemize}

<<genomic_initialize_ideogram_3, eval = FALSE, echo = c(1:2, 6:7, 11:12)>>=
cytoband = cytoband.df
circos.initializeWithIdeogram(cytoband, sort.chr = TRUE)
text(0, 0, "read from cytoband df\nsort.chr = TRUE", cex = 0.7)
text(-0.9, 0.9, "C", cex = 1.5)

cytoband = cytoband.df
circos.initializeWithIdeogram(cytoband, sort.chr = FALSE)
text(0, 0, "read from a data frame\nunique(cytoband[[1]])", cex = 0.7)
text(-0.9, 0.9, "D", cex = 1.5)

cytoband[[1]] = factor(cytoband[[1]], levels = paste0("chr", c(22:1, "X", "Y")))
circos.initializeWithIdeogram(cytoband, sort.chr = FALSE)
text(0, 0, "read from cytoband file\nfirst column converted to factor\nlevels = paste0('chr', c(22:1, 'X', 'Y'))", cex = 0.7)
text(-0.9, 0.9, "E", cex = 1.5)
@

\textbf{circlize} provides a function {\tt read.cytoband} which can read/download and process cytoband data.
In fact, {\tt circos.initializeWithIdeogram} calls {\tt read.cytoband} internally (actually, if there is no cytoband 
available, {\tt read.chromInfo} will be called later). Please refer to the
help page of the function for more details.

<<eval=FALSE>>=
cytoband = read.cytoband()
cytoband = read.cytoband(file)
cytoband = read.cytoband(df)
cytoband = read.cytoband(species)
@

\subsubsection{Pre-defined tracks}

After the initialization of the circos plot, the function will additionally create a track where there
are genomic axes and chromosome names, and create another track where there is an ideogram 
(depends on whether cytoband data is available).
{\tt plotType} can be used to control which kind of graphics need to be added 
(figure \ref{fig:genomic_initialize_ideogram} F, G).

<<genomic_initialize_ideogram_4, eval = FALSE, echo = c(1, 4, 8)>>=
circos.initializeWithIdeogram(plotType = c("axis", "labels"))
text(0, 0, "plotType = c('axis', 'labels')", cex = 0.7)
text(-0.9, 0.9, "F", cex = 1.5)
circos.initializeWithIdeogram(plotType = NULL)
text(0, 0, "plotType = NULL", cex = 0.7)
text(-0.9, 0.9, "G", cex = 1.5)

circos.clear()
@

\subsubsection{Other general settings}

Similar as general circos plot, the layout of circos plot can be controlled by {\tt circos.par}
(figure \ref{fig:genomic_initialize_ideogram} H, I).

<<genomic_initialize_ideogram_5, eval = FALSE, echo = c(1:3, 7:9)>>=
circos.par("start.degree" = 90)
circos.initializeWithIdeogram()
circos.clear()
text(0, 0, "'start.degree' = 90", cex = 0.7)
text(-0.9, 0.9, "H", cex = 1.5)

circos.par("gap.degree" = rep(c(2, 4), 12))
circos.initializeWithIdeogram()
circos.clear()
text(0, 0, "'gap.degree' = rep(c(2, 4), 12)", cex = 0.7)
text(-0.9, 0.9, "I", cex = 1.5)
@


\begin{figure}
\centering
\includegraphics[width=\textwidth]{figures/genomic_initialize_ideogram-1.pdf}
\caption{Different ways to initialize genomic circos plot. A) default; B) select subset of chromosomes; C) read from a data frame; D) the first column of the data frame is a factor; E) sort chromosomes; F) do not add ideogram; G) initialize the plot while plot nothing; H) set start degree of the plot; I) set gap degree.}
\label{fig:genomic_initialize_ideogram}
\end{figure}


\subsection{Customize ideogram}

The default behaviour of {\tt circos.initializeWithIdeogram} is both initializing the layout
and adding two tracks. But if {\tt plotType} is set to {\tt NULL}, the circular layout is only 
initialized but nothing is added. It provides possibility for users to completely design
their own style of plots. In the following example, we use different colors to represent 
chromosomes and change the style of chromosome names (figure \ref{fig:genomic_customize_ideogram}).

<<genomic_customize_ideogram, out.width = "0.8\\textwidth", fig.cap = "Customize ideogram.">>=
par(mar = c(1, 1, 1, 1))
circos.initializeWithIdeogram(plotType = NULL)
circos.trackPlotRegion(ylim = c(0, 1), panel.fun = function(x, y) {
    chr = get.cell.meta.data("sector.index")
    xlim = get.cell.meta.data("xlim")
    ylim = get.cell.meta.data("ylim")
    circos.rect(xlim[1], 0, xlim[2], 0.5, col = rand_color(1))
    circos.text(mean(xlim), 0.9, chr, cex = 0.5, facing = "clockwise", 
        niceFacing = TRUE)
}, bg.border = NA)
circos.clear()
@


\subsection{Initialize with general genomic category}

Cytoband data is just a special case of genomic category. {\tt circos.genomicInitialize}
can initialize circular layout with any kind of genomic categories. In fact, 
{\tt circos.initializeWithIdeogram} is implemented by {\tt circos.genomicInitialize}. 
The input data for the function is also a data frame with at least three columns. The first 
column is genomic category (for cytoband data, it is chromosome name), and the next two columns 
are genomic positions in each genomic category. The range of each category will be inferred 
from the minimum position and the maximum position in corresponding category. 
In the following example, a circos plot is initialized with three genes.

<<eval = FALSE>>=
df = data.frame(
    name  = c("TP53",  "TP63",    "TP73"),
    start = c(7565097, 189349205, 3569084),
    end   = c(7590856, 189615068, 3652765))
circos.genomicInitialize(df)
@

Note it is not necessary that the record for each gene is one row.

As explained in previous section, the order of genomic categories is controlled by the first column of {\tt df}
which depends on whether it is a factor or a simple vector. Alternative names can be assigned to each category
and the order of alternative names correspond to the order of genomic categories.

<<eval=FALSE>>=
circos.genomicInitialize(df)
circos.genomicInitialize(df, sector.names = c("tp53", "tp63", "tp73"))
circos.genomicInitialize(df, plotType)

circos.par(gap.degree = 2)
circos.genomicInitialize(df)
@

In following example, we plot the transcripts for TP53, TP63 and TP73 in a circular layout (figure \ref{fig:genomic_gene_model}).
The object {\tt tp\_family} is a list of data frames which contain positions of exons for
each transcript.

<<>>=
load(paste0(system.file(package = "circlize"), "/extdata/tp_family.RData"))
names(tp_family)
names(tp_family[["TP53"]])
head(tp_family[["TP53"]][[1]])

df = data.frame(gene = names(tp_family),
                start = sapply(tp_family, function(x) min(unlist(x))),
                end = sapply(tp_family, function(x) max(unlist(x))))
df
@

In the following code, we first create a track which identifies three genes.

<<genomic_gene_model_1, eval = FALSE>>=
circos.genomicInitialize(df)
circos.genomicTrackPlotRegion(ylim = c(0, 1), 
    bg.col = c("#FF000040", "#00FF0040", "#0000FF40"), 
    bg.border = NA, track.height = 0.05)
@

Next, we put each transcript line by line. It is simply adding lines and rectangles. 
The usage of {\tt circos.genomicTrackPlotRegion} will be introduced in later sections.

<<genomic_gene_model_2, eval = FALSE>>=
n = max(sapply(tp_family, length))
circos.genomicTrackPlotRegion(ylim = c(0.5, n + 0.5), 
    panel.fun = function(region, value, ...) {
        gn = get.cell.meta.data("sector.index")
        tr = tp_family[[gn]]  # all transcripts for this gene
        for(i in seq_along(tr)) {
            # for each transcript
            current_tr_start = min(tr[[i]]$start)
            current_tr_end = max(tr[[i]]$end)
            circos.lines(c(current_tr_start, current_tr_end), 
                c(n - i, n - i), col = "#CCCCCC")
            circos.genomicRect(tr[[i]], ytop = n - i + 0.4, 
                ybottom = n - i - 0.4, col = "orange", border = NA)
        }
}, bg.border = NA, track.height = 0.3)
circos.clear()
@

<<genomic_gene_model, echo = FALSE, out.width = '\\textwidth', fig.cap = "Alternative transcripts for genes.">>=
par(mar = c(1, 1, 1, 1))
<<genomic_gene_model_1>>
<<genomic_gene_model_2>>
@

\subsection{Zooming}

{\tt circos.genomicInitialize} is implemented by the general {\tt circos.initialize} function,
so same strategy can be applied to zoom chromosomes (figure \ref{fig:genomic_zoom}). Just be careful with the order of chromosomes.

<<>>=
cytoband = read.cytoband()
df = cytoband$df
chromosome = cytoband$chromosome

# copy regions for the two zoomed chromosomes
zoom_df = df[df[[1]] %in% chromosome[1:2], ]
zoom_df[[1]] = paste0("zoom_", zoom_df[[1]])
df2 = rbind(df, zoom_df)

# attach ranges for two zoomed chromosomes
xrange = c(cytoband$chr.len, cytoband$chr.len[1:2])
normal_sector_index = seq_along(chromosome)
zoomed_sector_index = length(chromosome) + 1:2

# normalize in normal chromsomes and zoomed chromosomes separately
sector.width = c(xrange[normal_sector_index] / sum(xrange[normal_sector_index]), 
                 xrange[zoomed_sector_index] / sum(xrange[zoomed_sector_index])) 
@

Here we can only use {\tt circos.genomicInitialize} to deal with zoomed chromosomes because chromosome 
with name "zoom\_chr1" is not a normal chromosome.

<<genomic_zoom_1, eval = FALSE>>=
par(mar = c(1, 1, 1, 1))
circos.par(start.degree = 90)
circos.genomicInitialize(df2, sector.width = sector.width)
@

In following tracks, data frames should also be extended with zoomed chromosomes.
We can write a simple function to do the extension.

<<>>=
extend_zoomed_chromosome_in_bed = function(bed, chromosome, prefix = "zoom_") {
    zoom_bed = bed[bed[[1]] %in% chromosome, , drop = FALSE]
    zoom_bed[[1]] = paste0(prefix, zoom_bed[[1]])
    rbind(bed, zoom_bed)
}
@

Then add a new track:

<<genomic_zoom_2, eval = FALSE>>=
bed = generateRandomBed(100)
circos.genomicTrackPlotRegion(extend_zoomed_chromosome_in_bed(bed, chromosome[1:2]),
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, pch = 16, cex = 0.5)
})
@

Also add link from original chromosome to the zoomed chromosome.

<<genomic_zoom_3, eval = FALSE>>=
circos.link("chr1", get.cell.meta.data("cell.xlim", sector.index = "chr1"),
    "zoom_chr1", get.cell.meta.data("cell.xlim", sector.index = "zoom_chr1"),
    col = "#00000020", border = NA)
circos.clear()
@

<<genomic_zoom, echo = FALSE, out.width = "0.8\\textwidth", fig.cap = "Zoom chromosomes.">>=
<<genomic_zoom_1>>
<<genomic_zoom_2>>
<<genomic_zoom_3>>
@

Same strategy can be applied if you only want to zoom sub regions inside a certain chromosome (just make them as pseudo chromosomes).

\section{Create plotting regions}

In following sections, chromosome will be used as the type of genomic category. 
In this section, we assume {\tt data} is simply a data frame. For more complex
situations and behaviour of the functions, we will introduce in the next section.

Similar as {\tt circos.trackPlotRegion}, {\tt circos.genomicTrackPlotRegion} also
accepts a self-defined function {\tt panel.fun} which is applied in every cell but with 
different form.

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(data, panel.fun = function(region, value, ...) {
    circos.genomicPoints(region, value, ...)
})
@

Inside {\tt panel.fun}, users can use low-level genomic graphic functions to add basic
graphics in each cell. {\tt panel.fun} expects two arguments {\tt region} and {\tt value}.
{\tt region} is a data frame containing start position and end position in the current chromosome which is extracted
from {\tt data}, and you can think it corresponds to values on x-axes. {\tt value} is also 
a data frame which contains other columns in {\tt data}, and you can think it corresponds
to values on y-axes. Besides, there should be a third arguments {\tt ...} which is 
mandatory and is used to pass user-invisible variables to inner functions (which we will
explain in next sections).


Following codes demonstrate values for {\tt region} and {\tt value} when used inside {\tt panel.fun}.

<<echo = 2:6>>=
pdf(NULL)
bed = generateRandomBed(nc = 1)
head(bed)
circos.initializeWithIdeogram(plotType = NULL)
circos.genomicTrackPlotRegion(bed, panel.fun = function(region, value, ...) {
    if(get.cell.meta.data("sector.index") == "chr1") {
        print(head(region))
        print(head(value))
    }
})
circos.clear()
invisible(dev.off())
@

Since {\tt circos.genomicTrackPlotRegion} will create a new track, it needs values to
calculate range of y-values to arrange data points. Users can either specify the index of 
numeric columns in {\tt data} by {\tt numeric.column} ({\bf named index or numeric index})
or set {\tt ylim}. If none of them are set, the function will try to look for all numeric columns
in {\tt data} (of course, excluding the first three columns), and set them as {\tt numeric.column}.

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(data, ylim = c(0, 1),
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, ...)
})
circos.genomicTrackPlotRegion(data, numeric.column, 
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, ...)
})
@

Since genomic functions are implemented by basic circos functions, you can use {\tt circos.info}
anywhere to get information of sectors and tracks.

\subsection{Points}

{\tt circos.genomicPoints} is similar as {\tt circos.points}. The difference is {\tt circos.genomicPoints}
expects a data frame containing genomic regions and a data frame containing values. The data column for 
plotting should be indicated by {\tt numeric.column}. If the function is called
inside {\tt circos.genomicTrackPlotRegion} and users have been already set {\tt numeric.column}
in {\tt circos.genomicTrackPlotRegion}, proper value of {\tt numeric.column} will be
passed to {\tt circos.genomicPoints} through {\tt ...} in {\tt panel.fun}. Which means, you need
to add {\tt ...} as the final argument in {\tt circos.genomicPoints} to get such information.
If {\tt numeric.column} is not set in both place, {\tt circos.genomicPoints} will use all numeric columns detected in {\tt value}.

Note here {\tt numeric.column} is measured in {\tt value} while {\tt numeric.column} in
{\tt circos.genomicTrackPlotRegion} is measured in the complete data frame. There is a
difference of 3! When {\tt numeric.column} is passed to {\tt circos.genomicPoints} internally, 3 is subtracted automatically. 
If you use character index instead of numeric index, you do not need to worry about it. 

<<eval=FALSE>>=
circos.genomicPoints(region, value, ...)
circos.genomicPoints(region, value, numeric.column = c(1, 2))
circos.genomicPoints(region, value, cex, pch)
circos.genomicPoints(region, value, sector.index, track.index)
@

If there is only one numeric column, graphical parameters such as {\tt pch}, {\tt cex}
can be of length one or number of rows of {\tt region}. If there are more than
one numeric columns specified, points for each numeric column will be added iteratively,
and the graphical parameters should be either length one or number of numeric columns specified.

{\tt circos.genomicPoints} is implemented by {\tt circos.points}. The basic idea of the 
implementation is like follows.

<<eval = FALSE>>=
circos.genomicPoints = function(region, value, numeric.column = 1, ...) {
    x = (region[[2]] + region[[1]])/2
    y = value[[numeric.column]]
    circos.points(x, y, ...)
}
@

As shown above, {\tt circos.genomicPoints} also provides arguments {\tt sector.index} and {\tt track.index}.
The default values are current sector and current track. With these two arguments, you can
add points by {\tt circos.genomicPoints} outside {\tt circos.genomicPlotTrackRegion}. 

\subsection{Lines}

{\tt circos.genomicLines} is similar as {\tt circos.lines}. The setting of graphical parameters
is similar as {\tt circos.genomicPoints}. 

<<eval=FALSE>>=
circos.genomicLines(region, value, ...)
circos.genomicLines(region, value, numeric.column = c(1, 2))
circos.genomicLines(region, value, lwd, lty = "segment")
circos.genomicLines(region, value, area, baseline, border)
circos.genomicLines(region, value, sector.index, track.index)
@

For {\tt lty}, we additionally provide a new option {\tt segment} by which each genomic
interval will represent as a 'horizontal' line.

\subsection{Text}

For {\tt circos.genomicText}, the position of text can be specified either by {\tt numeric.column} (index)
or a separated vector {\tt y}. The labels of text can be specified either by {\tt labels.column} (index)
or a vector {\tt labels}.

<<eval=FALSE>>=
circos.genomicText(region, value, ...)
circos.genomicText(region, value, y, labels)
circos.genomicText(region, value, numeric.column, labels.column)
circos.genomicText(region, value, facing, niceFacing, adj)
circos.genomicText(region, value, sector.index, track.index)
@


\subsection{Rectangle}

For {\tt circos.genomicRect}, the positions of top and bottom of the rectangles can be
specified by {\tt ytop}, {\tt ybottom} or {\tt ytop.column}, {\tt ybottom.column} (index).

<<eval=FALSE>>=
circos.genomicRect(region, value, ytop = 1, ybottom = 0)
circos.genomicRect(region, value, ytop.column = 2, ybottom = 0)
circos.genomicRect(region, value, col, border)
@

One of the usage of {\tt circos.genomicRect} is to plot heatmap on the circle.
\textbf{circlize} provides a simple function {\tt colorRamp2} to interpolate colors. 
The arguments of {\tt colorRamp2} are break points and colors
which correspond to the break points. {\tt colorRamp2} returns a new function
which can be used to generate new colors. 

<<>>=
col_fun = colorRamp2(breaks = c(-1, 0, 1), colors = c("green", "black", "red"))
col_fun(c(-2, -1, -0.5, 0, 0.5, 1, 2))
col_fun = colorRamp2(breaks = -log10(c(1, 0.05, 1e-4)), 
    colors = c("green", "black", "red"))
p_value = c(0.8, 0.5, 0.001)
col_fun(-log10(p_value))
@

\section{More details on {\tt circos.genomicTrackPlotRegion}}

The behaviour of {\tt circos.genomicTrackPlotRegion} and {\tt panel.fun} 
will be different according to different input data (e.g. is it a simple data frame or a list
of data frames? If it is a data frame, how many numeric columns it has?) and different settings.

\subsection{Normal mode}

\subsubsection{Input is a data frame}

If input data is a simple data frame, {\tt region} in {\tt panel.fun} would be a data frame
containing start position and end position in the current chromosome which is extracted
from input data. {\tt value} is also a data frame which contains columns in input data
excluding the first three columns. Index of proper numeric columns will be passed by {\tt ...}
if it is set in {\tt circos.genomicTrackPlotRegion}.
So if users want to use such information, they need to pass {\tt ...} to low-level
genomic function such as {\tt circos.genoimcPoints} as well.

If there are more than one numeric columns, all columns that are specified will be added to
the plot.

<<eval=FALSE>>=
bed = generateRandomBed(nc = 2)
# just note `numeric.column` is measured in `bed`
circos.genomicTrackPlotRegion(bed, numeric.column = 4, 
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, ...)
        circos.genomicPoints(region, value)
        # here `numeric.column` is measured in `value`
        circos.genomicPoints(region, value, numeric.column = 1)
})
@

\subsubsection{Input is a list of data frames}

If input data is a list of data frames, {\tt panel.fun} is applied on each data frame iteratively.
Under such situation, {\tt region} and {\tt value} will contain corresponding data
in the current data frame. The index for the current data frame can be get by {\tt getI(...)}.
For {\tt numeric.column} argument, if input data is a list of data frame, the length can only be one
or the number of data frames, which means, there is only one numeric column that will be used in each data frame.

<<echo=TRUE, eval=FALSE>>=
bedlist = list(generateRandomBed(), generateRandomBed())
circos.genomicTrackPlotRegion(bedlist,
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicPoints(region, value, col = i, ...)
})

# column 4 in the first bed and column 5 in the second bed
circos.genomicTrackPlotRegion(bedlist, numeric.column = c(4, 5),
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicPoints(region, value, col = i, ...)
})
@

\subsection{{\tt stack} mode}

{\tt circos.genomicTrackPlotRegion} also supports a {\tt stack} mode. Under {\tt stack} mode,
{\tt ylim} is re-defined inside the function. The y-axis will be splitted into several layers
with equal height and graphics are put on 'horizontal' layers (y = 1, 2, ...).

\subsubsection{Input is a data frame}

Under {\tt stack} mode, when input data is a single data frame containing one or more numeric columns, 
each numeric column defined in {\tt numeric.column} will be treated as a single unit. 
{\tt ylim} is re-defined to {\tt c(0.5, n+0.5)} in which {\tt n} is number of numeric columns. 
{\tt panel.fun} will be applied iteratively on each numeric column. In each iteration, 
in {\tt panel.fun}, {\tt region} is still the genomic regions in current genomic category, 
but {\tt value} only contains current numeric column plus all non-numeric columns.
All low-level genomic graphic functions will be applied on each `horizontal line' {\tt y = i} 
in which {\tt i} is the index of current numeric column. The value of {\tt i} can 
be obtained by {\tt getI(...)}.

<<eval=FALSE>>=
bed = generateRandomBed(nc = 2)
circos.genomicTrackPlotRegion(bed, stack = TRUE,
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicPoints(region, value, col = i, ...)
})
@

\subsubsection{Input is a list of data frame}

When input data is a list containing data frames, each data frame will be treated 
as a single unit. The situation is quite similar as described previously.
{\tt ylim} is re-defined to {\tt c(0.5, n+0.5)} in which {\tt n} is number of data frames. 
{\tt panel.fun} will be applied iteratively on each data frame. In each
iteration, in {\tt panel.fun}, {\tt region} is still the genomic regions in current 
chromosome, and {\tt value} contains columns in current data frame excluding the first three columns.
Still, graphics by low-level genomic functions will be added on the `horizontal' lines.

<<echo=TRUE, eval=FALSE>>=
bedlist = list(generateRandomBed(), generateRandomBed())
circos.genomicTrackPlotRegion(bedlist, stack = TRUE,
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicPoints(region, value, ...)
})
@

Under {\tt stack} mode, the differentce between using a data frame with multiple numeric columns
and using a list of data frames is if using a data frame, all layers share the same genomic positions
while if using a list of data frames, the genomic positions for each layer can be different.

\section{Applications}

In this section, we will show several real examples of making genomic plots under different modes.

\subsection{Points}

To make plots more clear to look at, we only add graphics in the 1/4 of the circle and
initialize the plot only with chromosome 1 (figure \ref{fig:genomic_application_points}).

<<genomic_application_points_0, eval = FALSE>>=
par(mar = c(1, 1, 1, 1))
set.seed(999)

circos.par("track.height" = 0.1, start.degree = 90,
    canvas.xlim = c(0, 1), canvas.ylim = c(0, 1), gap.degree = 270)
circos.initializeWithIdeogram(chromosome.index = "chr1", plotType = NULL)
@

In track A, it is the most simple way to add points.

<<genomic_application_points_A, eval = FALSE, echo = 1:5>>=
### track A
bed = generateRandomBed(nr = 300)
circos.genomicTrackPlotRegion(bed, panel.fun = function(region, value, ...) {
    circos.genomicPoints(region, value, pch = 16, cex = 0.5, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "A", adj = c(1.1, 0.5))
@

In track B, under {\tt stack} mode, points are added in one horizontal line (here we
additionally add the dashed line) because there is only one numeric column in {\tt bed}.

<<genomic_application_points_B, eval = FALSE, echo = 1:10>>=
### track B
bed = generateRandomBed(nr = 300)
circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, pch = 16, cex = 0.5, ...)
        
        i = getI(...)
        cell.xlim = get.cell.meta.data("cell.xlim")
        circos.lines(cell.xlim, c(i, i), lty = 2, col = "#00000040")
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "B", adj = c(1.1, 0.5))
@

In track C, the input data is a list of two data frames. In the plot, sizes of the points
correspond to the values of regions, and colors correspond to different data frames.

<<genomic_application_points_C, eval = FALSE, echo = 1:10>>=
### track C
bed1 = generateRandomBed(nr = 300)
bed2 = generateRandomBed(nr = 300)
bed_list = list(bed1, bed2)
circos.genomicTrackPlotRegion(bed_list, 
    panel.fun = function(region, value, ...) {
        cex = (value[[1]] - min(value[[1]]))/(max(value[[1]]) - min(value[[1]]))
        i = getI(...)
        circos.genomicPoints(region, value, cex = cex, pch = 16, col = i, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "C", adj = c(1.1, 0.5))
@

In track D, plot the data frame list under {\tt stack} mode.

<<genomic_application_points_D, eval = FALSE, echo = 1:9>>=
### track D
circos.genomicTrackPlotRegion(bed_list, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        cex = (value[[1]] - min(value[[1]]))/(max(value[[1]]) - min(value[[1]]))
        i = getI(...)
        circos.genomicPoints(region, value, cex = cex, pch = 16, col = i, ...)
        cell.xlim = get.cell.meta.data("cell.xlim")
        circos.lines(cell.xlim, c(i, i), lty = 2, col = "#00000040")
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "D", adj = c(1.1, 0.5))
@

In track E, the data frame has four numeric columns. Different colors are assigned to 
different columns. Note {\tt value} is a data frame with four columns.

<<genomic_application_points_E, eval = FALSE, echo = 1:6>>=
### track E
bed = generateRandomBed(nr = 300, nc = 4)
circos.genomicTrackPlotRegion(bed, 
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, cex = 0.5, pch = 16, col = 1:4, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "E", adj = c(1.1, 0.5))
@

In track F, the data frame has four columns but is plotted under {\tt stack} mode.
Note here {\tt value} is a data frame with only one column (but it will be exexuted
4 times in each cell).

<<genomic_application_points_F, eval = FALSE, echo = c(1:11, 14)>>=
### track F
bed = generateRandomBed(nr = 300, nc = 4)
circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        cex = (value[[1]] - min(value[[1]]))/(max(value[[1]]) - min(value[[1]]))
        i = getI(...)
        circos.genomicPoints(region, value, cex = cex, pch = 16, col = i, ...)
        
        cell.xlim = get.cell.meta.data("cell.xlim")
        circos.lines(cell.xlim, c(i, i), lty = 2, col = "#00000040")
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "F", adj = c(1.1, 0.5))
circos.clear()
@

<<genomic_application_points, echo = FALSE, out.width = "0.8\\textwidth", fig.cap = "Add points under different modes.">>=
<<genomic_application_points_0>>
<<genomic_application_points_A>>
<<genomic_application_points_B>>
<<genomic_application_points_C>>
<<genomic_application_points_D>>
<<genomic_application_points_E>>
<<genomic_application_points_F>>
@

\subsection{Lines}

Initialize the plot with one quarter of the circle (figure \ref{fig:genomic_application_lines}).

<<genomic_application_lines_0, eval = FALSE>>=
par(mar = c(1, 1, 1, 1))
circos.par("track.height" = 0.1, start.degree = 90,
    canvas.xlim = c(0, 1), canvas.ylim = c(0, 1), gap.degree = 270)
circos.initializeWithIdeogram(chromosome.index = "chr1", plotType = NULL)
@

In track A, it is the most simple way to add lines.

<<genomic_application_lines_A, eval = FALSE, echo = 1:6>>=
### track A
bed = generateRandomBed(nr = 500)
circos.genomicTrackPlotRegion(bed, 
    panel.fun = function(region, value, ...) {
        circos.genomicLines(region, value, type = "l", ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "A", adj = c(1.1, 0.5))
@

In track B, add lines for a list of data frames. Different colors refers to different
data frames.

<<genomic_application_lines_B, eval = FALSE, echo = 1:9>>=
### track B
bed1 = generateRandomBed(nr = 500)
bed2 = generateRandomBed(nr = 500)
bed_list = list(bed1, bed2)
circos.genomicTrackPlotRegion(bed_list, 
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicLines(region, value, col = i, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "B", adj = c(1.1, 0.5))
@

In track C, it is the {\tt stack} mode of adding lines. Here the width of each line corresponds
to the width of each genomic interval.

<<genomic_application_lines_C, eval = FALSE, echo = 1:6>>=
### track C
circos.genomicTrackPlotRegion(bed_list, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicLines(region, value, col = i, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "C", adj = c(1.1, 0.5))
@

In track D, lines for the four numeric columns are added with different colors.

<<genomic_application_lines_D, eval = FALSE, echo = 1:6>>=
### track D
bed = generateRandomBed(nr = 500, nc = 4)
circos.genomicTrackPlotRegion(bed, 
    panel.fun = function(region, value, ...) {
        circos.genomicLines(region, value, col = 1:4, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "D", adj = c(1.1, 0.5))
@

In track E, it is the {\tt stack} mode for a data frame with more than one numeric columns.

<<genomic_application_lines_E, eval = FALSE, echo = 1:7>>=
### track E
bed = generateRandomBed(nr = 500, nc = 4)
circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicLines(region, value, col = i, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "E", adj = c(1.1, 0.5))
@

In track F, we specify {\tt type} to {\tt segment} and use different colors for segments.
Note each segment is located at its corresponding y-value.

<<genomic_application_lines_F, eval = FALSE, echo = c(1:7, 10)>>=
### track F
bed = generateRandomBed(nr = 200)
circos.genomicTrackPlotRegion(bed, 
    panel.fun = function(region, value, ...) {
        circos.genomicLines(region, value, type = "segment", lwd = 2, 
            col = rand_color(nrow(region)), ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "F", adj = c(1.1, 0.5))
circos.clear()
@


<<genomic_application_lines, echo = FALSE, out.width = "0.8\\textwidth", fig.cap = "Add lines under different modes.">>=
<<genomic_application_lines_0>>
<<genomic_application_lines_A>>
<<genomic_application_lines_B>>
<<genomic_application_lines_C>>
<<genomic_application_lines_D>>
<<genomic_application_lines_E>>
<<genomic_application_lines_F>>
@

\subsection{Rectangles}

Again, initialize the plot with one quarter of the circle (figure \ref{fig:genomic_application_rect}).
Also, we want to map the values to colors by {\tt f}.

<<genomic_application_rect_0, eval = FALSE>>=
par(mar = c(1, 1, 1, 1))
circos.par("track.height" = 0.1, start.degree = 90,
    canvas.xlim = c(0, 1), canvas.ylim = c(0, 1), gap.degree = 270)
circos.initializeWithIdeogram(chromosome.index = "chr1", plotType = NULL)
f = colorRamp2(breaks = c(-1, 0, 1), colors = c("green", "black", "red"))
@

In track A, {\tt bed} has four numeric columns and {\tt stack} mode is used to arrange 
the heatmap.

<<genomic_application_rect_A, eval = FALSE, echo = 1:6>>=
### track A
bed = generateRandomBed(nr = 100, nc = 4)
circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        circos.genomicRect(region, value, col = f(value[[1]]), border = NA, ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "A", adj = c(1.1, 0.5))
@

In track B, add rectangles for a list of data frames. Comparing to track A, here genomic
regions for each data frame are different, so the positions of rectangles in the first layer
are different from that in the second layer.

Under {\tt stack} mode, {\tt ytop} and {\tt ybottom} can be set to adjust the height of 
rectangles. It would be straightforward because for each layer, the position is {\tt y = i}
and the height of each rectangle is 1.

<<genomic_application_rect_B, eval = FALSE, echo = 1:10>>=
### track B
bed1 = generateRandomBed(nr = 100)
bed2 = generateRandomBed(nr = 100)
bed_list = list(bed1, bed2)
circos.genomicTrackPlotRegion(bed_list, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicRect(region, value, ytop = i + 0.4, ybottom = i - 0.4,
            col = f(value[[1]]), ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "B", adj = c(1.1, 0.5))
@

In track C, the plot are the same as track B, but without {\tt stack} mode.
Note here we explicitly specify {\tt ylim}.

<<genomic_application_rect_C, eval = FALSE, echo = 1:7>>=
### track C
circos.genomicTrackPlotRegion(bed_list, ylim = c(0, 3), 
    panel.fun = function(region, value, ...) {
        i = getI(...)
        circos.genomicRect(region, value, ytop = i + 0.4, ybottom = i - 0.4, 
            col = f(value[[1]]), ...)
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "C", adj = c(1.1, 0.5))
@

In track D, bars are added to the base line ({\tt y = 0}).

<<genomic_application_rect_D, eval = FALSE, echo = c(1:10, 13)>>=
### track D
bed = generateRandomBed(nr = 200)
circos.genomicTrackPlotRegion(bed, 
    panel.fun = function(region, value, ...) {
        circos.genomicRect(region, value, ytop.column = 1, ybottom = 0, 
            col = ifelse(value[[1]] > 0, "red", "green"), ...)
            
        cell.xlim = get.cell.meta.data("cell.xlim")
        circos.lines(cell.xlim, c(0, 0), lty = 2, col = "#00000040")
})
pos = get.cell.meta.data("yplot")
text(0, mean(pos), "D", adj = c(1.1, 0.5))
circos.clear()
@

<<genomic_application_rect, echo = FALSE, out.width = "0.8\\textwidth", fig.cap = "Add rectangles under different modes.">>=
<<genomic_application_rect_0>>
<<genomic_application_rect_A>>
<<genomic_application_rect_B>>
<<genomic_application_rect_C>>
<<genomic_application_rect_D>>
@

\subsection{Mixed use of general circos functions}

{\tt panel.fun} is applied on each cell, which means, besides genomic functions, you can 
also use general circos functions to add more graphics. For example, some horizontal lines 
and texts are added to each cell and axes are put on top of each cell:

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(bed, ylim = c(-1, 1),
    panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, ...)

        cell.xlim = get.cell.meta.data("cell.xlim")
        for(h in c(-1, -0.5, 0, 0.5, 1)) {
            circos.lines(cell.xlim, c(0, 0), lty = 2, col = "grey")
        }
        circos.text(x, y, labels)
        circos.axis("top")
})
@

\section{Links}

{\tt circos.genomicLink} expects two data frames and it will add links from genomic
intervals in the first data frame to corresponding genomic intervals in the second
data frame (figure \ref{fig:genomic_links}). You can set graphical parameters as a single
scalar or a vector.

<<genomic_links, echo = 2:11, out.width = "\\textwidth", out.height = "0.5\\textwidth", fig.width = 12, fig.height = 6, fig.cap = "Add links from two sets of genomic regions.">>=
par(mfrow = c(1, 2), mar = c(1, 1, 1, 1))
bed1 = generateRandomBed(nr = 100)
bed1 = bed1[sample(nrow(bed1), 20), ]
bed2 = generateRandomBed(nr = 100)
bed2 = bed2[sample(nrow(bed2), 20), ]
circos.initializeWithIdeogram(plotType = c("axis", "labels"))
circos.genomicLink(bed1, bed2)
circos.clear()

circos.initializeWithIdeogram(plotType = c("axis", "labels"))
circos.genomicLink(bed1, bed2, col = rand_color(nrow(bed1), transparency = 0.5), 
    border = NA)
circos.clear()
par(mfrow = c(1, 1))
@

\section{Highlight chromosomes}

{\tt highlight.chromosome} provides a simple way to highlight chromosomes. Just remember
to use transparent filled colors. The position of the highlighted regions
can be fine-tuned by {\tt padding} argument which are percent of corresponding height 
and width in the highlighted regions.

In following examples, we first create five more tracks with random points.

<<genomic_highlight_1, eval = FALSE>>=
circos.par("track.height" = 0.1)
circos.initializeWithIdeogram(plotType = c("axis", "labels"))

for(i in 1:5) {
    bed = generateRandomBed(nr = 100)
    circos.genomicTrackPlotRegion(bed, panel.fun = function(region, value, ...) {
        circos.genomicPoints(region, value, pch = 16, cex = 0.5, ...)
    })
}
@

Following are several ways of highlighting.
we can highlight individual tracks by setting {\tt track.index}. If two tracks or
two chromosomes are neighbours, they will be taken as a union (figure \ref{fig:genomic_highlight}).

<<genomic_highlight_2, eval = FALSE>>=
highlight.chromosome(c("chrX", "chrY", "chr1"))
highlight.chromosome("chr3", col = "#00FF0040", padding = c(0.05, 0.05, 0.15, 0.05))
highlight.chromosome("chr5", col = NA, border = "red", lwd = 2, 
    padding = c(0.05, 0.05, 0.15, 0.05))
highlight.chromosome("chr7", col = "#0000FF40", track.index = c(2, 4, 5))
highlight.chromosome(c("chr9", "chr10", "chr11"), col = NA, border = "green", 
    lwd = 2, track.index = c(2, 4, 5))
highlight.chromosome(paste0("chr", c(1:22, "X", "Y")), col = "#FFFF0040", 
    track.index = 6)
circos.clear()
@

<<genomic_highlight, echo = FALSE, out.width = "0.8\\textwidth", fig.cap = "Highlight different chromosomes and tracks.">>=
par(mar = c(1, 1, 1, 1))
<<genomic_highlight_1>>
<<genomic_highlight_2>>
@

{\tt circos.info} would be helpful for you to find the correct track index. 
And if you only want to highlight some sub regions in one chromosome, use {\tt draw.sector} directly. 

\section{High-level genomic functions}

\textbf{circlize} implements several high-level functions which may help to visualize genomic data more effectively.

\subsection{Position transformation}

There is one representative situation when genomic position transformation needs to be applied. 
For example, there are regions which are quite densely located on the chromosome. 
If heatmap or text for those dense regions are to be plotted. They would be overlapped and 
hard to identify, also ugly to visualize. Thus, a way to transform original positions to 
new positions would help for the visualization.

Low-level genomic functions such as {\tt circos.genomicPoints} all accept an argument {\tt posTransform}
to apply user-defined position transformation. Value for {\tt posTransform} is a self-defined function which
only accepts at least one argument: a data frame with two columns (start position and end position). 
Because {\tt posTransform} is bound to low-level graphic functions, it is not necessary to add
chromosome information in the input data. There is only one requirement for position 
transformation: Number of rows of regions should be the same before and after the transformation. 

\subsubsection{The default position transformation function}

In \textbf{circlize}, there already provides a position transformation function {\tt posTransform.default} 
which distributes positions uniformly in the current chromosome.

Following code does the transformation. The points are plotted with the new transformed regions.

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(data, panel.fun = function(region, value, ...) {
    circos.genomicPoints(region, value, posTransform = posTransform.default, ...)
})
@

Of course we need identification lines to connect the untransformed regions to the transformed
regions. There is a function called {\tt circos.genomicPosTransformLines} which does the job.
Same value for {\tt posTransform} which was used before should be passed to {\tt circos.genomicPosTransformLines}.
Note {\tt circos.genomicPosTransformLines} will create a new track to put such identification lines.
In the function, {\tt direction} controls whether the position transformed track is inside or outside
the identification line track, and {\tt horizontalLine} controls whether add lines to identify the regions.


<<eval=FALSE>>=
circos.genomicPosTransformLines(data, posTransform = posTransform.default)
circos.genomicPosTransformLines(data, posTransform = posTransform.default,
    horizontalLine = "top")
circos.genomicPosTransformLines(data, posTransform = posTransform.default,
    direction = "outside")
@

In following example, we make heatmap for selected regions in the genome (figure \ref{fig:genomic_postransform} A).

<<genomic_postransform_1, eval = FALSE>>=
circos.par(cell.padding = c(0, 0, 0, 0))
circos.initializeWithIdeogram()
bed = generateRandomBed(nr = 100, nc = 4)

# note how 'horizontalLine' works
circos.genomicPosTransformLines(bed, posTransform = posTransform.default, 
    horizontalLine = "top", track.height = 0.1)

f = colorRamp2(breaks = c(-1, 0, 1), colors = c("green", "black", "red"))
circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        circos.genomicRect(region, value, col = f(value[[1]]), 
            border = f(value[[1]]), posTransform = posTransform.default, ...)
}, bg.border = NA)

circos.clear()
@

For the second example, the heatmap is plotted outside the identification lines, so here
{\tt direction} is set to {\tt outside} (figure \ref{fig:genomic_postransform} B). 
Also note how we add the ideogram inside the circle.

<<genomic_postransform_2, eval = FALSE>>=
circos.par(cell.padding = c(0, 0, 0, 0))
circos.initializeWithIdeogram(plotType = NULL)

circos.genomicTrackPlotRegion(bed, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        circos.genomicRect(region, value, col = f(value[[1]]), 
            border = f(value[[1]]), posTransform = posTransform.default, ...)
}, bg.border = NA)

circos.genomicPosTransformLines(bed, posTransform = posTransform.default, 
    direction = "outside", horizontalLine = "bottom", track.height = 0.1)

cytoband = read.cytoband()$df
circos.genomicTrackPlotRegion(cytoband, stack = TRUE, 
    panel.fun = function(region, value, ...) {
        circos.genomicRect(region, value, col = cytoband.col(value[[2]]), border = NA)
        xlim = get.cell.meta.data("xlim")
        ylim = get.cell.meta.data("ylim")
        circos.rect(xlim[1], ylim[1], xlim[2], ylim[2], border = "black")
}, track.height = 0.05, bg.border = NA)

circos.clear()
@

<<genomic_postransform, echo = FALSE, out.width = "0.6\\textwidth", out.height = "1.2\\textwidth", fig.width = 6, fig.height = 12, fig.cap = "Position transformation with {\\tt posTransform.default}. A) transformation is inside; B) transformation is outside.">>=
par(mfrow = c(2, 1), mar = c(1, 1, 1, 1))
<<genomic_postransform_1>>
text(-0.9, 0.9, "A", cex = 1.5)
<<genomic_postransform_2>>
text(-0.9, 0.9, "B", cex = 1.5)
par(mfrow = c(1, 1))
@

\subsubsection{Position transformation for text}

There is another position transformation function {\tt posTransform.text} provided in \textbf{circlize} that can smartly position
text on the circle. Normally, we don't want the text too away from the original position
and also we want to get avoid of text overlapping. {\tt posTransform.text} calculates the height of 
text and transform positions properly. Since such text position transformation relies on font size of text
and which track the text is in, the usage of {\tt posTransform.text} is a little bit complex. 
Currently, {\tt posTransform.text} only makes sense if it is called inside {\tt circos.genomicText} 
and only works when {\tt facing} is set to {\tt clockwise} or {\tt reverse.clockwise}.

In following example code, using {\tt posTransform.text} is quite similar as {\tt posTransform.default}
(but note this function is quit experimental):

<<eval=FALSE>>=
bed = generateRandomBed(nr = 400, fun = function(k) rep("text", k))
circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), 
    panel.fun = function(region, value, ...) {
        circos.genomicText(region, value, y = 0, labels.column = 1, 
            facing = "clockwise", adj = c(0, 0.5), 
            posTransform = posTransform.text, cex = 0.8)
}, track.height = 0.1, bg.border = NA)
@

For the next track where position transformation lines are plotted. Code will be more complex.
Since the calculation of the position of the line ends which corresponds to the transformed track
relies on text settings (e.g. font size, font family) and track for the text (i.e.
the track will affect the position of text), we need to go back to the track where text position transformation happened
and use the same settings for the text. Such information should be collected when plotting transformation lines.
The solution by \textbf{circlize} is to wrap {\tt posTransform.text} with such information,
then text information will be recovered and the transformation will be calculated in the correct track.

As we mentioned before, {\tt posTransform} expect at least one argument, so the argument {\tt value}
is completely fine here.

<<eval=FALSE>>=
i_track = get.cell.meta.data("track.index")  # the nearest track
# we put `y`, `labels`, ... into a self-defined function
# because these parameters will affect the text position
circos.genomicPosTransformLines(bed, direction = "outside",
    posTransform = function(region, value) 
        posTransform.text(region, y = 0, labels = value[[1]], 
            cex = 0.8, track.index = i_track)
)
@

If the transformation line track is plotted before the track where text position is transformed, 
things will be more complicated. Since {\tt circos.genomicPosTransformLines} 
needs information of the text, but at that moment, track which transforms text positions 
has not be created. The solution is first creating an empty track for position transformation lines, 
then adding the position transformation track. Finally go back to the 
transformation line track to add transformation lines.

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), track.height = 0.1, bg.border = NA)
i_track = get.cell.meta.data("track.index")  # remember this empty track, we'll come back

circos.genomicTrackPlotRegion(bed, ylim = c(0, 1),
    panel.fun = function(region, value, ...) {
        circos.genomicText(region, value, y = 1, labels.column = 1, 
            facing = "clockwise", adj = c(1, 0.5),
            posTransform = posTransform.text, cex = 0.8)
}, track.height = 0.1, bg.border = NA)
tr_track = get.cell.meta.data("track.index") # position transformation track

# because `circos.genomicPosTransformLines` is implemented by 
# `circos.trackPlotRegion`, it accepts `track.index` argument.
circos.genomicPosTransformLines(bed, 
    posTransform = function(region, value) 
        posTransform.text(region, y = 1, labels = value[[1]], 
            cex = 0.8, track.index = tr_track),
    direction = "inside", track.index = i_track
)
@

Padding of text after transformation can be set through {\tt padding} argument to adjust the space
between two neighbouring text. Of course, when add transformation lines, {\tt padding} should
also be wrapped in. 

<<eval=FALSE>>=
circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), 
    panel.fun = function(region, value, ...) {
        circos.genomicText(region, value, y = 0, labels.column = 1, 
            facing = "clockwise", adj = c(0, 0.5), posTransform = posTransform.text, 
            cex = 0.8, padding = 0.2)
}, track.height = 0.1, bg.border = NA)

i_track = get.cell.meta.data("track.index")  # previous track
circos.genomicPosTransformLines(bed, 
    posTransform = function(region, value) posTransform.text(region, y = 0, 
        labels = value[[1]], cex = 0.8, padding = 0.2, track.index = i_track),
    direction = "outside"
)
@

For examples and comparisons between {\tt posTransform.default} and {\tt posTransform.text} when
visualizing text, please refer to figure \ref{fig:genomic_text_pos_transformation}. Source code
is available in the help page of {\tt posTransform.text}.

<<genomic_text_pos_transformation, echo = FALSE, out.width = "\\textwidth", fig.cap = "Transformation of text positions.">>=
source("src/genomic-07-posTransformLinesText.R")
@

\subsubsection{How about if I have a lot of labels?}

Position transformation is applied inside each chromosome, which means, the transformed new positions
are still located in the original chromosome. This would cause a problem that if you have 
a lot of labels to put, they will overlap with each other (figure \ref{fig:genome_more_labels} B). Unfortunately there is no simple solution
to extend the positions outside of the chromosome.

Anyway, we provide a 'not-too-perfect' solution but the code is a little bit lengthy (figure \ref{fig:genome_more_labels} A). The idea
is to put labels and transformation lines in one same track. Here we use {\tt posTransform.text}
directly to calculate the transformed positions.

To make a clear demonstration, we only add graphics in the first quarter of the circle, only 
on one chromosome.

<<genome_more_labels_1, eval = FALSE>>=
set.seed(999)
bed = generateRandomBed(nr = 800, fun = function(k) rep("text", k))

par(mar = c(1, 1, 1, 1))
circos.par(start.degree = 75, canvas.xlim = c(0, 1), canvas.ylim = c(0, 1), 
    gap.degree = 300, cell.padding = c(0, 0, 0, 0), track.margin = c(0, 0))
circos.initializeWithIdeogram(plotType = NULL, chromosome.index = "chr1")
@

The steps inside {\tt panel.fun} are:

\begin{enumerate}
    \item Original positions are plotted as dots.
    \item In order to get rid of text overlapping, we separate the text into two groups.
    One group of text are plotted lower than the other group. In following code, the x-axis 
    is splitted into 6 intervals and {\tt region} and {\tt value} are separated into groups
    by looking at whether the middle points of corresponding regions are in the odd intervals (i.e., the first interval, the third interval, ...)
    or even intervals. Of course, in this step, users should determine their own rule to separate the groups.
    \item Text in different groups are assigned with different heights ({\tt y}).
    \item Do the position transformation on text and add labels.
    \item Add the position transformation lines manually. Note the line for lower text and higher text
    are different. We do something more to make sure lines which points to higher text do not overlap
    with lower text too much.
\end{enumerate}

<<genome_more_labels_2, eval = FALSE>>=
circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), 
    panel.fun = function(region, value, ...) {
    
    # original positions
    circos.genomicPoints(region, data.frame(rep(0, nrow(region))), pch = 16)
        
    xlim = get.cell.meta.data("xlim")
    breaks = seq(xlim[1], xlim[2], length.out = 7)
    midpoints = (region[[1]] + region[[2]])/2
    for(i in seq_along(breaks)[-1]) {
        # index for current interval
        l = midpoints >= breaks[i - 1] & midpoints < breaks[i]
        # if there is no data in this interval
        if(sum(l) == 0) next
        
        # sub-regions in this interval
        sub_region = region[l, , drop = FALSE]
        sub_value = value[l, , drop = FALSE]
        
        # note here i == 2, 4, 6, ... corresponds to odd intervals
        # text in odd intervals are in lower position
        if(i %% 2 == 0) {
            y = 0.2
        } else {
            y = 0.8
        }
        
        # get the transformed position and add text with new positions
        tr_region = posTransform.text(sub_region, y = y, labels = sub_value[[1]], 
            cex = 0.8, adj = c(0, 0.5))
        circos.genomicText(tr_region, sub_value, labels.column = 1, y = y, 
            adj = c(0, 0.5), facing = "clockwise", niceFacing = TRUE, cex = 0.8)
        
        # add position transformation lines for odd intervals
        if(i %% 2 == 0) {    
            for(i in seq_len(nrow(sub_region))) {
                x = c( (sub_region[i, 1] + sub_region[i, 2])/2,
                       (sub_region[i, 1] + sub_region[i, 2])/2,
                       (tr_region[i, 1] + tr_region[i, 2])/2,
                       (tr_region[i, 1] + tr_region[i, 2])/2)
                y = c(0, 0.2/3, 0.2/3*2, 0.2)
                circos.lines(x, y)
            }
        } else { # add position transformation lines for even intervals
            median_sub_region_midpoint = median(midpoints[l])
            sub_region_width = max(midpoints[l]) - min(midpoints[l])
            for(i in seq_len(nrow(sub_region))) {
                x = c( (sub_region[i, 1] + sub_region[i, 2])/2,
                       (sub_region[i, 1] + sub_region[i, 2])/2,
                       median_sub_region_midpoint + 
                           sub_region_width*(i - nrow(sub_region))/nrow(sub_region) * 0.2,
                       median_sub_region_midpoint + 
                           sub_region_width*(i - nrow(sub_region))/nrow(sub_region) * 0.2,
                       (tr_region[i, 1] + tr_region[i, 2])/2,
                       (tr_region[i, 1] + tr_region[i, 2])/2)
                y = c(0, 0.1, 0.2, 0.6, 0.7, 0.8)
                circos.lines(x, y)
            }
        }
    }
    
}, track.height = 0.2, bg.border = NA)

circos.clear()
@

You can self-define positions of the `necks' of the long transformation lines so that
they will not overlap with the lower text.

<<genome_more_labels, echo = FALSE, out.width = "\\textwidth", out.height = "0.5\\textwidth", fig.width = 14, fig.height = 7, fig.cap = "Position transformation for a lot of text. A) put text on two layers; B) put text on one layer.">>=
par(mfrow = c(1, 2), mar = c(1, 1, 1, 1))
<<genome_more_labels_1>>
text(0.1, 0.1, "A", cex = 1.5)

<<genome_more_labels_2>>
par(mar = c(1, 1, 1, 1))
circos.par(start.degree = 75, canvas.xlim = c(0, 1), canvas.ylim = c(0, 1), gap.degree = 300, cell.padding = c(0, 0, 0, 0), track.margin = c(0, 0))
circos.initializeWithIdeogram(plotType = NULL, chromosome.index = "chr1")
circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), panel.fun = function(region, value, ...) {
    circos.genomicText(region, value, y = 0, labels.column = 1, facing = "clockwise", adj = c(0, 0.5),
        posTransform = posTransform.text, cex = 0.8, niceFacing = F)
}, track.height = 0.1, bg.border = NA)
i_track = get.cell.meta.data("track.index")

circos.genomicPosTransformLines(bed, 
    posTransform = function(region, value) posTransform.text(region, y = 0, labels = value[[1]], cex = 0.8, track.index = i_track),
    direction = "outside"
)

circos.genomicTrackPlotRegion(bed, ylim = c(0, 1), panel.fun = function(region, value, ...) {
    circos.points( (region[[1]] + region[[2]])/2, rep(0.5, nrow(region)), pch = 16)
}, track.height = 0.02, bg.border = NA)

circos.clear()

text(0.1, 0.1, "B", cex = 1.5)
@


\subsection{Genomic density and Rainfall plot}


{\tt circos.genomicDensity} calculates how much a genomic window is covered by regions in {\tt bed}.
The input data can be a single data frame or a list of data frames.

<<eval = FALSE>>=
circos.genomicDensity(bed)
circos.genomicDensity(bed, baseline = 0)
circos.genomicDensity(bed, window.size = 1e6)
circos.genomicDensity(bedlist, col = c("#FF000080", "#0000FF80"))
@

Rainfall plot can be used to visualize distribution of regions. On the plot, y-axis
corresponds to the distance to neighbour regions (log10-based). So if there is a drop-down on
the plot, it means there is a cluster of regions at that area (figure \ref{fig:genomic_rainfall}).
The input data can be a single data frame or a list of data frames.

<<eval = FALSE>>=
circos.genoimcRainfall(bed)
circos.genoimcRainfall(bedlist, col = c("red", "green"))
@

Following example makes a rainfall plot for differentially methylated regions (DMR) and their
genomic densities. In the plot, red corresponds to hyper-methylated DMRs (gain of methylation) and blue corresponds
to hypo-methylated DMRs (loss of methylation) (figure \ref{fig:genomic_rainfall}).

<<genomic_rainfall, out.width = "0.8\\textwidth", fig.cap = "Genomic rainfall plot and densities.">>=
load(paste0(system.file(package = "circlize"), "/extdata/DMR.RData"))
par(mar = c(1, 1, 1, 1))
circos.initializeWithIdeogram(plotType = c("axis", "labels"))

bed_list = list(DMR_hyper, DMR_hypo)
circos.genomicRainfall(bed_list, pch = 16, cex = 0.4, col = c("#FF000080", "#0000FF80"))
circos.genomicDensity(DMR_hyper, col = c("#FF000080"), track.height = 0.1)
circos.genomicDensity(DMR_hypo, col = c("#0000FF80"), track.height = 0.1)
circos.clear()
@

\end{document}
