### Exercícios 8. Reamostragem e Simulação ### #Exercício 1. # 1.1) Ler os dados fornecidos no exercício: eutad <- read.table("palmadulto.txt", header = T, sep = "\t") head(eutad) tail(eutad) length(eutad$tag) # 1.2) Distância entre cada indivíduo: dist=matrix(NA, ncol=102, nrow=102) # 1.3) Distância observada de cada indivíduo: for(i in 1:101) { for(j in (i+1):102) { difx2=(eutad$gx[i]-eutad$gx[j])^2 dify2=(eutad$gy[i]-eutad$gy[j])^2 dist[i,j]<-sqrt(difx2 + dify2) dist[j,i]<-sqrt(difx2 + dify2) } } # 1.4) Distância média do vizinho próximo(MNN): (nn<-apply(dist, 1, min, na.rm=TRUE)) (mnn<-mean(nn)) # 1.5) Crie um vetor resultado, com 1000 NA's, para guardar os valores de cada simulação. resultado <- rep(NA,1000) # 1.6) Guarde o valor observado na primeira posição de resultado; resultado[1] <- mnn head(resultado) # 1.7) Crie um ciclo com contador (k) que vai de 2 a 1000; # 1.8) Dentro do ciclo: # 1.8.1) Crie o objetos xsim, um vetor com valores amostrados aleatoriamente de uma distribuição uniforme de 0 a 320 (o tamanho x da parcela), arredondando o valor para uma casa decimal2). Lembre-se de sortear o mesmo quantidade de valores que a população de palmito observado; # 1.8.2) Faça o mesmo que no passo anterior e guarde no objeto ysim; # 1.8.3) Como no tópico anterior, crie uma matriz para guardar as distâncias entre cada valor xy simulados; # 1.8.3) Crie os ciclos para o calculo das distâncias como no tópico anterior; # 1.8.4) Guarde o valor da distância média do vizinho mais próximo dos dados simulados na posição k do vetor resultado; # 1.8.5) Feche o ciclo. for (k in 2:1000) { xsim <- round(runif(102, min = 0, max = 320),1) ysim <- round(runif(102, min = 0, max = 320),1) sim.dist <- matrix(NA, ncol=102, nrow=102) for (i in 1:101) { for(j in (i+1):102) { difx2 <- (xsim[i]-xsim[j])^2 dify2 <- (ysim[i]-ysim[j])^2 dist[i,j] <- sqrt(difx2 + dify2) dist[j,i] <- sqrt(difx2 + dify2) } } resultado[k] <- mean(apply(dist, 1, min, na.rm=TRUE)) } # 1.9. Faça um histograma dos valores simulados e coloque uma linha vermelha # vertical na posição do valor observado; hist(resultado) abline(v = resultado[1], col = "red") # 1.10. Calcule a probabilidade de uma distribuição espacial aleatória gerar # valores iguais ou mais extremos do que o valor MNN observado. prob <- sum(resultado <= resultado[1])/(length(resultado)) #A probabilidade de gerar valores aleatórios iguais ou menores do que o MNN #observado(14.76015) neste caso foi de 0.02. ######################################################################### ### Exercício 2 ### # Simulando o teste de uma regressão linear # Conjunto de Dados: Massa do Corpo e do Cérebro de Vertebrados # Observações dos dados: Massas médias do corpo (kg) e do cérebro (g) de # 28 vertebrados terrestres (Rosseeuw & Leroy, 1987 1)). # Note que há valores faltantes. # Lendo os dados: animais <- read.table("animais.txt", header = T, sep = ";", dec = ",") str(animais) head(animais) # Retirando as linhas com valores faltantes: animais2 <- animais[(animais$brain!=""),] animais2 <- na.omit(animais) str(animais2) head(animais2) # 2.1. Calcule a inclinação da relação log(brain) ~ log(body) plot(brain ~ body, data = animais2) plot(log(brain) ~ log(body), data = animais2) ml <- lm(log(brain) ~ log(body), data = animais2) abline(ml, col = "red") summary(ml)$coefficients[2,1] # 0.7522607 med.x <- mean(log(animais2$body)) med.y <- mean(log(animais2$brain)) desv.x <- log(animais2$body) - med.x desv.y <- log(animais2$brain) - med.y desv.xy <- (desv.x)*(desv.y) soma.desv.xy <- sum(desv.xy) desv.x2 <- (desv.x)^2 sum(desv.x2) b <- soma.desv.xy / sum(desv.x2) b # b é 0.7522607, o mesmo apresentado pelo summary. # 2.2. Crie um vetor para guardar o resultado de simulações resultado <- rep(NA,1000) # 2.3. Guarde o valor observado no objeto criado em 2.2 resultado[1] <- b # 2.4. Abra um ciclo de 2 a 1000 # 2.5. Desordene o vetor brain e guarde no objeto sim_brain; # 2.6. Calcule a inclinação entre o log(sim_brain) ~ log(body) e guarde no vetor resultado; # 2.7. Feche o ciclo; for (i in 2:1000) { sim_brain <- sample(animais2$brain) med.x <- mean(log(animais2$body)) med.y <- mean(log(sim_brain)) desv.x <- log(animais2$body) - med.x desv.y <- log(sim_brain) - med.y desv.xy <- (desv.x)*(desv.y) soma.desv.xy <- sum(desv.xy) desv.x2 <- (desv.x)^2 sum(desv.x2) resultado[i] <- soma.desv.xy / sum(desv.x2) } # 2.8. Faça o histograma dos valores simulados e compare com o valor # observado da inclinação da relação hist(resultado) abline(v = resultado[1], col = "red") # 2.9. Calcule a probabilidade da inclinação observada ter sido gerada # por variáveis que são independentes prob <- sum(resultado >= resultado[1])/(length(resultado)) ## A probabilidade das variáveis serem independentes é de 0.001.