
library(MASS)
library(lattice)

round(dbinom(0:16,16,0.3),3)
plot(dbinom(0:16,16,0.3),type="h", xlab="Numero de mujeres", ylab="Probabilidad")
pbinom(2,16,0.3,lower.tail=T)
rbinom(100,16,0.3)

nn <- 20
muestra.1 <- rnorm(nn)
muestra.2 <- rnorm(nn)
cor.test(muestra.1,muestra.2)

plot(muestra.1,muestra.2)
abline(h=0, lty=2)
abline(v=0, lty=2)
title(paste("r =", round(cor(muestra.1, muestra.2), 3)))

cor.2 <- function(nn) {
   muestra.1 <- rnorm(nn)
   muestra.2 <- rnorm(nn)
   cor(muestra.1,muestra.2)
   }
resultados <- NULL
tam <- 1000
for (i in 1:tam)
   {
   resultados[i] <- cor.2(20)
   }
hist(resultados, xlab="r", main = "Coeficiente de Correlacion")

(paste("5% tiene r menor que ", round(sort(resultados)[tam/20], 3)))

(paste("5% tiene r mayor que ", round(sort(resultados)[tam - tam/20], 3)))

(paste("eta vale ", round(sort(abs(resultados))[tam-tam/10], 3)))


rango <- seq(-3,3,by =0.1)
muestra.1 <- dnorm(rango, 0, 1)
muestra.2 <- dnorm(rango, 0, 0.7)
mesh <-  outer(muestra.1, muestra.2, function(x,y) x*y)
persp(rango, rango, mesh, phi=20, theta=30, expand =0.8, xlab="N(0,1)",
    ylab="N(0,0.7)", zlab="FD", main="Distribucion Normal Bivariada")

ls()

rm(muestra.1,muestra.2,resultados,rango,rand.2,rand.1,nn,mesh,tam,i,cor.2)
ls()

data(iris)
help(iris)
str(iris)
summary(iris)

histogram(~ Petal.Length,data=iris)
histogram(~ Petal.Length|Species,data=iris)
densityplot(~ Petal.Length|Species,data=iris)
bwplot(~ Petal.Length|Species,data=iris)

attach(iris)
plot(Petal.Length,Petal.Width)
plot(Petal.Width ~ Petal.Length)

plot(Petal.Length,Petal.Width,pch=20, cex=1.2, xlab='Longitud del petalo (cm.)',
   ylab='Ancho del petalo (cm.)', main='Datos de Anderson sobre plantas Iris',
   col=c('slateblue', 'firebrick', 'darkolivegreen')[as.numeric(Species)])

abline(v=mean(Petal.Length), lty=2, col='red')
abline(h=mean(Petal.Width), lty=2, col='red')
abline(v=median(Petal.Length), lty=2, col='blue')
abline(h=median(Petal.Width), lty=2, col='blue')
grid()

points(mean(Petal.Length),mean(Petal.Width), cex=2, pch=23,
   col='black',bg='red')
points(median(Petal.Length),median(Petal.Width), cex=2, pch=23,
   col='black',bg='blue')
title(sub='Centroides: media (rojo) y mediana (azul)')

text(1,2.4, 'Tres especies de Iris', pos=4, col='navyblue')
legend(1,2.4,levels(Species), pch=20, bty='n',
   col=c('slateblue', 'firebrick', 'darkolivegreen'))

abline(lm(Petal.Width ~ Petal.Length), lty = 'longdash', col='red')
abline(lqs(Petal.Width ~ Petal.Length), lty = 2, col='blue')

boxplot(Petal.Length ~ Species, horizontal = T, col='lightblue',boxwex=.5,
   xlab='Longitud del petalo (cm)', ylab='Especie', main='Box plot agrupado')

coplot(Petal.Width ~ Petal.Length | Species, col=as.numeric(Species),
   pch=as.numeric(Species))

pairs(iris[,1:4], col=as.numeric(Species),main='Graficas a pares')

stars(iris[,1:4], key.loc=c(2,35), mar=c(2,2,10,2),
   main= 'Star plot de plantas individuales', frame=T)


