Hi Yan,
You're correct: multi2di() calls rtree(, equiprob = TRUE) which generates the
random topology. Actually, the unlabelled topologies are generated with equal
probabilities. For the labelled topologies, this is a bit more complicated. In
the case n = 4 tips, there are two unlabelled topologies:
(((,),),);
((,),(,));
with 12 possible arrangements of the tip labels for the first topology and 3
for the second one. A solution is to use different probabilities when splitting
the set of tips. Let's say the n tips are split into 2 subsets of sizes a and b
(with n > 2, a + b = n, and both a and b > 0). Then the weight of the pair
(a,b) is calculated with:
howmanytrees(a) * howmanytrees(b) * choose(n, a)
This seems to give correct frequencies for n = 4. It's more complicated to
check with larger values of n because of the large numbers of topologies, but
this looks OK.
Cheers,
Emmanuel
PS: thank you for spelling mistake on ape's homepage: this must have been there
since 2012.
----- Le 20 Oct 20, à 4:23, Yan Wong y...@pixie.org.uk a écrit :
> Thanks to Emmanuel for adding the equiprob=TRUE option to multi2di. However,
> I’m
> not convinced that it’s working correctly - I still get more balanced trees
> than unbalanced trees in ape 5.4.1
>
>> library(phytools)
>> library(phangorn)
>> library(ape)
>> packageVersion("ape")
> [1] '5.4.1'
>
>> reps <- do.call(c, lapply(1:100000, function(x)
> + multi2di(starTree(c('a','b','c','d')), equiprob=TRUE)))
>> balanced = read.tree(text="((a,c),(b,d));”) # test one of 15 4-tomy
>> resolutions
>> test <- sapply(reps, function(x) RF.dist(x, balanced, rooted=TRUE) == 0)
>> mean(test) # should be ~ 1/15 ~ 0.0667
> [1] 0.16625
>
> By the way, there is a spelling mistake on the ape homepage too ("ape souce
> code” - is missing an “r” in “source”).
>
> Cheers
>
> Yan
> [[alternative HTML version deleted]]
>
> _______________________________________________
> R-sig-phylo mailing list - R-sig-phylo@r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
> Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/
_______________________________________________
R-sig-phylo mailing list - R-sig-phylo@r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Searchable archive at http://www.mail-archive.com/r-sig-phylo@r-project.org/
