2  Introduction to ML

ML is about building an algorithm that can perform a certain tasks. The general steps in an ML are:

  1. Get some data
  2. Establish something that should be done with such data (= task)
  3. Find an algorithm that can perform the task
  4. Operationalize the goal by an error metric (=loss) that is to be minimized
  5. Train the algorithm to achieve a small loss on the data
  6. Test the performance of the algorithm
Example
  1. data = airquality
  2. predict Ozone
  3. algorithm = lm
  4. loss = residual sum of squares
  5. m = lm(Ozone ~ ., data = airquality)
  6. could predict(m) to data or hold-out

Our task in the course will now be to better understand

2.1 ML Tasks

There are three types of machine learning tasks:

  • Supervised learning
  • Unsupervised learning
  • Reinforcement learning

Supervised learning, you train algorithms to predict something (classes or values = regression) from some other data (= features), and you provide it with correct examples of the execution of the task (called training data). A linear regression is an example of supervised learning.

Unsupervised learning, on the other hand, is when you provide the features, but no examples of the correct execution of the task. Clustering techniques are examples of unsupervised learning

Reinforcement learning is a technique that mimics a game-like situation. The algorithm finds a solution through trial and error, receiving either rewards or penalties for each action. As in games, the goal is to maximize the rewards. We will talk more about this technique on the last day of the course.

For now, we will focus on the first two tasks, supervised and unsupervised learning (here a YouTube video explaining again the difference).

2.1.1 Test questions

In ML, predictors (or the explaining variables) are often called features:

In supervised learning the response (y) and the features (x) are known:

In unsupervised learning, only the features are known:

In reinforcement learning an agent (ML model) is trained by interacting with an environment:

Have a look at the two textbooks on ML (Elements of statistical learning and introduction to statistical learning) in our further readings at the end of the GRIPS course - which of the following statements is true?

2.2 Unsupervised Learning

In unsupervised learning, we want to identify patterns in data without having any examples (supervision) about what the correct patterns / classes are. As an example, consider the iris data set. Here, we have 150 observations of 4 floral traits:

iris = datasets::iris
colors = hcl.colors(3)
traits = as.matrix(iris[,1:4])
species = iris$Species
image(y = 1:4, x = 1:length(species) , z = traits,
      ylab = "Floral trait", xlab = "Individual")
segments(50.5, 0, 50.5, 5, col = "black", lwd = 2)
segments(100.5, 0, 100.5, 5, col = "black", lwd = 2)

Trait distributions of iris dataset

The observations are from 3 species and indeed those species tend to have different traits, meaning that the observations form 3 clusters.

pairs(traits, pch = as.integer(species), col = colors[as.integer(species)])

However, imagine we don’t know what species are, what is basically the situation in which people in the antique have been. The people just noted that some plants have different flowers than others, and decided to give them different names. This kind of process is what unsupervised learning does.

2.2.1 K-means Clustering

An example for an unsupervised learning algorithm is k-means clustering, one of the simplest and most popular unsupervised machine learning algorithms (see more on this in section “distance based algorithms”).

To start with the algorithm, you first have to specify the number of clusters (for our example the number of species). Each cluster has a centroid, which is the assumed or real location representing the center of the cluster (for our example this would be how an average plant of a specific species would look like). The algorithm starts by randomly putting centroids somewhere. Afterwards each data point is assigned to the respective cluster that raises the overall in-cluster sum of squares (variance) related to the distance to the centroid least of all. After the algorithm has placed all data points into a cluster the centroids get updated. By iterating this procedure until the assignment doesn’t change any longer, the algorithm can find the (locally) optimal centroids and the data points belonging to this cluster. Note that results might differ according to the initial positions of the centroids. Thus several (locally) optimal solutions might be found.

The “k” in K-means refers to the number of clusters and the ‘means’ refers to averaging the data-points to find the centroids.

