search for: pdident

...raction effects tim*trt*ctr. The response variable is ge. Here is my lme code: dat <- data.frame(ge=ge,trt=factor(trt),tim=factor(tim),ctr=factor(ctr)) dat$grp = as.factor(rep(1, nrow(dat))) # dim(dat) = 216x5 w <- lme(ge ~ 1,data=dat, random = list( grp = pdBlocked(list( pdIdent(~ trt - 1) ,pdIdent(~ tim - 1) ,pdIdent(~ ctr - 1) ,pdIdent(~ trt:tim - 1) ,pdIdent(~ trt:ctr - 1) ,pdIdent(~ tim:ctr - 1) ,pdIdent(~ trt:tim:ctr - 1...

...he code: library(nlme) fossil <- read.table("fossil.dat",header=T) x <- fossil$age y <- 100000*fossil$strontium.ratio knots <- seq(94,121,length=25) n <- length(x) X <- cbind(rep(1,n),x) Z <- outer(x,knots,"-") Z <- Z*(Z>0) fit <- lme(y~-1+X,random=pdIdent(~-1+Z)) When I ran the code fit <- lme(y~-1+X,random=pdIdent(~-1+Z)) I got an error message: Error in getGroups.data.frame(dataMix, groups) : Invalid formula for groups I was really puzzled. I asked Dr. Ngo and he said they did it in S-plus but not R. Does anyone knows how to do it in R?...

...rmula I used: fm$pmu <- getGroups(fm, ~1|pop/mum, level=2) fm$grp = as.factor(rep(1,nrow(fm))) fm$pl <- getGroups(fm, ~1|plot) fm$mo <- getGroups(fm, ~1|moist) fm$su <- getGroups(fm, ~1|sub) > fm1 <- lme(sqrt(mass) ~ iheight + moist*sub*pop, data=fm, random=list(grp=pdBlocked(list(pdIdent(~pl - 1), pdIdent(~pmu - 1), pdIdent(~pmu:su - 1), pdIdent(~pmu:mo - 1))))) Error in chol((value + t(value))/2) : non-positive definite matrix in chol I know the problem is with the random interaction terms, but I don't know how to overcome this. Any advice would be greatly appreciated. I...

...6]) gpdat <- groupedData(Y1 ~X1 + X2 + X3 + X4 + X5 + m_sum|.g, data = gpdat) the model fitted using: library(Matrix) library(nlme) m_sum <- rowSums(gpdat[,11:27]) mod1 <- lme(fixed = Y1 ~ X1 + X2 + X3 + X4 + X5 + m_sum, random = pdBlocked(list(pdIdent(~1), pdIdent (~ X6 - 1), pdIdent (~ X7 - 1), pdIdent(~ rand_mat -1))), data = gpdat) Which should recover the variance components: var_label var_est rand_mat_scalar 0.00021983 X6_scalar 0.62314002...

Dear R-help group, How can I define a lme with 3 factors(a,b,c), where c is nested in b, and a is crossed with b/c? I think that: lme(response ~ ..., data = Data, random = pdBlocked(list(pdIdent(~ a - 1), pdIdent(~ b - 1)))) is one part of the answer and: lme(response~..., data=Data, random=~1|b/c) is the other part of the answer but how can I combine them?? Could anybody please help me? Many greetings, Thomas PS: Is there good literatur available which could help working with mixe...

I've posted the following to R-help on May 15. It has reproducible R code for real data -- and a real (academic, i.e unpaid) consultion background. I'd be glad for some insight here, mainly not for myself. In the mean time, we've learned that it is to be expected for anova(*, "marginal") to be contrast dependent, but still are glad for advice if you have experience. Thank

...ews/2002-05/msg00095.html http://134.148.236.121/R/help/06/08/32763.html that this is the expected behavior because of the way grouping works with lme(). My questions are: 1. is this the correct random argument to lmer: anova(lme(Acc~A*B*C,random=list(Sub=pdBlocked(list( pdIdent(~1), pdIdent(~A-1), pdIdent(~B-1), pdIdent(~C-1)))),data=data)) 2. How much do the multiple comparisons depend on the random statement? 3. I'm also playing with lmer: Acc~A*B*C+(1|Sub) Is this the correct lmer call for the cr...

Douglas Bates wrote: > > The GLMM function in the lme4 package allows you to specify crossed > random effects within the random argument without the need for the > pdBlocked and pdIdent constructions. Simply ensure that your grouping > factors are defined in such a way that each distinct group has a > different level in the grouping factor (this is usually not a problem > for crossed grouping factors but can be a problem with nested factors) > and list them. For...

...ollapse="+"),"-1")) > fmla1<-as.formula(paste("~",paste(colnames(ldata[,14:18]),collapse="+"),"-1")) > fmla2<-as.formula(paste("~",paste(colnames(ldata[,19:23]),collapse="+"),"-1")) > Block=pdBlocked(list(pdIdent(fmla0),pdIdent(fmla1),pdIdent(fmla2))) > lme(fixed=Score ~ factor(time)-1,data=ldata,random=list(Block), + weights=varIdent(form= ~ 1|time), + correlation=corSymm(form= ~ 1|Block/ID)) and the Error message is Error in parse(file, n, text, prompt) : syntax error in "~ " I repe...

Try this data(Assay) as1 <- lme(logDens~sample*dilut, data=Assay, random=pdBlocked(list( pdIdent(~1), pdIdent(~sample-1), pdIdent(~dilut-1)))) update(as1,random=pdCompSymm(~sample-1)) update(as1,random=pdCompSymm(~sample-1)) update(as1,random=pdCompSymm(~sample-1)) update(as1,random=pdCompSymm(~sample-1)) I'm getting different results on differen...

