[R-sig-phylo] Error estimating RMA in phytools
Liam J. Revell
liamjrevell at gmail.com
Sat Jun 23 03:30:57 CEST 2012
Hi Oscar.
My educated guess after some investigating is this is because the
objects Chem2mean and Morpho2mean are data frames not matrices; and thus
Chem2mean[,1] does not contain names. Try first doing:
Chem2mean<-as.matrix(Chem2mean), and the same with Morpho2mean, and then
repeat the analysis.
You might also get this error (or one quite similar) if you tree has
zero-length terminal edges or y=k*x for non-zero scalar constant k
(i.e., r(x,y)==1 or -1).
All the best, Liam
--
Liam J. Revell
University of Massachusetts Boston
web: http://faculty.umb.edu/liam.revell/
email: liam.revell at umb.edu
blog: http://phytools.blogspot.com
On 6/22/2012 12:17 PM, Oscar Valverde wrote:
> Hi
>
> I am working on the correlation among traits of 36 plants. I am using
> the function phyl.RMA to estimate relationship between pairs of traits
> but I keep having a error message
>
> CNvrSRLphylo<-phyl.RMA(Chem2mean[,1],Morpho2mean[,1],treeARP,method="lambda")
> Error in solve.default(kronecker(R, C)) :
> system is computationally singular: reciprocal condition number = 0
>
> I tried running the code provided by phytools for this function and
> works fine. So it must be something about my tree that is not working.
> Thanks in advance for the help
>
>
>
> On Fri, Apr 20, 2012 at 12:19 PM, Jarrod Hadfield<j.hadfield at ed.ac.uk> wrote:
>> Hi,
>>
>> The warning:
>>
>> 1: glm.fit: algorithm did not converge
>>
>> is from a standard glm that MCMCglmm uses to obtain semi-reasonable starting
>> values. For a three (I think this is correct for your data?) category
>> response the starting values are obtained from 2 binomial glms (the presence
>> of category 2 and the presence of category 3). The fact that a simple glm
>> did not converge indicates a problem. You seem to suggest the warning only
>> appears with phylogenetic models, which seems odd. Is this correct? If so,
>> try the argument start=list(QUASI=FALSE) in MCMCglmm and report back.
>>
>> If not, I would try and diagnose the problem using a simple binomial glm for
>> each of the two categories in order to pinpoint the problem. The most likely
>> problem is complete separation whereby some categories are not represented
>> in some environments. Is this possible?
>>
>> With completely separated data MCMCglmm (and other software) is likely to
>> give spurious results unless something is done. One possibility is to fit a
>> more informative prior on the fixed effects which is approximately flat on
>> the probability scale (See the last half of section 2.6 in the CourseNotes).
>>
>> Cheers,
>>
>> Jarrod
>>
>>
>>
>>
>>
>>
>> Quoting Rafael Rubio de Casas<rafa at nescent.org> on Tue, 17 Apr 2012
>> 17:15:17 -0400:
>>
>>> Dear R-Sig-Phylo users,
>>> I am trying to use MCMCglmm to analyze the evolution of a discrete
>>> character and the influence of the environment on it, but I keep on
>>> having problems. All my mixed models give crazy results, with chains
>>> that do not converge. I always get teh following error message, which I
>>> believe is probably at the root of the problem
>>> Warning messages:
>>> 1: glm.fit: algorithm did not converge
>>>
>>> However, when I fit a fixed effects model (no phylogenetic correction),
>>> everything seems to go fine, and it is when I incorporate the phylogeny
>>> that things fall apart. I keep on getting the warning message,
>>> regardless of the priors (although there might be something to improve
>>> there). I was hoping somebody could give some advice.
>>> So here is a detailed description of what I am doing:
>>> My data consists of a trait with three levels, a set of environmental
>>> parameters and a phylogeny that describes the evolutionary relationships
>>> among the different taxa.
>>> Let the data be stored in a dataframe "df", including taxa names as
>>> "taxa", the trait of interest as a factor with three states "char", and
>>> the environmental parameter a continuous variable "env". The phylogeny
>>> is an ultrametric tree "phy". Below I include a simulation that creates
>>> a similar dataset but that for some reason is not problematic. My actual
>>> dataset is very big and somewhat messy. If someone is interested I can
>>> send it.
>>> The code I'm using is as follows
>>> #Fixed effects model:
>>> #Prior definition
>>> k<- length(levels(data$char))
>>> I<- diag( k-1 )
>>> J<- matrix( rep(1, (k-1)^2), c(k-1, k-1))
>>> IJ<- (1/3) * (diag(k-1) + J)
>>> prior = list(R = list(V = IJ, fix = 1))
>>> myMCMC<- MCMCglmm( char ~ trait + env:trait -1,
>>> rcov = ~us(trait):units,
>>> data = df, family="categorical",
>>> prior=prior )
>>>
>>> #Model using the phylogeny as a random factor
>>> #Prior definition . (There might be some problems here...)
>>> Prior.phyl = list(R = list(V = IJ, fix = 1),G = list( G1 = list(V = IJ,
>>> n = k-1 ) ) )
>>> Ainv<-inverseA(phy)$Ainv
>>> df$species=df$taxa
>>> myMCMC.phyl<- MCMCglmm( ch ~ trait + env:trait -1,
>>> random=~us(trait):species,
>>> rcov = ~us(trait):units,
>>> ginverse = list(species=Ainv),
>>> data = df, family="categorical",
>>> prior=prior.phyl)
>>>
>>> Thanks in advance,
>>> Rafa
>>>
>>> ####################
>>> ### Simulated data ###
>>> ####################
>>> ##NOTE: this data runs seamlessly
>>> a=rnorm(100,5,3)
>>> b=rnorm(100,10,3)
>>> c=rnorm(100,20,3)
>>> env=c(a,b,c)
>>> char=rep(c(1,2,3),c(100,100,100))
>>> df=as.data.frame(cbind(char,env))
>>> taxa.labels =paste(rep("sp",30),c(1:30),sep="")
>>>
>>> trait.data=rbind(cbind(taxa.labels[1:10],1),cbind(taxa.labels[11:20],2),cbind(taxa.labels[21:30],3))
>>> df$taxa<- trait.data[,1][match(df$char, df[,2])]
>>> df$char=as.factor(df$char)
>>> phy<- rcoal(30, tip.label=taxa.labels)
>>>
>>> --
>>> National Evolutionary Synthesis Center
>>> *NESCent<http://www.nescent.org/>*
>>> 2024 W. Main Street, Suite A200
>>> Durham, NC27705
>>> rafa at nescent.org<mailto:rfr6 at duke.edu>
>>> 919.668.9107
>>>
>>> [[alternative HTML version deleted]]
>>>
>>> _______________________________________________
>>> R-sig-phylo mailing list
>>> R-sig-phylo at r-project.org
>>> https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
>>>
>>>
>>
>>
>>
>> --
>> The University of Edinburgh is a charitable body, registered in
>> Scotland, with registration number SC005336.
>>
>> _______________________________________________
>> R-sig-phylo mailing list
>> R-sig-phylo at r-project.org
>> https://stat.ethz.ch/mailman/listinfo/r-sig-phylo
>
>
>
More information about the R-sig-phylo
mailing list