[R-sig-ME] Adding additional factor to previously working MCMCglmm model kills it. Problem with priors? Help?

[Tricia Markle]

Hi Jarrod,

Thanks for the response.

If I am understanding correctly, would I then end up with the code below?
If species is incorporated in the random part, would I then not include it
as a fixed factor?

tree<-read.tree("Plethodontidae_comb61_PL.phy")
species<-c("D._carolinensis_KHK103", "D._fuscus_KHK142",
"D._ochrophaeus_WKS05", "D._ocoee_B_KHK62", "D._orestes_KHK129",
"D._monticola_A",  "D._santeetlah_11775", "P_cinereus", "P_cylindraceus",
"P_glutinosus", "P_hubrichti", "P_montanus", "P_punctatus", "P_richmondi",
"P_teyahalee", "P_virginia", "P_wehrlei")
pruned.tree<-drop.tip(tree,tree$tip.label[-match(species,
tree$tip.label)])# Prune tree to just include species of interest
sptree<-makeNodeLabel(pruned.tree, method="number", prefix="node") #rename
nodes to be unique
treeAinv<-inverseA(sptree, nodes="TIPS")$Ainv

dataset$species.ide<-dataset$species
random=~us(1+Temp):species.ide

prior<-list(G=list(G1=list(V=diag(2), nu=2, alpha.mu=c(0,0),
alpha.V=diag(2)*1000)), R=list(V=diag(1), nu=0.002))

model6<-MCMCglmm(LVO2~1+Temp+Acclm+Range+Mass+Sex+Acclm*Range,
random=random, data=dataset, family="gaussian",
ginverse=list(species=treeAinv), prior=prior, nitt=300000, burnin=25000,
thin = 1000, verbose=FALSE)

Alternatively, someone mentioned to me that to reduce model complexity I
could group the handful of species I have with smaller samples sizes and
label them all as "other" in a second species list ("New_spp"). When I
added this new species assortment as a fixed factor in my original model it
runs! Is this a reasonable thing to do?

Thanks, Tricia

On Wed, Apr 29, 2015 at 12:46 AM, Jarrod Hadfield <j.hadfield at ed.ac.uk>
wrote:

> Hi,
>
> If you want iid species effects then put them in the random part. You will
> need to rename them so that they are not associated with the ginverse. For
> example,
>
> dataset$species.ide<-dataset$species
>
> and fit species.ide as a random effect.
>
> You only have 17 species. This is not enough to get precise estimates of
> the variance components and so you should expect prior sensitivity,
> particularly with respect to the phylogenetic part. I would try and
> simplify your model. I'm not sure how many of the fixed effects are
> species-level or individual-level but I would try and reduce the complexity
> of the fixed part of the model too (for example is the 4-way interaction
> needed?). The warning usually indicates that the model is overparameterised.
>
> Cheers,
>
> Jarrod
>
>
>
>
>  Quoting Tricia Markle <markl033 at umn.edu> on Wed, 29 Apr 2015 00:24:53
> -0500:
>
>  Hello,
>>
>>
>>
>> I have a working MCMCglmm model with phylogenetic consideration and repeat
>> measures. I realized after the fact that ?species? wasn?t properly
>> included
>> in the model. When I added this additional factor (of 17 levels), however,
>> I received the following error message:
>>
>> --------
>>
>>
>>
>> Error in MCMCglmm(LVO2 ~ 1 + Temp + Acclm + Lat_Ext + LMass + Sex +
>> species
>> +  :  ill-conditioned G/R structure: use proper priors if you haven't or
>> rescale data if you have
>>
>> In addition: Warning message:
>>
>> In MCMCglmm(LVO2 ~ 1 + Temp + Acclm + Lat_Ext + LMass + Sex + species +  :
>>
>>   some fixed effects are not estimable and have been removed. Use
>> singular.ok=TRUE to sample these effects, but use an informative prior!
>>
>>
>>
>> ----------
>>
>>
>>
>> I had some help with my priors originally so I am not sure the best way to
>> now tweak them to make them work? Could different priors help or is there
>> something ?wrong? with adding species as another factor?
>>
>>
>>
>>
>> My hypothesis centers around the remaining significant interaction term
>> Acclm:Lat_Ext ? whether the relationship of acclimation on oxygen
>> consumption (VO2) is differs depending on the latitudinal extent of a
>> species (while also considering a number of covariates).
>>
>>
>>
>> Here is my code:
>>
>>
>>
>> library(ape)
>>
>> library(MCMCglmm)
>>
>> dataset<-read.csv(file="RespData.csv", head=TRUE)
>>
>> attach(dataset)
>>
>> str(dataset) # confirming that sex, range, species, and ID are all factors
>>
>>
>>
>> #Phylogeny Component
>>
>> tree<-read.tree("Plethodontidae_comb61_PL.phy")
>>
>> species<-c("D._carolinensis_KHK103", "D._fuscus_KHK142",
>> "D._ochrophaeus_WKS05", "D._ocoee_B_KHK62", "D._orestes_KHK129",
>> "D._monticola_A",  "D._santeetlah_11775", "P_cinereus", "P_cylindraceus",
>> "P_glutinosus", "P_hubrichti", "P_montanus", "P_punctatus", "P_richmondi",
>> "P_teyahalee", "P_virginia", "P_wehrlei")
>>
>> pruned.tree<-drop.tip(tree,tree$tip.label[-match(species,
>> tree$tip.label)])# Prune tree to just include species of interest
>>
>> sptree<-makeNodeLabel(pruned.tree, method="number", prefix="node")
>>
>>
>>
>> treeAinv<-inverseA(sptree, nodes="TIPS")$Ainv
>>
>> random=~us(1+Temp):species
>>
>> prior<-list(G=list(G1=list(V=diag(2), nu=2, alpha.mu=c(0,0),
>> alpha.V=diag(2)*1000)), R=list(V=diag(1), nu=0.002))
>>
>>
>>
>> #Model 1
>>
>>
>> model1<MCMCglmm(LVO2~1+Temp+Acclm+Lat_Ext+LMass+Sex+species+Temp*Acclm+Temp*Lat_Ext+Acclm*Lat_Ext+Temp*Acclm*Lat*Ext,
>> random=random, data=dataset, family="gaussian",
>> ginverse=list(species=treeAinv), prior=prior, nitt=300000, burnin=25000,
>> thin = 1000, verbose=FALSE)
>>
>>
>>
>> Thank you kindly for your help!
>>
>>
>>
>> Tricia
>>
>>         [[alternative HTML version deleted]]
>>
>> _______________________________________________
>> R-sig-mixed-models at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>>
>>
>
> --
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336.
>
>
>

	[[alternative HTML version deleted]]