Friday, November 23, 2012

Bug fix for phyl.RMA

A user pointed out that phyl.RMA did not work (specifically, returned the incorrect signed RMA slope) when the relationship between x & y is negative. This is, in fact, correct. I have now fixed this problem by changing the following line of code:
beta1<-sqrt(temp$R[2,2]/temp$R[1,1])
to:
beta1<-sign(temp$R[1,2])*sqrt(temp$R[2,2]/temp$R[1,1])
the difference being that in the latter case I pre-multiply the evolutionary variance ratio by the sign of the interspecific covariance. Updated code is here and will be in the next version of phytools.

This update introduces the problem that our statistic for hypothesis testing the RMA regression slope, based on equation [8] of McArdle (1988), is only defined for positive RMA regression slopes. To address this, I have boldly assumed that for h0: β < 0 the null distribution is symmetric about 0 & consequently if b & β have the same sign, I use equation [8] on the absolute values of b & β - otherwise I set T to zero & return a warning. This (at least) has the desirable effect of using McArdle equation [8] without modification with b and β are both > 0 (typical of hypothesis tests about allometry).

Note that, unlike in LS regression, for RMA we cannot normally test h0: β = 0. This is because β = 0 implies that σy = 0.

No comments:

Post a Comment