[R] lme and SAS Proc mixed

Spencer Graves spencer.graves at pdf.com
Mon Jul 3 20:31:06 CEST 2006

	  Your example is entirely too complicated for me to parse in the time 
available, but I have a few questions that I hope might help:

	  First, have you examined str(fit.lme) plus all the other help pages 
listed under "See Also" in the "lme" help page, especially "lmeObject"? 
  With luck, this may answer your questions.

	  Second, are all your random effects nested, or are some crossed?  If 
the model is all standard nesting with no complicated inhomogeneity of 
variance or correlation structure, you may be working too hard:  The 
'nlme' package is designed to make standard nested mixed-model analysis 
as easily as Doug Bates could make it with a general tool in the S 
language paradigm.  Functions like 'pdBlocked' are relatively 
inefficient attempts to adapt the 'nlme' paradigm to support crossed 
random effects, etc.  If you have both crossed and nested random 
effects, you might try Bates' current project, 'lmer' associated with 
the 'lme4' package.  Unfortunately, that package currently has fewer 
helper functions for things like confidence intervals.  (It's Bates' 
leading edge development effort.)

	  Third, can generate a much simpler, self-contained problem that 
exhibits the same features as your more complicated example?  If your 
email had contained a few lines of R code that someone like me could 
copy and paste into R and get what you got, it would increase your 
chances of getting a more informative reply quicker.  Without it, I'm 
left to guessing.

	  Fourth, if you want to know, "which is more correct", have you 
considered Monte Carlo?  The nlme package includes a function 
simulate.lme, which might help you.  The lme4 package includes mcmcsamp.

	  Hope this helps.
	  Spencer Graves

Kellie J. Archer, Ph.D. wrote:
> I am trying to use lme to fit a mixed effects model to get the same
> results as when using the following SAS code:
> proc mixed;
> class refseqid probeid probeno end;
> model expression=end logpgc / ddfm=satterth;
> random probeno probeid / subject=refseqid type=cs;
> lsmeans end / diff cl; run;
> There are 3 genes (refseqid) which is the large grouping factor, with
> 2 probeids nested within each refseqid, and 16 probenos nested within
> each of the probeids.
> I have specified in the SAS Proc Mixed procedure that the
> variance-covariance structure is to be compound symmetric. Therefore,
> the variance-covariance matrix is a block diagonal matrix of the form,
> V_1  0   0
> 0   V_2  0
> 0    0   V3
> where each V_i represents a RefSeqID. Moreover, for each V_i the
> structure within the block is 
> v_{11}   v{12}
> v_{21}   v{22}
> where v_{11} and v_{22} are different probeids nested within the
> refseqid, and so are correlated. The structure of
> both v_{11} and v_{22} are compound symmetric, and v_{12} and v{21}
> contain a constant for all elements of the matrix.
> I have tried to reproduce this using lme, but it is unclear from the
> documentation (and Pinheiro & Bates text) how the pdBlocked and
> compound symmetric structure can be combined.
> fit.lme<-lme(expression~End+logpgc,random=list(RefSeqID=pdBlocked(list
> (~1,~ProbeID-1),pdClass="pdSymm")),data=dataset,correlation=corCompSym
> m(form=~1|RefSeqID/ProbeID/ProbeNo))
> The point estimates are essentially the same comparing R and SAS for
> the fixed effects, but the 95% confidence intervals are much shorter
> using lme(). In order to find the difference in the algorithms used by
> SAS and R I tried to extract the variance-covariance matrix to look at
> its structure. I used the getVarCov() command, but it tells me that
> this function is not available for nested structures. Is there another
> way to extract the variance-covariance structure for nested models?
> Does anyone know how I could get the var-cov structure above using
> lme?
> Kellie J. Archer, Ph.D.
> Assistant Professor, Department of Biostatistics
> Fellow, Center for the Study of Biological Complexity
> Virginia Commonwealth University
> 1101 East Marshall St., B1-066
> Richmond, VA 23298-0032
> phone: (804) 827-2039
> fax: (804) 828-8900
> e-mail: kjarcher at vcu.edu
> website: www.people.vcu.edu/~kjarcher
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide! http://www.R-project.org/posting-guide.html

More information about the R-help mailing list