multcomp/doBy comparação multipla com interação
Caros colegas Estou fazendo uma análise de um experimento DIC com arranjo fatorial (Gest e Manej). A resposta medida é o peso vivo. A análise de variância mostrou que há efeito de interação entre os dois fatores avaliados. Eu gostaria de saber onde está esta diferenças. A minha dúvida é o seguinte como fazer a matrix de contrastes para a interação? estou utilizando o pacote doBy. A outra dúvida é a seguinte: o R utiliza como padrão o contraste tratamento o qual compara todos os tratamento como o grupo controle( em meu caso Gest=0 e Manej=1) como faço para alterar esta forma de contraste, para um que compare todos os tratamentos entre si ?
dput(abate[,1:5]) structure(list(Animal = c(1L, 2L, 3L, 4L, 5L, 7L, 8L, 9L, 12L, 14L, 15L, 17L, 18L, 19L, 21L, 22L, 23L, 25L, 26L, 28L, 30L, 32L, 34L, 35L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 46L, 47L, 48L, 49L, 50L, 53L, NA, NA), Gest = structure(c(4L, 2L, 3L, 2L, 2L, 4L, 3L, 4L, 2L, 2L, 3L, 4L, 1L, 4L, 2L, 1L, 4L, 1L, 2L, 4L, 3L, 3L, 1L, 4L, 3L, 1L, 2L, 3L, 1L, 1L, 3L, 1L, 2L, 3L, 1L, 2L, 1L, 1L, NA, NA), .Label = c("0", "100", "130", "140"), class = "factor"), Manej = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, NA, NA ), .Label = c("1", "2"), class = "factor"), rep = c(2L, 3L, 1L, 2L, 4L, 3L, 1L, 2L, 2L, 5L, 4L, 3L, 1L, 4L, 1L, 2L, 4L, 4L, 3L, 1L, 3L, 5L, 6L, 1L, 2L, 5L, 1L, 2L, 4L, 3L, 3L, 3L, 4L, 4L, 5L, 5L, 1L, 2L, NA, NA), PV = c(49000L, 40800L, 51800L, 41900L, 35000L, 40500L, 35400L, 42300L, 36600L, 36600L, 38400L, 46400L, 26900L, 50200L, 35500L, 35500L, 48400L, 40500L, 32500L, 42800L, 36500L, 50000L, 31200L, 33500L, 46500L, 40900L, 36500L, 31200L, 31500L, 35500L, 38700L, 30000L, 37900L, 42700L, 35500L, 42000L, 34500L, 36700L, NA, NA)), .Names = c("Animal", "Gest", "Manej", "rep", "PV"), class = "data.frame", row.names = c(NA, -40L))
abate$Gest<-factor(abate$Gest) abate$Manej<-factor(abate$Manej) C<-lm(PV~Gest*Manej,abate) Anova(C,type="III") -- Fernando Antônio de Souza Zootecnista, Msc. Nutricão Animal, UFMG lattes: http://lattes.cnpq.br/6519538815038307
Fernando, Parabéns, seu código é totalmente reproduzível (com excessão da função Anova() que tá ali sem carregar o pacote car, mas eu conhecia). CMR é meio raro na nessa lista. Para resolver seu problema você fazer assim. abate <- structure(list(Animal = c(1L, 2L, 3L, 4L, 5L, 7L, 8L, 9L, 12L, 14L, 15L, 17L, 18L, 19L, 21L, 22L, 23L, 25L, 26L, 28L, 30L, 32L, 34L, 35L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 46L, 47L, 48L, 49L, 50L, 53L, NA, NA), Gest = structure(c(4L, 2L, 3L, 2L, 2L, 4L, 3L, 4L, 2L, 2L, 3L, 4L, 1L, 4L, 2L, 1L, 4L, 1L, 2L, 4L, 3L, 3L, 1L, 4L, 3L, 1L, 2L, 3L, 1L, 1L, 3L, 1L, 2L, 3L, 1L, 2L, 1L, 1L, NA, NA), .Label = c("0", "100", "130", "140"), class = "factor"), Manej = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, NA, NA ), .Label = c("1", "2"), class = "factor"), rep = c(2L, 3L, 1L, 2L, 4L, 3L, 1L, 2L, 2L, 5L, 4L, 3L, 1L, 4L, 1L, 2L, 4L, 4L, 3L, 1L, 3L, 5L, 6L, 1L, 2L, 5L, 1L, 2L, 4L, 3L, 3L, 3L, 4L, 4L, 5L, 5L, 1L, 2L, NA, NA), PV = c(49000L, 40800L, 51800L, 41900L, 35000L, 40500L, 35400L, 42300L, 36600L, 36600L, 38400L, 46400L, 26900L, 50200L, 35500L, 35500L, 48400L, 40500L, 32500L, 42800L, 36500L, 50000L, 31200L, 33500L, 46500L, 40900L, 36500L, 31200L, 31500L, 35500L, 38700L, 30000L, 37900L, 42700L, 35500L, 42000L, 34500L, 36700L, NA, NA)), .Names = c("Animal", "Gest", "Manej", "rep", "PV"), class = "data.frame", row.names = c(NA, -40L)) abate$Gest <- factor(abate$Gest) abate$Manej <- factor(abate$Manej) abate <- na.omit(abate) xtabs(~Gest+Manej, abate) C <- lm(PV~Gest*Manej, abate) anova(C) car::Anova(C,type="III") # interação não significativa à 5%! D <- lm(PV~Gest+Manej, abate) X <- doBy::popMatrix(D, effect="Gest") choose(nlevels(abate$Gest), 2) cb <- combn(nrow(X), 2) Xc <- X[cb[1,],]-X[cb[2,],] Xc rownames(Xc) <- apply(combn(levels(abate$Gest), 2), 2, paste, collapse="-") require(multcomp) summary(glht(D, linfct=Xc)) # em caso de interação significativa fazer # isso aqui para cada nível que deseja fixar levels(abate$Manej) X <- doBy::popMatrix(C, effect="Gest", at=list(Manej="1")) X # o resto segue normal... # para automatizar os desdobramentos você pode usar a lapply() ou # similares. A opção de constrate nada mais é que uma restrição paramétrica no vetor de efeitos que te leva para uma particular solução ou conjunto de estimativas. Como você falou, o treatment anula (impõe/assume) que o efeito do primeiro nível de cada fator é zero, o soma zero impõe que a soma deles é zero. Em todas elas, o que você precisa estimar ao final são p-1 parâmetros/efeitos pois a restrição cuidou de um. Todos contra todos na saída do summary() não tem como especificar porque isso requer mais que p-1 hipóteses você só pode ir com p-1. Mas isso não é problema porque uma vez estimado os efeitos, qualquer função e quantidade de funções estimáveis a partir deles pode ser obtida. À 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 ==========================================================================
Caro Walmes, muito obrigado pela atenção. Foi muito útil suas sugestões. entretanto não estou entendendo algo e gostaria que vc me esclarecesse Coloquei desta vez um CMR com um banco de dados em que a interação deu significativa. Seguindo o comando que você colocou para no caso de interação entretanto não apresenta nenhum efeito significante.. Não entendo o que aconteceu. Você pode me explicar? dput(abate[,1:8])structure(list(Animal = c(1L, 2L, 3L, 4L, 5L, 7L, 8L, 9L, 12L, 14L, 15L, 17L, 18L, 19L, 21L, 22L, 23L, 25L, 26L, 28L, 30L, 32L, 34L, 35L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 46L, 47L, 48L, 49L, 50L, 53L, NA, NA), Gest = structure(c(4L, 2L, 3L, 2L, 2L, 4L, 3L, 4L, 2L, 2L, 3L, 4L, 1L, 4L, 2L, 1L, 4L, 1L, 2L, 4L, 3L, 3L, 1L, 4L, 3L, 1L, 2L, 3L, 1L, 1L, 3L, 1L, 2L, 3L, 1L, 2L, 1L, 1L, NA, NA), .Label = c("0", "100", "130", "140"), class = "factor"), Manej = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, NA, NA ), .Label = c("1", "2"), class = "factor"), rep = c(2L, 3L, 1L, 2L, 4L, 3L, 1L, 2L, 2L, 5L, 4L, 3L, 1L, 4L, 1L, 2L, 4L, 4L, 3L, 1L, 3L, 5L, 6L, 1L, 2L, 5L, 1L, 2L, 4L, 3L, 3L, 3L, 4L, 4L, 5L, 5L, 1L, 2L, NA, NA), PV = c(49000L, 40800L, 51800L, 41900L, 35000L, 40500L, 35400L, 42300L, 36600L, 36600L, 38400L, 46400L, 26900L, 50200L, 35500L, 35500L, 48400L, 40500L, 32500L, 42800L, 36500L, 50000L, 31200L, 33500L, 46500L, 40900L, 36500L, 31200L, 31500L, 35500L, 38700L, 30000L, 37900L, 42700L, 35500L, 42000L, 34500L, 36700L, NA, NA), PCVZ = c(39841.2, 35457.1, 42343.2, 32989.7, 26566.5, 30339.9, 25965.5, 28475.2, 30815.7, 28162.6, 28103.5, 37905.8, 19174.5, 34319.4, 27963.2, 26609.4, 40852.9, 34781.1, 23411.8, 43490.8, 28900.1, 38240.8, 22463.9, 24631.9, 36483.9, 33688.2, 31838.7, 26048.1, 25646.6, 22815.6, 31369.6, 22173.7, 30371.2, 36740.8, 25773, 38171.6, 26401.4, 26711.2, NA, NA), CAR.PVZ = c(0.474383302, 0.518936969, 0.493585747, 0.512281106, 0.508158771, 0.451550598, 0.427490324, 0.431954824, 0.490009962, 0.511316427, 0.462575836, 0.485413842, 0.547602284, 0.416673951, 0.507810265, 0.541162146, 0.474874489, 0.554899069, 0.478391239, 0.494357427, 0.46020602, 0.489006506, 0.534190412, 0.418156943, 0.482404567, 0.558058905, 0.515096408, 0.45684714, 0.538083021, 0.547870755, 0.471794349, 0.518632434, 0.520229691, 0.435483169, 0.52768401, 0.458456025, 0.537850266, 0.539099703, NA, NA), PCARC = c(18900L, 18400L, 20900L, 16900L, 13500L, 13700L, 11100L, 12300L, 15100L, 14400L, 13000L, 18400L, 10500L, 14300L, 14200L, 14400L, 19400L, 19300L, 11200L, 21500L, 13300L, 18700L, 12000L, 10300L, 17600L, 18800L, 16400L, 11900L, 13800L, 12500L, 14800L, 11500L, 15800L, 16000L, 13600L, 17500L, 14200L, 14400L, NA, NA)), .Names = c("Animal", "Gest", "Manej", "rep", "PV", "PCVZ", "CAR.PVZ", "PCARC"), class = "data.frame", row.names = c(NA, -40L)) library(car) library(doBy) library(multcomp) abate$Gest<-factor(abate$Gest) abate$Manej<-factor(abate$Manej) C<-lm(PCARC~Gest*Manej,abate) Anova(C,type="III") levels(abate$Manej) X<-popMatrix(C, effect="Gest", at=list(Manej="2")) cb <- combn(nrow(X), 2) Xc <- X[cb[1,],]-X[cb[2,],] summary(glht(C, linfct=Xc)) popMeans(C, effect="Gest", at=list(Manej="2")) Em 3 de dezembro de 2012 16:39, Walmes Zeviani <walmeszeviani@gmail.com>escreveu:
Fernando,
Parabéns, seu código é totalmente reproduzível (com excessão da função Anova() que tá ali sem carregar o pacote car, mas eu conhecia). CMR é meio raro na nessa lista. Para resolver seu problema você fazer assim.
abate <- structure(list(Animal = c(1L, 2L, 3L, 4L, 5L, 7L, 8L, 9L, 12L, 14L, 15L, 17L, 18L, 19L, 21L, 22L, 23L, 25L, 26L, 28L, 30L, 32L, 34L, 35L, 37L, 38L, 39L, 40L, 41L, 42L, 43L, 44L, 46L, 47L, 48L, 49L, 50L, 53L, NA, NA), Gest = structure(c(4L, 2L, 3L, 2L, 2L, 4L, 3L, 4L, 2L, 2L, 3L, 4L, 1L, 4L, 2L, 1L, 4L, 1L, 2L, 4L, 3L, 3L, 1L, 4L, 3L, 1L, 2L, 3L, 1L, 1L, 3L, 1L, 2L, 3L, 1L, 2L, 1L, 1L, NA, NA), .Label = c("0", "100", "130", "140"), class = "factor"), Manej = structure(c(1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L, 2L, 1L, 2L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, NA, NA ), .Label = c("1", "2"), class = "factor"), rep = c(2L, 3L, 1L, 2L, 4L, 3L, 1L, 2L, 2L, 5L, 4L, 3L, 1L, 4L, 1L, 2L, 4L, 4L, 3L, 1L, 3L, 5L, 6L, 1L, 2L, 5L, 1L, 2L, 4L, 3L, 3L, 3L, 4L, 4L, 5L, 5L, 1L, 2L, NA, NA), PV = c(49000L, 40800L, 51800L, 41900L, 35000L, 40500L, 35400L, 42300L, 36600L, 36600L, 38400L, 46400L, 26900L, 50200L, 35500L, 35500L, 48400L, 40500L, 32500L, 42800L, 36500L, 50000L, 31200L, 33500L, 46500L, 40900L, 36500L, 31200L, 31500L, 35500L, 38700L, 30000L, 37900L, 42700L, 35500L, 42000L, 34500L, 36700L, NA, NA)), .Names = c("Animal", "Gest", "Manej", "rep", "PV"), class = "data.frame", row.names = c(NA, -40L))
abate$Gest <- factor(abate$Gest) abate$Manej <- factor(abate$Manej) abate <- na.omit(abate)
xtabs(~Gest+Manej, abate)
C <- lm(PV~Gest*Manej, abate) anova(C) car::Anova(C,type="III")
# interação não significativa à 5%! D <- lm(PV~Gest+Manej, abate)
X <- doBy::popMatrix(D, effect="Gest")
choose(nlevels(abate$Gest), 2) cb <- combn(nrow(X), 2) Xc <- X[cb[1,],]-X[cb[2,],] Xc rownames(Xc) <- apply(combn(levels(abate$Gest), 2), 2, paste, collapse="-")
require(multcomp)
summary(glht(D, linfct=Xc))
# em caso de interação significativa fazer # isso aqui para cada nível que deseja fixar
levels(abate$Manej) X <- doBy::popMatrix(C, effect="Gest", at=list(Manej="1")) X
# o resto segue normal... # para automatizar os desdobramentos você pode usar a lapply() ou # similares.
A opção de constrate nada mais é que uma restrição paramétrica no vetor de efeitos que te leva para uma particular solução ou conjunto de estimativas. Como você falou, o treatment anula (impõe/assume) que o efeito do primeiro nível de cada fator é zero, o soma zero impõe que a soma deles é zero. Em todas elas, o que você precisa estimar ao final são p-1 parâmetros/efeitos pois a restrição cuidou de um. Todos contra todos na saída do summary() não tem como especificar porque isso requer mais que p-1 hipóteses você só pode ir com p-1. Mas isso não é problema porque uma vez estimado os efeitos, qualquer função e quantidade de funções estimáveis a partir deles pode ser obtida.
À 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 ==========================================================================
_______________________________________________ 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.
Isso que aconteceu com você acontece com frequência. Veja, sua interação foi significativa à 4% Sum Sq Df F value Pr(>F) (Intercept) 1254528000 1 301.3933 < 2e-16 *** Gest 39011053 3 3.1241 0.04044 * Manej 28043758 1 6.7374 0.01448 * Gest:Manej 39571188 3 3.1689 *0.03859 ** Residuals 124872833 30 Bem, agora você vai fazer comparações duas à duas. O seu número de hipóteses avaliadas é grande pois você fará k(k-1)/2 contrastes, onde k é o número de níveis. Quando eleva-se o número de hipóteses eleva-se também a chance de erro tipo I. Se cada hipótese tem 95% de confiança (ou 5% de significância), ao testar 10 hipóteses *independentes*, a chance de *pelo menos uma* ser significativa por mero acaso é 1-0.95^10=0.40, ou seja, bem distante do 5%. Ou seja, o nível de significância global fica longe do nominal. Para contornar isso faça a engenharia inversa, quanto deve ser nível individual para o global ser 5%? Diversos tipos de testes de comparação multipla existem justamente para fazer essa engenharia inversa (Tukey, SNK, t protegido por Bonferroni, etc). Cada um deles procura resolver o mesmo problema: a manutenção do nível global de significância para comparações múltiplas. A glht() faz correção nos p-valores dos testes t. Vários métodos estão disponíveis e o padrão é o single-step. Mudar o método muda o apelo e consequentemente o rigor, uns serão mais liberais (tenderão a dar mais diferenças) que outros. Quando você pega p-valores do F perto do 5%, como no seu caso, se o número de comparações for grande (6 no seu caso), existe chance de não haver diferenças. Por exemplo, para ser significativo pelo critério de Bonferroni, o nível individual de cada teste t teria que ser inferior à 5%/6=0.833, ou seja, uma hipótese individual tem que ser significativa à 1% para representar 5% no global. É natural você pagar esse preço, afinal, você tá testando várias hipóteses. Sou dono de uma barraca de tiro ao alvo numa feira. O prêmio é R$ 10 se acertar o alvo com 3 tentativas. Se você diz "quanto eu ganho se eu acertar com 10 tentativas?", lógico que não vou te pagar R$ 10. Pagaria bem menos, você terá muito mais chances. Tenho que manter a esperança matemática no mesmo lugar. Você pode trocar a opção de correção de p-valor para que veja diferenças mas não é recomendado fazer. Os testes de hipótese bem como as hipóteses devem ser definidos antes do experimento/análise dos dados e não decididos durante análise para favorecer certos resultados. Em caráter exploratório, você pode fazer o teste à 10% para recomendar estratégias/selecionar níveis para realizar um experimento futuro. X1 <- popMatrix(C, effect="Gest", at=list(Manej="1")) X2 <- popMatrix(C, effect="Gest", at=list(Manej="2")) cb <- combn(nrow(X), 2) Xc1 <- X1[cb[1,],]-X1[cb[2,],] Xc2 <- X2[cb[1,],]-X2[cb[2,],] summary(glht(C, linfct=Xc1)) summary(glht(C, linfct=Xc2)) summary(glht(C, linfct=Xc1), test=adjusted(*type="none"*)) summary(glht(C, linfct=Xc2), test=adjusted(*type="none"*))
summary(glht(C, linfct=Xc1)) Linear Hypotheses: Estimate Std. Error t value Pr(>|t|) 1 == 0 -180 1290 -0.139 0.9990 2 == 0 -1760 1290 -1.364 0.5307 3 == 0 -3710 1369 -2.711 0.0508 . 4 == 0 -1580 1290 -1.224 0.6163 5 == 0 -3530 1369 -2.579 0.0680 . 6 == 0 -1950 1369 -1.425 0.4940 (Adjusted p values reported -- *single-step method*)
summary(glht(C, linfct=Xc1), test=adjusted(type="none")) Linear Hypotheses: Estimate Std. Error t value Pr(>|t|) 1 == 0 -180 1290 -0.139 0.890 2 == 0 -1760 1290 -1.364 0.183 3 == 0 -3710 1369 -2.711 0.011 * 4 == 0 -1580 1290 -1.224 0.230 5 == 0 -3530 1369 -2.579 0.015 * 6 == 0 -1950 1369 -1.425 0.165 (Adjusted p values reported -- *none method*)
À 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 ==========================================================================
Muito obrigado Walmes! Mais uma vez foi muito esclarecedor sua resposta. Penso que é mais sensato manter a analis com o default da função glht mesmo. grato. Em 3 de dezembro de 2012 22:14, Walmes Zeviani <walmeszeviani@gmail.com>escreveu:
Isso que aconteceu com você acontece com frequência. Veja, sua interação foi significativa à 4%
Sum Sq Df F value Pr(>F) (Intercept) 1254528000 1 301.3933 < 2e-16 *** Gest 39011053 3 3.1241 0.04044 * Manej 28043758 1 6.7374 0.01448 * Gest:Manej 39571188 3 3.1689 *0.03859 ** Residuals 124872833 30
Bem, agora você vai fazer comparações duas à duas. O seu número de hipóteses avaliadas é grande pois você fará k(k-1)/2 contrastes, onde k é o número de níveis. Quando eleva-se o número de hipóteses eleva-se também a chance de erro tipo I. Se cada hipótese tem 95% de confiança (ou 5% de significância), ao testar 10 hipóteses *independentes*, a chance de *pelo menos uma* ser significativa por mero acaso é 1-0.95^10=0.40, ou seja, bem distante do 5%. Ou seja, o nível de significância global fica longe do nominal. Para contornar isso faça a engenharia inversa, quanto deve ser nível individual para o global ser 5%?
Diversos tipos de testes de comparação multipla existem justamente para fazer essa engenharia inversa (Tukey, SNK, t protegido por Bonferroni, etc). Cada um deles procura resolver o mesmo problema: a manutenção do nível global de significância para comparações múltiplas. A glht() faz correção nos p-valores dos testes t. Vários métodos estão disponíveis e o padrão é o single-step. Mudar o método muda o apelo e consequentemente o rigor, uns serão mais liberais (tenderão a dar mais diferenças) que outros. Quando você pega p-valores do F perto do 5%, como no seu caso, se o número de comparações for grande (6 no seu caso), existe chance de não haver diferenças. Por exemplo, para ser significativo pelo critério de Bonferroni, o nível individual de cada teste t teria que ser inferior à 5%/6=0.833, ou seja, uma hipótese individual tem que ser significativa à 1% para representar 5% no global.
É natural você pagar esse preço, afinal, você tá testando várias hipóteses. Sou dono de uma barraca de tiro ao alvo numa feira. O prêmio é R$ 10 se acertar o alvo com 3 tentativas. Se você diz "quanto eu ganho se eu acertar com 10 tentativas?", lógico que não vou te pagar R$ 10. Pagaria bem menos, você terá muito mais chances. Tenho que manter a esperança matemática no mesmo lugar.
Você pode trocar a opção de correção de p-valor para que veja diferenças mas não é recomendado fazer. Os testes de hipótese bem como as hipóteses devem ser definidos antes do experimento/análise dos dados e não decididos durante análise para favorecer certos resultados. Em caráter exploratório, você pode fazer o teste à 10% para recomendar estratégias/selecionar níveis para realizar um experimento futuro.
X1 <- popMatrix(C, effect="Gest", at=list(Manej="1")) X2 <- popMatrix(C, effect="Gest", at=list(Manej="2")) cb <- combn(nrow(X), 2) Xc1 <- X1[cb[1,],]-X1[cb[2,],] Xc2 <- X2[cb[1,],]-X2[cb[2,],] summary(glht(C, linfct=Xc1)) summary(glht(C, linfct=Xc2)) summary(glht(C, linfct=Xc1), test=adjusted(*type="none"*)) summary(glht(C, linfct=Xc2), test=adjusted(*type="none"*))
summary(glht(C, linfct=Xc1)) Linear Hypotheses: Estimate Std. Error t value Pr(>|t|) 1 == 0 -180 1290 -0.139 0.9990 2 == 0 -1760 1290 -1.364 0.5307 3 == 0 -3710 1369 -2.711 0.0508 . 4 == 0 -1580 1290 -1.224 0.6163 5 == 0 -3530 1369 -2.579 0.0680 . 6 == 0 -1950 1369 -1.425 0.4940 (Adjusted p values reported -- *single-step method*)
summary(glht(C, linfct=Xc1), test=adjusted(type="none")) Linear Hypotheses: Estimate Std. Error t value Pr(>|t|) 1 == 0 -180 1290 -0.139 0.890 2 == 0 -1760 1290 -1.364 0.183 3 == 0 -3710 1369 -2.711 0.011 * 4 == 0 -1580 1290 -1.224 0.230 5 == 0 -3530 1369 -2.579 0.015 * 6 == 0 -1950 1369 -1.425 0.165 (Adjusted p values reported -- *none method*)
À 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 ==========================================================================
_______________________________________________ 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.
participantes (3)
-
Fernando A. Souza -
Fernando Antonio de souza
-
Walmes Zeviani