Description Usage Arguments Details Value Note See Also Examples
View source: R/kpPlotTranscripts.R
Plot gene transcripts on the genome, with options to add strand markers and to differentiate between coding and non-coding exons.
1 2 3 4 5 6 7 8 9 | kpPlotTranscripts(karyoplot, data, y0=NULL, y1=NULL, non.coding.exons.height=0.5,
detail.level=2,
add.strand.marks=TRUE, mark.height=0.20, mark.width=1, mark.distance=4,
add.transcript.names=TRUE, transcript.names=NULL, transcript.name.position="left", transcript.name.cex=1,
col="black", border=NULL, coding.exons.col=NULL, coding.exons.border.col=NULL,
non.coding.exons.col=NULL, non.coding.exons.border.col=NULL,
introns.col=NULL, marks.col=NULL, transcript.name.col=NULL,
ymax=NULL, ymin=NULL, r0=NULL, r1=NULL,
data.panel=1, clipping=TRUE, ...)
|
karyoplot |
(a |
data |
(a |
y0 |
(numeric) The bottom of the transcripts in the y axis. It can have a different value for each transcript and values will be recycled if needed. If null, it will be set to the minimum y value in the data.panel, usually 0. (Defaults to NULL) |
y1 |
(numeric) The top of the transcripts in the y axis. It can have a different value for each transcript and values will be recycled if needed. If null, it will be set to the maximum y value in the data.panel, usually 1. (Defaults to NULL) |
non.coding.exons.height |
(numeric) The height of the non.coding exons relative to the transcript height. For example, if 0.5, non-coding exons will have a height half the size of the coding ones. (default 0.5) |
detail.level |
(numeric: 1 or 2) The detail level of the transcript representation: 1 will plot only boxes representing the transcripts, 2 will plot detailed structure of the trasncripts (coding and non-coding exons and introns). (Defaults to 2) |
add.strand.marks |
(boolean) Whether strand marks should be plotted or not. Strand marks are small arrows along the introns (or whole transcripts if detail level=1). (defaults to TRUE) |
mark.height |
(numeric) The height of the strand marks in "coding exons heights", that is, if mark.height is 0.5, the mark will have a height of half the height of an exon. (defaults to 0.2) |
mark.width |
(numeric) The width of the strand marks, in mark heights. mark.width=1 will produce arrow heads with a slope pf 45 degrees. A value higher than 1 will produce smaller angles and a value below 1 larger angles with more vertical lines. (defaults to 1, 45 degrees) |
mark.distance |
(numeric) The distance between marks, in mark widths. A distance of 2, will add a space of 2*mark.width between consecutive marks. (defaults to 4) |
add.transcript.names |
(boolean) Whether to add transcript names to the plot (defailts to TRUE) |
transcript.names |
(named character) A named character vector with the labels of the transcripts. If not null, it will be used as a dictionary, so transcript ids should be names and desired labels the values. If NULL, the transcript ids will be used as labels. (defaults to null) |
transcript.name.position |
(character) The position of the text relative to the rectangle. Can be "left", "right", "top", "bottom" or "center". (Defaults to "left") |
transcript.name.cex |
(numeric) The cex value to plot the transcript labels. (defaults to 1) |
col |
(color) The color of the transcripts and labels. It is possible to specify different colors for each element class (transcript names, exons, strand marks...). All elements with no explicit color will be plotted using col. (defaults to "black) |
border |
(color) The color of the border of rectangles representing genes, transcripts and exons. Every element class may have its own specific color using the appropiate parameters. The ones with no explicit color will use border. At the same time, if border is NULL, it will default to col. (fesults to NULL) |
coding.exons.col |
(color) The fill color of the rectangles representing the coding exons. If NULL, it will use col. (defaults to NULL) |
coding.exons.border.col |
(color) The color of the border of the coding exons. If NULL, it will use border. (defaults to NULL) |
non.coding.exons.col |
(color) The fill color of the rectangles representing the non-coding exons. If NULL, it will use col. (defaults to NULL) |
non.coding.exons.border.col |
(color) The color of the border of the non-coding exons. If NULL, it will use border. (defaults to NULL) |
introns.col |
(color) The color of the lines representing the introns. If NULL, it will use col. (defaults to NULL) |
marks.col |
(color) The color of the arrows representing the strand. If NULL, it will use col. (defaults to NULL) |
transcript.name.col |
(color) The color of the transcript labels. If NULL, it will use col. (defaults to NULL) |
ymax |
(numeric) The maximum value of |
ymin |
(numeric) The minimum value of |
r0 |
(numeric) r0 and r1 define the vertical range of the data panel to be used to draw this plot. They can be used to split the data panel in different vertical ranges (similar to tracks in a genome browser) to plot differents data. If NULL, they are set to the min and max of the data panel, it is, to use all the available space. (defaults to NULL) |
r1 |
(numeric) r0 and r1 define the vertical range of the data panel to be used to draw this plot. They can be used to split the data panel in different vertical ranges (similar to tracks in a genome browser) to plot differents data. If NULL, they are set to the min and max of the data panel, it is, to use all the available space. (defaults to NULL) |
data.panel |
(numeric) The identifier of the data panel where the data is to be plotted. The available data panels depend on the plot type selected in the call to |
clipping |
(boolean) Only used if zooming is active. If TRUE, the data representation will be not drawn out of the drawing area (i.e. in margins, etc) even if the data overflows the drawing area. If FALSE, the data representation may overflow into the margins of the plot. (defaults to TRUE) |
... |
The ellipsis operator can be used to specify any additional graphical parameters. Any additional parameter will be passed to the internal calls to the R base plotting functions. |
This is one of the high-level, or specialized, plotting functions of
karyoploteR. It takes a list with the transcripts, coding and non-coding
exons and creates a traditional boxes and line representation of the
transcripts. With y0 and y1, it is possible to specify a different vertical
position and different height for each transcript.
It can add little arrows,
strand marks, along the introns to show the transcript strand. The marks
appearance can be customized using 4 different parameters specifying the
height (relative to the height of the transcript), width of each mark
(relative to the height), distance between marks (relative to the width) and
color of the marks. Marks are centered on the space they have available and
if the available space to too tight for a single mark, no mark will be
plotted. The direction of the marks is based on the transcript strand as
specified in data$transcripts
object.
Two detail levels are available: detail.level=1
will represent
transcripts as solid boxes (optionally with strand marks along the whole
transcript); detail.level=2
will represent the internal structure of
the transcripts -coding and non coding exons, introns, and optionally the
strand marks only in the introns.
Returns the original karyoplot object, unchanged.
IMPORTANT: The direction of the strand marks is taken from the strand
information in the data$transcripts
object. If transcripts have no
strand information there (they have strand="*"
), no marks will be
drawn.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 |
#Build the data objet expected by transcripts
transcripts <- c(toGRanges("chr1", 100, 1000),
toGRanges("chr1", 1500, 3000))
names(transcripts) <- c("T1", "T2")
strand(transcripts) <- c("+", "-")
coding.exons <- list("T1"=c(toGRanges("chr1", 200, 300),
toGRanges("chr1", 500, 800),
toGRanges("chr1", 900, 950)),
"T2"=c(toGRanges("chr1", 2200, 2300),
toGRanges("chr1", 2500, 2510),
toGRanges("chr1", 2700, 2800)))
non.coding.exons <- list("T1"=c(toGRanges("chr1", 100, 199),
toGRanges("chr1", 951, 1000)),
"T2"=c(toGRanges("chr1", 1500, 1700),
toGRanges("chr1", 1900, 1950),
toGRanges("chr1", 2100, 2199),
toGRanges("chr1", 2801, 3000)))
data <- list(transcripts=transcripts, coding.exons=coding.exons, non.coding.exons=non.coding.exons)
#Create a simple example plot
karyoplot <- plotKaryotype(zoom=toGRanges("chr1", 0, 3200))
kpAddBaseNumbers(karyoplot, tick.dist = 400)
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0, r1=0.1)
#Create a plot with different variants of the transcripts
karyoplot <- plotKaryotype(zoom=toGRanges("chr1", 0, 3200))
kpAddBaseNumbers(karyoplot, tick.dist = 400)
#Standard
kpPlotTranscripts(karyoplot, data=data, r0=0, r1=0.1)
#Customize colors
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.11, r1=0.21, transcript.name.position = "right", non.coding.exons.col = "#888888", non.coding.exons.border.col = "#888888", marks.col="#CC6666", transcript.name.col="#CC6666")
#Change vertical position and transcript height
kpPlotTranscripts(karyoplot, data=data, y0=c(0, 0.4), y1=c(0.2, 1), r0=0.25, r1=0.8, add.transcript.names = TRUE, transcript.name.position = "left", add.strand.marks = TRUE, clipping=FALSE)
#Change detail level, colors and transcript names
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.9, r1=1, detail.level = 1, add.transcript.names = TRUE, transcript.name.position = "top", add.strand.marks = TRUE, clipping=FALSE, col="blue", marks.col="white", transcript.names=list(T1="Transcript 1", T2="Transcript 2"))
#Create a plot with different variants of the strand marks
karyoplot <- plotKaryotype(zoom=toGRanges("chr1", 0, 3200))
kpAddBaseNumbers(karyoplot, tick.dist = 400)
#Standard
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0, r1=0.1)
#No marks
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.15, r1=0.25, add.strand.marks=FALSE)
#Change the strand marks height
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.3, r1=0.4, mark.height=1)
#Change the mark width
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.45, r1=0.55, mark.width=2)
#Change the mark distance
kpPlotTranscripts(karyoplot, data=data, y0=0, y1=1, r0=0.6, r1=0.7, mark.distance=1.5)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.