setwd("/Users/subasishdas1/Dropbox/---- TRB_2017 ----/Raul Tire Deb/Resources")
td <- read.csv("TD Data for TRB Paper 20160610.csv")
head(td)
##   AvgOfNUM_LANES AveTrkADT AveNonTrkADT AvgOfSPD_MAX AvgOfSPD_MIN
## 1           2.22    694.60     14327.62        55.00         4.44
## 2           4.00   9028.06     27510.82        75.00         0.00
## 3           4.00   5166.92     42560.84        75.00         1.00
## 4           4.00    534.93     25287.64        53.37         0.05
## 5           8.89  14360.77    191908.35        60.20         0.00
## 6           4.00   4674.39     35438.10        64.54         0.00
##   Total_TD_ctL Total_TD_wL Total_TD_ct Total_TD_w Runs Cond_P Round FLen
## 1         9.56       24.18          27      33.11    2   1.00     2 5.24
## 2       140.00      284.06        1096     687.41    4   0.23     2 4.96
## 3         5.00        8.47          56      24.79    4   1.00     2 5.66
## 4         1.00        1.50           2       1.82    4   1.00     2 4.33
## 5        66.61      175.65         174     248.71    5   0.20     2 7.80
## 6         4.00        6.61         141      31.24    4   0.93     1 6.62
##                Type_GR                       FUNC_SYS_DESC
## 1         Only Spatial                URBAN MINOR ARTERIAL
## 2         Only Spatial            URBAN PRIN ARTERIAL (IH)
## 3         Only Spatial URBAN PRIN ARTERIAL (OTHER FREEWAY)
## 4         Only Spatial         URBAN PRIN ARTERIAL (OTHER)
## 5              Imputed            URBAN PRIN ARTERIAL (IH)
## 6 Temporal And Spatial                    RURAL INTERSTATE
names(td)
##  [1] "AvgOfNUM_LANES" "AveTrkADT"      "AveNonTrkADT"   "AvgOfSPD_MAX"  
##  [5] "AvgOfSPD_MIN"   "Total_TD_ctL"   "Total_TD_wL"    "Total_TD_ct"   
##  [9] "Total_TD_w"     "Runs"           "Cond_P"         "Round"         
## [13] "FLen"           "Type_GR"        "FUNC_SYS_DESC"
td1 <- td[c(1:5, 8)]
names(td1)
## [1] "AvgOfNUM_LANES" "AveTrkADT"      "AveNonTrkADT"   "AvgOfSPD_MAX"  
## [5] "AvgOfSPD_MIN"   "Total_TD_ct"

Correlation Plot

library(corrplot)
corrp <- cor(td1)
corrplot(corrp, method = "circle")

SVR Model

library(e1071)

### Total_TD_ct
dim(td1)
## [1] 83  6
svm.model <- svm(Total_TD_ct~., data = td1, cost = 100, gamma =1)
svm.pred  <- predict(svm.model, td1[c(1:5)])

svm_err_all <- td1$Total_TD_ct- svm.pred
(svm_rmse1 <- sqrt(mean(svm_err_all^2)))
## [1] 32.60432
pl <- cbind(obs= td1$Total_TD_ct, pred=abs(svm.pred))
pl <- data.frame(pl)

library(ggplot2)
p1 <- ggplot(pl, aes(obs, pred))
p1a <- p1 + geom_point(colour = "#00A08A", size = 1.5)+ theme_bw()+ geom_abline(intercept=0, slope=1, color="#F98400",size = 1.5,linetype=4)+labs(x = "Observed \n (a)", y="Predicted")+labs(title="Total Count of Tire Debris")
p1a

### Total_TD_w
td2 <- td[c(1:5, 9)]
svm.model <- svm(Total_TD_w~., data = td2, cost = 100, gamma =1)
svm.pred  <- predict(svm.model, td2[c(1:5)])

svm_err_all <- td2$Total_TD_w- svm.pred
(svm_rmse1 <- sqrt(mean(svm_err_all^2)))
## [1] 34.91371
pl <- cbind(obs= td2$Total_TD_w, pred=abs(svm.pred))
pl <- data.frame(pl)

