turn on suggestions

Auto-suggest helps you quickly narrow down your search results by suggesting possible matches as you type.

Showing results for

Find a Community

- Home
- /
- Analytics
- /
- Stat Procs
- /
- proc mixed estimated V matrix for different subjec...

Topic Options

- Subscribe to RSS Feed
- Mark Topic as New
- Mark Topic as Read
- Float this Topic for Current User
- Bookmark
- Subscribe
- Printer Friendly Page

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

06-24-2015 04:34 PM

Hello SAS community,

Following is the model that I'm interested to extract blocks of the V matrix using hold option with known or predefined variance parameters. I would like to get V matrix that contains G-side and R-side structure at plot level, i.e.a matrix of 240 x 240 in the following example.

I was able to get V matrix for location as subject, but I'm not sure how to change the subject to a plot level to get V with both G and R structures.

*Generate numbers;

data one;

do har_y = 1 to 2; * two harvest years;

do loc = 1 to 5; * five locations ;

do blk = 1 to 4; * four blocks;

do cult = 1 to 6;*six cultivars;

yield= 1; * dummy value

output;

end;

end;

end;

end;

*Sort

proc print; run;

proc sort data=one out=two;by har_y loc blk cult;

run;

proc mixed lognote data=two;

class cult loc har_y blk;

model yield=cult har_y cult*har_y;

random har_y/subject=loc type=csh;

random har_y/subject=loc*cult type=csh;

random har_y/subject=loc*blk type=csh v ;

repeated har_y/subject=loc*blk*cult type=csh r;

parms(6.1387) (3.9508) ( 0.7087) (0.09453) (0.1454) (0.8500) (0.01970) (0.04729) (0) (0.8706) (0.7041) (0.5151)/noiter hold= 1 2 3 4 5 6 7 8 9 10 11 12; *parameters used to get V matrix

ods output Tests3 = vcfixed CovParms=vcparms FitStatistics= vcfitstat v=v1 r=r1 ;

run;

Thank you in advance for your time and consideration.

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

06-26-2015 01:11 PM

This isn't going to help at all, but I am surprised that you don't get a message in the output that the G matrix is not positive definite.

It appears that the repeated statement has a subject of loc*blk*cult, which I assume is at the plot level. With what you have, I guess I would have approached this somewhat differently:

proc mixed lognote data=two;

class cult loc har_y blk;

model yield=cult har_y cult*har_y;

random intercept cult blk/subject=loc type=csh v;

repeated har_y/subject=loc*blk*cult type=csh r;

Note that this fits 3 random effects (loc, loc*cult and loc*blk) each with different values and constant covariances between all effects. Putting harvest year in as both a random and repeated effect seems hard to justify in my opinion--I would see it as random if measurements were on different plots in the two years, and repeated if on the same plot.

See if this helps any.

Steve Denham

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

06-26-2015 03:38 PM

Thanks, Steve for your suggestion on the model fitting. I was trying different variance-covariance structures on the field trials data set that consists of multiple harvest.Yes,I had a problem with not positive definite G matrix (only for few trials), but I resolved that by providing different initial parms. The three random statements in my original question was to model har_y at three random subject levels, i.e. loc, loc*cult and loc*blk.

RVS

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

06-26-2015 04:39 PM

Hi Steve,

I did replicate proc mixed code into proc glimmix, this provided a V matrix consists of 240 x 240 to the plot level. Do you think it is appropriate?

proc glimmix data=two;

class cult loc har_y blk;

model yield=cult har_y cult*har_y;

random har_y/subject=loc type=un;

random har_y/subject=loc*cult type=un;

random har_y/subject=loc*blk type=un;

random har_y/ residual subject=loc*blk*cult type=un v;

ods output Tests3 = vcfixed CovParms=vcparms FitStatistics= vcfitstat v=glimmix ;

run;

Thank you for your time and consideration.

RVS

- Mark as New
- Bookmark
- Subscribe
- Subscribe to RSS Feed
- Highlight
- Email to a Friend
- Report Inappropriate Content

06-29-2015 09:04 AM

If you are looking to fit a random slopes model (in har_y), this looks good. I still think of it as a random intercepts model, with a single G side random statement:

random intercept cult blk/subject=loc type=un;

I end up thinking this because of the R side factors including all of these variance components in an unstructured covariance structure matrix. My experience would be that the output would throw a message that the G matrix was not positive definite--basically saying that the R part was explaining all of the variability. However, with sufficient data, it could be done.

Steve Denham