Deliver4.Rmd

---
title: "Forecasting modeling of numeric target"
author: "Katerina Dimitrova, Jose Romero, Sergi Munoz"
date: "May 11, 2018"
output: pdf_document
editor_options: 
  chunk_output_type: console
chunk_output_type: console
---

```{r,include=FALSE}
#knitr::opts_chunk$set(echo = TRUE)

```

```{r,include=FALSE}
setwd("C://Users/JOSE CAMILO/Desktop/ADEI-masterrrr") #Change 
#setwd("D:/ADEI/ADEI.git/trunk") #Change

```

```{r,include=FALSE}
rm(list=ls())
requiredPackages <- c("effects","FactoMineR","car", "factoextra","RColorBrewer","ggplot2","missMDA","mvoutlier","dplyr","ggmap","ggthemes","knitr","MVA")
missingPackages <- requiredPackages[!(requiredPackages %in% installed.packages()[,"Package"])]

if(length(missingPackages)) install.packages(missingPackages)
lapply(requiredPackages, require, character.only = TRUE)

# Useful function
calcQ <- function(x) {
  s.x <- summary(x)
  iqr<-s.x[5]-s.x[2]
  list(souti=s.x[2]-3*iqr, mouti=s.x[2]-1.5*iqr, min=s.x[1], q1=s.x[2], q2=s.x[3],        q3=s.x[5], max=s.x[6], mouts=s.x[5]+1.5*iqr, souts=s.x[5]+3*iqr ) }

```

#MODELING
A statistical model is an expression that attempts to explain patterns in the observed values of a response variable by relating the response variable to a set of predictor variables and parameters.
Consider the following familiar statistical model:

y = mx+c

This simple statistical model relates a response variable (y) to a single predictor variable (x) as a straight line according to the values of two constant parameters:
m - the degree to which y changes per unit of change in x (gradient of line) c - the value of y when x = 0 (yintercept).


In complex systems, variables are typically the result of many influential and interacting factors and therefore simple models usually fail to fully explain a response variable.
Consequently, the statistical model also has an error component that represents the portion of the
response variable that the model fails to explain. Hence, statistical models are of the form:

response variable = model + error




#Linear models in R
Les variables que utilitzarem per a començar a fer el model són les variables contínues més significatives per a nosaltres.
Per començar hem de veure quines són les variables més relacionades amb la nostra variable resposta Total_amount.
```{r} 
load("Taxi5000_raw_DataDefinitivev1.RData")

# Quines variables són les més adientes per crear el nostre model? La resposta és 
#sencilla, les variables més realcionades amb la nostra variable resposta són les 
#seleccionades com a possibles candidates per a formar part del nostre model
condes(df[,vars_con],12)
# Les variables més relacionades amb Total_amount (variable resposta) són:
# Trip_distance,Fare_amount,trip_distance_km,Trip_distance,travel_time,trip_length

# Feature Selection 
cor(df[,c("Total_amount",vars_con)],method="spearman") #coeficient no parametric

# Quantes variables agafem? Hem d'intentar arribar al millor model que pugem representar.
```
Les variables més correlacionades amb la nostra variable resposta Total_amount són les que tenen una correlació propera a 1. En el nostre cas les variables més correlacionades són: Trip_distance, Fare_amount, trip_distance_km, Trip_distance, travel_time, trip_length, però no totes sòn útils.

# Multiple Linear Regression issues
*Target numeric Total_amount*