plot(Petal.Length, Petal.Width)
(p <- identify(Petal.Length, Petal.Width))
iris[p,]


xyplot(Petal.Width ~ Petal.Length, data=iris, groups=Species, auto.key=T)
xyplot(Petal.Width ~ Petal.Length|Species, data=iris, groups=Species, auto.key=T)

pl1 <- cloud(Sepal.Length ~ Petal.Length * Petal.Width, group=Species,
   data=iris, pch=20, main='Datos de Anderson sobre plantas Iris, Todas las Especies',
   screen=list(z=30, x=-60))

data(volcano)
pl2 <- wireframe(volcano, shade=TRUE, aspect=c(61/87, 0.4),
   ligth.source = c(10,0,10), zoom=1.1, box=F, scales=list(draw=F),
   xlab='',ylab='', zlab='', main='Volcan Maunga Whau, Auckland')

pl3 <- levelplot(volcano, col.regions=gray(0:16/16),
   main='Grafica de Nivel, Volcan Maunga Whau, Auckland')

pl4 <- contourplot(volcano, at=seq(floor(min(volcano)/10)*10,
   ceiling(max(volcano)/10)*10, by=10), 
   main='Curvas de Nivel, Volcan Maunga Whau, Auckland',
   sub='Intervalo entre curvas 10 m.', region=T,
   col.regions=terrain.colors(100))

print(pl1, split=c(1,1,2,2), more=T)
print(pl2, split=c(2,1,2,2), more=T)
print(pl3, split=c(1,2,2,2), more=T)
print(pl4, split=c(2,2,2,2), more=F)
rm(pl1, pl2, pl3, pl4)


model <- lm(formula = Petal.Length ~ Species, data = iris)
summary(model)

anova(model)

par(mfrow=c(1,1))
plot(Sepal.Length ~ Petal.Length, data = iris, col=as.numeric(Species),pch=20)

model <- lm(Sepal.Length ~ Petal.Length, data=iris)
summary(model)

plot(residuals(model) ~ fitted(model),col=as.numeric(iris$Species),pch=20)
abline(h=0)

model <- lm(Sepal.Length ~ Petal.Length * Species, data=iris)
summary(model)

plot(residuals(model) ~ fitted(model),col=as.numeric(iris$Species),pch=20)
abline(h=0)

model.2 <- lm(Sepal.Length ~ Petal.Length / Species, data=iris)
summary(model.2)
anova(model,model.2)

pca <- prcomp(iris[,1:4],scale=T)
summary(pca)
screeplot(pca)

pca$rotation
biplot(pca)


petal.v <- iris$Petal.Length[iris$Species == "versicolor"]
par(mfrow=c(1,2))
hist(petal.v, breaks=12, xlab="Longitud de los petalos", main="Iris versicolor")
qqnorm(petal.v)
qqline(petal.v)
shapiro.test(petal.v)

par(mfrow=c(1,1))
mu <- mean(petal.v)
sigma <- sd(petal.v)
puntos <- seq(3,5.5,length=100)
puntos.normal <- dnorm(puntos,mean=mu,sd=sigma)
hist(petal.v, breaks=12, xlab="Longitud de los petalos", main="Iris versicolor",probability=T)
lines(puntos, puntos.normal)

petal.v.2 <- sort(rnorm(length(petal.v),mu,sigma))
summary(petal.v.2)


sepals.2 <- predict(model, data.frame(Petal.Length = petal.v.2,Species = "versicolor"),
   interval = 'prediction',level=0.9)
plot(petal.v.2,sepals.2[,'fit'], type='l',col='blue',lwd=1.5,
   xlab='Longitud de petalos simulada', ylab='Prediccion de longitud de sepalos',
   main='Iris versicolor')
lines(petal.v.2,sepals.2[,'lwr'], col='red', lty=2)
lines(petal.v.2,sepals.2[,'upr'], col='red', lty=2)
points(petal.v, iris$Sepal.Length[iris$Species == 'versicolor'],pch=20)
rug(petal.v.2)
rug(sepals.2[,1],side=2)


ls()
rm(iris, model, model.2, mu, p, pca, petal.v, petal.v.2,
    puntos, puntos.normal, sigma)
ls()
