[R] Accounting for correlated random effects in coxme with matrix from a phylogeny rather than pedigree
Hello All, I have a problem with running the mixed effects Cox regression model using a distance matrix from a phylogeny rather than a pedigree. I searched previous posts and didn't find any directly relevant previous posts. I am interested in using a mixed effects Cox regression model to determine the best predictors of time to recruitment in 80 different reintroduced plant populations representing a total of 31 species. I will like to account for correlated random effects that result from phylogenetic relationships amongst species. Dr. Therneau's 2015 article on Mixed Effects Cox Models provide a very helpful template for me to do this with the coxme function in R. In this article, the correlation structure due to genetic relationships amongst individuals was defined using a kinship matrix derived from a pedigree. Instead of a pedigree, I have a phylogeny for these 31 species. Hence, I used the inverseA function in the MCMCglmm package to generate an inverse additive genetic relatedness matrix from the phylogeny for these 31 species. And then fed it in as input to the varlist argument in my mixed effects cox regression model (using function coxme). I got an error message (please see below). Based on the error, one thought I had was to convert the inverseA matrix from a “dgCMatrix” to “bdsmatrix” but this was not successful either. I have also unsuccessfully tried to use a pairwise phylogenetic distance matrix. Is there a better way to do this? I basically just want to account for the correlated random effects due to phylogenetic relatedness amongst the 31 species represented in the dataset for the Cox regression model. Please see my code below and I welcome suggestions on how best to make this work. Thank you. #Load packages library(MCMCglmm) library(asremlPLUS) library(ape) source("read.newick.R") mytree <- read.newick(file="Phylo_2015Sept15.txt") mytree6 <- makeNodeLabel(mytree, method="number", prefix = "node")#Make sure each node is uniquely labeled IA <- inverseA(mytree6, scale=TRUE) #generate inverse of the additive genetic relatedness matrix (A) from phylogeny #Use IA as input in correlated random effects model. Doesn't work. fit2 <- coxme(Surv(surv.time, recruitment) ~ pred1 + pred2 + sixcatD1 + sixcatD2 + sixcatD3 + (1|species), data = traitcox, varlist=coxmeMlist(IA, rescale=F)) Error in as(x, "bdsmatrix") : no method or default for coercing “dgCMatrix” to “bdsmatrix” [[alternative HTML version deleted]] __ R-help@r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
[R] Help!
Please I need some help using R to analyze my data. What I would like to do is to repeat the same basic process (e.g. linear regression between wood density and distance from pith) for at least 240 data subsets within the main data-frame. Within the main data-frame, these data subsets will be defined by three variables namely,  species, individual and core (i.e. 20 species, at least 6 individuals of each species, and 2 cores from each individual).  Whereas I can write the code to carry out this process for each subset, I am unable to successfully instruct R to automatically carry out the process for each of these subsets (perhaps using loops). So to illustrate what I have done so far with the codes below I was able to run a regression analysis for core âaâ of individual 1 in the species âApeimeâ. But rather than do this 240 times, I would like to tell R to repeat the process automatically using loops or any method that works.  Code:  RG2-BCI[BCI$Species == APEIME BCI$Individual == 1 BCI$Core == a, ] plot(x=RG2$DP..cm., y=RG2$WD..g.cm3, xlab=Distance from pith cm, main=APEIME1a, ylab=Wood density g/cm3) RG2lm-lm(RG2$WD..g.cm3~RG2$DP..cm.) summary(RG2lm)  Thanks  Oyomoare [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
Re: [R] Help!
From: Oyomoare Osazuwa-Peters oyomo...@yahoo.com Subject: Re: [R] Help! To: Erich Neuwirth erich.neuwi...@univie.ac.at Date: Monday, September 20, 2010, 5:16 PM Thanks for responding to my request for help. I understand what you mean about the repeated measurements methods for the two cores. The thing though is to answer my research question, the data I really need is the radial gradient (equals the slope from a regression of the response variable (WD) on the predictor (DP)) for each core. Then, I can be begin to test for the effects of species, individuals and core using an appropriate test (likely nested anova). For now I am in the initial process of getting radial gradients and having problems with the code that would instruct R to do it all at once. My main problem is when I define the subsetting indices to be species, individual and core at the same time for the whole data frame, so that R performs the operation for each of the 240 data subsets automatically, it doesn't work. But it works when I define only a single subset of the data like I showed in my first mail. Oyomoare --- On Mon, 9/20/10, Erich Neuwirth erich.neuwi...@univie.ac.at wrote: From: Erich Neuwirth erich.neuwi...@univie.ac.at Subject: Re: [R] Help! To: r-help@r-project.org Date: Monday, September 20, 2010, 5:02 PM If you assume that the variance is the same in all your subsets, you can do an lm analysis with your subset classification as a factor. You could also analyze the interaction between factors and between factors and your numeric independent variable. You also should consider repeated measurement methods since you are taking 2 cores from the same individuals. On 9/20/2010 11:46 PM, Oyomoare Osazuwa-Peters wrote: Please I need some help using R to analyze my data. What I would like to do is to repeat the same basic process (e.g. linear regression between wood density and distance from pith) for at least 240 data subsets within the main data-frame. Within the main data-frame, these data subsets will be defined by three variables namely, à species, individual and core (i.e. 20 species, at least 6 individuals of each species, and 2 cores from each individual). à Whereas I can write the code to carry out this process for each subset, I am unable to successfully instruct R to automatically carry out the process for each of these subsets (perhaps using loops). So to illustrate what I have done so far with the codes below I was able to run a regression  analysis for core ââ¬Ëaââ¬â¢ of individual 1 in the species ââ¬ÅApeimeââ¬ï¿½. But rather than do this 240 times, I would like to  tell R to repeat the process automatically using loops or any method that works.  Ã Code:  Ã RG2-BCI[BCI$Species == APEIME BCI$Individual == 1 BCI$Core == a, ] plot(x=RG2$DP..cm., y=RG2$WD..g.cm3, xlab=Distance from pith cm, main=APEIME1a, ylab=Wood density g/cm3) RG2lm-lm(RG2$WD..g.cm3~RG2$DP..cm.) summary(RG2lm)  Ã Thanks  Ã Oyomoare         [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. -Inline Attachment Follows- __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. [[alternative HTML version deleted]] __ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.