cboettig/geiger — Issue #1: Bug in geiger::fitContinuous with meserr 
Author: cboettig
Created: 2012-02-02 20:24:03+00:00
State: open

I'm not sure if this is the best place to make a bug report, but I thought others might be interested.

I was trying to replicate, with geiger, a demonstration of estimating lambda with measurement error from Liam's blog:

http://phytools.blogspot.com/2011/11/including-measurement-error-in.html

Replicating this with fitContinuous, I received the following error:

> fitContinuous(tree,xe,model="lambda",meserr=se)
Error in meserr[o, ] : incorrect number of dimensions

The problem is related to this bit of code in fitContinuous, used to determine whether meserr is a single value, a vector, or a matrix:

else if (is.vector(meserr)) {
if (!is.null(names(meserr))) {
o <- match(rownames(td$data), names(meserr))
if (length(o) != ntax)
stop("meserr is missing some taxa from the tree")
meserr <- as.matrix(meserr[o,])

Removing the comma from the final line fixed the problem. The values I receive from geiger after removing the comma are identical to the values I get from phytools::phylosig for the same data, so I am assuming I haven't broken anything important.

The full code I used is at the end of this message, for those interested.

Best,

Matt Johnson

> library(geiger)
Loading required package: ape
Loading required package: MASS
Loading required package: mvtnorm
Loading required package: msm
Loading required package: ouch
Loading required package: subplex
> library(phytools)
Loading required package: mnormt
Loading required package: phangorn
Loading required package: quadprog
Loading required package: igraph
Loading required package: rgl
> source("mjContinuous.R") #fitContinuous with offending comma removed
> set.seed(2)
> gen.lambda=0.7
> tree = drop.tip(birthdeath.tree(b=1,d=0,taxa.stop=501),"501")
> se = sqrt(rchisq(n=500,df=2))
> names(se) = tree$tip.label
>
> x= fastBM(lambdaTree(tree,gen.lambda),sig2=1)
> e = rnorm(n=500,sd=se)
> xe=x+e
>
> fitContinuous(tree,x,model="lambda")
Fitting lambda model:
$Trait1
$Trait1$lnl
[1] -1031.783

$Trait1$beta
[1] 1.221674

$Trait1$lambda
[1] 0.7613218

$Trait1$aic
[1] 2069.567

$Trait1$aicc
[1] 2069.615

$Trait1$k
[1] 3

> fitContinuous(tree,xe,model="lambda")
Fitting lambda model:
$Trait1
$Trait1$lnl
[1] -1154.987

$Trait1$beta
[1] 1.412728

$Trait1$lambda
[1] 0.5495113

$Trait1$aic
[1] 2315.975

$Trait1$aicc
[1] 2316.023

$Trait1$k
[1] 3

> fitContinuous(tree,xe,model="lambda",meserr=se)
Error in meserr[o, ] : incorrect number of dimensions
>
> mjContinuous(tree,xe,model="lambda",meserr=se) #with problematic comma removed
Fitting lambda model:
$Trait1
$Trait1$lnl
[1] -1141.14

$Trait1$beta
[1] 1.170196

$Trait1$lambda
[1] 0.7484273

$Trait1$aic
[1] 2288.281

$Trait1$aicc
[1] 2288.329

$Trait1$k
[1] 3