search for: varcorr

If you are indeed using lme and not lmer then the needed function is VarCorr(). However, 2 recommendations. First, this is a busy list and better emails subject headers get better attention. Second, I would recommend using lmer as it is much faster. However, VarCorr seems to be incompatible with lmer and I do not know of another function to work with lmer. Hence, a better...

Folks In trying to use lmer for a hierarchical model, I encountered the following message: Error in UseMethod("VarCorr") : no applicable method for 'VarCorr' applied to an object of class "mer" foo.mer <- lmer(y ~ TP + (TP|M),data=joe.q) > head(joe.q[,1:5]) TP M AB Trt y 1 1 Jan A NN 19.20002 2 1 Jan A NN 19.06378 3 1 Jan A NN 21.36292 4 1 Jan A NN 16.85191...

I have a loop that I would like to use to extract the "stddev" for each itteration so I can average the "stddev" for all the runs. It would be helpful to know how to extract the "stddev" for each run from the VarCorr. Thanks MCruns<-1000 sighatlvec<-rep(NA,MCruns) sighatbvec<-rep(NA,MCruns) sighatevec<-rep(NA,MCruns) for(mc in 1:MCruns) { samples<-(sqrt(sigbsq)*rnorm(6)) labs<-(sqrt(siglsq)*rnorm(2)) errors<-(sqrt(sigesq)*rnorm(24)) y<-errors+samples[sampfac]+labs[labfac] lmer<-lme...

...) where E is the individual level for consistency analysis, A-D are other fixed and random effects that I have to control for. Following Nakagawa and Scheilzeth I can work out the repeatability estimate using the following (as it is a binomial and the residual variance is fixed at 1). attr(lme4::VarCorr(model)$E, "stddev")^2 / (1*(pi^2)/3 + attr(lme4::VarCorr(model)$E, "stddev")^2 + attr(lme4::VarCorr(model)$C, "stddev")^2 + attr(lme4::VarCorr(model)$D, "stddev")^2 ) My question is can I use stimulate(model) to generate values that I can then use to do p...

Dear All, I would like to ask you to help me with my memeory. I remember using some function that would list all the possible methods I could apply to an object. Say, if I had an object of class=lme, it would tell me that that I could do stuff like qqnorm(myobjct), or VarCorr(myobject). In general, a very complete list. I though this list of all possible methods would pop out by typing methods(class = lme) or, methods(class = whatever class my object is), but, if I type: methods(class = lme), I get: "anova.lme" "plot.lme" "simulate.l...

...s namespace is loaded but not attached. I now have problems to fix this issue. It is easy to get rid of two of the three package warnings by using import(lme4) import(survival) in my NAMESPACE. Yet, I cannot import nlme as I then get the following warnings: Warning: replacing previous import ?VarCorr? when loading ?lme4? Warning: replacing previous import ?lmList? when loading ?lme4? I know that nlme is imported in lme4 but I cannot import both at the same time. Can anyone help me here? I do *not* want or can use importFrom as I rely on multiple functions and also rely on non-exported methods...

...el-2 variance component may be obtained via subscripting, but what about the level-1 variance, viz., the 3.215072 term? (actually this term squared) Didn't see anything in the archives on this. Cheers, David > fm <- lmer( dv ~ time.num*drug + (1 | Patient.new), data=dat.new ) > VarCorr(fm) $Patient.new 1 x 1 Matrix of class "dpoMatrix" (Intercept) (Intercept) 8.519916 attr(,"sc") scale 3.215072 > VarCorr(fm)[[1]][1] [1] 8.519916 > VarCorr(fm)[[2]][1] Error in VarCorr(fm)[[2]] : subscript out of bounds ########################...

