3
# if(!is.null(phy$maps))
4
## FIXME write the characters/states block (but not matrix block) as well.
6
## FIXME support writing multiphylos, list of multiphylos to nexml
12
#' @param phy a phy object containing simmap phy$maps element,
13
#' from the phytools pacakge
14
#' @param state_ids a named character vector giving the state
15
#' names corresponding to the ids used to refer to each state
16
#' in nexml. If null ids will be generated and states taken from
17
#' the phy$states names.
18
#' @return a nexml representation of the simmap
23
#' phy <- nexml_to_simmap(simmap_ex)
24
#' nex <- simmap_to_nexml(phy)
25
simmap_to_nexml <- function(phy, state_ids = NULL){
28
## Hack to deal with S3 class issues when coercing to S4
29
if(class(phy) == c("simmap", "phylo"))
32
## Create the NeXML object
33
nexml <- as(phy, "nexml")
35
if(!is.null(phy$states)){
36
nexml <- add_characters(data.frame(states = as.integer(as.factor(phy$states))), nexml)
38
## FIXME doesn't have states
39
chars_ids <- get_state_maps(nexml)[[1]]
40
char_id <- names(chars_ids)
41
state_ids <- reverse_map(chars_ids[[1]])
42
# can assume no other states added yet since works on a phy, not nexml
46
if(!is.null(phy$maps))
47
nexml <- simmap_edge_annotations(maps = phy$maps,
49
state_ids = state_ids,
57
simmap_edge_annotations <- function(maps, nexml, state_ids = NULL, char_id = "simmapped_trait"){
60
## if state ids are not given
61
if(is.null(state_ids)){
62
state_ids <- levels(as.factor(names(unlist(maps))))
63
names(state_ids) <- state_ids
67
# Loop over all edges, adding the simmap annotation to each:
68
for(i in 1:length(nexml@trees[[1]]@tree[[1]]@edge)){
69
## FIXME check to assure this is always the correct order??
71
## Read the mapping of the current edge
74
## Generate the list of XML "stateChange" nodes
75
mapping <- lapply(1:length(edge_map), function(j){
76
## A node has an id, a length and a state
77
meta(property = "simmap:stateChange",
78
children = list(meta(property = "simmap:order",
80
meta(property = "simmap:length",
81
content = edge_map[[j]]),
82
meta(property = "simmap:state",
83
content = state_ids[[names(edge_map[j])]] )
88
reconstruction <-meta(property = "simmap:reconstruction",
89
children =c(list(meta(property="simmap:char",
93
## Insert the reconstructions into a <meta> element in each nexml edge
94
nexml@trees[[1]]@tree[[1]]@edge[[i]]@meta <-
95
c(meta(type = "LiteralMeta",
96
property = "simmap:reconstructions",
97
children = list(reconstruction)))
100
## Return the entire nexml object
101
nexml <- add_namespaces(c(simmap = "https://github.com/ropensci/RNeXML/tree/master/inst/simmap.md"),
110
## Returns list of multiPhylo ...
115
#' @param nexml a nexml object
116
#' @return a simmap object (phylo object with a $maps element
117
#' for use in phytools functions).
121
#' phy <- nexml_to_simmap(simmap_ex)
122
#' nex <- simmap_to_nexml(phy)
123
nexml_to_simmap <- function(nexml){
125
## Get the statemap, if available
127
characters <- get_characters(nexml)
129
## loop over trees blocks
130
out <- lapply(nexml@trees, function(trees){
131
phys <- lapply(trees@tree,
133
get_otu_maps(nexml)[[trees@otus]],
134
get_state_maps(nexml)[[1]][[1]]
136
phys <- lapply(phys, characters_to_simmap, characters)
139
class(phys) <- "multiPhylo"
144
if(length(out) == 1){
145
if(length(out[[1]]) > 1){
146
flatten_multiphylo(out)
154
characters_to_simmap <- function(phy, characters){
155
out <- as.character(characters[[1]]) # coerce factor to string
156
names(out) <- rownames(characters)
162
# given the nexml tree:
163
tree_to_simmap <- function(tree, otus, state_maps = NULL){
164
maps <- lapply(tree@edge, function(edge){
166
reconstructions <- sapply(edge@meta, function(x) x@property == "simmap:reconstructions")
167
if(any(reconstructions))
168
reconstruction <- edge@meta[[which(reconstructions)]]@children
170
else { # handle exceptions
171
warning("no simmap data found")
172
return(toPhylo(tree, otus)) ## no simmap found
176
# lapply(reconstruction, function(reconstruction){ # for each reconstruction
177
stateChange <- sapply(reconstruction[[1]]@children, function(x) x@property == "simmap:stateChange")
179
values <- sapply(reconstruction[[1]]@children[which(stateChange)], function(stateChange){
180
# phytools only supports one reconstruction of one character per phy object
181
property <- sapply(stateChange@children, function(x) x@property)
182
names(stateChange@children) <- property # clean labels
183
sapply(stateChange@children, function(x) x@content)
185
out <- as.numeric(values["simmap:length", ])
186
ordering <- as.numeric(values["simmap:order",])
187
if(!is.null(state_maps))
188
states <- state_maps[values["simmap:state", ]]
190
states <- values["simmap:state", ]
193
out <- out[ordering] # sort according to explicit order
201
phy <- toPhylo(tree, otus)
204
## create the rest of the maps elements
207
## Return phylo object
214
#' @title A nexml class R object that includes simmap annotations
215
#' @description A nexml object with simmap stochastic character mapping
216
#' annotations added to the edges, for use with the RNeXML package
217
#' parsing and serializing NeXML into formats that work with the ape and
218
#' phytools packages.
221
#' @format a \code{nexml} instance
222
#' @source Simulated tree and stochastic character mapping based on
223
#' Revell 2011 (doi:10.1111/j.2041-210X.2011.00169.x)
224
#' @author Carl Boettiger
229
## Extend directly with XML representation instead of S4
231
#setClass("simmap:reconstructions",
232
# slots = c(reconstruction = "ListOfreconstruction"))
233
#setClass("ListOfreconstruction", contains = "list")
234
#setClass("simmap:stateChange", contains = "IDTagged")
added
removed
R/simmap.R
