tests/testthat/test_picx.R

# Per http://r-pkgs.had.co.nz/tests.html

library(hutan)
library(ape)
context("picx")

# Prepare test data
## Make the tree
phy_text = "(((ENSGALP00000006308:296,ENSOANP00000015173:296)Amniota:296,(ENSMODP00000020505:296,ENSMODP00000024085:296)Euteleostomi:296)Euteleostomi:296,(((ENSGGOP00000016325:490,(((((ENSMUSP00000093719:30.66666667,ENSMUSP00000100528:30.66666667)Mus_musculus:30.66666664,ENSMUSP00000097030:61.3333333)Mus_musculus:30.6666667,((ENSPTRP00000018140:7,ENSP00000389103:7)Homininae:22,ENSMMUP00000039495:29)Catarrhini:63)Euarchontoglires:132.6666667,(ENSMMUP00000035786:92,ENSMUSP00000099541:92)Euarchontoglires:132.6666667)Eutheria:132.6666667,(((ENSMMUP00000010107:89.33333333,ENSMMUP00000009270:89.33333333)Macaca_mulatta:89.33333334,(ENSMMUP00000034357:89.33333333,ENSMMUP00000036086:89.33333333)Macaca_mulatta:89.33333334)Eutheria:89.33333333,ENSMMUP00000009522:268)Eutheria:89.33333333)Eutheria:132.6666667)Eutheria:132.6666667,((ENSMMUP00000040814:207.5555556,ENSMMUP00000034588:207.5555556)Eutheria:207.5555556,(ENSMUSP00000100709:207.5555556,ENSMUSP00000100711:207.5555556)Eutheria:207.5555556)Eutheria:207.5555556)Eutheria:132.6666667,((((ENSMMUP00000038489:128.3055556,ENSMMUP00000000543:128.3055556)Eutheria:128.3055556,ENSGGOP00000022232:256.6111111)Eutheria:128.3055556,ENSMMUP00000022209:384.9166667)Theria:128.3055556,(ENSGGOP00000025227:271.1111111,(ENSPTRP00000056551:29,ENSMMUP00000035858:29)Catarrhini:242.1111111)Eutheria:242.1111111)Theria:242.1111111)Theria:132.6666667)Euteleostomi;"

phy = read.tree( text=phy_text )

## make the tip values
x = c(0.0653820383843777, 0.0581881115273513, 0.146262584851995, 0.146165720301541, 0.122051423671061, 0.061922747428527, 0.0720437904893511, 0.071819687505455, 0.116688417148195, 0.0449859561687792, 0.128921905870532, 0.125443194067582, 0.08165720212685811, 0.134465324836396, 0.341592962693034, 0.161406059961964, 0.173168792917929, 0.138992117848763, 0.190340915795125, 0.145133341227623, 1, 0.217477395678853, 0.995751317106596, 0.854659489405832, 1, 0.769083978252658, 1, 0.619784102971771, 0.882213750595532)

#names(x) = c("ENSGALP00000006308", "ENSOANP00000015173", "ENSMODP00000020505", "ENSMODP00000024085", "ENSGGOP00000016325", "ENSMUSP00000093719", "ENSMUSP00000100528", "ENSMUSP00000097030", "ENSPTRP00000018140", "ENSP00000389103", "ENSMMUP00000039495", "ENSMMUP00000035786", "ENSMUSP00000099541", "ENSMMUP00000010107", "ENSMMUP00000009270", "ENSMMUP00000034357", "ENSMMUP00000036086", "ENSMMUP00000009522", "ENSMMUP00000040814", "ENSMMUP00000034588", "ENSMUSP00000100709", "ENSMUSP00000100711", "ENSMMUP00000038489", "ENSMMUP00000000543", "ENSGGOP00000022232", "ENSMMUP00000022209", "ENSGGOP00000025227", "ENSPTRP00000056551", "ENSMMUP00000035858")

# Calculate the contrasts under different methods
pics = data.frame(
	ape_pic=ape::pic( x, phy, scaled=TRUE, var.contrasts=FALSE ),
	brownian_picx=picx( x, phy, model_method="BM" ) ,
	ou_picx=picx( x, phy, model_method="OU" ) 
)

r2 = summary(lm(brownian_picx~ape_pic, pics))$r.squared
test_that("hutan::picx and ape::pic give similar result",{
	expect_true( r2 > 0.95 )
})
caseywdunn/hutan documentation built on June 12, 2020, 3:12 a.m.