[R-sig-ME] Fitting a fixed effect model while modelling correlation in the residuals in 2 directions
Thomas Bishop
t.bishop at usyd.edu.au
Tue Sep 15 05:57:12 CEST 2009
Dear R users,
I have recently started to use the nlme and lme4 packages for fitting
mixed effects models. Previously I used Genstat for doing this.
Current I am struggling with analysing an experiment while modelling the
correlation in the residual term in two directions, along rows and down
columns.
The data is below. It is an agricultural experiment with variety as a
fixed effect and yield as the response. The experiment was laid out as
a lattice square design but I wish to ignore this blocking structure and
model the spatial correlation in the residuals. The layout is 10 rows
by 15 columns of plots, each planted with a different variety. I wish
to model the correlation with an AR(1) model.
I have below example code based on the nlme package.
I can model correlation in 1 direction and get the same results as I do
with Genstat.
slate.ar1 <- gls(Yield~as.factor(Variety),corr = corAR1(form = ~
1|FieldRow),data=Slate)
summary(slate.ar1)
#Generalized least squares fit by REML
# Model: Yield ~ as.factor(Variety)
# Data: Slate
# AIC BIC logLik
# 552.6045 628.969 -249.3022
#
#Correlation Structure: AR(1)
# Formula: ~1 | FieldRow
# Parameter estimate(s):
# Phi
#0.7325946
But to model spatial correlation in 2 directions along rows (FieldRow)
and down columns (Field Column) I get an error message.
slate.ar1*ar1 <- gls(Yield~as.factor(Variety),corr = corAR1(form = ~
1|(FieldRow+FieldColumn),data=Slate))
#Error in corAR1(form = ~1 | (FieldRow + FieldColumn), data = Slate) :
# unused argument(s) (data = Slate)
Thank you for any advice you can give.
Tom
Data:
Variety Yield FieldRow FieldColumn
1 10.03 1 1
2 13.56 1 2
4 14.12 1 3
3 12.39 1 4
5 15.08 1 5
19 19.67 1 6
23 15.72 1 7
2 19.69 1 8
6 17.47 1 9
15 15.98 1 10
18 16.3 1 11
25 16.33 1 12
9 12.55 1 13
11 12.77 1 14
2 15.72 1 15
6 15.31 2 1
7 15.4 2 2
9 12.5 2 3
8 16.58 2 4
10 11.85 2 5
8 16.05 2 6
12 15.5 2 7
16 15 2 8
25 16.42 2 9
4 15.04 2 10
5 16.8 2 11
7 15.26 2 12
16 14.52 2 13
23 14.8 2 14
14 14.82 2 15
21 11.26 3 1
22 14 3 2
24 13.29 3 3
23 12.87 3 4
25 15.55 3 5
11 13.95 3 6
20 16.96 3 7
24 15.7 3 8
3 14.04 3 9
7 12.85 3 10
6 14.73 3 11
13 17.61 3 12
22 16.95 3 13
4 13.64 3 14
20 17.9 3 15
11 12.61 4 1
12 14.23 4 2
14 11.1 4 3
13 17.35 4 4
15 16.17 4 5
22 18.2 4 6
1 13.51 4 7
10 12.97 4 8
14 14.12 4 9
18 15.06 4 10
24 15.12 4 11
1 13.55 4 12
15 15.24 4 13
17 14.78 4 14
8 13.71 4 15
16 14.58 5 1
17 20.36 5 2
19 21.19 5 3
18 19.12 5 4
20 18.93 5 5
5 17.48 5 6
9 14.5 5 7
13 17.4 5 8
17 14.5 5 9
21 15.23 5 10
12 13.64 5 11
19 16.9 5 12
3 13.34 5 13
10 12.39 5 14
21 15.57 5 15
3 16.23 6 1
18 18.62 6 2
8 16.45 6 3
13 18.88 6 4
23 15.27 6 5
16 16.06 6 6
24 18.42 6 7
10 11.86 6 8
13 14.62 6 9
2 12.42 6 10
10 10.82 6 11
4 13.04 6 12
17 12.67 6 13
11 12.66 6 14
23 12 6 15
1 13.31 7 1
16 14.17 7 2
6 16.11 7 3
11 14.54 7 4
21 17.9 7 5
12 17.67 7 6
20 19.17 7 7
1 12.64 7 8
9 10.6 7 9
23 9.51 7 10
12 11.3 7 11
6 12.66 7 12
24 12.89 7 13
18 12.6 7 14
5 11.74 7 15
5 12.11 8 1
20 14.11 8 2
10 11.83 8 3
15 15.5 8 4
25 16.6 8 5
4 15.26 8 6
7 16.81 8 7
18 15.45 8 8
21 12.9 8 9
15 9.76 8 10
19 12.4 8 11
13 11.81 8 12
1 9.17 8 13
25 12.87 8 14
7 9.75 8 15
2 13.88 9 1
17 14.53 9 2
7 13.84 9 3
12 16.69 9 4
22 17.38 9 5
25 18.45 9 6
3 17 9 7
14 15.28 9 8
17 13.73 9 9
6 12.4 9 10
21 12.52 9 11
20 15.91 9 12
8 14.28 9 13
2 15.09 9 14
14 12.73 9 15
4 14.43 10 1
19 16.67 10 2
9 15.49 10 3
14 14.59 10 4
24 17.22 10 5
8 15.83 10 6
11 14.9 10 7
22 16.07 10 8
5 13.15 10 9
19 11.74 10 10
3 14.43 10 11
22 16.49 10 12
15 14.07 10 13
9 13.15 10 14
16 13.18 10 15
More information about the R-sig-mixed-models
mailing list