Creem un nou model m10, amb les variables que més es relacionen amb la nostra variable resposta. 
```{r}
# Creem un nou model amb les variables Fare_amount, 
# trip_distance_km, travel_time, Tip_amount, espeed.
m10<-lm(Total_amount~Fare_amount+trip_distance_km
        +travel_time+Tip_amount+espeed,data=df)
summary(m10)
#Explicació:
# 1. L'estimació de B0 és 1.1175243 amb un error estandar 0.0407993.
# El contrast de H0: B0 = 0 davant H1: B0 != 0 llença un valor estadistic
# de 27.391 i un p-valor inferior a 2e-16.

# 2. L'estimació de B1 és 0.9864011 amb un error estandar 0.0053612.
# El contrast de H0: B1 = 0 davant H1: B1 != 0 llença un valor estadistic
# de 183.989 i un p-valor inferior a 2e-16.

# 3. L'estimació de B2 és 0.0585972 amb un error estandar 0.0102656. 
# El contrast de H0: B2 = 0 davant H1: B2 != 0 llença un valor estadistic
# de 5.708 i un p-valor de 1.21e-08.

# 4. L'estimació de B3 és -0.0041999 amb un error estandar 0.0019113.
# El contrast de H0: B3 = 0 davant H1: B3 != 0 llença un valor estadistic
# de -2.197 i un p-valor de 0.028.

# 5. L'estimació de B4 és 1.0491685 amb un error estandar 0.0066837. 
#. El contrast de H0: B4 = 0 davant H1: B4 != 0 llença un valor estadistic
# de 156.974 i un p-valor inferior a 2e-16.

# 6. L'estimació de B5 és 0.0003075 amb un error estandar 0.0014119.
# El contrast de H0: B5 = 0 davant H1: B5 != 0 llença un valor estadistic
# de 0.218 i un p-valor de 0.828.

# La recta ajustada apareix, per tant, especificada a través dels cinc 
#coeficients:
# Total_amount = 1.1175243+0.9864011*Fare_amount+0.0585972*trip_distance_km
# -0.0041999*travel_time+1.0491685*Tip_amount+0.0003075*espeed

# 7. L'error estandar de l'ajust te un valor de 0.7687.

# 8. El coeficient R^2 te el valor 0.9907, que indica que el 99.07% 
# de tota la variabilitat que te el fenomen relatiu al preu pagat.

# Les variables més significatives són les que tenen un p-valor menor a 0.05,
# ja que no podem rebutjar la hipotesi nul·la. En aquest cas les més
# significatives són: Fare_amount, trip_distance_km, travel_time,Tip_amount.
# Observem que podem descartar la variable espeed, almenys de moment.
```

```{r}
# Test net effects
# All net effects are significant?
Anova(m10)
# No tots els efectes son significatius, ja que obtenim variables que
# tenen un p-valor major al nostre llindar 0.05
```



```{r}

# Apliquem la comanda step al model anterior per determinar
# quina combinació es la més adient i poder treure variables
# que no són necessaries
m11<-step(m10)
m11<-step(m10,k=log(nrow(df)))  # BIC
vif(m11) 
# Observem que l'unica variable que compleix la nostra cota
# es Tip_amount, per tant podem agafar aquesta variable per 
# al nostre model, i hem de continuar treballant per tdeterminar
# quina es la millor entre Fare_amount trip_distance_km
```

Apliquem la comanda AIC per determinar quin és el millor model fins al moment
``` {r}
AIC(m10,m11) 
# Comparem tos el models que tenim fins el momnet per 
# determinar quin és el millor.
# Segons el dit a classe el millor model es el que te el Valor
# de AIC més petit, en el nostre cas és m10, però aquest mostra colinealitat.
```

Creem el model m20 i m21 amb totes les variables contínues.
```{r}
vars_con_model <- vars_con[c(1:5,7:8,10:11,13:17)]
m20<-lm(Total_amount~.,data= df[,c("Total_amount",vars_con_model)])
summary(m20)
vif(m20)

m21<-step(m20,k=log(nrow(df)))  # BIC
summary(m21)
vif(m21)

# Test Fisher: m20 - m21 H0 Equivalents
anova(m21,m20)

# All net effects are significant?
Anova(m21)

vif(m21)  # Check association between explanatory vars
round(cor(df[,vars_con],use="pairwise.complete.obs"),dig=2)

```
En el model m20 observem que hi ha variables que no són necessàries, ja que hi ha moltes correlacions implícites entre les mateixes variables del model. El mateix passa amb el model m21. A l'hora de comparar l'equivalència entre els dos model observem que el p-valor és de 0.2552 per tant podem dir que no són equivalents.


 
```{r}
# En aquest moment ja comencem a tenir una ida de quines variables
# continues necessitem per al nostre model
AIC(m20,m21,m10,m11)
```
Observem que de tots els models que hem vist fins al moment el millor és el m21.

