I just tried this using the code below, and I can't see a trend. I calculated
the mean tip value across the species in each simulation, and made a box plot
over 10000 simulations. Although sim.char() doesn't do this, one can simply
draw these numbers from a multivariate normal distribution - and so their
expected mean is zero.
rr<-c(0.001, 0.01, 0.1, 1, 10, 100, 1000)
phy<-birthdeath.tree(b=1, d=0, taxa.stop=20)
allMeans<-matrix(nrow=10000, ncol=length(rr))
for(i in 1:length(rr)) {
mm<-rr[i]
ss<-sim.char(phy, as.matrix(rr[i]), nsims=10000)
allMeans[,i]<-apply(ss, 3, FUN=mean)
}
meanmeans<-apply(allMeans, 2, FUN=mean)
cc<-as.numeric(allMeans)
xx<-rep(rr, each=10000)
plot(as.factor(xx), cc)
On May 20, 2011, at 12:32 PM, Dean Adams wrote:
> Hi all,
>
> I'm observing a curious pattern in continuous trait data simulated on a
> phylogeny that I cannot explain. I thought I'd throw it out to the group for
> ideas. My apologies if this was addressed in an earlier thread.
>
> I'm exploring four different R functions for simulating continuous data under
> Brownian Motion: 'rTraitCont' in Ape, 'sim.char' in Geiger,
> 'transformPhylo.sim' from MotMot, and 'fastBM' from Liam Revell's
> PhyTools-beta. For these simulations I generate a tree (in this case a
> perfectly-balanced tree) and simulate 100 data sets on the same phylogeny
> using a particular initial BM rate parameter (sigma). For each simulation, I
> calculate the mean & variance of the simulated tips data, and then summarize
> these across 100 simulations to understand some general properties of the 4
> data-generating functions. I'm using a wide range of initial sigmas to see
> how things fall out (sigma = 0.0001, 0.001. 0.01 ... 1000).
>
> At sigma = 1.0, the mean of the variance in tips data across simulations is
> similar for all 4 methods. And as expected, as sigma increases or decreases,
> the mean variation among tips values also increases or decreases. For each
> sigma value, mean levels of variation are also similar for 3 of the 4
> methods: with rTraitCont having much larger relative levels of variation as
> sigma >>1 , and much smaller relative levels of variation as sigma << 1. (I
> suspect this might be due to a scaling difference in the functions: sim.char,
> fastBM, and transformPhylo.sim all use ~sqrt(sigma*branchlengths) when
> generating random normal data, but I'm not certain what rTraitCont uses.)
>
> However, the most curious finding is that for all methods, as sigma
> increases, so too does the mean trait value across the tips (and the converse
> occurs as sigma decreases). This observation is curious to me, as one should
> not see a predictable shift in the mean under Brownian motion. I thought
> this might be due to simulating too few data sets and taking the mean of the
> mean, but the pattern remained when 10,000 simulated data sets were generated
> (though obviously, the effect was smaller, with mean values closer to zero).
> I then repeated the simulations for different numbers of taxa (N=16, 32, 64,
> 128) and the pattern was still present. So I think this observation is
> robust.
>
> Presently, I am at a loss to explain this observation, though I have a
> nagging suspicion that I'm missing something obvious. For instance, it may be
> related to the Central limit theorem, as increasing the number of simulated
> data sets decreased the magnitude of the mean effect across them. But that
> doesn't fully explain why the mean deviation from zero increases
> systematically with an increasing rate parameter (sigma).
>
> Any thoughts on this would be greatly appreciated.
>
> Best,
>
> Dean
>
> --
> Dr. Dean C. Adams
> Associate Professor
> Department of Ecology, Evolution, and Organismal Biology
> Department of Statistics
> Iowa State University
> Ames, Iowa
> 50011
> www.public.iastate.edu/~dcadams/
> phone: 515-294-3834
>
> _______________________________________________
> R-sig-phylo mailing list
> [email protected]
> https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
Luke Harmon
Assistant Professor
Biological Sciences
University of Idaho
208-885-0346
[email protected]
[[alternative HTML version deleted]]
_______________________________________________
R-sig-phylo mailing list
[email protected]
https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
