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 <walmeszeviani@gmail.com>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: walmes@ufpr.br skype: walmeszeviani twitter: @walmeszeviani homepage: http://www.leg.ufpr.br/~walmes linux user number: 531218 ==========================================================================
2012/10/20 tiago souza marçal <tiagosouzamarcal@hotmail.com>
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: eder@leg.ufpr.br To: R-br@listas.c3sl.ufpr.br 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 R-br@listas.c3sl.ufpr.br 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 R-br@listas.c3sl.ufpr.br 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 R-br@listas.c3sl.ufpr.br 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.