...trap functions in the library "bootstrap" to estimate confidence intervals of ICC values calculated in lme. In lme, the ICC is calculated as tau/(tau+sigma-squared). So, for instance the ICC in the following example is 0.116: > tmod<-lme(CINISMO~1,random=~1|IDGRUP,data=TDAT) > VarCorr(tmod) IDGRUP = pdLogChol(1) Variance StdDev (Intercept) 0.1829931 0.427777 Residual 1.3907732 1.179310 > 0.18299/(0.18299+1.39077) [1] 0.1162757 Using the bootstrap library, I can set up theta to do the ICC as follows: >theta<-function(x,DATA){tmod<-lme(CINISMO~1,ra...

Hello, has someone written by chance a function to extract the variance-covariance matrix from a lmer-object? I've noticed the VarCorr function, but it gives unhandy output. Regards, Roel de Jong

...a simple random-effects model using lme: model<-lme(y~1,random=~1|groupA/groupB) which returns the output for the StdDev of the Random effects, and model AIC etc as expected. Until yesterday I was using R v. 2.0, and had no problem in calling the variance components of the above model using VarCorr(model), together with their 95% confidence intervals using intervals() - although for some response variables a call to intervals() returns the error: Cannot get confidence intervals on var-cov components: Non-positive definite approximate variance-covariance. I have now installed R v. 2.3.1 and...

...4 loaded at same time!) # OrthoFem<-Orthodont[Orthodont$Sex=="Female",] # http://tolstoy.newcastle.edu.au/R/help/06/03/23758.html library(lme4) fm1OrthF. <- lmer(distance~age+(age|Subject), data=OrthoFem) lmer(distance~age+(age|Subject), data=OrthoFem)@bVar$Subject[2,2,]* (attr(VarCorr(lmer(distance~age+(age| Subject),data=OrthoFem)),"sc")^2)[1] (attr(VarCorr(fm1OrthF.),"sc")^2)[1] fm1.s <- coef(fm1OrthF.)$Subject fm1.s.var <- fm1OrthF. at bVar$Subject*(attr(VarCorr(fm1OrthF.),"sc")^2)[1] fm1.s0.s <- sqrt(fm1.s.var[1,1,]) fm1.s0.a &l...

...*beta + error and error = U + W(t) + Z where U is the random intercept (normally distributed), W(t) the stationary Gaussian process and Z also a normally distributed (the residual) rv. Each of these three random variables have a variance which I am not sure to which output in lme they belong to. VarCorr gives the intercept and residual variance which I assume belong to U and Z respectively. The output of lme gives another estimate called "range" which I assume belongs to the parameter estimate needed for the Gaussian correlation. Are my assumptions correct? And where can I get the varia...

...le.com/ANOVA-problem-td4609062.html, but now focussed on the problem of extracting variance components. I have got my mixed effects model now /narrow$line<-as.factor(narrow$line) rg.lmer <- lapply(split(narrow, narrow$gene),function(x)lmer(value~sex+(1|line:sex)+(1|line),data=x)) str(lme4::VarCorr(rg.lmer[[1]])) (rg.lmer[[1]])/ and want the variance components from each model (one per gene) put in to a data frame, this works for one gene: /## SUCCESS! for one gene... alpha<-as.matrix(lme4::VarCorr(rg.lmer[[1]])) beta<-data.frame(c(alpha,attr(alpha,"sc")^2)) colnames(beta)&l...

Here is a piece of code fitting a model to a (part) of a dataset, just for illustration. I can extract the random interaction and the residual variance in group meth==1 using VarCorr, but how do I get the other residual variance? Is there any way to get the other variances in numerical form directly - it seems a litte contraintuitive to use "as.numeric" when extracting estimates, it's a bit like good old days writing SAS-programs that reads the SAS output files...

...thought it would be appropriate to follow this procedure: > library(nlme) > fm1Rail.lme <- lme(travel ~ 1, data = Rail, random = ~ 1 | Rail) > summary(fm1Rail.lme) ... Random effects: Formula: ~1 | Rail (Intercept) Residual StdDev: 24.80547 4.020779 ... > as.numeric(VarCorr(fm1Rail.lme)[1]) / (as.numeric(VarCorr(fm1Rail.lme)[1]) + as.numeric(VarCorr(fm1Rail.lme)[2] )) [1] 0.9743987 This gives a reasonable, high repeatability. First, I would like to know, whether this computation of repeatability is correct, and if so, how I would have to proceed if there was also a v...

Dear all, Here is the question: For example, using the "petrol" data offered with R. pet3.lme<-lme(Y~SG+VP+V10+EP,random=~1|No,data=petrol) pet3.lme$sigma gives the residual StdDev. But I can't figure out how to extract the "(intercept) StdDev", although it is in the print out if I do "summary(pet3.lme)". In

Dear all, I'm working with a model to estimate components of variance by using the lme() function. The model is as the following: model=lme(fixed=X~1+as.factor(station),data=myData,na.action=na.exclude,rand om=~N+1|spliceOrHoming) Where X is the response measured variable, station is treated as fixed effects factor, N is a continuous variable, and spliceOrHoming is a (ordered)

...rawley's book (page 361 onwards) so I can check what I am doing. I am tryg to reproduce the nested analysis on page 368: model<-aov(Glycogen~Treatment/Rat/Liver + Error(Treatment/Rat/Liver), rats) using lme. The code: model1<- lme(Glycogen~Treatment, random = ~1|Rat/Liver, data=rats) VarCorr(model1) Variance StdDev Rat = pdLogChol(1) (Intercept) 20.6019981 4.538942 Liver = pdLogChol(1) (Intercept) 0.0540623 0.232513 Residual 42.4362241 6.514309 Does NOT give me the same variance componets I find in Crawley's book (page 371 o...

Hi All. I need a bit of help extracting the residual error variance from the VarCorr structure from lmer. #Here's a 2-way random effects model lmer.1 <- lmer(rating ~ (1|person)+(1|rater), data = dat) #Get the structure vc.fit <- VarCorr(lmer.1) #results in..... $person 1 x 1 Matrix of class "dpoMatrix" (Intercept) (Intercept) 0.7755392 $rat...

...<2.2 coxme 2.2 lmekin 2.2 ------------------------------------------------------ beta fixef fixef fixef fixef fixef b ranef ranef ---- ranef ranef var(beta) vcov vcov ---- vcov vcov var(b) VarCorr VarCorr ranef VarCorr VarCorr Notice that "ranef" in the prior release was out of step with everyone else, which became obvious when working on lmekin. Also, the constructs "fit$coef$random" and "fit$coef$fixed" no longer work due to some internal rea...