Random effects for subsampling in RCBD effects model dimensions

[mthorne]

I am working with a dataset from an agricultural weed management study. The data are counts per sq. meter. There are 2 study sites, 6 treatments, 4 replications per treatment, and 2 subsamples per plot. Plot is the experimental unit and I am considering study site as a random effect. This is a 6-year trial but I have already determined that I need to evaluate each year separately because of an interaction with year sequence and treatment. Using a negative binomial distribution doesn't always work because of zeros resulting from 100% effectiveness of the treatments. 

This is the model that seems to work the best and produce dimensions that seem realistic.

ods graphics on;
proc glimmix data=ltf method=laplace plots=studentpanel; by seq;
class site seq trt sub rep;
tctm2 = sqrt(ctm2);
model tctm2 = trt / dist=normal link=identity;
random intercept / subject=rep(sub);
lsmeans trt / ilink diff lines plot=meanplot;
run; quit;
ods graphics off;

It produces this structure. (site=field site, seq=year, trt=treatment, sub=subplot, rep=replication)

Class Level Information
Class	Levels	Values
Site	2	ED ST
Seq	1	6
Trt	6	1 2 3 4 5 6
sub	2	1 2
Rep	4	1 2 3 4

Number of Observations Read	96
Number of Observations Used	96

Dimensions
G-side Cov. Parameters	1
R-side Cov. Parameters	1
Columns in X	7
Columns in Z per Subject	1
Subjects (Blocks in V)	8
Max Obs per Subject	12

If I try to add "site" as a random effect, the model structure seems all wrong, as follows.

ods graphics on;
proc glimmix data=ltf method=laplace plots=studentpanel; by seq;
class site seq trt sub rep;
tctm2 = sqrt(ctm2);
model tctm2 = trt / dist=normal link=identity;
random intercept / subject=site;
random intercept / subject=rep(sub);
lsmeans trt / ilink diff lines plot=meanplot;
run; quit;
ods graphics off;

Class Level Information
Class	Levels	Values
Site	2	ED ST
Seq	1	6
Trt	6	1 2 3 4 5 6
sub	2	1 2
Rep	4	1 2 3 4

Number of Observations Read	96
Number of Observations Used	96

Dimensions
G-side Cov. Parameters	2
R-side Cov. Parameters	1
Columns in X	7
Columns in Z	10
Subjects (Blocks in V)	1
Max Obs per Subject	96

The number of subjects is wrong, and the number of obs per subject is the total number of obs.

Any thoughts on this analysis?

Thanks!

[alexben · Accepted solution]

---

@mthorne wrote:

I am working with a dataset from an agricultural weed management study. The data are counts per sq. meter. There are 2 study sites, 6 treatments, 4 replications per treatment, and 2 subsamples per plot. Plot is the experimental unit and I am considering study site as a random effect. This is a 6-year trial but I have already determined that I need to evaluate each year separately because of an interaction with year sequence and treatment. Using a negative binomial distribution doesn't always work because of zeros resulting from 100% effectiveness of the treatments. 

This is the model that seems to work the best and produce dimensions that seem realistic.

ods graphics on;
proc glimmix data=ltf method=laplace plots=studentpanel; by seq;
class site seq trt sub rep;
tctm2 = sqrt(ctm2);
model tctm2 = trt / dist=normal link=identity;
random intercept / subject=rep(sub);
lsmeans trt / ilink diff lines plot=meanplot;
run; quit;
ods graphics off;

It produces this structure. (site=field site, seq=year, trt=treatment, sub=subplot, rep=replication)

Class Level InformationClass Levels ValuesSiteSeqTrtsubRep

2	ED ST
1	6
6	1 2 3 4 5 6
2	1 2
4	1 2 3 4

Number of Observations ReadNumber of Observations Used

96

96

DimensionsG-side Cov. ParametersR-side Cov. ParametersColumns in XColumns in Z per SubjectSubjects (Blocks in V)Max Obs per Subject

1

1

7

1

8

12

If I try to add "site" as a random effect, the model structure seems all wrong, as follows.

ods graphics on;
proc glimmix data=ltf method=laplace plots=studentpanel; by seq;
class site seq trt sub rep;
tctm2 = sqrt(ctm2);
model tctm2 = trt / dist=normal link=identity;
random intercept / subject=site;
random intercept / subject=rep(sub);
lsmeans trt / ilink diff lines plot=meanplot;
run; quit;
ods graphics off;

Class Level InformationClass Levels ValuesSiteSeqTrtsubRep

2	ED ST
1	6
6	1 2 3 4 5 6
2	1 2
4	1 2 3 4

Number of Observations ReadNumber of Observations Used

96

96

DimensionsG-side Cov. ParametersR-side Cov. ParametersColumns in XColumns in ZSubjects (Blocks in V)Max Obs per Subject

2

1

7

10

1

96

The number of subjects is wrong, and the number of obs per subject is the total number of obs.

Any thoughts on this analysis?

Thanks!

---

To properly include site as a random effect in your GLIMMIX model while keeping the structure correct, you can try nesting rep within site, rather than adding a second intercept statement. Here’s a revised model:

sas

Copy code

ods graphics on; proc glimmix data=ltf method=laplace plots=studentpanel; by seq; class site seq trt sub rep; tctm2 = sqrt(ctm2); model tctm2 = trt / dist=normal link=identity; random intercept / subject=site; random intercept / subject=rep(site*sub); lsmeans trt / ilink diff lines plot=meanplot; run; quit; ods graphics off;

This way, you should see more realistic subject structure and parameter estimates. Let me know if this helps!

[mthorne]

Thank you so much! This looks good!