Walmes, Obrigado pela explicação, mais creio que possa não ter sido claro em minha descrição, não tenho bloco dentro dos locais, são os locais o meu "Bloco". Se eu não repetir essas cultivar dentro dos locais eu não tenho a possibilidade de estimar a variância dentro do local. Minha maior duvida e sobre a codificação das repetições que se colocando 2 e 1, como se fosse rep 1 e 2 isso estaria correto.
dados <- expand.grid(Variedade=factor(c(letters[1:3],letters[1:5])),Local=factor(1:NL)) dados$Rep <-(rep(c(2,2,2,1,1,1,1,1),t=NL)) dados$Resp <- sort(rnorm(nrow(dados)))+10 desta forma minha duvida vem mais em saber se o modelo declarado desta forma: m1 <- lmer(Resp~1+(Rep|Local)+(1|Variedade),dados) esta estimando os componentes de forma correta. Então a grande questão é como declarar a Rep. Obrigado Em 20 de outubro de 2012 20:59, Walmes Zeviani <[email protected]>escreveu: > Éder, > > Na simulação faltou o efeito de bloco. Bem, nas formulas para efeito > aleatório da nlme e lme4 nos temos que declarar usando termos que vão antes > de depois da barra em pé | (pipe). Antes da barra nos colocamos um fator de > efeito fixo e após um de efeito aleatório. Assim, essas expecificações tem > a seguinte interpretação > > * (1|local) : 1 representa o parâmetro beta_0 ou intercepto e com isso > dizemos que os níveis de local são de efeito aleatório e incidem no > intercepto, ou seja, esses desvios somam-se ao intercepto e existe um para > cada nível de local; > * (1|local)+(x|local) : considerando x numérico (ex: dose, tempo) estamos > dizendo que o efeito aleatório de local incide em b_0 e b_1 que o > coeficiente angular, taxa de incremento por unidade de x, e isso é o modelo > de intercepto e inclinação aleatório, note que x tem efeito fixo. > * (1|local)+(trat|local): considerando trat como fator, isso significa que > o efeito aleatório de local causa desvios em cada um dos parâmetros > associados aos níveis de trat (efeito fixo). Dificilmente você ira declarar > uma formular com um fator antes do pipe, aqui será estimado uma variância > para bloco dentro de cada nível de trat, pouco usado. Na prática poderia-se > ter (ano|local), que significa que para ano o efeito de local tem uma > variância diferente (demora bastante para estimar). > * (1|local:bloco): com isso você declara que existe uma fonte de variação > cujos níveis são as combinações (produto cartesiano) dos níveis dos fatores > envolvidos e que o efeito deste incide no intercepto. > > Bem, com isso o termo que você incluiu (0+Rep|Local) não tem sentido > olhando para sua descrição do experimento, embora coisas sejam estimadas. > Na realidade, se você repete níveis de um fator dentro do bloco você > aumenta a precisão na estimativa daqueles repetidos e aumenta a informação > para estimação do erro resídual. Você não precisa declarar essa fonte de > variação, ela "sai por diferença". Como é o ultimo termo de efeito > aleatório (que todo modelo tem) não precisa ser declarado. Então acredito > que você não esteja definindo corretamente o modelo que você descreveu. Se > você muda a estrutura do efeito aleatório algum impacto terá nas > estimativas dos efeitos fixos, por isso que m0 e m1 diferem. No CMR abaixo > eu simulo do modelo que entendi da sua descrição e faço ajustes de alguns > modelos alternativos. > > > #------------------------------------------------------------------------------------------ > # Experimento com efeito de locais, locais/bloco e variedade > # dados artificiais (usando muitos níveis para estimativas ficarem > próximas dos valores > # usados na simulação, facilita depuração) > nl <- 20; nb <- 5; nc <- 12; nr <- 3; nu <- nc+nr > # efeito de local ~N(0, 2), 20 locais > l <- rnorm(nl, 0, 2) > # efeito de local/bloco ~N(0, 1), 5 blocos por local > lb <- rnorm(nl*nb, 0, 1) > # efeito de local/bloco/unidade ~N(0, 0.5), 8 unidades por bloco > lbu <- rnorm(nl*nb*nu, 0, 0.5) > # efeito de cultivar ~N(0, 1.5), 5 cultivares > c <- rnorm(nc, 0, 1.5) > da1 <- expand.grid(ano=gl(2,1), local=gl(nl,1), bloco=gl(nb,1), > cultivar=gl(nc,1), r=1) > da2 <- expand.grid(ano=gl(2,1), local=gl(nl,1), bloco=gl(nb,1), > cultivar=gl(nr,1), r=2) > da <- rbind(da1, da2) > Z <- with(da, cbind(model.matrix(~-1+local), > model.matrix(~-1+local:bloco), > model.matrix(~-1+cultivar))) > dim(Z) > ranef <- c(l, lb, c) > length(ranef) > X <- model.matrix(~ano, da) > y <- X%*%c(1,0)+Z%*%ranef+lbu > > require(lme4) > > # modelo declarado é o mesmo usado para simular > m0 <- lmer(y~ano+(1|local)+(1|local:bloco)+(1|cultivar), data=da) > summary(m0) > > # local fixo, bloco dentro de local e cultivar aleatório > m1 <- lmer(y~ano+local+(1|local:bloco)+(1|cultivar), data=da) > summary(m1) > > # cultivar fixo, bloco dentro de local aleatório > m2 <- lmer(y~ano+cultivar+(1|local/bloco), data=da) # local/bloco = > local+local:bloco > summary(m2) # cultivares mais repetidas, 1:3, tem menor erro padrão > > # efeito de local/bloco com variâncias diferentes nos anos (omiti o outros > efeitos) > m3 <- lmer(y~(ano|local/local), data=da) > summary(m3) # cultivares mais repetidas, 1:3, tem menor erro padrão > > > #------------------------------------------------------------------------------------------ > > À disposição. > Walmes. > > ========================================================================== > Walmes Marques Zeviani > LEG (Laboratório de Estatística e Geoinformação, 25.450418 S, 49.231759 W) > Departamento de Estatística - Universidade Federal do Paraná > fone: (+55) 41 3361 3573 > VoIP: (3361 3600) 1053 1173 > e-mail: [email protected] > skype: walmeszeviani > twitter: @walmeszeviani > homepage: http://www.leg.ufpr.br/~walmes > linux user number: 531218 > ========================================================================== > > > 2012/10/20 tiago souza marçal <[email protected]> > >> Não sei se eu entendi bem mas a sua duvida é saber qual dos modelos >> abaixo usar? >> >> Yijk = m + Gi + B/Ajk + Aj + GAij + Eijk >> >> Yijk = m + Gi + Bk + Aj + GAij + Eijk >> >> Att. >> >> Tiago. >> >> ------------------------------ >> Date: Sat, 20 Oct 2012 15:36:31 -0300 >> From: [email protected] >> To: [email protected] >> Subject: Re: [R-br] Blocos com algumas repetições do material >> >> >> >> Boa tarde Pessoal, >> >> Tenho experimentos em blocos(Locais) e dentro deles gostaria de por >> algumas de repetições de alguns tratamentos para estimar a variância dentro >> do bloco, que desconfio ser grande. Por exemplo, tenho 5 Cultivares >> gostaria de repetir 3 delas dentro do bloco, tendo ao todo 8 parcelas, o >> ideal seria ter uma repetição de tudo, porem na pratica isso não >> é possível. >> Minhas duvidas são: >> 1) Estou definindo de maneira correta essas repetições dentro do modelo? >> 2) Como o método de calculo dos dois modelos abaixo é igual para os >> efeitos fixos, por que o intercepto esta com valores diferentes? >> 3) Neste dois casos(m0 e m1, abaixo) o BLUP seria a média mais o efeito >> da cultivar? >> Obrigado pela ajuda >> >> require(lme4) >> ### Simulando dados >> dados <- >> expand.grid(Variedade=factor(c(letters[1:3],letters[1:5])),Local=factor(1:5)) >> dados$Rep <- factor(rep(c(2,2,2,1,1,1,1,1),t=5)) >> dados$Resp <- sort(rnorm(nrow(dados))) >> dados >> >> ### Modelos >> # Sem efeito dentro do local >> m0 <- lmer(Resp~1+(1|Local)+(1|Variedade),dados) >> summary(m0) >> ranef(m0) >> # Com efeito dentro do local >> m1 <- lmer(Resp~1+(0+Rep|Local)+(1|Variedade),dados) >> summary(m1) >> ranef(m1) >> >> >> _______________________________________________ R-br mailing list >> [email protected] >> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br Leia o guia de >> postagem (http://www.leg.ufpr.br/r-br-guia) e forne�a c�digo m�nimo >> reproduz�vel. >> >> _______________________________________________ >> R-br mailing list >> [email protected] >> https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br >> Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça >> código mínimo reproduzível. >> > > > _______________________________________________ > R-br mailing list > [email protected] > https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br > Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça > código mínimo reproduzível. >
_______________________________________________ R-br mailing list [email protected] https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br Leia o guia de postagem (http://www.leg.ufpr.br/r-br-guia) e forneça código mínimo reproduzível.
