---
title: "Trying to be rich"
output: rmdformats::readthedown
---

We will work on the EuStockMarkets dataset which has 4 financial indexes for various European countries. Our goal is to predict if the CAC40 will increase next day.

## Data pre-processing

Let's have a look at the data.
```{r}
plot(EuStockMarkets)
```
As expected we can see that all indexed increase with time and there seems to have no anomalies in our data. To be rich we need to predict if the CAC40 will increase next day. Hence we will work with increments and values of yesterday for the other indexes. 

```{r}
inc_cac <- diff(EuStockMarkets[,"CAC"])## increments for CAC40
```
We need to do that for each index and the way we code it above is not practical. It's better to do this way
```{r}
increments <- apply(EuStockMarkets, 2, diff)##2 means that we apply the function diff on the columns of EuStockMarkets
```

For real life modelling, we cannot use future observations so we have to work with the following dataset
```{r}
nobs <- nrow(increments)
outcome <- factor(ifelse(increments[-1,"CAC"] > 0, "rich", "not rich"))
df <- data.frame(outcome = outcome, increments[-nobs,])
```

Now we can have a look at the data.
```{r}
plot(df[,-1], col = df$outcome)
```
Clearly and as expected, variables are (strongly) correlated. This is unfortunate. Why? Because it means that although we have 4 features, we don't have so many information.

## Modelling

We have a binary outcome so the use of the logistic regression model is sensible. Let's try!
```{r}
first_fit <- glm(outcome ~ ., data = df, family = binomial)
summary(first_fit)
```


The model is very poor, no *statistically* significant variable. Maybe we can complexify our model using squared values
```{r}
fit_second <- glm(outcome ~ .^2, data = df, family = binomial)
summary(fit_second)
```
Not a good improvement. Let's work with the sign of the increments as well.
```{r}
pos_or_neg <- function(x)
  ifelse(x > 0, "positive", "negative")

new_feat <- apply(df[,-1], 2, pos_or_neg)
colnames(new_feat) <- paste(colnames(new_feat), "sign", sep = "_")
df <- data.frame(df, new_feat)
```

Last attempt for feature engineering.
```{r}
fit <- glm(outcome ~ .^2, data = df, family = binomial)
summary(fit)
```

Ok not so bad. Now we do model selection using the *stepAIC* function (a bit harsh since I don't like automatic model selection but this is live coding so...)
```{r}
library(MASS)
fit <- stepAIC(fit, trace = 0)
summary(fit)
```

## Analyzing the odds

Although it is not mandatory for this application to analyze the odds, we want to earn mony but, better, learn a statistical analysis. I will do it only intrepret for variable *DAX_signpositive*. As $\exp(-4.239e-01) \approx 0.65$, having a positive increment for DAX as opposed to a negative one, decrease the odds by amount of $35%$.

## Model performance: Will we be rich?

Trying to assess the performance of our binary classifier using the *caret* package that we need to install first
```{r eval = FALSE}
install.packages("caret")
```
and load
```{r}
library(caret)
```


```{r}
pred <- factor(ifelse(predict(fit, type = "response") > 0.5, "rich", "not rich"))
confusionMatrix(pred, df$outcome)
```

Now to get the ROC curve from the *pROC* package we install it
```{r eval = FALSE}
install.packages("pROC")
```
adn load it as before
```{r}
library(pROC)
```

We are now all set to plot the ROC curve
```{r}
test_prob = predict(fit, newdata = df[,-1], type = "response")
test_roc = roc(df$outcome ~ test_prob, plot = TRUE, print.auc = TRUE)
```
We can see that our logistic classifier is slightly better than the random coin classifier. Unfortunately we won't be millionaires :-(