Creem un nou model m23.
```{r}
# Creem una nou model amb les variables trip_distance_km, 
# Fare_amount, Extra, MTA_tax, Tip_amount, Tolls_amount
m23<-lm(Total_amount~trip_distance_km+Fare_amount+Extra+MTA_tax+Tip_amount+Tolls_amount,data=df)
summary(m23)

Anova(m23)
# Observem que totes les variables son significatives

vif(m23)
# Encara observem que hi ha variables que estan relacionades entre elles,
# ja que superen la cota de 3 en el valor resultant 

AIC(m23)
```
Observem que aquest nou model m23 té un AIC més petit que el model m22 però no volem un model amb correlacions entre les pròpies variables.

Creem el model m24.
```{r}
# Creem un nou model amb les variables que creiem més significatives
# i que no estan ralacionades entre elles
m24<-lm(Total_amount~trip_distance_km
        +Extra+Tip_amount+Tolls_amount,data=df)
summary(m24)

Anova(m24) # All Net effects are significant
# condicionat a tenir 3 variables en el model la quarta 
# es significativa 

vif(m24)
# Observem que no hi ha colinealitat entre les variables
# utilitzades per aquet model

AIC(m23,m24)
# Agafem el model amb un AIC més petit
```

# Diagnostics
L'homoscedasticitat és una propietat fonamental del model de regressió lineal general i està dins dels seus supòsits clàssics bàsics. Es diu que hi ha homoscedasticitat quan la variància dels errors de la regressió són els mateixos per a cada observació i (d'1 a n observacions).

Creem el model m25.
```{r}
# Heterocedastricitat (variança no constant)
# Creem un nou model amb les variables que hens han semblat 
# les millors per a fer el model segons l'apartat anterior 
# i apliquem el logaridme a la variable resposta.
m25 <-lm(log(Total_amount)~trip_distance_km
        +Extra+Tip_amount+Tolls_amount,data=df)
summary(m25)
# Explicació:
# 1. L'estimació de B0 és 1.839703 amb un error estandar 0.006091.

# 2. L'estimació de B1 és 0.123742 amb un error estandar 0.001059. 

# 3. L'estimació de B2 és 0.110499 amb un error estandar 0.008970. 

# 4. L'estimació de B3 és 0.065952 amb un error estandar 0.001967. 

# 5. L'estimació de B4 és 0.015307 amb un error estandar 0.004860. 

# La recta ajustada apareix, per tant, especificada a través 
# dels quatre coeficients:
# log(Total_amount) = 1.839703+0.123742*trip_distance_km+
# 0.110499*Extra+0.065952*Tip_amount+0.015307*Tolls_amount

# 6. L'error estandar de l'ajust te un valor de 0.2264.


Anova(m25)
# Observem que totes les variables són significant ja que tenen el p-valor inferir a 0.05

# apliquem l'exponencial per poder treure el logaritme aplicat anteriorment
# i així poder fer el predict de forma correcta.
vif(m25)
# No hi ha colinealitat entre les variables

par(mfrow=c(2,2))
plot(m25)
par(mfrow=c(1,1))
```

Creem el model m26.
```{r}
# No cal transformar, pero cal termes de ordre superior
m26 <-lm(Total_amount~poly(trip_distance_km,2)
        +Extra+poly(Tip_amount,2)+Tolls_amount,data=df)
summary(m26)

Anova(m26)

anova(m26,m23)

marginalModelPlot(m26) #EXTRA Y Tolls_amount com a variables factors

Boxplot(rstudent(m26))
llout<-which(abs(rstudent(m26))>5); length(llout)

```