...ing lme and from that, calculate the intraclass correlation as a measure of rater reliability. So I analysed it using a two-factor, unreplicated random effects model, with sitenumber and rater as crossed random effects. The call to lme was: fit <- lme(score ~ 1, random=list(level=pdBlocked(list(pdIdent(~sitenumber-1), pdIdent(~rater-1)))),data=scores) which makes sense, and I then calculated the intraclass correlation using the within sitenumber (ie among raters) variance, the within rater variance and the residual variance (using VarCorr(fit)): rho = var_among/(var_among + var_within + var_res...

...I am using lme function in nlme package. Here is the R code I have 1. using lme: First I created a grouped data using gx <- groupedData(y ~ 1 | Sample, data=x) gx$dummy <- factor(rep(1,nrow(gx))) then I run the lme fm<- lme(y ~ 1, data=gx, random=list(dummy=pdBlocked(list(pdIdent(~Sample-1), pdIdent(~Operator-1), pdIdent(~Operator:Run-1))))) finally, I use VarCorr to extract the variance components vc <- VarCorr(fm) Variance StdDev Operator:Run 1.595713e-10:20 1.263215e-05:20 Sam...

...ome difficulty bringing commands over from S-PLUS to R. I have the following setup (modified from Ngo and Wand, 2004): set.seed(1) x <- runif(200) y <- sin(3*pi*x) + rnorm(200)*.4 ## library(splines) z <- ns(x, 4) The following runs without error on S-PLUS f <- lme(y ~ 1, random = pdIdent(~ -1 + z)) But in R I get library(nlme) > f <- lme(y ~ 1, random = pdIdent(~ -1 + z)) Error in getGroups.data.frame(dataMix, groups) : Invalid formula for groups Does the S-PLUS lme have some default setting that R doesn't? Any help would be appreciated. Thanks, -roger

...block, sample within block and dilut within block is random. This model can be fitted with (where I define Assay2 to get an ordinary data frame rather than a grouped data object): Assay2 <- as.data.frame(Assay) fm2<-lme(logDens~sample*dilut, data=Assay2, random=list(Block = pdBlocked(list(pdIdent(~1), pdIdent(~sample-1),pdIdent(~dilut-1))) )) Now, block has only 2 levels so I prefer to treat it as fixed: fm3<-lme(logDens~Block+sample*dilut, data=Assay2, random=list(Block = pdBlocked(list(pdIdent(~sample-1),pdIdent(~dilut-1))) )) This works fine until I try a summary() or anova(),...

Hello all: I am interested in computing what the multilevel modeling literature calls a multiple membership model. More specifically, I am working with a data set involving clients and providers. The clients are the lower-level units who are nested within providers (higher-level). However, this is not nesting in the usual sense, as clients can belong to multple providers, which I understand

...2010 1.137097586 10.9 3 2010:191 I want to include both nesting attempt (hack) and respirometer channel (ch) as random factors in a model trying to explain variation in oxygen consumption. From Pinheiro & Bates (2000), I've gathered that this model could be fit making use of pdBlocked and pdIdent, so I've tried fitting the below model: m1.bmr<-with(bmred.df,lme(mlmin~treat*year+massout,random=pdBlocked(list(pdIdent(~hack-1),pdIdent(~ch-1))) )) However, my model fails with the following error message: Error in getGroups.data.frame(dataMix, groups) : Invalid formula for groups I...

...to fit a fixed effect of Cross, and estimate random effects for females and males. (data available as an *.rda for anyone who wants it...) Following the advice from various postings from R-help and from Pinhero and Bates, I can successfully (I think) code an lme() model, using pdBlocked, pdIdent and a dummy grouping variable in a grouped data object. > load("~/RbinaryData.rda") > dd<-data.frame(RbinaryData,dum=factor(rep(1,120))) > names(dd) [1] "Cross" "Female" "Male" "Mate1" "dum" > summary(dd) Cross...

hello! this is a question, how can i specify the random part in the GLMM-call (of the lme4 library) for compound matrices just in the the same way as they defined in the lme-Call (of the nlme library). For example i would just need random=list(my.Subject=pdBlocked(list(pdIdent(~... , ...),pdIdent(~... , ...)))) this specification , if i also attach library(nlme) , is not accepted in the GLMM-call, though the simple form random=list(my.Subject=pdIdent(~...,...)) is accepted. what is the analogous of pdBlocke...

...andom blocked structure (i.e. ~var.a1+var.a2+var.a3, ~var.b1+var.b2,~var.c1+var.c2+var.c3+var.c4), and each one part of it i would like to model as Identity-structure matrix. So i had, in symbols of nlme-package, and for only one cluster-variable my.Subject: random=list(my.Subject=pdBlocked(list(pdIdent(~var.a1+var.a2+var.a3,...),pdIdent(~var.b1+var.b2,...),pdIdent(~var.c1+var.c2+var.c3+var.c4),...))) As the lme4-package doesn't use the pdMat-classes for specification of the random-part in GLMM, i used, with advice of Douglas Bates, the explicit specification in the GLMM-call (this call...

...plement a NON NESTED random effects structure in lme() In particular I've tried the following translation from lmer to lme, as suggested from some web example mod1<-lmer(y~x*z+(x*z|factorA1/factorB)+(x*z|factorA2/factorB)) # y,x,z continuous mod2<-lme(y~x*z, random= pdBlocked(list(pdIdent(~1|factorA1/factorB ),pdIdent(~1|factorA2/factorB)))) In detail check how I've tried to state in mod1 that Iwant to evaluate randomness in the interaction x*z (i.e intercept, slope, interaction) grouped by by a general nesting structure that sets factorA1 and factorA2 as same level effects (h...