A typical pipeline for using k-means clustering looks the same as for other algorithms. After having visualized the data, we fit a model, visualize the results and have a look at the performance by use of the confusion matrix. By setting a fixed seed, we can ensure that results are reproducible.

set.seed(123)

#Reminder: traits = as.matrix(iris[,1:4]).

kc = kmeans(traits, 3)
print(kc)
K-means clustering with 3 clusters of sizes 50, 62, 38

Cluster means:
  Sepal.Length Sepal.Width Petal.Length Petal.Width
1     5.006000    3.428000     1.462000    0.246000
2     5.901613    2.748387     4.393548    1.433871
3     6.850000    3.073684     5.742105    2.071053

Clustering vector:
  [1] 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
 [38] 1 1 1 1 1 1 1 1 1 1 1 1 1 2 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2
 [75] 2 2 2 3 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 2 3 2 3 3 3 3 2 3 3 3 3
[112] 3 3 2 2 3 3 3 3 2 3 2 3 2 3 3 2 2 3 3 3 3 3 2 3 3 3 3 2 3 3 3 2 3 3 3 2 3
[149] 3 2

Within cluster sum of squares by cluster:
[1] 15.15100 39.82097 23.87947
 (between_SS / total_SS =  88.4 %)

Available components:

[1] "cluster"      "centers"      "totss"        "withinss"     "tot.withinss"
[6] "betweenss"    "size"         "iter"         "ifault"      

Visualizing the results. Color codes true species identity, symbol shows cluster result.

plot(iris[c("Sepal.Length", "Sepal.Width")],
     col =  colors[as.integer(species)], pch = kc$cluster)
points(kc$centers[, c("Sepal.Length", "Sepal.Width")],
       col = colors, pch = 1:3, cex = 3)

We see that there are are some discrepancies. Confusion matrix:

table(iris$Species, kc$cluster)
            
              1  2  3
  setosa     50  0  0
  versicolor  0 48  2
  virginica   0 14 36

If you want to animate the clustering process, you could run

library(animation)

saveGIF(kmeans.ani(x = traits[,1:2], col = colors),
        interval = 1, ani.width = 800, ani.height = 800)

Elbow technique to determine the probably best suited number of clusters:

set.seed(123)

getSumSq = function(k){ kmeans(traits, k, nstart = 25)$tot.withinss }

#Perform algorithm for different cluster sizes and retrieve variance.
iris.kmeans1to10 = sapply(1:10, getSumSq)
plot(1:10, iris.kmeans1to10, type = "b", pch = 19, frame = FALSE,
     xlab = "Number of clusters K",
     ylab = "Total within-clusters sum of squares",
     col = c("black", "red", rep("black", 8)))

Often, one is interested in sparse models. Furthermore, higher k than necessary tends to overfitting. At the kink in the picture, the sum of squares dropped enough and k is still low enough. But keep in mind, this is only a rule of thumb and might be wrong in some special cases.

2.3 Supervised Learning

The two most prominent branches of supervised learning are regression and classification. The basic distinction between the two is that classification is about predicting a categorical variable, and regression is about predicting a continuous variable.

2.3.1 Regression

The random forest (RF) algorithm is possibly the most widely used machine learning algorithm and can be used for regression and classification. We will talk more about the algorithm later.

For the moment, we want to go through a typical workflow for a supervised regression: First, we visualize the data. Next, we fit the model and lastly we visualize the results. We will again use the iris data set that we used before. The goal is now to predict Sepal.Length based on the information about the other variables (including species).

Fitting the model:

library(randomForest)
set.seed(123)

Sepal.Length is a numerical variable:

str(iris)
'data.frame':   150 obs. of  5 variables:
 $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
 $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
 $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
 $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
 $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...
hist(iris$Sepal.Length)

The randomForest can be used similar to a linear regression model, we can specify the features using the formula syntax (~. means that all other variables should be used as features):

m1 = randomForest(Sepal.Length ~ ., data = iris)   # ~.: Against all others.
print(m1)

