rm(list=ls())

##eolica <- data.frame(
##    vento = c(5.00, 6.00, 3.40, 2.70, 10.0, 9.70, 9.55, 3.05, 8.15, 6.20, 2.90, 6.35, 4.60, 5.80, 7.40, 3.60, 7.85, 8.80, 7.00, 5.45, 9.10, 10.2, 4.10, 3.95, 2.45),
##    energia = c(1.582, 1.822, 1.057, 0.500, 2.236, 2.386, 2.294, 0.558, 2.166, 1.866, 0.653, 1.930, 1.562, 1.737, 2.088, 1.137, 2.179, 2.112, 1.800, 1.501, 2.303, 2.310, 1.194, 1.144, 0.123)
)
##eolica <- eolica[order(eolica$vento), ]

eolica <- read.csv2("http://www.leg.ufpr.br/~paulojus/CE092/dados/eolica.csv")

with(eolica, plot(energia ~ vento))
with(eolica, plot(energia ~ vento, 
                  ylim=c(0, 2.5), xlim = c(0, 10.5)))

with(eolica, abline(lm(energia~vento)))

##
## Regressão (Grau 1)
with(eolica, plot(energia ~ vento))
preds <- data.frame(vento = seq(2,11, length=101))

(seg1.0 <- lm(energia ~ vento, data = eolica))
preds$pseg1.0 <- predict(seg1.0, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg1.0, lty = 3))
##
## Regressão segmentada
##
## um nó
eolica$xc5 <- with(eolica, ifelse(vento<5, 0, vento-5))
preds$xc5 <- with(preds, ifelse(vento<5, 0, vento-5))

(seg1.1 <- lm(energia ~ vento + xc5, data = eolica))
preds$pseg1.1 <- predict(seg1.1, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg1.1))
abline(v = 5, lty = 3, col = "gray")

## dois nós
eolica$xc7.5 <- with(eolica, ifelse(vento<7.5, 0, vento - 7.5))
preds$xc7.5 <- with(preds, ifelse(vento<7.5, 0, vento - 7.5))

(seg1.2 <- lm(energia ~ vento + xc5 + xc7.5, data = eolica))
preds$pseg1.2 <- predict(seg1.2, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg1.2))
abline(v = c(5, 7.5), lty = 3, col = "gray")

##
## Regressão segmentada (grau 2)
## - restrição de continuidade da função
## - restrição de continuidade na primeira derivada da função
##
(seg2.0 <- lm(energia ~ vento + I(vento^2), data = eolica))
preds$pseg2.0 <- predict(seg2.0, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg2.0, lty = 3))

eolica$xc <- with(eolica, ifelse(vento<5, 0, I(vento-5)))
preds$xc <- with(preds, ifelse(vento<5, 0, I(vento-5)))
eolica$x2c2 <- with(eolica, ifelse(vento<5, 0, I(vento^2)-5^2))
preds$x2c2 <- with(preds, ifelse(vento<5, 0, I(vento^2)-5^2))

(seg2.1 <- lm(energia ~ vento + I(vento^2) + xc + x2c2, data = eolica))
preds$pseg2.1 <- predict(seg2.1, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg2.1))
abline(v = 5, lty = 3, col = "gray")

eolica$xc2 <- with(eolica, ifelse(vento<5, 0, I((vento-5)^2)))
preds$xc2 <- with(preds, ifelse(vento<5, 0, I((vento-5)^2)))

(seg2.2 <- lm(energia ~ vento + I(vento^2) + xc2, data = eolica))
preds$pseg2.2 <- predict(seg2.2, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg2.2))


with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg2.0, lty = 3))
with(preds, lines(vento, pseg2.1))
with(preds, lines(vento, pseg2.2, col = 4))
abline(v = 5, lty = 3, col = "gray")

##
## Regressão segmentada (grau 3)
## - restrição de continuidade da função
## - restrição de continuidade na primeira derivada da função
## - restrição de continuidade na segunda derivada da função
##
(seg3.0 <- lm(energia ~ vento + I(vento^2) + I(vento^3), data = eolica))
preds$pseg3.0 <- predict(seg3.0, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg3.0, lty = 3))

eolica$xc3 <- with(eolica, ifelse(vento<5, 0, I((vento-5)^3)))
preds$xc3 <- with(preds, ifelse(vento<5, 0, I((vento-5)^3)))

(seg3.3 <- lm(energia ~ vento + I(vento^2) + I(vento^3) + xc3, data = eolica))
preds$pseg3.3 <- predict(seg3.3, newdata = preds)
with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg3.3))
abline(v = 5, lty = 3, col = "gray")

with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg3.0, lty = 3))
with(preds, lines(vento, pseg3.3, col = 4))
abline(v = 5, lty = 3, col = "gray")


(seg3.bs <- lm(energia ~ splines::bs(vento, knot = 5), data = eolica))
preds$pseg3.bs <- predict(seg3.bs, newdata = preds)
##with(eolica, plot(energia ~ vento))
with(preds, lines(vento, pseg3.bs, col = 2))
abline(v = 5, lty = 3, col = "gray")
