0
$\begingroup$

My question is whether it’s possible to compute lsmeans defined in this SAS algorithm if the design matrix is not in GLM form. In particular, in R, if one feeds that design to model.matrix(), then Reference coding is used, so A will occupy 2 columns, B will occupy 1 column, A*B will occupy 2 columns, and C will take 1 column.

The meaning of columns in GLM and Ref coding is quite different. For instance, if we have just one factor with two levels, A1 and A2, then in GLM format the intercept corresponds to the response averaged over A1 and A2, but in Ref the intercept corresponds to response at A1.

Apparently, the column space of GLM and Ref is the same, so I wonder if there is a way to represent each GLM column as a linear combination of Ref columns. In that case, one could take the LSM vector defined in GLM terms and apply it to a model fitted in Ref coding. Likewise, a contrast vector can be set as a difference of two LSM vectors and then estimated in Ref coding.

I know there are some R packages, like "contrast" that can take a fitted model in Ref coding and a contrast or lsmeans specification in string format, like "LSM(A1, B1)" or "A1 vs A2". However, here I am interested in a more specific solution that assumes we have already implemented the SAS algorithm for LSM vector. I am looking for a GLM->Ref "adapter" of sorts.

Thanks,

James

$\endgroup$
  • 1
    $\begingroup$ I'm pretty sure the emmeans package has extensive support for contrasts. Without having checked, I'd expect that you can pass whichever contrasts you desire. The only question should be how difficult it is and if possibly the contrasts you need can be created automatically or if you need to pass the contrast matrix. $\endgroup$ – Roland Feb 18 at 9:04
  • $\begingroup$ Have you looked at the CONTRAST and ESTIMATE procedures in GLM: support.sas.com/documentation/cdl/en/statug/63962/HTML/default/… and support.sas.com/documentation/cdl/en/statug/63962/HTML/default/… $\endgroup$ – StatsStudent Feb 18 at 15:39
  • $\begingroup$ Those two links describe how things work in GLM coding. I myself posted a similar link for LSMEANS. Again, the question is how to represent each column of GLM design as a linear combination of columns from Ref design. $\endgroup$ – James Feb 18 at 16:21
  • 1
    $\begingroup$ You seem to be asking about R and SAS at the same time, so I’m confused. Anyway, if it’s about R and you use the emmeans package, you can estimate any contrast and it doesn’t matter what coding was used in fitting the model. $\endgroup$ – Russ Lenth Feb 18 at 16:35
  • 1
    $\begingroup$ The joint_tests() function does type III tests. For understanding what emmeans does, I suggest reading the “basics” vignette. $\endgroup$ – Russ Lenth Feb 19 at 14:25
0
$\begingroup$

It turns out that lsmeans coefficients computed according to SAS guidelines are exactly the same for GLM and Ref coding. Consider crossed factors FactorA (A1, A2, A3) and FactorB (B1, B2). Suppose there is also a continuous covariate, Cov. If one runs:

model.matrix(~ FactorA + FactorB + FactorA:FactorB + Cov)

in R, the design matrix is exactly the same as that used in PROC GLM where ref="A1" and ref="B1" are specified in CLASS statement. Even though Ref design matrix has seven columns, PROC GLM will report 13 lsmeans coefficients. They are identical to lsmeans coefficients reported when PROC GLM or MIXED are run without the "ref" specification which amounts to using GLM coding for the design matrix.

Therefore, if one knows how to obtain the lsmeans coefficients with SAS guidelines, they can be applied directly to a model fitted in R, if one just matches the labels between the two. In this example, SAS will generate (13 - 7) = 6 lsmeans coefficients whose labels are not found in the model.matrix() column names, but those can be ignored because Ref coding assumes the corresponding regression coefficients are zeros.

| cite | improve this answer | |
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.