Qual é a alturas média das årvores em uma floresta?
Departamento de Ecologia IB-USP
Qual é a alturas média das årvores em uma floresta?
Fazer afirmaçÔes sobre as caracterĂsticas de uma população, com base em informaçÔes dadas por amostras.
População Ă© qualquer conjunto de elementos que vocĂȘ estĂĄ observando
Amostra é qualquer subconjunto desta população.
## [1] 22.72789
## [1] 32.3065
## [1] 5.683881
A variĂĄvel aleatĂłria X Ă© uma variĂĄvel que tem um valor Ășnico (determinado aleatoriamente) para cada resultado de um experimento (lembre-se esse termo veio da teoria de probabilidade).
Exemplos:
Altura das ĂĄrvores de uma floresta (!!)
NĂșmero de presas capturadas por um predador em um determinado dia;
Comprimento de um peixe adulto selecionado aleatoriamente.
As variĂĄveis aleatĂłrias podem ser discretas ou contĂnuas.
Função que define as probabilidades P(X) para cada valor possĂvel da variĂĄvel X.
Distribuição de Probabilidades.
Discreta: nĂșmero de sucessos em n tentativas . p - prob sucesso.
NĂșmero sementes predadas em cada experimento contendo 20 semente.
NĂșmero de animais infectados em relação ao nĂșmero de capturados em cada local
Discreta: Probabilidade de uma sĂ©rie de eventos ocorrem em um perĂodo fixo de tempo, ĂĄrea, volume, quadrante, etc.
Dados de contagem
Abundùncia de espécies em um local
NĂșmero de indivĂduos capturados por tempo
dpois(5,10)
SimĂ©trica, contĂnua
Metade (50%) dos dados estå acima (e abaixo) da média
68% dados estå dentro de 1 desvio padrão da média
95% estå dentro de 2 desvios padrÔes da média
Virtualmente todos os valores estão dentro de 3 desvios padrÔes da média
"Toda soma de variĂĄveis aleatĂłrias independentes de mĂ©dia finita e variĂąncia limitada Ă© aproximadamente Normal, desde que o nĂșmero de termos da soma seja suficientemente grande."
Independentemente do tipo de distribuição da população, na medida em que o tamanho da amostra aumenta, a distribuição das médias amostrais tende a uma distribuição Normal.
levels(arv$sp)
## [1] "sp1" "sp2"
Serå que a altura média de uma espécie nessa floresta é diferente da altura média da outra espécie?
Ou seja, serå que as espécies possuem uma mesma distribuição de probabilidades e as diferenças encontradas são devido ao acaso ou serå que cada espécie é uma variåvel aleatória com médias diferentes?
mean(arv$alt[arv$sp=="sp2"]) - mean(arv$alt[arv$sp=="sp1"])
## [1] 5.109466
Posso afirmar que as médias são diferentes?
H0 - não hå diferença entre as alturas das duas espécies.
A média da altura da espécie 1 é igual à média da altura da espécie 2 e a diferença observada se deve ao acaso.
H1 - hå diferença entre as alturas das duas espécies.
A média da altura da espécie 1 é diferente da média da altura da espécie 2, elas vem de distribuiçÔes diferentes.
AFIRMAĂĂO: nĂŁo hĂĄ diferença entre as alturas das espĂ©cies,
VERDADE: hå diferenças.
AFIRMAĂĂO: hĂĄ diferença entre as alturas das espĂ©cies,
VERDADE: não hå diferenças.
ERRO DO TIPO I (α): rejeitar a hipótese nula dada que ela é verdadeira.
ERRO DO TIPO II (B): aceitar a hipĂłtese nula dado que ela Ă© falsa.
Contrastar a probabilidade de que a diferença que vocĂȘ encontrou nas mĂ©dias das alturas das duas espĂ©cies pode se considerada devido ao acaso ou nĂŁo.
Qual é a probabilidade de que as duas amostras pertencem à mesma população de medidas?
Podemos ou não rejeitar a hipótese nula (de que a diferença é ao acaso), baseado em alguma probabilidade de estarmos cometendo algum erro (α).
"O nĂșmero de observaçÔes independentes menos o nĂșmero de parĂąmetros estimados."
Quantas observaçÔes independentes foram usadas para calcular a média e a variùncia dos dados sob a hipótese nula?
Nosso exemplo: 100 ĂĄrvores - 2 parĂąmetros = 98 graus de liberdade
à uma medida que reflete a diferença esperada entre a a situação quando H0 é correta e quando a H1 é a correta.
Tendo o t e os graus de liberdade recorremos à distribuição t para saber a probabilidade de obter o valor observado dado que a H0 é verdade.
Se esta probabilidade for muito pequena, a chance de estarmos caindo no erro do tipo 1 Ă© pequena.
t.test(arv$alt~arv$sp, var.equal=T)
## ## Two Sample t-test ## ## data: arv$alt by arv$sp ## t = -5.0125, df = 98, p-value = 2.387e-06 ## alternative hypothesis: true difference in means is not equal to 0 ## 95 percent confidence interval: ## -7.132311 -3.086621 ## sample estimates: ## mean in group sp1 mean in group sp2 ## 20.17316 25.28262
Grande famĂlia de modelos estatĂsticos que modelam uma relação linear entre a variĂĄvel dependente Y e a(s) variĂĄvel(is) independente(s) X.
teste t, ANOVA, ANCOVA, regressĂŁo, etc
LMs: y tem distribuição normal
GLMs: y tem outra distribuição (binomial, poisson) que é linearizada através de uma funçÔa de ligação
Comparando alturas de 3 espécies de plantas
H0: que as médias das alturas das 3 espécies não diferem.
SutraĂmos cada altura das ĂĄrvores pela mĂ©dia geral e elevamos ao quadrado.
(media.geral <- mean(arv3$alt))
## [1] 23.28284
(medias.sps <- apply(arv4, 2, mean))
## sp1 sp2 sp3 ## 20.44245 21.78334 27.62274
dif.geral <- arv4 - media.geral (ss.especies <- dif.geral^2)
## sp1 sp2 sp3 ## [1,] 4.00337202 17.121094 185.447876 ## [2,] 11.60801175 21.347420 30.089610 ## [3,] 0.07824755 3.323227 106.599051 ## [4,] 15.75374565 26.834280 24.948725 ## [5,] 9.15530094 3.247214 77.710675 ## [6,] 224.43179074 23.066347 4.313177 ## [7,] 15.28543337 7.866884 1.539904 ## [8,] 53.88994194 35.359990 7.051024 ## [9,] 7.21577099 15.631997 26.264064 ## [10,] 1.52250306 25.688498 46.856173
(ss.total <- sum(ss.especies))
## [1] 1033.251
(ss.sp1=sum((arv4$sp1-medias.sps["sp1"])^2))
## [1] 262.2655
(ss.sp2=sum((arv4$sp2-medias.sps["sp2"])^2))
## [1] 157.0018
(ss.sp3=sum((arv4$sp3-medias.sps["sp3"])^2))
## [1] 322.4728
(ss.intra=ss.sp1+ss.sp2+ss.sp3)
## [1] 741.7401
medias.sps
## sp1 sp2 sp3 ## 20.44245 21.78334 27.62274
media.geral
## [1] 23.28284
ss.entre=10*sum((medias.sps-media.geral)^2) ss.entre
## [1] 291.5112
#conferindo os cĂĄlculos ss.intra+ss.entre
## [1] 1033.251
ss.total
## [1] 1033.251
Soma dos quadrados / graus de liberdade
graus de liberdade entre grupos = 3 - 1
graus de liberdade dentro dos brupos = 30 - 3
ms.entre=ss.entre/2 ms.intra=ss.intra/27 ms.entre
## [1] 145.7556
ms.intra
## [1] 27.47186
F = desvio médio entre grupos / desvio médio intra grupos
#F - razĂŁo das variĂąncias F.sps=ms.entre/ms.intra F.sps
## [1] 5.305634
Vamos ver na distribuição F qual a probabilidade de termos econtrado o valor F.sps
ao acaso:
Colocando os dados calculados em nossa tabela de anova
A função que constrĂłi esta tabela de ANOVA e faz o teste estatĂstico Ă© a aov
. Vamos comparar:
anova.sps <- aov(alt~sp, data=arv3) summary(anova.sps)
## Df Sum Sq Mean Sq F value Pr(>F) ## sp 2 291.5 145.76 5.306 0.0114 * ## Residuals 27 741.7 27.47 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Agora nĂłs sabemos de ondem vem os nĂșmeros nesta tabela e como foi calculado o P, certo?
veja os grĂĄficos
faça teste à posteriori - Tukey HSD
TukeyHSD(anova.sps)
## Tukey multiple comparisons of means ## 95% family-wise confidence level ## ## Fit: aov(formula = alt ~ sp, data = arv3) ## ## $sp ## diff lwr upr p adj ## sp2-sp1 1.340893 -4.4708805 7.152667 0.8360276 ## sp3-sp1 7.180300 1.3685259 12.992073 0.0132022 ## sp3-sp2 5.839406 0.0276327 11.651180 0.0487468
lm.anova.sp <- lm(alt~sp, data=arv3) summary.aov(lm.anova.sp) #mesmo que o summary do aov
## Df Sum Sq Mean Sq F value Pr(>F) ## sp 2 291.5 145.76 5.306 0.0114 * ## Residuals 27 741.7 27.47 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(lm.anova.sp)
## ## Call: ## lm(formula = alt ~ sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -12.1407 -3.0002 -0.0742 3.2719 9.2780 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 20.442 1.657 12.334 1.32e-12 *** ## spsp2 1.341 2.344 0.572 0.57202 ## spsp3 7.180 2.344 3.063 0.00492 ** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 5.241 on 27 degrees of freedom ## Multiple R-squared: 0.2821, Adjusted R-squared: 0.229 ## F-statistic: 5.306 on 2 and 27 DF, p-value: 0.01139
alt.nutr <- lm(alt ~ nutri, data = arv3); summary(alt.nutr)
## ## Call: ## lm(formula = alt ~ nutri, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -6.7176 -1.7171 0.0088 1.4124 9.2072 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 14.38711 1.05832 13.594 7.44e-14 *** ## nutri 0.58006 0.05971 9.714 1.82e-10 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.906 on 28 degrees of freedom ## Multiple R-squared: 0.7712, Adjusted R-squared: 0.763 ## F-statistic: 94.37 on 1 and 28 DF, p-value: 1.818e-10
Y continuo
alt.sp.nutr <- lm(alt ~ nutri + sp, data = arv3); summary(alt.sp.nutr)
## ## Call: ## lm(formula = alt ~ nutri + sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.6682 -0.9554 -0.1585 1.3581 7.6134 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 13.3971 1.1034 12.142 3.23e-12 *** ## nutri 0.5256 0.0561 9.369 8.09e-10 *** ## spsp2 1.6421 1.1423 1.437 0.16251 ## spsp3 3.8326 1.1965 3.203 0.00357 ** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.553 on 26 degrees of freedom ## Multiple R-squared: 0.836, Adjusted R-squared: 0.817 ## F-statistic: 44.16 on 3 and 26 DF, p-value: 2.401e-10
alt.sp.nutr3 <- lm(alt ~ nutri*sp, data=arv3); summary(alt.sp.nutr3)
## ## Call: ## lm(formula = alt ~ nutri * sp, data = arv3) ## ## Residuals: ## Min 1Q Median 3Q Max ## -5.6724 -0.9444 -0.3313 1.0416 7.0702 ## ## Coefficients: ## Estimate Std. Error t value Pr(>|t|) ## (Intercept) 13.40167 1.42727 9.390 1.66e-09 *** ## nutri 0.52527 0.08906 5.898 4.38e-06 *** ## spsp2 2.93717 1.97080 1.490 0.149 ## spsp3 0.43636 2.75666 0.158 0.876 ## nutri:spsp2 -0.10095 0.12423 -0.813 0.424 ## nutri:spsp3 0.17187 0.14350 1.198 0.243 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## Residual standard error: 2.474 on 24 degrees of freedom ## Multiple R-squared: 0.8578, Adjusted R-squared: 0.8282 ## F-statistic: 28.96 on 5 and 24 DF, p-value: 2.003e-09
summary.aov(alt.sp.nutr3)
## Df Sum Sq Mean Sq F value Pr(>F) ## nutri 1 796.8 796.8 130.179 3.53e-11 *** ## sp 2 66.9 33.5 5.467 0.0111 * ## nutri:sp 2 22.6 11.3 1.846 0.1796 ## Residuals 24 146.9 6.1 ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Outras distribuiçÔes de probabilidades além da Normal
VariĂąncias nĂŁo homogĂȘneas
Mais comuns:
Poisson â Ășteis para dados de contagem
Binomial â Ășteis para dados com proporçÔes, ou presença/ausĂȘncia
Uma hipótese sobre a distribuição da variåvel resposta Yi.
Especificação do modelo (variåveis X).
3.Função de ligação entre a média Yi e as variåveis X. Função que lineariza a regressão
Distribuição | link function |
---|---|
normal | identity |
poisson | log |
binomial | logit |
Gamma | reciprocal |
Y - NĂșmero de atropelamentos de anuros em uma estrada
X - distĂąncia a um Parque Nacional.
mod <- glm(TOT.N ~ D.PARK, data = atrop, family = poisson) summary(mod)
## ## Call: ## glm(formula = TOT.N ~ D.PARK, family = poisson, data = atrop) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -8.1100 -1.6950 -0.4708 1.4206 7.3337 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 4.316e+00 4.322e-02 99.87 <2e-16 *** ## D.PARK -1.059e-04 4.387e-06 -24.13 <2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for poisson family taken to be 1) ## ## Null deviance: 1071.4 on 51 degrees of freedom ## Residual deviance: 390.9 on 50 degrees of freedom ## AIC: 634.29 ## ## Number of Fisher Scoring iterations: 4
anova(mod, test = "Chisq")
## Analysis of Deviance Table ## ## Model: poisson, link: log ## ## Response: TOT.N ## ## Terms added sequentially (first to last) ## ## ## Df Deviance Resid. Df Resid. Dev Pr(>Chi) ## NULL 51 1071.4 ## D.PARK 1 680.55 50 390.9 < 2.2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Y - Quantidade de veados infectados dentre os capturados
X - quantidade de vegetação aberta
y <- cbind(veados$infectado, veados$amostrado - veados$infectado) mod.veado <- glm(y ~ openland, data=veados, family="binomial") summary(mod.veado)
## ## Call: ## glm(formula = y ~ openland, family = "binomial", data = veados) ## ## Deviance Residuals: ## Min 1Q Median 3Q Max ## -3.8052 -1.1463 0.6707 1.3416 3.5483 ## ## Coefficients: ## Estimate Std. Error z value Pr(>|z|) ## (Intercept) 2.4820 0.1736 14.29 <2e-16 *** ## openland -6.0336 0.4432 -13.61 <2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 ## ## (Dispersion parameter for binomial family taken to be 1) ## ## Null deviance: 376.18 on 31 degrees of freedom ## Residual deviance: 121.30 on 30 degrees of freedom ## AIC: 204.65 ## ## Number of Fisher Scoring iterations: 5
anova(mod.veado, test = "Chisq")
## Analysis of Deviance Table ## ## Model: binomial, link: logit ## ## Response: y ## ## Terms added sequentially (first to last) ## ## ## Df Deviance Resid. Df Resid. Dev Pr(>Chi) ## NULL 31 376.18 ## openland 1 254.88 30 121.30 < 2.2e-16 *** ## --- ## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Teste de Shapiro-Wilks
shapiro.test(arv3$alt)
## ## Shapiro-Wilk normality test ## ## data: arv3$alt ## W = 0.96262, p-value = 0.3609
shapiro.test(arv$alt[arv$sp == "sp1"])
## ## Shapiro-Wilk normality test ## ## data: arv$alt[arv$sp == "sp1"] ## W = 0.9906, p-value = 0.9591
shapiro.test(arv$alt[arv$sp == "sp2"])
## ## Shapiro-Wilk normality test ## ## data: arv$alt[arv$sp == "sp2"] ## W = 0.97453, p-value = 0.3502
Teste de Bartlett:
bartlett.test(arv3$alt~arv3$sp)
## ## Bartlett test of homogeneity of variances ## ## data: arv3$alt by arv3$sp ## Bartlett's K-squared = 1.1109, df = 2, p-value = 0.5738
Teste de Levene:
leveneTest(arv3$alt~arv3$sp)
## Levene's Test for Homogeneity of Variance (center = median) ## Df F value Pr(>F) ## group 2 0.512 0.605 ## 27
RESĂDUO = Y(observado) - Y(dados ajustados)
A premissa de normalidade na verdade recai sobre os resĂduos do modelo.
resid(anova.sps)
## 1 6 11 16 21 26 ## 0.8395548 -0.5666556 3.1201253 -1.1287012 5.8661704 -12.1406501 ## 31 36 41 46 51 56 ## 6.7500566 -4.5005793 0.1541789 1.6065001 -2.6382600 -3.1208224 ## 61 66 71 76 81 86 ## 3.3224765 -3.6806771 -0.3024985 6.3022481 -1.3052922 7.4459310 ## 91 96 101 106 111 116 ## -2.4542276 -3.5688778 9.2780228 1.1454976 5.9847796 -9.3347719 ## 121 126 131 136 141 146 ## 4.4754638 -2.2630830 -5.5808309 -6.9952785 0.7849457 2.5052549
#residuals(anova.sps) # o mesmo que a função de cima