-
Notifications
You must be signed in to change notification settings - Fork 6
/
plot_mtdna_tree.R
117 lines (108 loc) · 4.68 KB
/
plot_mtdna_tree.R
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
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
library(tidyverse)
library(ape)
library(treeio)
library(ggtree)
# set working directory
setwd("")
# import tree
mt.tree <- read.iqtree("partitions.nex.treefile")
# remove outgroup
mt.tree.pruned <- drop.tip(mt.tree, "JN632674")
# check the node labels
#plot.phylo(as.phylo(mt.tree.pruned), no.margin = TRUE, cex = 0.75)
#nodelabels(cex = 0.75, frame = "circle")
# create the basic plot
mt.t1 <- ggtree(mt.tree.pruned,
layout = "roundrect",
ladderize = TRUE,
right = TRUE) +
# add scale bar
geom_treescale(x = 0, y = 0, fontsize = 2.3) +
# add tip labels
geom_tiplab(size = 2.3, align = TRUE, linesize = 0.1) +
# add node labels (support values)
geom_nodepoint(aes(color = UFboot, subset = !is.na(as.numeric(UFboot))),
size = 1.5) +
# annotate clades
# note: Some clade annotations are knowngly wrong due to the placement of a
# few individuals outside its assigned (sub)species. To get the tree
# annotations as shown in the paper requires image editing. This was
# done just to facilitate the image editing process later on.
geom_cladelabel(node = 54, label = "Northern",
color = c("#000000", "#000000"),
offset = 0.02, offset.text = 0.003,
barsize = 0, angle = -90, hjust = 0.5,
fontsize = 3.2, align = TRUE) +
geom_cladelabel(node = 61, label = "West African",
color = c("#F0E442", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 55, label = "Kordofan",
color = c("#E69F00", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 65, label = "Nubian",
color = c("#D55E00", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 72, label = "Reticulated",
color = c("#CC79A7", "#000000"),
offset = 0.02, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 3.2, align = TRUE) +
geom_cladelabel(node = 79,
label = expression(paste("Masai ", italic("s. l."))),
color = c("#000000", "#000000"),
offset = 0.02, offset.text = 0.003,
barsize = 0, angle = -90, hjust = 0.5,
fontsize = 3.2, align = TRUE) +
geom_cladelabel(node = 80, label = "Luangwa",
color = c("#006046", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 87,
label = expression(paste("Masai ", italic("s. str."))),
color = c("#009E73", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 96, label = "Southern",
color = c("#000000", "#000000"),
offset = 0.02, offset.text = 0.003,
barsize = 0, angle = -90, hjust = 0.5,
fontsize = 3.2, align = TRUE) +
geom_cladelabel(node = 91, label = "South African",
color = c("#56B4E9", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
geom_cladelabel(node = 96, label = "Angolan",
color = c("#0072B2", "#000000"),
offset = 0.0125, offset.text = 0.003,
barsize = 1, angle = -90, hjust = 0.5,
fontsize = 2.3, align = TRUE) +
# add legend
scale_colour_distiller("UFBoot support",
palette = "Spectral",
direction = 1) +
theme(legend.title = element_text(face = "bold", size = 7),
legend.text = element_text(size = 6),
legend.key.height = unit(0.75, "line"),
legend.position = "bottom")
# flip clades
mt.t2 <- flip(mt.t1, 52, 78) %>%
flip(79, 96) %>%
flip(34, 81) %>%
flip(27, 53) %>%
flip(54, 72) %>%
flip(69, 71)
# check node numbers interactively
#identify(mt.t2)
# show plot
mt.t2
# save plot in '.tif' format
ggsave("mtdna_tree.tiff", width = 87, height = 170, units = "mm", dpi = 300)