Call:
 randomForest(formula = Sepal.Length ~ ., data = iris) 
               Type of random forest: regression
                     Number of trees: 500
No. of variables tried at each split: 1

          Mean of squared residuals: 0.1364625
                    % Var explained: 79.97

In statistics we would use a linear regression model:

mLM = lm(Sepal.Length~., data = iris)

As many other ML algorithms, the RF is not interpretable, so we don’t get coefficients that connect the variables to the response. But, at least we get the variable importance which is similar to an anova, telling us which variables were the most important ones:

varImpPlot(m1)

Our liner model would report linear effects, however, the lm cannot keep up with the flexibility of a random forest!

coef(mLM)
      (Intercept)       Sepal.Width      Petal.Length       Petal.Width 
        2.1712663         0.4958889         0.8292439        -0.3151552 
Speciesversicolor  Speciesvirginica 
       -0.7235620        -1.0234978 

And the finally, we can use the model to make predictions using the predict method:

plot(predict(m1), iris$Sepal.Length, xlab = "Predicted", ylab = "Observed")
abline(0, 1)

To understand the structure of a random forest in more detail, we can use a package from GitHub.

reprtree:::plot.getTree(m1, iris)

Here, one of the regression trees is shown.

2.3.2 Classification

With the random forest, we can also do classification. The steps are the same as for regression tasks, but we can additionally see how well it performed by looking at the confusion matrix. Each row of this matrix contains the instances in a predicted class and each column represents the instances in the actual class. Thus the diagonals are the correctly predicted classes and the off-diagonal elements are the falsely classified elements.

Species is a factor with three levels:

str(iris)
'data.frame':   150 obs. of  5 variables:
 $ Sepal.Length: num  5.1 4.9 4.7 4.6 5 5.4 4.6 5 4.4 4.9 ...
 $ Sepal.Width : num  3.5 3 3.2 3.1 3.6 3.9 3.4 3.4 2.9 3.1 ...
 $ Petal.Length: num  1.4 1.4 1.3 1.5 1.4 1.7 1.4 1.5 1.4 1.5 ...
 $ Petal.Width : num  0.2 0.2 0.2 0.2 0.2 0.4 0.3 0.2 0.2 0.1 ...
 $ Species     : Factor w/ 3 levels "setosa","versicolor",..: 1 1 1 1 1 1 1 1 1 1 ...

Fitting the model (syntax is the same as for the regression task):

set.seed(123)
library(randomForest)
m1 = randomForest(Species ~ ., data = iris)
print(m1)

Call:
 randomForest(formula = Species ~ ., data = iris) 
               Type of random forest: classification
                     Number of trees: 500
No. of variables tried at each split: 2

        OOB estimate of  error rate: 4.67%
Confusion matrix:
           setosa versicolor virginica class.error
setosa         50          0         0        0.00
versicolor      0         47         3        0.06
virginica       0          4        46        0.08
varImpPlot(m1)

Predictions:

head(predict(m1))
     1      2      3      4      5      6 
setosa setosa setosa setosa setosa setosa 
Levels: setosa versicolor virginica

Confusion matrix:

table(predict(m1), as.integer(iris$Species))
            
              1  2  3
  setosa     50  0  0
  versicolor  0 47  4
  virginica   0  3 46

Our model made a few errors.

Visualizing results ecologically:

plot(iris$Petal.Width, iris$Petal.Length, col = iris$Species, main = "Observed")

plot(iris$Petal.Width, iris$Petal.Length, col = predict(m1), main = "Predicted")

Visualizing one of the fitted models:

reprtree:::plot.getTree(m1, iris)

Confusion matrix:

knitr::kable(table(predict(m1), iris$Species))

2.4 Exercise - Supervised Learning

Using a random forest on the iris dataset, which parameter would be more important (remember there is a function to check this) to predict Petal.Width?
Task: First deep neural network

