10 Plotting alongside phylogeny
Jupyter interactive version:
This guide shows the files to plot idiograms alongside a phylogeny
10.1 Load package
visit gitlab for installation instructions https://gitlab.com/ferroao/idiogramFISH
#load package
library(idiogramFISH) 10.2 ggtree of iqtree and monocentrics
idiogramFISH comes with two trees and data.frames with chr. and marks’ data for the correspondent OTUs, first we will plot a tree produced with iqtree (Nguyen et al., 2015)
Load the iqtree:
We will use phytools for that (Revell, 2012)
require(ggplot2)
require(phytools)
require(ggpubr)
require(grid) #pushViewport
require(ggtree)
# list.files(system.file('extdata', package = 'my_package') )
# find path of iqtree file
iqtreeFile <- system.file("extdata", "eightSpIqtree.treefile", package = "idiogramFISH")
# load file as phylo object
iqtreephylo <- read.newick(iqtreeFile) # phytools
# transform tree
iqtreephyloUM <- force.ultrametric(iqtreephylo, method= "extend") # phytoolsMake a ggtree (Yu et al., 2018)
ggtreeOf8 <- ggtree(iqtreephyloUM) + geom_tiplab(size=6)Modify optionally graphical parameters with ggplot and ggpubr: (Wickham, 2016; Kassambara, 2019)
gbuil2 <- ggplot_build(ggtreeOf8) # get ggplot_built
gtgbuild <- ggplot_gtable(gbuil2) # get gtable from ggplot_built
gtgbuild$layout$clip[gtgbuild$layout$name == "panel"] <- "off" # modify gtable
ggtreeOf8b <- as_ggplot(gtgbuild) # back to ggplot
gtgbuildgg2 <- ggtreeOf8b + theme(plot.margin = unit(c(1,9.5,3,1.5), "cm") ) # top right bottom left - modify marginsOrder OTUs in data.frame of chr. data
Apply order of phylogeny to the data.frame
# Let's get the order of species in tree
ggtreeOf8TIPS <- ggtreeOf8$data[which(ggtreeOf8$data$isTip),]
desiredOrder <- rev(ggtreeOf8TIPS[order(ggtreeOf8TIPS$y),]$label)
# make a vector without missing OTUs
desiredFiltered <- intersect(desiredOrder,allChrSizeSample$OTU)
# establish desired order
allChrSizeSample$OTU <- factor(allChrSizeSample$OTU, levels = desiredFiltered)
# order
allChrSizeSample <- allChrSizeSample[order(allChrSizeSample$OTU),]Now we have to establish where are the OTUs in the tree, that don’t have chr. data
# Establish position of OTUs before missing data OTUs
matchres <- match(desiredOrder,desiredFiltered)
matchres[is.na(matchres)] <- "R"
reps <- rle(matchres)
posOTUsBeforeMissing <- as.numeric(matchres[which(matchres=="R")-1][which(matchres[which(matchres=="R")-1]!="R")] )
# This are the OTUs that come before missing chr. data OTUs
BeforeMissing <- desiredFiltered[posOTUsBeforeMissing]
# This is the amount of missing OTUs, spaces to add (ghost karyotypes)
valuesOfMissRepsBefore <- reps$lengths[which(reps$values=="R")]Plotting
Now we are ready to plot adding those arguments for addMissingOTUAfter and missOTUspacings
# plot to png file
png(file="firstplot.png" ,width=962,height=962 )
par(omi=rep(0,4),
mar=c(0,1,2,1),
mfrow=c(1,2) ) # one row two columns
par(fig=c(0,.3,0,1)) # location of left ghost plot
plot.new() # ghost plot to the left
par(fig=c(.3,1,0,1)) # location of right plot
plotIdiograms(dfChrSize = allChrSizeSample, # data.frame of Chr. Sizes
dfMarkPos = allMarksSample, # d.f. of Marks (inc. cen. marks)
dfMarkColor = mydfMaColor, # d.f. of mark characteristics
squareness = 4, # squareness of vertices
lwd.chr=.5, # width of lines
orderChr = "name", # order chr. by name
centromereSize = 1.3, # apparent cen. size
chrWidth =.75, # width of chr.
chrSpacing = .25, # horizontal spacing of chr.
indexIdTextSize=.4, # font size of indices and chr. names
karHeight = 4.8, # karyotype vertical relative size without spacing
karHeiSpace = 6.5, # karyotype vertical relative size with spacing
nameChrIndexPos=4, # move the name of chr. indexes to left
morpho="both", # add chr. morphology
chrIndex = "both", # add chr. indices
karIndex = TRUE, # add karyotype indices
yTitle= "", # remove units title of ruler
markLabelSpacer = 0 # spaces from rightmost chr. to legend
,ylimTopMod = -.1 # modify ylim top margin
,ylimBotMod=2.6 # modify ylim bottom margin
,rulerPos = -0.5 # position of rulers
,rulerNumberSize = .35 # font size of ruler number
,rulerNumberPos = .4 # position of ruler numbers
,ruler.tck = -.004 # tick size and orient.
,addMissingOTUAfter = BeforeMissing # OTUs after which there are ghost karyotypes - empty spaces
,missOTUspacings = valuesOfMissRepsBefore # number of ghost karyotypes
)
# plot to the left the ggtree
pushViewport(viewport(layout = grid.layout(1, 2)))
pushViewport(viewport(layout.pos.col = 1, layout.pos.row = 1))
print(gtgbuildgg2,newpage=F)
# close png
dev.off()
10.3 plot of revBayes tree and holocentrics
Now we are going to plot a tree from revBayes (Höhna et al., 2017)
First, load the revBayes tree:
require(ggplot2)
require(phytools)
require(ggpubr)
require(grid) #pushViewport
require(ggtree)
require(treeio)
# find path of iqtree file
revBayesFile <- system.file("extdata", "revBayesTutorial.tree", package = "idiogramFISH")
# load file as phylo object
revBayesPhylo <- read.beast(revBayesFile) # ggtree or treeio
# transform tree
revBayesPhyloUM <- force.ultrametric(revBayesPhylo@phylo, method= "extend") # phytoolsOrder OTUs in data.frame of chr. data
# Get order of OTUs in tree
is_tip <- revBayesPhyloUM$edge[,2] <= length(revBayesPhyloUM$tip.label)
ordered_tips <- revBayesPhyloUM$edge[is_tip, 2]
desiredorderRevB <- rev(revBayesPhyloUM$tip.label[ordered_tips])
# Ceate some holocentrics' data
allChrSizeSampleHolo <- allChrSizeSample
allChrSizeSampleHolo <- allChrSizeSampleHolo[,c("OTU","chrName","longArmSize")]
colnames(allChrSizeSampleHolo)[which(names(allChrSizeSampleHolo)=="longArmSize")]<-"chrSize"
allMarksSampleHolo <- allMarksSample
allMarksSampleHolo <- allMarksSampleHolo[which(allMarksSampleHolo$chrRegion!="cen"),]
allMarksSampleHolo <- allMarksSampleHolo[c("OTU","chrName","markName","markDistCen","markSize")]
colnames(allMarksSampleHolo)[which(names(allMarksSampleHolo)=="markDistCen")] <- "markPos"
allMarksSampleHolo[which(allMarksSampleHolo$markName=="5S"),]$markSize <- .5
# Apply order of phylogeny to data.frame
# make a vector without missing OTUs
desiredFiltered <- intersect(desiredorderRevB,allChrSizeSampleHolo$OTU)
# establish desired order
allChrSizeSampleHolo$OTU <- factor(allChrSizeSampleHolo$OTU, levels = desiredFiltered)
# order
allChrSizeSampleHolo <- allChrSizeSampleHolo[order(allChrSizeSampleHolo$OTU),]
Now we have to establish where are the OTUs in the tree, that don’t have chr. data
# Establish position of OTUs before missing data OTUs
matchres <- match(desiredorderRevB,desiredFiltered)
matchres[is.na(matchres)] <- "R"
reps <- rle(matchres)
posOTUsBeforeMissing <- as.numeric(matchres[which(matchres=="R")-1][which(matchres[which(matchres=="R")-1]!="R")] )
# This are the OTUs that come before missing chr. data OTUs
BeforeMissingPlot2 <- desiredFiltered[posOTUsBeforeMissing]
# This is the amount of missing OTUs, spaces to add (ghost karyotypes)
valuesOfMissRepsBeforePlot2 <- reps$lengths[which(reps$values=="R")]Plotting
Now we are ready to plot adding those arguments for addMissingOTUAfter and missOTUspacings
# plot to png file
png(file=paste0("secondplot.png" ),width=962,height=700)
{
par(omi=rep(0,4) , mar=c(0,0,0,0), mfrow=c(1,2))
par(fig=c(0,.27,0,1))
par(mar=c(2,0,2,0)) # b l t r
plot(revBayesPhyloUM)
par(fig=c(0.27,1,0,1), new=TRUE)
par(mar=c(0,0,0,0)) # b l t r
plotIdiograms(allChrSizeSampleHolo, # chr. size data.frame
dfMarkPos = allMarksSampleHolo, # data.frame of marks' positions
dfMarkColor = mydfMaColor, # d.f. of mark characteristics
squareness = 4, # vertices squareness
karHeight = 2.8, # karyotype height
karHeiSpace = 4.5, # vertical size of kar. including spacing
yTitle="",
karIndex = TRUE, # add karyotype index
indexIdTextSize=.4 # font size of indices and chr. names
,addMissingOTUAfter = BeforeMissingPlot2 # add ghost OTUs after these names
,missOTUspacings = valuesOfMissRepsBeforePlot2 # how many ghosts, respectively
,lwd.chr=.5 # line width
,markLabelSpacer = 0 # dist. of legend to rightmost chr.
,legendWidth = 2.3 # width of square or dots of legend
,rulerPos = -0.5 # position of ruler
,rulerNumberSize = .35 # font size of number of ruler
,rulerNumberPos = .4 # position of ruler number
,ruler.tck=-.004 # tick of ruler size and orient.
,ylimTopMod = -4 # modify ylim of top
,ylimBotMod = -4 # modify ylim of bottom
,xlimRightMod = 3 # modify xlim right argument
,xModifier = 100 # separ. among chromatids
)
}
# close png
dev.off()
10.4 plot of revBayes tree and holocentrics and monocentrics
Create data.frames with both types of karyotypes (Wickham, 2011)
require(plyr)
# Select this OTU from the monocen.
monosel<-c("Species_F","Species_C","Species_A")
# chr.
allChrSizeSampleSel <- allChrSizeSample [which(allChrSizeSample$OTU %in% monosel ),]
# marks
allMarksSampleSel <- allMarksSample [which(allMarksSample$OTU %in% monosel ),]
# Select the others from the holocen.
holosel <- setdiff(unique(allChrSizeSampleHolo$OTU),monosel)
# chr.
allChrSizeSampleHoloSel <- allChrSizeSampleHolo[which(allChrSizeSampleHolo$OTU %in% holosel ),]
# marks
allMarksSampleHoloSel <- allMarksSampleHolo [which(allMarksSampleHolo$OTU %in% holosel ),]
# merge chr d.fs
mixChrSize <- plyr::rbind.fill(allChrSizeSampleSel,allChrSizeSampleHoloSel)
# merge marks' d.fs
mixMarks <- plyr::rbind.fill(allMarksSampleSel,allMarksSampleHoloSel)
Order data.frame and determine missing karyotypes
# make a vector without missing OTUs
desiredFiltered <- intersect(desiredorderRevB, mixChrSize$OTU)
# establish desired order
mixChrSize$OTU <- factor(mixChrSize$OTU, levels = desiredFiltered)
# order data.frame
mixChrSize <- mixChrSize[order(mixChrSize$OTU),]
# Establish position of OTUs before missing data OTUs
matchres <- match(desiredorderRevB,desiredFiltered)
matchres[is.na(matchres)] <- "R"
reps <- rle(matchres)
posOTUsBeforeMissing <- as.numeric(matchres[which(matchres=="R")-1][which(matchres[which(matchres=="R")-1]!="R")] )
# This are the OTUs that come before missing chr. data OTUs
BeforeMissingPlot2 <- desiredFiltered[posOTUsBeforeMissing]
# This is the amount of missing OTUs, spaces to add (ghost karyotypes)
valuesOfMissRepsBeforePlot2 <- reps$lengths[which(reps$values=="R")]Plotting
Now we are ready to plot adding those arguments for addMissingOTUAfter and missOTUspacings
# plot to png file
png(file=paste0("thirdplot.png" ),width=1100,height=1000)
{
par(omi=rep(0,4) , mar=c(0,0,0,0), mfrow=c(1,2))
par(fig=c(0,.25,0,1))
par(mar=c(1,0,0,0))
plot(revBayesPhyloUM)
par(fig=c(0.25,1,0,1), new=TRUE)
par(mar=c(0,0,0,0))
plotIdiograms(mixChrSize, # chr. size data.frame
dfMarkPos = mixMarks, # data.frame of marks' positions (inc. cen. marks)
dfMarkColor = mydfMaColor, # d.f. of mark characteristics
origin="b", # position measured from bottom of chr.
karHeight = 2.8, # vertical size of kar. including spacing
karHeiSpace = 4.5, # vertical size of kar. including spacing
squareness = 5, # vertices squareness
chrSpacing = .25, # horizontal spacing among chr.
yTitle = "",
karIndex = TRUE # add karyotype index
,indexIdTextSize=.4 # font size of indices and chr. names
,distTextChr=0.7 # dist. among chr. and chr. name
,addMissingOTUAfter = BeforeMissingPlot2 # add ghost OTUs after these names
,missOTUspacings = valuesOfMissRepsBeforePlot2 # how many ghosts, respectively
,lwd.chr=.5 # line width
,markLabelSpacer = 0 # dist. of legend to rightmost chr.
,legendWidth = 2 # width of square or dots of legend
,ylimTopMod = -2 # modify ylim of top
,ylimBotMod = -2 # modify ylim of bottom
,rulerPos = -0.5 # position of ruler
,rulerNumberSize = .35 # font size of number of ruler
,rulerNumberPos = .4 # position of ruler number
,ruler.tck=-.004 # ruler tick size and orient.
,OTUfont=3 # italics
,OTUfamily="Courier New" # for OTU name
,xModMonoHoloRate = 5 # factor (quotient) to shrink separation of chromatids of holocen.
)
}
# close png
dev.off()