Dear colleagues,

I am trying to analyze the data of a longitudinal experiment in which we 
measured gonadal size of birds each month for a 9 months period, over 
the whole gonadal cycle (from growth until regression). Birds were 
divided in two treatments, control and experiment, and in each treatment 
we had two groups of birds depending on their origin (urban vs. rural). 
The experiment was perfectly balanced, with 10 birds per origin per 
treatment (hence a total of 40 individuals). Data is attached at the 
bottom of this e-mail.
I analyzed this data with nlme, including julian date, quadratic julian 
date, origin, treatment and all relevant interactions in the model. I 
accounted for autocorrelation and heteroscedasticity using CAR1 and 
varIdent arguments. I also included random intercepts and slopes in the 
model to account for non independence of data in the repeated measure 
design. The model looks like this:

 > mod <- 
lme(log(gonad)~julian+julian2+treatment+origin+julian*treatment*origin+julian2*treatment*origin,random=~julian|id,
corr = 
corCAR1(form=~julian),weights=varComb(varIdent(form=~julian|treatment),varIdent(form=~julian|origin)),
control = list(niterEM = 100000),data=gonads)

During model selection I excluded all non relevant-non significant 
interactions. Eventually the final model includes only julian:treatment 
and julian:origin interactions, which are significant. Output is also at 
the bottom of this e-mail.

Now my question: I am interested to know at which julian dates the 
treatment and origin effects are significant. I would avoid post-hoc 
tests since I think that they are not really needed and anyway 
complicated for such models. Could I use contrasts to get p-values of 
effects at different julian dates? If so, what would be the best way to 
write down the contrasts matrix? Alternatively, what would you suggest I 
should do to answer my question?

I would appreciate if you could help help me.
Thank you for your time,

Davide Dominoni




Model output:

Linear mixed-effects model fit by REML
  Data: gonads
        AIC      BIC    logLik
   401.1266 458.9096 -185.5633

Random effects:
  Formula: ~julian | id
  Structure: General positive-definite, Log-Cholesky parametrization
             StdDev      Corr
(Intercept) 0.012697642 (Intr)
julian      0.001408309 -0.998
Residual    0.347322171

Correlation Structure: Continuous AR(1)
  Formula: ~julian | id
  Parameter estimate(s):
Phi
0.2

Combination of variance functions:

  Structure: Different standard deviations per stratum
  Formula: ~julian | treatment
  Parameter estimates:
experimental      control
    1.0000000    0.9907898

  Structure: Different standard deviations per stratum
  Formula: ~julian | origin
  Parameter estimates:
    rural    urban
1.000000 1.031409

Fixed effects: log(gonad) ~ julian + julian2 + treatment + origin + 
julian:treatment + julian:origin + treatment:origin

                                                                 
Value                  Std.Error          DF             t-value       
p-value
(Intercept)                                              
0.3779310        0.06163952   312        6.131309      0.0000
julian                                                       
0.0247056        0.00108908  312       22.684912      0.0000
julian2                                                   
-0.0001274        0.00000497   312     -25.602676      0.0000
treatmentexperimental                        0.3653295        
0.07748821     36        4.714647      0.0000
originurban                                            0.2343046        
0.07797856     36        3.004731      0.0048
julian:treatmentexperimental            -0.0027794        0.00074267   
312      -3.742475       0.0002
julian:originurban                                -0.0017778        
0.00074313   312      -2.392261       0.0173
treatmentexperimental:originurban -0.1040218        0.09261304     
36      -1.123187       0.2688

  Correlation:
                                                                 
(Intr)    julian  julian2 trtmnt orgnrb jln:tr jln:rg
julian                                                       -0.608
julian2                                                     0.323 -0.813
treatmentexperimental                       -0.621  0.182  0.005
originurban                                           -0.612  0.181 
-0.004  0.337
julian:treatmentexperimental             0.342 -0.331 -0.009 -0.554  0.007
julian:originurban                                 0.333 -0.333  0.011  
0.010 -0.548 -0.015
treatmentexperimental:originurban  0.358  0.006 -0.004 -0.577 -0.584 
-0.005 -0.010

Standardized Within-Group Residuals:
        Min         Q1        Med         Q3        Max
-2.6727564 -0.8101380  0.1647109  0.7735809  1.9034670

Number of Observations: 356
Number of Groups: 40


head(gonads)

id 	origin 	treatment 	     julian 	   julian2 	    gonad
30945 	rural 	control 	-15 	-225 	1.9
30945 	rural 	control 	21 	441 	1.9
30945 	rural 	control 	42 	1764 	2.43
30945 	rural 	control 	61 	3721 	4.775
30945 	rural 	control 	84 	7056 	6.95
30945 	rural 	control 	105 	11025 	7
30945 	rural 	control 	133 	17689 	6.9
30945 	rural 	control 	161 	25921 	2
30945 	rural 	control 	189 	35721 	1.6
30991 	urban 	control 	-15 	-225 	1.7
30991 	urban 	control 	21 	441 	1.8
30991 	urban 	control 	42 	1764 	2.3
30991 	urban 	control 	61 	3721 	4.1
30991 	urban 	control 	84 	7056 	6.4
30991 	urban 	control 	105 	11025 	7.7
30991 	urban 	control 	133 	17689 	6.8
30991 	urban 	control 	161 	25921 	2.4
30991 	urban 	control 	189 	35721 	2



-- 
Davide M. Dominoni - PhD student
Schlossalle 2
78315 Radolfzell
Germany

ddominoni@orn.mpg.de
mobile: 0049 (0) 17662376221
office: 0049 (0) 7732150124


	[[alternative HTML version deleted]]


