Lab 4 - Logistic Regression, LDA, QDA, KNN

We will now try to train models; the target :Direction has two classes: Up and Down; it needs to be interpreted as a categorical object, and we will mark it as a ordered factor to specify that 'Up' is positive and 'Down' negative (for the confusion matrix later):

y = coerce(y, OrderedFactor)
classes(y[1])

2-element CategoricalArrays.CategoricalArray{String,1,UInt8}:
 "Down"
 "Up"

Note that in this case the default order comes from the lexicographic order which happens to map to our intuition since D comes before U.

cm = countmap(y)
categories, vals = collect(keys(cm)), collect(values(cm))
Plots.bar(categories, vals, title="Bar Chart Example", legend=false)
ylabel!("Number of occurrences")

Seems pretty balanced.

Let's now try fitting a simple logistic classifier (aka logistic regression) not using :Year and :Today:

LogisticClassifier = @load LogisticClassifier pkg=MLJLinearModels
X2 = select(X, Not([:Year, :Today]))
classif = machine(LogisticClassifier(), X2, y)

import MLJLinearModels ✔
untrained Machine; caches model-specific representations of data
  model: LogisticClassifier(lambda = 2.220446049250313e-16, …)
  args: 
    1:	Source @283 ⏎ ScientificTypesBase.Table{AbstractVector{ScientificTypesBase.Continuous}}
    2:	Source @894 ⏎ AbstractVector{ScientificTypesBase.OrderedFactor{2}}

Let's fit it to the data and try to reproduce the output:

fit!(classif)
ŷ = MLJ.predict(classif, X2)
ŷ[1:3]

3-element CategoricalDistributions.UnivariateFiniteVector{ScientificTypesBase.OrderedFactor{2}, String, UInt8, Float64}:
 UnivariateFinite{ScientificTypesBase.OrderedFactor{2}}(Down=>0.493, Up=>0.507)
 UnivariateFinite{ScientificTypesBase.OrderedFactor{2}}(Down=>0.519, Up=>0.481)
 UnivariateFinite{ScientificTypesBase.OrderedFactor{2}}(Down=>0.519, Up=>0.481)

Note that here the ŷ are scores. We can recover the average cross-entropy loss:

cross_entropy(ŷ, y) |> mean |> r3

0.691

in order to recover the class, we could use the mode and compare the misclassification rate:

ŷ = predict_mode(classif, X2)
misclassification_rate(ŷ, y) |> r3

0.478

Well that's not fantastic...

Let's visualise how we're doing building a confusion matrix, first is predicted, second is truth:

cm = confusion_matrix(ŷ, y)

          ┌─────────────┐
          │Ground Truth │
┌─────────┼──────┬──────┤
│Predicted│ Down │  Up  │
├─────────┼──────┼──────┤
│  Down   │ 145  │ 141  │
├─────────┼──────┼──────┤
│   Up    │ 457  │ 507  │
└─────────┴──────┴──────┘

We can then compute the accuracy or precision, etc. easily for instance:

@show false_positive(cm)
@show accuracy(ŷ, y)  |> r3
@show accuracy(cm)    |> r3  # same thing
@show positive_predictive_value(ŷ, y) |> r3   # a.k.a. precision
@show recall(ŷ, y)    |> r3
@show f1score(ŷ, y)   |> r3

false_positive(cm) = 457
accuracy(ŷ, y) |> r3 = 0.522
accuracy(cm) |> r3 = 0.522
positive_predictive_value(ŷ, y) |> r3 = 0.526
recall(ŷ, y) |> r3 = 0.782
f1score(ŷ, y) |> r3 = 0.629

Let's now train on the data before 2005 and use it to predict on the rest. Let's find the row indices for which the condition holds

train = 1:findlast(X.Year .< 2005)
test = last(train)+1:length(y);

We can now just re-fit the machine that we've already defined just on those rows and predict on the test:

fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)
accuracy(ŷ, y[test]) |> r3

0.48

Well, that's not very good... Let's retrain a machine using only :Lag1 and :Lag2:

X3 = select(X2, [:Lag1, :Lag2])
classif = machine(LogisticClassifier(), X3, y)
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)
accuracy(ŷ, y[test]) |> r3

0.56

Interesting... it has higher accuracy than the model with more features! This could be investigated further by increasing the regularisation parameter but we'll leave that aside for now.

We can use a trained machine to predict on new data:

Xnew = (Lag1 = [1.2, 1.5], Lag2 = [1.1, -0.8])
ŷ = MLJ.predict(classif, Xnew)
ŷ |> pprint

CategoricalDistributions.UnivariateFinite{ScientificTypesBase.OrderedFactor{2}, String, UInt8, Float64}[UnivariateFinite{ScientificTypesBase.OrderedFactor{2}}(Down=>0.521, Up=>0.479), UnivariateFinite{ScientificTypesBase.OrderedFactor{2}}(Down=>0.504, Up=>0.496)]

Note: when specifying data, we used a simple NamedTuple; we could also have defined a dataframe or any other compatible tabular container. Note also that we retrieved the raw predictions here i.e.: a score for each class; we could have used predict_mode or indeed

mode.(ŷ)

2-element CategoricalArrays.CategoricalArray{String,1,UInt8}:
 "Down"
 "Down"

‎

Let's do a similar thing but with a LDA model this time:

BayesianLDA = @load BayesianLDA pkg=MultivariateStats

classif = machine(BayesianLDA(), X3, y)
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)

accuracy(ŷ, y[test]) |> r3

import MLJMultivariateStatsInterface ✔
0.56

Note: BayesianLDA is LDA using a multivariate normal model for each class with a default prior inferred from the proportions for each class in the training data. You can also use the bare LDA model which does not make these assumptions and allows using a different metric in the transformed space, see the docs for details.

LDA = @load LDA pkg=MultivariateStats
using Distances

classif = machine(LDA(dist=CosineDist()), X3, y)
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)

accuracy(ŷ, y[test]) |> r3

import MLJMultivariateStatsInterface ✔
0.548

‎

Bayesian QDA is available via ScikitLearn:

BayesianQDA = @load BayesianQDA pkg=MLJScikitLearnInterface

import MLJScikitLearnInterface ✔
MLJScikitLearnInterface.BayesianQDA

Using it is done in much the same way as before:

classif = machine(BayesianQDA(), X3, y)
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)

accuracy(ŷ, y[test]) |> r3

0.599

‎

We can use K-Nearest Neighbors models via the NearestNeighbors package:

KNNClassifier = @load KNNClassifier

knnc = KNNClassifier(K=1)
classif = machine(knnc, X3, y)
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)
accuracy(ŷ, y[test]) |> r3

import NearestNeighborModels ✔
0.5

Pretty bad... let's try with three neighbors

knnc.K = 3
fit!(classif, rows=train)
ŷ = predict_mode(classif, rows=test)
accuracy(ŷ, y[test]) |> r3

0.532

A bit better but not hugely so.

‎

Since we have a probabilistic classifier, we can also check metrics that take scores into account such as the area under the ROC curve (AUC):

ŷ = MLJ.predict(classif, rows=test)

auc(ŷ, y[test])

0.7388551666996884

We can also display the curve itself

fprs, tprs, thresholds = roc_curve(ŷ, y[test])

plot(fprs, tprs, linewidth=2, size=(800,600))
xlabel!("False Positive Rate")
ylabel!("True Positive Rate")

‎