From: Ignacio Colonna <iacolonn_at_uiuc.edu>

Date: Sun 20 Mar 2005 - 15:14:59 EST

m1<-corMatrix(ar1)

plot(m1$"1/name of first treatment"[,1])

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide!

http://www.R-project.org/posting-guide.html

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sun Mar 20 15:20:39 2005

Date: Sun 20 Mar 2005 - 15:14:59 EST

Emma,

I am not an expert, but I have been trying to fit similar models. Adding to Keith's reply to your question, I can suggest what I concluded was the most reasonable model for my case. Based on Keith's Model1, you might also want to allow for a correlation among years within each experimental unit (I am assuming the experiment was conducted at the exact same location over the 3 years).

Say you want to impose an autoregressive, order 1 structure (you can change this to any other structure you may consider appropriate)

To do this at the block*treatment unit (the smallest experimental unit size in your experiment), you can add to keith's code:

correlation=corAR1(form=~1|block/treatment)

thus the entire code would be

Model1<-lme(mg~treatment + year + treatment:year, random=~1|block,
correlation=corAR1(form=~1|block/treatment),data=magnesium)

This results in a model with a certain covariance among all exp.units within the same block, plus another covariance between pairs of years within the same exp.unit, and this covariance decreases as the difference in time increases.

You can see graphically the structure of this covariance by doing:

rho<-0.3 ar1<-corAR1(value=rho,form=~1|block/treatment) ar1<-Initialize(ar1,data=yourdata)

m1<-corMatrix(ar1)

plot(m1$"1/name of first treatment"[,1])

Now, I really hope someone more knowledgeable is checking this out there and will point out whether this is incorrect, as I have used it for my analysis assuming was the correct approach.

Ignacio

-----Original Message-----

From: r-help-bounces@stat.math.ethz.ch

[mailto:r-help-bounces@stat.math.ethz.ch] On Behalf Of Keith Wong
Sent: Friday, March 18, 2005 5:41 PM

To: r-help@stat.math.ethz.ch

Subject: [R] Re: Repeated Measures, groupedData and lme

Hello,

I'm an R-newbie, but I've been learning to use lme for repeated measures experiments as well.

If I understand correctly:

Outcome variable: Mg (Kg/ha)

Subject/grouping variable: block

Condition/treatment: treatment (19 levels) Repeated factor: time (3 levels: 99, 02, 04)

I think if you specify the model formula in the lme call, then the formula structure specified in the groupedData object is ignored.

One suggestion for the model:

Model1<-lme(mg~treatment + year + treatment:year, random=~1|block, data=magnesium)

If the question of interest is the treatment:year interaction

Or

Model2 <- lme(mg~treatment, random=~1|block, data=magnesium)

Hope this helps ... and hope the experts chime in at this point to give more guidance.

Keith

------quoting original post---

Hello

I am trying to fit a REML to some soil mineral data which has been collected over the time period 1999 - 2004. I want to know if the 19 different treatments imposed, differ in terms of their soil mineral content. A tree model of the data has shown differences between the treatments can be attributed to the Magnesium, Potassium and organic matter content of the soil, with Magnesium being the primary separating variable.

I am looking at soil mineral data were collected : 99, 02, 04.

In the experiment, there are 19 different treatments (treatmentcontrol, treatment6TFYM, treatment 12TFYM etc), which are replicated in 3 blocks.

For the magnesium soil data, I have created the following groupedData object:

magnesium<-groupedData(Mg~year|treatment, inner=~block) Where mg=magnesium Kg/ha

If so is the following command correct:

Model1<-lme(mg~treatment, random=block|year, data=magnesium)?

Should I be fitting the model using the whole data set rather than the groupedData object?

Thank you very much for your help

Emma Pilgrim

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide!

http://www.R-project.org/posting-guide.html

R-help@stat.math.ethz.ch mailing list

https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html Received on Sun Mar 20 15:20:39 2005

*
This archive was generated by hypermail 2.1.8
: Fri 03 Mar 2006 - 03:30:51 EST
*