Deep neural networks are currently the state of the art in unsupervised learning. Their ability to model different types of data (e.g. graphs, images) is one of the reasons for their rise in recent years. However, their use beyond tabular data (tabular data == features have specific meanings) requires extensive (programming) knowledge of the underlying deep learning frameworks (e.g. TensorFlow or PyTorch), which we will teach you in two days. For tabular data, we can use packages like cito, which work similarly to regression functions like lm and allow us to train deep neural networks in one line of code.

A demonstration with the iris dataset:

library(cito)

# always scale your features when using DNNs
iris_scaled = iris
iris_scaled[,1:4] = scale(iris_scaled[,1:4])

# the default architecture is 3 hidden layers, each with 10 hidden nodes (we will talk on Wednesday more about the architecture)
# Similar to a lm/glm we have to specify the response/loss family, for multi-target (3 species) we use the softmax loss function
model = dnn(Species~., lr = 0.1,data = iris_scaled, loss = "softmax", verbose = FALSE)

DNNs are not interpretable, i.e. no coefficients (slopes) that tell us how the features affect the response, however, similar to the RF, we can calculate a ‘variable importance’ which is similar to an anova:

summary(model)
Summary of Deep Neural Network Model

Feature Importance:
      variable importance_1
1 Sepal.Length     2.290604
2  Sepal.Width     5.880450
3 Petal.Length    78.079413
4  Petal.Width    65.322820

Average Conditional Effects:
                Response_1  Response_2  Response_3
Sepal.Length -0.0008044059  0.02891530 -0.02811095
Sepal.Width   0.0018329182  0.02917479 -0.03100772
Petal.Length -0.0018289219 -0.17218131  0.17401022
Petal.Width  -0.0022610827 -0.13096507  0.13322612

Standard Deviation of Conditional Effects:
              Response_1 Response_2 Response_3
Sepal.Length 0.002954661 0.09159783 0.09178804
Sepal.Width  0.006163875 0.11736476 0.11675831
Petal.Length 0.007047333 0.63252216 0.63205226
Petal.Width  0.008666864 0.47562443 0.47497322

Predictions

head(predict(model, type = "response"))
        setosa   versicolor    virginica
[1,] 0.9999356 6.438165e-05 2.101880e-16
[2,] 0.9993612 6.388194e-04 3.643038e-15
[3,] 0.9998255 1.745036e-04 9.812130e-16
[4,] 0.9996601 3.399171e-04 3.991938e-15
[5,] 0.9999580 4.191478e-05 1.591412e-16
[6,] 0.9999479 5.205601e-05 4.544236e-16

We get three columns, one for each species, and they are probabilities.

plot(iris$Sepal.Length, iris$Sepal.Width, col = apply(predict(model, type = "response"), 1, which.max))

Performance:

table(apply(predict(model), 1, which.max), as.integer(iris$Species))
   
     1  2  3
  1 50  0  0
  2  0 49  1
  3  0  1 49

Task:

  • predict Sepal.Length instead of Species (classification -> regression)
  • Use the ‘mse’ loss function
  • Plot predicted vs observed

Regression:

losses such as “mse” (mean squared error) or the “msa” (mean absolute error) are used for regression tasks

model = dnn(Sepal.Length~., lr = 0.1,data = iris_scaled, loss = "mse")

summary(model)
Summary of Deep Neural Network Model

Feature Importance:
      variable importance_1
1  Sepal.Width     1.877960
2 Petal.Length    28.846397
3  Petal.Width     1.545056
4      Species     1.715402

Average Conditional Effects:
             Response_1
Sepal.Width   0.2516230
Petal.Length  1.4459693
Petal.Width  -0.2616228

Standard Deviation of Conditional Effects:
             Response_1
Sepal.Width  0.09583777
Petal.Length 0.50299187
Petal.Width  0.10573136
plot(iris_scaled$Sepal.Length, predict(model))

Calculate \(R^2\):

cor(iris_scaled$Sepal.Length, predict(model))**2
          [,1]
[1,] 0.8819304