# dane BostonHousing
library(mlbench)
data(BostonHousing)
      # funkcje obliczająca błędy MAE, MSE, RMSE
mae <- function(y.true, y.pred) { sum(abs(y.pred-y.true))/length(y.true) }
mse <- function(y.true, y.pred) { sum((y.pred-y.true)^2)/length(y.true) }
rmse <- function(y.true, y.pred) { sqrt(mse(y.true, y.pred)) }

      # podzielmy dane na trenujące i testowe losowo w stosunku ok. 2:1
rbh <- runif(nrow(BostonHousing))
bh.train <- BostonHousing[rbh>=0.33,]
bh.test <- BostonHousing[rbh<0.33,]

      # model liniowy
bh.lm <- lm(medv ~ ., bh.train)
summary(bh.lm)
mae(bh.test$medv, predict(bh.lm, bh.test))
rmse(bh.test$medv, predict(bh.lm, bh.test))

      # model pamięciowy kNN
      # użyjemy implementacji z pakietu kknn

if (! "kknn" %in% row.names(installed.packages()))
  install.packages("kknn")
library(kknn)
bh.knn1 <- kknn(medv ~ ., bh.train, bh.test, k=1)
bh.knn3 <- kknn(medv ~ ., bh.train, bh.test, k=3)
bh.knn5 <- kknn(medv ~ ., bh.train, bh.test, k=5)
bh.knn7 <- kknn(medv ~ ., bh.train, bh.test, k=7)
mae(bh.test$medv, bh.knn1$fitted.values)
rmse(bh.test$medv, bh.knn1$fitted.values)
mae(bh.test$medv, bh.knn3$fitted.values)
rmse(bh.test$medv, bh.knn3$fitted.values)
mae(bh.test$medv, bh.knn5$fitted.values)
rmse(bh.test$medv, bh.knn5$fitted.values)
mae(bh.test$medv, bh.knn7$fitted.values)
rmse(bh.test$medv, bh.knn7$fitted.values)

      # to samo z uprzednim przeskalowaniem atrybutów numerycznych
      # średnie na zbiorze trenującym
bh.m <- apply(bh.train[,c(-4,-14)], 2, mean)
      # odchylenia standardowe na zbiorze trenującym
bh.sd <- apply(bh.train[,c(-4,-14)], 2, sd)
      # skalowanie atrybutów ze zbiorów trenującego i testowego
      # atrybuty nieskalowane kopiujemy bez zmian
bh.train.s <- as.data.frame(scale(bh.train[,c(-4,-14)], bh.m, bh.sd))
bh.train.s$chas <- bh.train$chas
bh.train.s$medv <- bh.train$medv
bh.test.s <- as.data.frame(scale(bh.test[,c(-4,-14)], bh.m, bh.sd))
bh.test.s$chas <- bh.test$chas
bh.test.s$medv <- bh.test$medv

bh.knn1s <- kknn(medv ~ ., bh.train.s, bh.test.s, k=1)
bh.knn3s <- kknn(medv ~ ., bh.train.s, bh.test.s, k=3)
bh.knn5s <- kknn(medv ~ ., bh.train.s, bh.test.s, k=5)
bh.knn7s <- kknn(medv ~ ., bh.train.s, bh.test.s, k=7)
mae(bh.test.s$medv, bh.knn1s$fitted.values)
rmse(bh.test.s$medv, bh.knn1s$fitted.values)
mae(bh.test.s$medv, bh.knn3s$fitted.values)
rmse(bh.test.s$medv, bh.knn3s$fitted.values)
mae(bh.test.s$medv, bh.knn5s$fitted.values)
rmse(bh.test.s$medv, bh.knn5s$fitted.values)
mae(bh.test.s$medv, bh.knn7s$fitted.values)
rmse(bh.test.s$medv, bh.knn7s$fitted.values)
      # czy wygląda na to, że kknn ma wbudowane skalowanie?

      # model SVR (support vector regression)
library(e1071)
bh.svm <- svm(medv ~ ., bh.train)
mae(bh.test$medv, predict(bh.svm, bh.test))
rmse(bh.test$medv, predict(bh.svm, bh.test))

      # strojenie parametru C
par.cost=2^(-2:6)
bh.tuneC <- tune(svm, medv ~ ., data=bh.train,
                 ranges=list(cost=par.cost),
                 tunecontrol=tune.control(nrepeat=10))
      # uzyskana najlepsza wartość
bh.tuneC$best.parameters
      # błędy walidacji krzyżowej dla poszczególnych wartości
bh.tuneC$performances[order(bh.tuneC$performances$error),]
      # błąd modelu dla najlepszej wartości C na zbiorze testowym
mae(bh.test$medv, predict(bh.tuneC$best.model, bh.test))
rmse(bh.test$medv, predict(bh.tuneC$best.model, bh.test))