library(ggplot2)
p1 <- ggplot(pl, aes(obs, pred))
p2 <- p1 + geom_point(colour = "#00A08A", size = 1.5)+ theme_bw()+ geom_abline(intercept=0, slope=1, color="#F98400",size = 1.5,linetype=4)+labs(x = "Observed \n (b)", y="Predicted")+labs(title="Total Weight of Tire Debris")
p2

### Total_TD_ctL
td3 <- td[c(1:5, 6)]
svm.model <- svm(Total_TD_ctL~., data = td3, cost = 100, gamma =1)
svm.pred  <- predict(svm.model, td3[c(1:5)])

svm_err_all <- td3$Total_TD_ctL- svm.pred
(svm_rmse1 <- sqrt(mean(svm_err_all^2)))
## [1] 10.74034
pl <- cbind(obs= td3$Total_TD_ctL, pred=abs(svm.pred))
pl <- data.frame(pl)

library(ggplot2)
p1 <- ggplot(pl, aes(obs, pred))
p3 <- p1 + geom_point(colour = "#00A08A", size = 1.5)+ theme_bw()+ geom_abline(intercept=0, slope=1, color="#F98400",size = 1.5,linetype=4)+labs(x = "Observed \n (c)", y="Predicted")+labs(title="Total Count of Large Tire Debris")
p3

### Total_TD_wL
td4 <- td[c(1:5, 7)]
svm.model <- svm(Total_TD_wL~., data = td4, cost = 100, gamma =1)
svm.pred  <- predict(svm.model, td4[c(1:5)])

svm_err_all <- td4$Total_TD_wL- svm.pred
(svm_rmse1 <- sqrt(mean(svm_err_all^2)))
## [1] 23.32193
pl <- cbind(obs= td4$Total_TD_wL, pred=abs(svm.pred))
pl <- data.frame(pl)

library(ggplot2)
p1 <- ggplot(pl, aes(obs, pred))
p4 <- p1 + geom_point(colour = "#00A08A", size = 1.5)+ theme_bw()+ geom_abline(intercept=0, slope=1, color="#F98400",size = 1.5,linetype=4)+labs(x = "Observed \n (d)", y="Predicted")+labs(title="Total Weight of Large Tire Debris ")
p4

library(grid)
grid.newpage()
pushViewport(viewport(layout = grid.layout(2, 2)))
vplayout <- function(x, y)
  viewport(layout.pos.row = x, layout.pos.col = y)
print(p1a, vp = vplayout(1, 1))
print(p2, vp = vplayout(1, 2))
print(p3, vp = vplayout(2, 1))
print(p4, vp = vplayout(2, 2))

######## TUNING
tuneResult <- tune(svm, Total_TD_ct~., data = td1,
              ranges = list(epsilon = seq(0,1,0.1), cost =1.77^(2:8))
)
print(tuneResult)
## 
## Parameter tuning of 'svm':
## 
## - sampling method: 10-fold cross validation 
## 
## - best parameters:
##  epsilon   cost
##      0.3 3.1329
## 
## - best performance: 62533.22
# best performance: MSE = 8.371412, RMSE = 2.89 epsilon 1e-04 cost 4
# Draw the tuning graph
plot(tuneResult)

svm.model <- svm(Total_TD_ct~., data = td1, cost = 50, epsilon =0.3, gamma=1)
svm.pred  <- predict(svm.model, td1[c(1:5)])

pl <- cbind(obs= td1$Total_TD_ct, pred=abs(svm.pred))
pl <- data.frame(pl)
pl1 <- sum(pl$pred-pl$obs)/sum(pl$obs)
head(pl1)
## [1] 0.1643002
library(ggplot2)
p1 <- ggplot(pl, aes(obs, pred))
p1a <- p1 + geom_point(colour = "#00A08A", size = 1.5)+ theme_bw()+ geom_abline(intercept=0, slope=1, color="#F98400",size = 1.5,linetype=4)+labs(x = "Observed", y="Predicted")
p1a