That LMatrix!
Locked
 
|
I need some help with the Lmatrix. Here’s the syntax. glm yc by mc3r with c1 c2 x/design c1 c2 x mc3r x*mc3r/print parameter/ lmatrix 'mm2' c1 -.016 c2 .016 x 2 mc3r 1 -1 0 x*mc3r 2 -2 0. Mc3r is three categories. The rest of the variables continuous covariates. I want to compare the first mc3r category to the last (the reference) category for different values of x while holding c1 and c2 at their mean values. The error
is that L is “not estimable”. I’ve read about estimableness but I don’t understand how to apply that knowledge.
(I don’t think that emmeans can be used here but if I am wrong, I’d like to know that fact and, more importantly, how to structure that subcommand.) Thanks, Gene Maguin |
|
|
I was too limited in my imagination. I had thought that emmeans could not be used when there was covariate by factor interaction. I was wrong so that part of my post is retracted. I’d still like to understand how to code the lmatrix. Gene Maguin |
|
|
Only Ryan Black knows for sure...
-----Original Message-----
=====================
To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the
command. To leave the list, send the command
SIGNOFF SPSSX-L
For a list of commands to manage subscriptions, send the command
INFO REFCARD
From: Maguin, Eugene <[hidden email]> To: SPSSX-L <[hidden email]> Sent: Fri, Apr 22, 2016 4:31 pm Subject: That LMatrix! I was too limited in my imagination. I had thought that emmeans could not be used when there was covariate by factor interaction. I was wrong so that part of my post is retracted. I’d still like to understand how to code the lmatrix.
Gene Maguin
|
|
|
In reply to this post by Maguin, Eugene
Gene, If you want to the compare the first category of mc3r to the last category of mc3r with x=2, then the code should be: /lmatrix 'mm2' mc3r 1 0 -1 x*mc3r 2 0 -2 Given the specified design matrix, the coefficients associated with terms c1 and c2 are 0 in the contrast of interest; therefore these terms can be excluded from the lmatrix code. Ryan On Fri, Apr 22, 2016 at 3:21 PM, Maguin, Eugene <[hidden email]> wrote:
|
Administrator
|
Personally, I find it easier to understand LMATRIX sub-commands when one uses the ALL keyword, and explicitly lists the values of all variables. For Gene's model, suppose the means of c1 and c2 are both = 50, and the value of X = 40 for the first contrast of interest. At least in the first draft of the syntax, I would include 3 LMATRIX commands as follows:
/LMATRIX 'M1: mc3r=1, x = 40' ALL 1 1 0 0 50 50 40 40 0 0 /LMATRIX 'M3: mc3r=3, x = 40' ALL 1 0 0 1 50 50 40 0 0 40 /LMATRIX 'M1-M2' ALL 0 1 0 -1 0 0 0 40 0 -40 For the first one (M1), the "estimate" = the fitted value of Y (or yc in Gene's model) for the first category of mc3r with x = 40 and c1 and c2 both = 50. For the second one (M2), the "estimate" = the fitted value of Y (or yc in Gene's model) for the third category of mc3r with x = 40 and c1 and c2 both = 50. And for the third one, the "estimate" is the M1-M3 contrast. You may find it convenient to use Excel to generate the coefficients--those reading via Nabble can see the attached example. Finally, one reason to include the Y1 and Y2 contrasts in the first draft of your syntax is that you can check the values of their estimates against values obtained via EMMEANS. If they don't match, you know your LMATRIX coefficients aren't correct. This is untested, but I think the following EMMEANS sub-command should produce a table that shows the same values as my M1 and M3 contrasts. /EMMEANS = TABLES(mc3r) WITH(c1=MEAN c2=MEAN X=40) Ironically, having done all the work to get to this point, I now suspect that you could get the contrast you want by just tacking COMPARE on to that EMMEANS line: /EMMEANS = TABLES(mc3r) WITH(c1=MEAN c2=MEAN X=40) COMPARE(mc3r) Gene, all of this is untested (in this particular instance), but give it a try, and let us know how it pans out! ;-) LMATRIX_coefficients_for_Gene.xlsx
--
Bruce Weaver bweaver@lakeheadu.ca http://sites.google.com/a/lakeheadu.ca/bweaver/ "When all else fails, RTFM." PLEASE NOTE THE FOLLOWING: 1. My Hotmail account is not monitored regularly. To send me an e-mail, please use the address shown above. 2. The SPSSX Discussion forum on Nabble is no longer linked to the SPSSX-L listserv administered by UGA (https://listserv.uga.edu/). |
|
|
Bruce,
I agree with your conceptual approach, and I have advocated for estimating the means first and then working one's way up to the contrast of interest. A couple of follow-up points, however. You set the values of c1 and c2 equal to their respective means using EMMEANS but the LMATRIX tells us the coefficients for those terms cancel out (equal 0) after subtracting the coefficients of like terms for this contrast. In other words, one can set the values of c1 and c2 to whatever value because the terms cancel out for this contrast anway. Put another way, understanding mathematically and conceptually why the c1 and c2 terms cancel out would allow one to achieve a deeper understanding of the traditional ANCOVA type model Gene is fitting *with respect to* c1 and c2. Finally, why are we still using GLM when the MIXED procedure is available, which does everything GLM and more? :-) Let's now generate some data to support my points. See below my name. Ryan -- *Generate Data. set seed 98765432. new file. input program. loop ID= 1 to 3000. compute x = rv.normal(0,1). compute c1 = rv.normal(0,1). compute c2 = rv.normal(0,1). IF (ID<= 1000) mc3r=1. IF (ID> 1000 and ID <= 2000) mc3r=2. IF (ID> 2000) mc3r=3. compute yc = -1.5 + .25*c1 + .50*c2 + 0.9*x + 2.2*(mc3r=1) + 2.8*(mc3r=2) + 0.5*x + 0.2*(mc3r=1)*x + 1.3*(mc3r=2)*x + rv.normal(0,1). end case. end loop. end file. end input program. execute. glm yc by mc3r with c1 c2 x /design c1 c2 x mc3r x*mc3r /print parameter /emmeans = TABLES(mc3r) WITH(x=2, c1=mean, c2=mean) COMPARE /emmeans = TABLES(mc3r) WITH(x=2, c1=4, c2=4) COMPARE /lmatrix 'mm2' mc3r 1 0 -1 c1 0 c2 0 x*mc3r 2 0 -2. mixed yc by mc3r with c1 c2 x /fixed=c1 c2 x mc3r x*mc3r /print=solution /test 'xbar(mc3r=1) at c1=-.016, c2=.016, x=2)' intercept 1 c1 -.016 c2 .016 x 2 mc3r 1 0 0 x*mc3r 2 0 0 /test 'xbar(mc3r=3) at c1=-.016, c2=.016, x=2)' intercept 1 c1 -.016 c2 .016 x 2 mc3r 0 0 1 x*mc3r 0 0 2 /test 'xbar(mc3r=3) - xbar(mc3r=1) at x=2)' intercept 0 c1 0 c2 0 x 0 mc3r 1 0 -1 x*mc3r 2 0 -2 /test 'xbar(mc3r=3) - xbar(mc3r=1) at x=2)' mc3r 1 0 -1 x*mc3r 2 0 -2. > On Apr 23, 2016, at 8:12 AM, Bruce Weaver <[hidden email]> wrote: > > Personally, I find it easier to understand LMATRIX sub-commands when one uses > the ALL keyword, and explicitly lists the values of all variables. For > Gene's model, suppose the means of c1 and c2 are both = 50, and the value of > X = 40 for the first contrast of interest. At least in the first draft of > the syntax, I would include 3 LMATRIX commands as follows: > > /LMATRIX 'M1: mc3r=1, x = 40' ALL 1 1 0 0 50 50 40 40 0 0 > /LMATRIX 'M3: mc3r=3, x = 40' ALL 1 0 0 1 50 50 40 0 0 40 > /LMATRIX 'M1-M2' ALL 0 1 0 -1 0 0 0 40 0 -40 > > For the first one (M1), the "estimate" = the fitted value of Y (or yc in > Gene's model) for the first category of mc3r with x = 40 and c1 and c2 both > = 50. > > For the second one (M2), the "estimate" = the fitted value of Y (or yc in > Gene's model) for the third category of mc3r with x = 40 and c1 and c2 both > = 50. > > And for the third one, the "estimate" is the M1-M3 contrast. > > You may find it convenient to use Excel to generate the coefficients--those > reading via Nabble can see the attached example. > > Finally, one reason to include the Y1 and Y2 contrasts in the first draft of > your syntax is that you can check the values of their estimates against > values obtained via EMMEANS. If they don't match, you know your LMATRIX > coefficients aren't correct. This is untested, but I think the following > EMMEANS sub-command should produce a table that shows the same values as my > M1 and M3 contrasts. > > /EMMEANS = TABLES(mc3r) WITH(c1=MEAN c2=MEAN X=40) > > Ironically, having done all the work to get to this point, I now suspect > that you could get the contrast you want by just tacking COMPARE on to that > EMMEANS line: > > /EMMEANS = TABLES(mc3r) WITH(c1=MEAN c2=MEAN X=40) COMPARE(mc3r) > > Gene, all of this is untested (in this particular instance), but give it a > try, and let us know how it pans out! ;-) > LMATRIX_coefficients_for_Gene.xlsx > <http://spssx-discussion.1045642.n5.nabble.com/file/n5732009/LMATRIX_coefficients_for_Gene.xlsx> > > > Ryan Black wrote >> Gene, >> >> If you want to the compare the first category of mc3r to the last category >> of mc3r with x=2, then the code should be: >> >> /lmatrix 'mm2' mc3r 1 0 -1 x*mc3r 2 0 -2 >> >> Given the specified design matrix, the coefficients associated with terms >> c1 and c2 are 0 in the contrast of interest; therefore these terms can be >> excluded from the lmatrix code. >> >> Ryan >> >> >> On Fri, Apr 22, 2016 at 3:21 PM, Maguin, Eugene < > >> emaguin@ > >> > wrote: >> >>> I need some help with the Lmatrix. Here’s the syntax. >>> >>> >>> >>> glm yc by mc3r with c1 c2 x/design c1 c2 x mc3r x*mc3r/print parameter/ >>> >>> lmatrix 'mm2' c1 -.016 c2 .016 x 2 mc3r 1 -1 0 x*mc3r 2 -2 0. >>> >>> >>> >>> Mc3r is three categories. The rest of the variables continuous >>> covariates. >>> I want to compare the first mc3r category to the last (the reference) >>> category for different values of x while holding c1 and c2 at their mean >>> values. The error is that L is “not estimable”. I’ve read about >>> estimableness but I don’t understand how to apply that knowledge. >>> >>> >>> >>> (I don’t think that emmeans can be used here but if I am wrong, I’d like >>> to know that fact and, more importantly, how to structure that >>> subcommand.) >>> >>> Thanks, Gene Maguin >>> ===================== To manage your subscription to SPSSX-L, send a >>> message to > >> LISTSERV@.UGA > >> (not to SPSSX-L), with no body text >>> except the command. To leave the list, send the command SIGNOFF SPSSX-L >>> For >>> a list of commands to manage subscriptions, send the command INFO REFCARD >> >> ===================== >> To manage your subscription to SPSSX-L, send a message to > >> LISTSERV@.UGA > >> (not to SPSSX-L), with no body text except the >> command. To leave the list, send the command >> SIGNOFF SPSSX-L >> For a list of commands to manage subscriptions, send the command >> INFO REFCARD > > > > > > ----- > -- > Bruce Weaver > [hidden email] > http://sites.google.com/a/lakeheadu.ca/bweaver/ > > "When all else fails, RTFM." > > NOTE: My Hotmail account is not monitored regularly. > To send me an e-mail, please use the address shown above. > > -- > View this message in context: http://spssx-discussion.1045642.n5.nabble.com/That-LMatrix-tp5732002p5732009.html > Sent from the SPSSX Discussion mailing list archive at Nabble.com. > > ===================== > To manage your subscription to SPSSX-L, send a message to > [hidden email] (not to SPSSX-L), with no body text except the > command. To leave the list, send the command > SIGNOFF SPSSX-L > For a list of commands to manage subscriptions, send the command > INFO REFCARD ===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD |
|
Administrator
|
I should hope you do like this approach, Ryan--I think I learned it from your posts! ;-)
Good question about why we aren't estimating models like this via MIXED. But bear in mind that other procedures like GLM or REGRESSION have some features that one might want to take advantage of. E.g., if one wants to enter variables in several blocks, the REGRESSION procedure is very nice, because it allows you to do this via one command with multiple /ENTER sub-commands, and it gives you the change in R^2 and associated F-test. To do the same thing via MIXED is a bit of a pain: You need a separate command for each block; you need to "roll your own" to do the likelihood ratio test comparing each model to the previous one; and you need to be very careful to ensure that you have the same cases used for every model in the series. If anyone from SPSS is paying attention, I think that it would be a huge improvement to MIXED it had the equivalent of REGRESSION's multiple /ENTER sub-commands and LR tests comparing each model to the previous one. Cheers, Bruce
--
Bruce Weaver bweaver@lakeheadu.ca http://sites.google.com/a/lakeheadu.ca/bweaver/ "When all else fails, RTFM." PLEASE NOTE THE FOLLOWING: 1. My Hotmail account is not monitored regularly. To send me an e-mail, please use the address shown above. 2. The SPSSX Discussion forum on Nabble is no longer linked to the SPSSX-L listserv administered by UGA (https://listserv.uga.edu/). |
|
|
Bruce, You mentioned multiple ENTER statements in the REGRESSION procedure, so I assume you are referring to fixed effects only models. Further, I assume you are referring to a multiple partial F test whereby one tests whether the effect of 2 or more variables significantly contributes above and beyond the variable(s) already in the model. This can be easily accomplished in the MIXED procedure in a single run using the TEST sub-command. Suppose we have three continuous predictors (x1, x2, x3) and we want to test whether adding variables x2 and x3 significantly predicts y above and beyond the model with only variable x1. The TEST statement can easily perform this test as follows: /test 'effect of x2 and x3 above and beyond x1' x2 1; x3 1. Below my name is the simulation and MIXED procedure code that demonstrates how the TEST sub-command produces identical multiple partial F test results (F, df, and p-value) to the REGRESSION procedure with multiple ENTER statements and the CHANGE keyword. Please write back if I misunderstood your question. Ryan -- *Generate Data. set seed 98765432. new file. input program. loop ID= 1 to 50. compute x1 = rv.normal(0,1). compute x2 = rv.normal(0,1). compute x3 = rv.normal(0,1). compute y = -1.5 + .40*x1 + .05*x2 + 0.10*x3 + rv.normal(0,1). end case. end loop. end file. end input program. execute. REGRESSION /STATISTICS COEFF OUTS R ANOVA CHANGE /DEPENDENT y /METHOD=ENTER x1 /METHOD=ENTER x1 x2 x3. mixed y with x1 x2 x3 /fixed=x1 x2 x3 /print=solution /test 'effect of x2 and x3 above and beyond x1' x2 1; x3 1. On Sat, Apr 23, 2016 at 1:22 PM, Bruce Weaver <[hidden email]> wrote: I should hope you do like this approach, Ryan--I think I learned it from your |
|
Administrator
|
Actually, I was thinking about comparing a fixed effects only model to a model with a random intercept, or random intercept model vs random intercept and slope. AFAIK, in order to carry out those likelihood ratio tests, one must use a separate MIXED command for each model, and then extract the -2LL values and number of model parameters from the output, and compute their own test statistic and p-value. To avoid doing all this "manually", one could use OMS to store the needed bits of output. But this is very cumbersome compared to the multiple /ENTER sub-commands for REGRESSION.
I'll give your code for fixed effects only a try when I get back to the office--no SPSS on this computer. Cheers, Bruce
--
Bruce Weaver bweaver@lakeheadu.ca http://sites.google.com/a/lakeheadu.ca/bweaver/ "When all else fails, RTFM." PLEASE NOTE THE FOLLOWING: 1. My Hotmail account is not monitored regularly. To send me an e-mail, please use the address shown above. 2. The SPSSX Discussion forum on Nabble is no longer linked to the SPSSX-L listserv administered by UGA (https://listserv.uga.edu/). |
|
|
Hi Bruce,
I thought you were implying that you were using regression over the mixed on occasion because of the multiple enter statement, which is what led me to my last post. Apologies for my misunderstanding your comment. If you want to test, for example, the subject variance estimate differs significantly from zero, there are three approaches that immediately come to mind: 1. Test the H0 that the difference in subject means = 0 using the TEST sub-command. 2. Test the H0 that the subject variance = 0 using the WALD test with asymptotic errors offered in the PRINT sub-command with the specified keyword TESTCOV. 3. Test the H0 that the subject variance = 0 using a Likelihood Ratio test (LRT) which as you correctly pointed out requires multiple runs. SAS offers LRTs as a sub-command for mixed models. I hope SPSS catches up. Frankly, there are many additional features I'd like to see for MIXED in SPSS. Ryan > On Apr 24, 2016, at 11:49 AM, Bruce Weaver <[hidden email]> wrote: > > Actually, I was thinking about comparing a fixed effects only model to a > model with a random intercept, or random intercept model vs random intercept > and slope. AFAIK, in order to carry out those likelihood ratio tests, one > must use a separate MIXED command for each model, and then extract the -2LL > values and number of model parameters from the output, and compute their own > test statistic and p-value. To avoid doing all this "manually", one could > use OMS to store the needed bits of output. But this is very cumbersome > compared to the multiple /ENTER sub-commands for REGRESSION. > > I'll give your code for fixed effects only a try when I get back to the > office--no SPSS on this computer. > > Cheers, > Bruce > > > Ryan Black wrote >> Bruce, >> >> You mentioned multiple ENTER statements in the REGRESSION procedure, so I >> assume you are referring to fixed effects only models. Further, I assume >> you are referring to a multiple partial F test whereby one tests whether >> the effect of 2 or more variables significantly contributes above and >> beyond the variable(s) already in the model. >> >> This can be easily accomplished in the MIXED procedure in a single run >> using the TEST sub-command. >> >> Suppose we have three continuous predictors (x1, x2, x3) and we want to >> test whether adding variables x2 and x3 significantly predicts y above and >> beyond the model with only variable x1. >> >> The TEST statement can easily perform this test as follows: >> >> /test 'effect of x2 and x3 above and beyond x1' x2 1; x3 1. >> >> Below my name is the simulation and MIXED procedure code that demonstrates >> how the TEST sub-command produces identical multiple partial F test >> results >> (F, df, and p-value) to the REGRESSION procedure with multiple ENTER >> statements and the CHANGE keyword. >> >> Please write back if I misunderstood your question. >> >> Ryan >> >> -- >> >> *Generate Data. >> set seed 98765432. >> new file. >> input program. >> loop ID= 1 to 50. >> compute x1 = rv.normal(0,1). >> compute x2 = rv.normal(0,1). >> compute x3 = rv.normal(0,1). >> compute y = -1.5 + .40*x1 + .05*x2 + 0.10*x3 + rv.normal(0,1). >> end case. >> end loop. >> end file. >> end input program. >> execute. >> >> REGRESSION >> /STATISTICS COEFF OUTS R ANOVA CHANGE >> /DEPENDENT y >> /METHOD=ENTER x1 >> /METHOD=ENTER x1 x2 x3. >> >> mixed y with x1 x2 x3 >> /fixed=x1 x2 x3 >> /print=solution >> /test 'effect of x2 and x3 above and beyond x1' x2 1; x3 1. >> >> --- snip --- > > > > > > ----- > -- > Bruce Weaver > [hidden email] > http://sites.google.com/a/lakeheadu.ca/bweaver/ > > "When all else fails, RTFM." > > NOTE: My Hotmail account is not monitored regularly. > To send me an e-mail, please use the address shown above. > > -- > View this message in context: http://spssx-discussion.1045642.n5.nabble.com/That-LMatrix-tp5732002p5732015.html > Sent from the SPSSX Discussion mailing list archive at Nabble.com. > > ===================== > To manage your subscription to SPSSX-L, send a message to > [hidden email] (not to SPSSX-L), with no body text except the > command. To leave the list, send the command > SIGNOFF SPSSX-L > For a list of commands to manage subscriptions, send the command > INFO REFCARD ===================== To manage your subscription to SPSSX-L, send a message to [hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD |
|
|
In reply to this post by Ryan
Ryan, Thank you.
From: SPSSX(r) Discussion [mailto:[hidden email]]
On Behalf Of Ryan Black Gene, If you want to the compare the first category of mc3r to the last category of mc3r with x=2, then the code should be: /lmatrix 'mm2' mc3r 1 0 -1 x*mc3r 2 0 -2 Given the specified design matrix, the coefficients associated with terms c1 and c2 are 0 in the contrast of interest; therefore these terms can be excluded from the lmatrix code. Ryan On Fri, Apr 22, 2016 at 3:21 PM, Maguin, Eugene <[hidden email]> wrote:
===================== To manage your subscription to SPSSX-L, send a message to
[hidden email] (not to SPSSX-L), with no body text except the command. To leave the list, send the command SIGNOFF SPSSX-L For a list of commands to manage subscriptions, send the command INFO REFCARD
|
«
Return to SPSSX Discussion
|
1 view|%1 views
| Free forum by Nabble | Edit this page |