#Francisco M. Ulloa 8090851

#Exercício 8, Reamostragem e Simulação

#1

setwd("F:\\Libros R\\1 R Panch\\8")

eutad<-read.table("palmadulto.txt", header=TRUE, sep="")
str(eutad)

plot(eutad$gx, eutad$gy)
dist<-matrix(NA, ncol=102, nrow=102)

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)
            }
          
    }

(nn<-apply(dist, 1, min, na.rm=TRUE))

(mnn<-mean(nn))#Distancia media do vizinho mais proximo

#Simulando

resultado <- rep(NA,1000)

## Guardo o primeiro valor observado

resultado[1] <- mnn
head(resultado)

## E agora o loop: diferencas sob a hipotese nula

for(i in 2:1000){ 
  xsim<-runif(nrow(eutad),0,320)
  ysim<-runif(nrow(eutad),0,320)

  dist <- matrix(NA, ncol = 102, nrow = 102)

	for(k in 1:101){
		for(j in (k+1):102){
			difx2 <- (xsim[k] - xsim[j])^2
			dify2 <- (ysim[k] - ysim[j])^2
			dist[k,j] <- sqrt(difx2 + dify2)
			dist[j,k] <- sqrt(difx2 + dify2)
		}
	}

	nn_random <- apply(dist, 1, min, na.rm = TRUE)
	mnn_random <- mean(nn_random)
	resultado[i] <- mnn_random 
}	


hist(resultado)
abline(v = mnn, col = "red", lty = 2, lwd = 2)

prob <- sum(resultado <= mnn) / length(resultado)
prob 
## p é menor que 0,05%, a distribuição é agregada no espaço

#########################################################

#2, desafio

animais <- read.table("animais.txt", header=T, as.is=T, sep=";", dec="," )
str(animais)
head(animais)
summary(animais)
animais <- animais[-which(is.na(animais$brain)),]
summary(animais)

##inclinação da relação log(brain) ~ log(body)

inclinacao <- coef(lm(log(brain) ~ log(body), data = animais))[2]

## vetor para guardar o resultado de simulações;

resultado <- rep(NA, 1000)
resultado[1] <- inclinacao 

for(i in 2:1000){
	sim_brain <- sample(animais$brain)
	resultado[i] <- coef(lm(log(sim_brain) ~ log(animais$body)))[2]
}

#Fazer histograma e compare com o valor observado da inclinação da relação

hist(resultado)
abline(v = inclinacao, col = "red", lty = 2, lwd = 2)


## Probabilidade da inclinação observada ter sido gerada 
## por variáveis que são independentes

p <- sum(resultado >= inclinacao) / length(resultado)
p


## Muito pequena a probabilidade dessa inclinação ser ao acaso