# Using factors as explanatory variables
Com hem vist un comportament extrany en les variables Extra i Tolls_amount probarem a utilitzarel com a factors.

Extra
```{r}
df$f.extra<-0
df$f.extra[df$Extra>0]<-1
df$f.extra[df$Extra>0.5]<-2
m27<-lm(Total_amount~poly(trip_distance_km,2)
        +f.extra+poly(Tip_amount,2)+Tolls_amount,data=df)
summary(m27)
Anova(m27)
```

Tolls_amount com a factor no ens dóna una millora, ja que en el moment de fer el BIC amb el model m27 observem que el model m37 és més gran.
```{r}
m37<-lm(Total_amount~poly(trip_distance_km,2)
        +f.extra+f.toll+poly(Tip_amount,2),data=df)
summary(m37)
par(mfrow=c(2,2))
plot(m37)
par(mfrow=c(1,1))
BIC(m27,m37)

```

Creem el model  m30
```{r}
m30<-lm(log(Total_amount)~trip_distance_km+Tip_amount+(f.extra+f.passenger+Payment_type+Tolls_amount)*(f.espeed)+f.ttime,data=df)
m30<-step(m30,k=log(nrow(df))) 
vif(m30) 
summary(m30) 
Anova(m30) 
BIC(m30)
```

# RESIDUAL ANALYSIS
Tenim el model m30 com a millor model. 

````{r}
# Residual Analysis
par(mfrow=c(2,2))
plot(m30)
par(mfrow=c(1,1))
```
RESIDUAL vs FITTED: En el gràfic no podem identificar una tendència de la variància dels residus a reduir-se ni augmentar a mesura que l'adherència augmenta. Observem que la línia vermella és pràcticament una línia recta. També observem uns punts que tenen algun valor estrany.


NORMAL Q-Q: Aquest grafic ens mostra les diferencies entre la distribució de probabilitat de la nostra població de la qual hem tret la mostra i la nostra distribució usada per la comparació. En el nostre grafic observem que en les puntes hi ha una separació marcada respecte a la distribució de la població, deguda a outliers.


SCALE-LOCATION: Aquest gràfic ens mostra si els residus es distribueixen per igual al llarg del rang dels predictors. Així és com podem verificar la suposició d'igual variança (homocedasticitat). Observem una línia horitzontal que ens diu que quan més augmenta Fitted values els residus augmenten.

RESIDUALS VS LEVERAGE: Aquest gràfic ens ajuda a trobar casos influents si n'hi ha. No tots els outliers són influents en l'anàlisi de regressió lineal. Encara que algunes dades tenen valors extrems, és possible que no siguin influents per determinar una línia de regressió. Això significa que els resultats no serien molt diferents si els incloem o els exclouen de l'anàlisi. Segueixen la tendència en la majoria dels casos i en realitat no els importa; no són influents. D'altra banda, alguns casos podrien ser molt influents, fins i tot si semblen estar dins d'un rang raonable dels valors. Poden ser casos extrems contra una línia de regressió i poden alterar els resultats si els exclouen de l'anàlisi. Una altra forma de plantejar-la és que no s'acompanyen amb la tendència en la majoria dels casos.


```{r}
marginalModelPlots(m30)

```

# INFLUENT DATA
```{r}
Boxplot(cooks.distance(m30),id.n=3)
abline(h=0.08,col="red",lwd=2)
llcoo<-which(cooks.distance(m30)>0.08);length(llcoo)

df[llcoo,]   
df<-df[-llrem]

```
Report:
En les observacions que tenen una distància de Cook més gran a 0.8 observem que tenen una relació molt peculiar entre el total pagat i la distància del viatge. Observem que hi ha distancies molt petites amb un preu pagat total molt elevat. (fixar-se en Trip_distance i Total_amount)