Fix to biplot S3 method for object class "phyl.pca"

I just pushed a fix to the S3 biplot method for objects of class "phyl.pca".

The source of the bug is closely related to my recent post on the discrepancy between what are called 'loadings' in princomp vs. phyl.pca. A few years the gentleman who kindly supplied the S3 method code used by phyl.pca made exactly the inverse mistake of princomp.

First, here is our bi-plot with the existing method:

library(phytools)
tree

## 
## Phylogenetic tree with 26 tips and 25 internal nodes.
## 
## Tip labels:
##  A, B, C, D, E, F, ...
## 
## Rooted; includes branch lengths.

X

##           [,1]        [,2]        [,3]       [,4]       [,5]
## A  0.419211564 -0.07072903 -0.67729412 -1.4275547 -0.8956189
## B  0.423491819 -0.07230281  0.04353057 -1.8236536 -1.5941026
## C  0.674572520 -0.97438998 -0.44582192 -1.2064721 -1.1968349
## D  0.144589136 -0.78968124  0.23939251 -0.1604712  1.2397129
## E  0.037492764 -0.33373486 -1.21897435  0.1363836  1.2665389
## F -0.049170578 -0.24004146 -1.25930976 -0.4204259  1.4010466
## G  0.270753737 -0.16233552 -1.44227716 -0.6901865  1.1583580
## H  0.307947940 -0.17681354  0.80860301  0.4294457  0.3390226
## I  0.121981063 -1.54597206 -0.60546618 -0.6466339  0.9344554
## J -0.420320783  0.21857516 -0.14983841  0.3039046  1.1689236
## K -0.250202949 -0.15091852 -0.26843245  1.0483096  1.4964830
## L  0.470691820 -1.64133733  0.01950930  0.9664008 -1.6606181
## M -0.859279349 -0.79613244  0.21614994 -0.6497791  0.9240529
## N -0.531125300 -0.83481129  0.84552758 -0.8766522  1.6286559
## O -0.362032884 -0.48093998  0.35350781  0.1044871  1.0196886
## P -0.271775777  0.11416429 -0.04792096  0.5975930  1.5266678
## Q -0.241609452  0.15047775  0.10962659  0.4692869  1.4156464
## R  0.250219313 -0.19964343  0.70567779 -0.8044968  1.4153869
## S  0.652031092 -0.57999496  0.41286288  0.2249271  1.5921390
## T  0.259519957 -0.52699138  0.47710290  0.4022593  1.5436005
## U  1.355538164 -0.97106032  0.90405187  0.9945186  2.4721401
## V  1.490160952  0.73566840  0.01132505 -0.6147261  2.8499558
## W  0.778375058  0.59260693 -0.13285445 -0.7414536  3.3113367
## X  0.542971500  0.73316942 -1.09947325 -2.6681790  3.1622164
## Y -0.259842325  0.62505613 -1.01765120 -1.3903978  3.3034612
## Z  0.001442499  0.56993438 -0.88314656 -1.4145177  2.7669227

pca<-phyl.pca(tree,X,method="lambda")
pca

## Phylogenetic pca
## Standard deviations:
##       PC1       PC2       PC3       PC4       PC5 
## 1.3889973 1.2627293 0.9716614 0.7283548 0.6591273 
## Loads:
##              PC1          PC2         PC3         PC4         PC5
## [1,]  0.27723112 -0.288933764  0.88668656  0.12016754  0.19750162
## [2,]  0.29061568  0.369744598 -0.25981002 -0.35072584  0.76702125
## [3,] -0.21999919 -0.266820125  0.29765195 -0.86158989 -0.22242651
## [4,] -0.98858119 -0.005420764  0.08204672  0.04923774  0.11628342
## [5,] -0.01938295  0.965462070  0.24065480 -0.02126167 -0.09560596
## lambda:
## [1] 1.002837

biplot(pca)

Here is our fixed biplot method code:

biplot.phyl.pca<-function(x,...){
    to.do<-list(...)
    if(hasArg(choices)){ 
        choices<-list(...)$choices
        to.do$choices<-NULL
    } else choices<-c(1,2)
    to.do$x<-x$S[,choices]
    to.do$y<-x$Evec[,choices]
    do.call(biplot,to.do)
}

and here is the result we should have obtained:

biplot(pca)

That's it!