Walmes, o mestre da didática. Parabéns!!<div><br></div><div>Walmes não deu para assistir a sua comunicação. Esepro que tenha transcorrido tudo bem.</div><div><br></div><div>Abs<br><br><div class="gmail_quote">Em 9 de maio de 2012 16:05, Walmes Zeviani <span dir="ltr"><<a href="mailto:walmeszeviani@gmail.com" target="_blank">walmeszeviani@gmail.com</a>></span> escreveu:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><font face="trebuchet ms,sans-serif">Ivan,<br><br>Vamos primeiro primeiro apodar as arestas a fim de fazer a coisa rolar. Você tem um experimento com um fator de 5 níveis nominais, denominado de F. Cada animal é uma unidade experimental do qual se observa uma resposta y e uma covariável x. Acredita-se que x afeta y e portanto quer-se remover esse efeito para poder estudar a variação explicada por F. Pode-se supor que F e x tenham efeitos aditivos (sem interação) ou pode-se assumir que existe interação. Na ausência de interação as diferenças entre níveis de F são as mesmas sem importar o nível de x, mas o mesmo não vale para presença de interação, de forma que, você deve definir antecipadamente em que nível(is) de x quer comparar níveis de F. Não sei o CMR é exatamente o que você tem, pois o CMR não parece ser um caso de modelo misto porque animal é a unidade experimental e não um fator de agrupamento (embora para Poisson este seja um modelo estimável), enfim... vou considerar um glm fixo, daí é só fazer as devidas extensões para o glmm.<br>

<br><span style="font-family:courier new,monospace">da <- data.frame(F=gl(5,20), x=runif(5*20))</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">da$Ey <- with(da, exp(1+as.numeric(F)/2+(as.numeric(F)-2)/3*x))</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">da$y <- with(da, rpois(F, lambda=exp(1+as.numeric(F)/2+(as.numeric(F)-2)/3*x)))</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">with(da, tapply(y, F, mean))</span><br style="font-family:courier new,monospace">

<br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">require(lattice)</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">xyplot(y~x, groups=F, data=da, type="a")</span><br style="font-family:courier new,monospace">

<br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">ac0 <- glm(y~F*x, data=da, family=poisson)</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">anova(ac0, test="Chisq")</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">betas <- coef(ac0)</span><br style="font-family:courier new,monospace"><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">require(doBy)</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace"># médias populacionais marginais em x=0.5</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">mpm <- popMatrix(ac0, effect="F", at=list(x=0.5))</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">mpm%*%betas</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">exp(mpm%*%betas)</span><br style="font-family:courier new,monospace">

<br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace"># diferença nas médias populacionais marginais duas à duas</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">nlevels(da$F)</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">comp <- t(combn(1:nlevels(da$F), 2))</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">rownames(comp) <- apply(comp, 1, paste, collapse="-")</span><br style="font-family:courier new,monospace">

<br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">dmpm <- t(apply(comp, 1, function(vs) mpm[vs[1],]-mpm[vs[2],]))</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">rownames(dmpm) <- rownames(comp)</span><br style="font-family:courier new,monospace">

<br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">apply(dmpm, 1, function(x) x%*%betas)</span><br style="font-family:courier new,monospace"><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">require(multcomp)</span><br style="font-family:courier new,monospace"><span style="font-family:courier new,monospace">summary(glht(ac0, linfct=mpm))  # médias populacionais marginais</span><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace">summary(glht(ac0, linfct=dmpm)) # contraste entre as médias populacionais marginais</span><br style="font-family:courier new,monospace"><br style="font-family:courier new,monospace">

<span style="font-family:courier new,monospace"># na verdade isso não é a média, mas sim o exp() disso</span><br style="font-family:courier new,monospace"><br>E é isso aí Ivan. Esse procedimento aí que fiz foi o que deixou um determinado software famoso por dar esses resultados aos usuários. Trata-se de algo que, depois de fazer pela primeira vez e entender, passa a ser trivial.<br>

<br>À disposição.<br>Walmes.<br><br clear="all"></font><span style="font-family:trebuchet ms,sans-serif">==========================================================================</span><br style="font-family:trebuchet ms,sans-serif">

<span style="font-family:trebuchet ms,sans-serif">Walmes Marques Zeviani</span><br style="font-family:trebuchet ms,sans-serif"><span style="font-family:trebuchet ms,sans-serif">LEG (Laboratório de Estatística e Geoinformação, 25.450418 S, 49.231759 W)</span><br style="font-family:trebuchet ms,sans-serif">

<span style="font-family:trebuchet ms,sans-serif">Departamento de Estatística - Universidade Federal do Paraná</span><br style="font-family:trebuchet ms,sans-serif"><span style="font-family:trebuchet ms,sans-serif">fone: <a href="tel:%28%2B55%29%2041%203361%203573" value="+554133613573" target="_blank">(+55) 41 3361 3573</a></span><br style="font-family:trebuchet ms,sans-serif">

<span style="font-family:trebuchet ms,sans-serif">VoIP: (3361 3600) 1053 1173</span><br style="font-family:trebuchet ms,sans-serif"><span style="font-family:trebuchet ms,sans-serif">e-mail: <a href="mailto:walmes@ufpr.br" target="_blank">walmes@ufpr.br</a></span><br style="font-family:trebuchet ms,sans-serif">

<span style="font-family:trebuchet ms,sans-serif">twitter: @walmeszeviani</span><br style="font-family:trebuchet ms,sans-serif"><span style="font-family:trebuchet ms,sans-serif">homepage: <a href="http://www.leg.ufpr.br/%7Ewalmes" target="_blank">http://www.leg.ufpr.br/~walmes</a></span><br style="font-family:trebuchet ms,sans-serif">

<span style="font-family:trebuchet ms,sans-serif">linux user number: 531218</span><br style="font-family:trebuchet ms,sans-serif"><span style="font-family:trebuchet ms,sans-serif">==========================================================================</span><br>

<br>_______________________________________________<br>
R-br mailing list<br>
<a href="mailto:R-br@listas.c3sl.ufpr.br">R-br@listas.c3sl.ufpr.br</a><br>
<a href="https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br" target="_blank">https://listas.inf.ufpr.br/cgi-bin/mailman/listinfo/r-br</a><br>
Leia o guia de postagem (<a href="http://www.leg.ufpr.br/r-br-guia" target="_blank">http://www.leg.ufpr.br/r-br-guia</a>) e forneça código mínimo reproduzível.<br></blockquote></div><br><br clear="all"><div><br></div>-- <br>
<div>Guilherme Moraes Ferraudo</div><div>Jaboticabal/Campinas - SP</div><div><a href="http://lattes.cnpq.br/2096118558794430" target="_blank">http://lattes.cnpq.br/2096118558794430</a></div><br>
</div>