# Matriz de datos
load(url("http://pagines.uab.cat/plopez/sites/pagines.uab.cat.plopez/files/IDH2014.rda"))

# install.packages("Hmisc")
library(Hmisc)
library(Deducer)

# Asignación de nombres a las variables
Y=IDH2014$Lifeexpectancy
X1=IDH2014$GNIpercapita2011
X2=IDH2014$Schooling

########## Análisis de regresión 1: regresión simple con las variables originales

# Diagrama de dispersión
plot(Y ~ X1)
ggplot() + geom_point(aes(x=X1, y=Y)) + geom_smooth(aes(x=X1, y=Y),method = 'lm') 
ggplot() + geom_point(aes(x=X1, y=Y)) + geom_smooth(aes(x=X1, y=Y),method = 'lm') +
  geom_text(aes(x=X1, y=Y,label=IDH2014$HDIrank),parse = FALSE)
ggplot() + geom_point(aes(x=X1, y=Y)) + geom_smooth(aes(x=X1, y=Y),method = 'lm') +
  geom_text(aes(x=X1, y=Y,label=IDH2014$Country),parse = FALSE,
            size=4,angle=35,hjust=0,vjust=1,alpha=0.5)

# Correlación (diferentes alternativas)
cor.matrix(variables=d(Y,X1), test=cor.test, method='pearson', aternative="two.sided") # Con Deducer
cor(Y,X1,use="complete.obs")
cor.test(Y,X1)
rcorr(Y,X1) # Con Hmisc

# Análisis de regresión simple
AR=lm(Y ~ X1)
summary(AR)
anova(AR,type='III')
plot(AR) # Gráficos para diagnóstico

# Para guardar como variables: residuos, residuos estandarizados, resiudos studentizados, 
#                              valores predichos y distancia de Cook 
tmp=cooks.distance(AR)
IDH2014[names(tmp),"cooks"]=tmp
tmp=residuals(AR)
IDH2014[names(tmp),"Residuals"]=tmp
tmp=rstandard(AR)
IDH2014[names(tmp),"resid.standardized"]=tmp
tmp=rstudent(AR)
IDH2014[names(tmp),"resid.studentized"]=tmp
IDH2014[["predicted.link"]] = predict(AR, newdata=IDH2014)
# Para editar IDH2014
edit(IDH2014)
# Casos extremos con residuos estandarizados por encima de 3 unidades de desviación
IDH2014.sub<-subset(IDH2014,resid.standardized < -3 | resid.standardized >3)
edit(IDH2014.sub)

########## Análisis de regresión 2: regresión simple con la transformación logarítmica

LX1=log(X1)

ggplot() + geom_point(aes(x=LX1, y=Y)) + geom_smooth(aes(x=LX1, y=Y),method = 'lm') 
ggplot() + geom_point(aes(x=LX1, y=Y)) + geom_smooth(aes(x=LX1, y=Y),method = 'lm') +
  geom_text(aes(x=LX1, y=Y,label=IDH2014$HDIrank),parse = FALSE)
ggplot() + geom_point(aes(x=LX1, y=Y)) + geom_smooth(aes(x=LX1, y=Y),method = 'lm') +
  geom_text(aes(x=LX1, y=Y,label=IDH2014$Country),parse = FALSE,
            size=4,angle=35,hjust=0,vjust=1,alpha=0.5)

cor.matrix(variables=d(Y,LX1), test=cor.test, method='pearson', aternative="two.sided") # Con Deducer
AR=lm(Y ~ LX1)
summary(AR)
anova(AR,type='III')
plot(AR)

tmp=cooks.distance(AR)
IDH2014[names(tmp),"cooks"]=tmp
tmp=residuals(AR)
IDH2014[names(tmp),"Residuals"]=tmp
tmp=rstandard(AR)
IDH2014[names(tmp),"resid.standardized"]=tmp
tmp=rstudent(AR)
IDH2014[names(tmp),"resid.studentized"]=tmp
IDH2014[["predicted.link"]] = predict(AR, newdata=IDH2014)
edit(IDH2014) # Para editar IDH2014
IDH2014.sub<-subset(IDH2014,resid.standardized > 3) # Casos extremos 
edit(IDH2014.sub)

########## Análisis de regresión 3: regresión múltiple

# Diagrama de dispersión en 3D
# install.packages("scatterplot3d")
library(scatterplot3d)
G3D=scatterplot3d(Y ~ LX1+X2,pch=16)

# Ajuste del plano de regresión
AR=lm(Y ~ LX1+X2)
G3D$plane3d(AR)

# Gráfico con Rotación
# install.packages("rgl")
library(rgl)
plot3d(Y,LX1,X2)
# install.packages("Rcmdr")
library(Rcmdr)
library(nlme, pos = 38)
library(mgcv, pos = 38)
scatter3d(Y ~ LX1 + X2, surface=FALSE, residuals=TRUE, bg="white", axis.scales=TRUE, grid=TRUE, ellipsoid=FALSE)
scatter3d(Y ~ LX1+X2)

# Análisis de regresión múltiple
AR=lm(Y ~ LX1+X2)
summary(AR)
anova(AR,type='III')
plot(AR)

# Regresión con variables estandarizadas
YS=scale(Y)
LX1S=scale(LX1)
X2S=scale(X2)
ARS=lm(YS ~ LX1S+X2S)
summary(ARS)

# Cálculo de la Importancia Relativa de cada predictor
# install.packages("relaimpo")
library(relaimpo)
calc.relimp(ARS)

# Diagnósticos: colinealidad (FIV: factor de inflación de la varianza) y residuos
library(car)
vif(ARS)
outlierTest(ARS)
summary(influence.measures(ARS))
