[R-sig-ME] R-sig-mixed-models Digest, Vol 105, Issue 8
John Morrongiello
john.morrongiello at unimelb.edu.au
Fri Sep 11 01:28:20 CEST 2015
Hi David
Just thinking, given you say month could have a seasonally cyclical effect (i.e. December should be very similar to January) another, option to model the mean month effect with a cyclic cubic regression spline smoother in gamm4. A smoother in this sense deals a bit with Thierry's point of not assuming a symmetrical pattern for month, you are estimating a lot less parameters in the model, and importantly you a not left with potential hard to explain large 'steps' in predicted monthly values at the end of the year.
In your case, this might look like:
library(gamm4)
gamm4(trait ~s(month, bs='cc',k=4),random=~1|ID,data)## play around with k yourself to set the desired level of smoother 'wriggliness'
Where this leaves the random regression part of your question, I'm not sure. You can fit a random month slope for month in the above model, but I'm yet to find a good interpretation of its meaning.
##example data (sorry just had poisson code at hand, concept similar with normal)
x <- runif(100)#imagine these are your months
fac <- sample(1:20,100,replace=TRUE)
eta <- x^2*3 + fac/20; fac <- as.factor(fac)
y <- rpois(100,exp(eta))
## fit model with just cyclic smoother for x
M1 <- gamm4(y~s(x,bs='cc'),random=~(1|fac))
plot(M1$gam)## ensures max(x) cycles to min(x) (i.e. December into January)
## fit model with cyclic smoother for x and potentially random slope for x?
M2 <- gamm4(y~s(x,bs='cc'),random=~(x|fac))
plot(M2$gam)
In this fictious case, the random x slopes result in a better model fit
AIC(M1$mer,M2$mer)
I guess my question would be: how do you interpret the x's (potential random slopes) below?
ranef(M2$mer)$fac
Can they be interpreted as reaction norms, and ultimately are they performing a similar role to a polynomial random slope in David's original example? Perhaps naievely, if x is positive, is it an overall increase in smoother 'wriggliness' i.e. increase in edf (opposite for negative)? Or is it something else? Would someone know how to visualise these random slopes dependant on the mean smoother so I can get a feel for what is going on?
Cheers
John
--
Dr. John R. Morrongiello
School of BioSciences
University of Melbourne
Victoria 3010, Australia
T: +61 3 8344 8929
M: +61 403 338 554
E: john.morrongiello at unimelb.edu.au
W: morrongiellolab.com
----------------------------------------------------------------------
Message: 1
Date: Thu, 10 Sep 2015 09:51:37 +0000
From: "Wilson, Alastair" <A.Wilson at exeter.ac.uk>
To: David Villegas R?os <chirleu at gmail.com>, "Thierry Onkelinx"
<thierry.onkelinx at inbo.be>
Cc: r-sig-mixed-models <r-sig-mixed-models at r-project.org>
Subject: Re: [R-sig-ME] non-linear slope in random regression model
Message-ID:
<AMSPR03MB3926F7C1117A865B3258B38B2510 at AMSPR03MB392.eurprd03.prod.outlook.com>
Content-Type: text/plain; charset="UTF-8"
Hi David,
Replying off list as last message has been held for moderation (I suspect uni has changed e-mail aliases or something) and this is not really an e-mail about R anyway!0. Fitting month as a factor makes a lot of sense to me for reasons Thierry gave. To reiterate though your random effect model is something you need to decide. Simplistically, imagine you fixed effects has ?taken out? all monthly variation in the mean. The question now is, do you want to model individual reaction norms as constant, a straight line, or something curvy (in which case polynomials are not a bad place to start).
Dingemanse et al have suggested using higher order polynomials to look at behavioural reaction norms ? something that has long been done in quant gen for longitudinal data (e.g. growth, milk yields with age etc). In evo ecol there is a long tradition of stickling with straight line reaction norms, but I think this is largely for semantic convenience, not necessarily because they fit people?s data brilliantly. If you centre the X axis then you can interpret the elevation of an individuals reaction norm as its phenotype in an average environment, while the slope can be used as its ?plasticity?. If you add a quadratic term people don?t know what to call it! Of course what it really means is that the ?plasticity? slope of the reaction norm, itself changes according to where you are on the X axis.
Not sure if that helps. Whatever you decide you should make sure you include covariance between slopes and elevations in your model (and with high order terms if you include them), to do otherwise forces very strong (and unrealistic) assumptions on the reaction norms that can distort things a lot. You should also consider whether it is appropriate to assume a homogeneous residual variance across all your months as violations of this can generate spurious support for IxE (i.e. among-individual variance in reaction norm slopes).
Best
A
From: David Villegas R?os [mailto:chirleu at gmail.com]
Sent: 10 September 2015 10:33
To: Thierry Onkelinx; Wilson, Alastair
Cc: r-sig-mixed-models
Subject: Re: [R-sig-ME] non-linear slope in random regression model
Thanks Thierry and Alastair.
I do have a lot of IDs (~290) but not so many month per ID (a mean of 10).
So from your comments I conclude that if I run model 2, or as suggested by Thierry this equivalent model, in which month is a factor in the fixed part...
trait~as.factor(month),random=~month|ID
...and I get significant between individual differences in the slope of the random effect, the conclusion is that individuals differ in their non-linear (as described by factor(month) or poly(month,3) in the fixed part) relationship between trait and month. Right?
If this is correct then I would not upgrade to a higher order random slope to keep things simple!
Best,
David
2015-09-10 11:08 GMT+02:00 Thierry Onkelinx <thierry.onkelinx at inbo.be<mailto:thierry.onkelinx at inbo.be>>:
Dear David,
Polynomials are technically still a linear model since they result is in linear combination of variable.
The only straightforward answer to you question is: it depends on the data. If you have plenty data, then I'd go for model 1. If the data doesn't require a second and 3th order polynomial, then their random effect variance will be very small. And the model will be reduced to model 2.
Note than plenty data means a lot of different id AND a lot of months per id. If you have only one year of data then random=~poly(month, 1)|ID is about as complex as you can go.
A 3th order polynomial seems a bit odd to me to model seasonality. I'd rather expect an even order. You might want to consider fitting month as a factor in the fixed effects. Then you model the seasonality without assumptions on the pattern.
Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature and Forest team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance Kliniekstraat 25
1070 Anderlecht
Belgium
To call in the statistician after the experiment is done may be no more than asking him to perform a post-mortem examination: he may be able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher The plural of anecdote is not data. ~ Roger Brinner The combination of some data and an aching desire for an answer does not ensure that a reasonable answer can be extracted from a given body of data. ~ John Tukey
2015-09-10 10:47 GMT+02:00 David Villegas R?os <chirleu at gmail.com<mailto:chirleu at gmail.com>>:
Hello,
I'm fitting some random regression models to investigate variation over time of a response trait.
The time variable is "month", and by fitting a random regression model I want to investigate variation in plasticity across individuals, i.e., differences in the "slope" between trait and time across individuals.
The relationship between the response trait and month is non-linear.
Basically, it describes a seasonal cycle.
I have considered two model candidates:
*Model 1*: fitting a polynomial of month in the fixed effects part to describe the non-linearity and get the population mean effect of month, and then a random slope using again the polynomial for month, to ge the individual differences.
trait~poly(month,3),random=~poly(month,3)|ID
*Model 2*: fitting a polynomial of month in the fixed effect part to describe the non-linearity, and then the individual-specific deviation from the fixed-effect means is modelled as a funtion of month (linear), assuming that the non-linearity is already accounted for with the fixed effect.
trait~poly(month,3),random=~month|ID
*My question is*:
Is it neccesary to include the non-linearity in the random part if it was already included in the fixed-effects part?
The idea of fitting model 2 comes from the following reference (page 488):
Dingemanse, N. J., Barber, I., Wright, J., & Brommer, J. E. (2012).
Quantitative genetics of behavioural reaction norms: genetic correlations between personality and behavioural plasticity vary across stickleback populations.*Journal of evolutionary biology*, *25*(3), 485-496.
Thank you.
David
[[alternative HTML version deleted]]
_______________________________________________
R-sig-mixed-models at r-project.org<mailto:R-sig-mixed-models at r-project.org> mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
[[alternative HTML version deleted]]
------------------------------
Message: 2
Date: Thu, 10 Sep 2015 09:05:03 +0000
From: "Wilson, Alastair" <A.Wilson at exeter.ac.uk>
To: David Villegas R?os <chirleu at gmail.com>, r-sig-mixed-models
<r-sig-mixed-models at r-project.org>
Subject: Re: [R-sig-ME] non-linear slope in random regression model
Message-ID:
<AMSPR03MB392F0A3C8B41CD9623630DBB2510 at AMSPR03MB392.eurprd03.prod.outlook.com>
Content-Type: text/plain; charset="iso-8859-1"
Hi David,
It is not necessary - it may or may not be useful/interesting depending on your questions and structure of your data. In this sort of random regression model you are using your fixed effect structure to account for the relationship between your X variable (here time) and the response mean. If a 3rd order polynomial does a good job then great! In the random effect structure you are modelling the individual level effects (i.e. deviations from the fixed effect mean attributable to individual identity). These could be constant (0 order), a linear (1st order) function of X or something else (e.g. a quadratic, or indeed some other type of function of X).
There is no necessary connection between the functional form you use in the fixed effects and that in the random effects. For plasticity studies a linear regression is common in the random effect but higher orders could be used (and may prove better fits).
For more on these sorts of models (they have some issues as well as some strength!) I would suggest looking in the animal breeding/quant gen literature, for instance work by Karin Meyer and Mark Kirkpatrick here
http://rstb.royalsocietypublishing.org/content/360/1459/1443
Alastair Wilson
-----Original Message-----
From: R-sig-mixed-models [mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of David Villegas R?os
Sent: 10 September 2015 09:47
To: r-sig-mixed-models
Subject: [R-sig-ME] non-linear slope in random regression model
Hello,
I'm fitting some random regression models to investigate variation over time of a response trait.
The time variable is "month", and by fitting a random regression model I want to investigate variation in plasticity across individuals, i.e., differences in the "slope" between trait and time across individuals.
The relationship between the response trait and month is non-linear.
Basically, it describes a seasonal cycle.
I have considered two model candidates:
*Model 1*: fitting a polynomial of month in the fixed effects part to describe the non-linearity and get the population mean effect of month, and then a random slope using again the polynomial for month, to ge the individual differences.
trait~poly(month,3),random=~poly(month,3)|ID
*Model 2*: fitting a polynomial of month in the fixed effect part to describe the non-linearity, and then the individual-specific deviation from the fixed-effect means is modelled as a funtion of month (linear), assuming that the non-linearity is already accounted for with the fixed effect.
trait~poly(month,3),random=~month|ID
*My question is*:
Is it neccesary to include the non-linearity in the random part if it was already included in the fixed-effects part?
The idea of fitting model 2 comes from the following reference (page 488):
Dingemanse, N. J., Barber, I., Wright, J., & Brommer, J. E. (2012).
Quantitative genetics of behavioural reaction norms: genetic correlations between personality and behavioural plasticity vary across stickleback populations.*Journal of evolutionary biology*, *25*(3), 485-496.
Thank you.
David
[[alternative HTML version deleted]]
_______________________________________________
R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
More information about the R-sig-mixed-models
mailing list