Wednesday, October 23, 2024

An Ntip S3 method for the "cophylo" object

A few days ago a phytools user wrote to me in frustration because the S3 method Ntip (exported by the ape package) didn’t work when applied to a "cophylo" object from phytools.

The reason this is the case is because a "cophylo" object doesn’t inherit the class "phylo" – that is, it doesn’t have "phylo" (or "multiPhylo") among its class attributes… which is good, in fact, but it is also not conformed like a "phylo" or "multiPhylo" object.

Just for fun, let’s see how easy it is to create a new S3 method for a different class type, when the generic already exists.

library(phytools)
data("wasp.trees")
wasp.trees
## 2 phylogenetic trees
data("wasp.data")
wasp.data
##        Pleistodontes                Sycoscapter
## 1      P._greenwoodi             S._4_{obliqua}
## 2  P._xanthocephalus        S._2_{aff._obliqua}
## 3        P._plebejus   S._14_{hesperidiiformis}
## 4           P._rieki          S._12_{xylosycia}
## 5         P._blandus         S._11_{glandifera}
## 6         P._regalis         S._8_{pleurocarpa}
## 7          P._nitens          S._10_{crassipes}
## 8    P._schizodontes          S._9_{triradiata}
## 9      P._imperialis         S._15_{rubiginosa}
## 10      P._athysanus          S._1_{brachypoda}
## 11 P._astrabocheilus         S._7_{subpuberula}
## 12       P._proximus        S._6_{lilliputiana}
## 13    P._macrocainus         S._5_{cerasicarpa}
## 14       P._cuneatus            S._3_{playpoda}
## 15      P._froggatti S._australis_{macrophylla}

These are the trees and data that we’re going to use to create our "cophylo" object.

First, let’s see what Ntip does with a "multiPhylo" object like wasp.trees.

Ntip(wasp.trees)
## [1] 19 15

This is just the number of tips in our first & second trees of the object, respectively.

Now to create our "cophylo" object, we can proceed to run phytools::cophylo as follows.

obj<-cophylo(wasp.trees[[1]],wasp.trees[[2]],
  wasp.data,print=TRUE)
## Rotating nodes to optimize matching...
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 359.577777777778
## objective: 357.777777777778
## objective: 357.777777777778
## objective: 345.111111111111
## objective: 345.111111111111
## objective: 345.111111111111
## objective: 322.777777777778
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 322.777777777778
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## objective: 201.17728531856
## Done.
plot(obj,link.type="curved",type="cladogram",lwd=3,
  fsize=0.7,part=0.45)

plot of chunk unnamed-chunk-32

Next, let’s see what happens when we don’t have an method for the object class & try to use Ntip.

Ntip(obj)
# Error in UseMethod("Ntip") : 
#   no applicable method for 'Ntip' applied to an object of class "cophylo"

Not applicable method? No problem. Let’s create one.

First, check the function definition of the Ntip. Ours should match:

args(Ntip)
## function (phy) 
## NULL

Now, just make it.

Ntip.cophylo<-function(phy) sapply(phy$trees,Ntip)

Does it work? Sure!

Ntip(obj)
## [1] 19 15

That’s all folks!

No comments:

Post a Comment

Note: due to the very large amount of spam, all comments are now automatically submitted for moderation.