Assignment 9
EMATM0061: Statistical Computing and Empirical Methods, TB1, 2024
Introduction
This is the 9th assignment for Statistical Computing and Empirical Methods. This assignment is mainly based on Lectures 25, 26, 27 (see the Blackboards). Please note that you don’t need to submit this assignment.
Load the tidyverse package:
library(tidyverse)
Some questions require the “glmnet” and “QSARdata” packages. If they are not already installed on your computer, please install them using “install.packages()” first.
1. Linear discriminant analysis
(Q1)
Describe the probabilistic model that underpins linear discriminant analysis.
(Q2)
In this question, we will train a linear discriminant analysis model to carry out the classification task to predict whether a hawk belongs to either the “Sharp-shinned” (SS) or the “Cooper’s” (CH) species of hawks, based on a four dimensional feature vector containing the weight, and the lengths of the wing, the tail and the hallux, which are generated by the following code.
library(Stat2Data)
data(Hawks)
hawks_total <- Hawks %>% select( Weight, Wing, Hallux, Tail, Species)
%>%
filter(Species=='SS' | Species =='CH') %>% drop_na() %>%
mutate(Species=as.numeric(Species=='SS'))
Assume that we have the following train-test split of our dataset.
num_total <- hawks_total %>% nrow() # number of penguin data
num_train <- floor(num_total*0.6) # number of train examples
num_test <- num_total-num_train # number of test samples
set.seed(0) # set random seed for reproducibility
test_inds <- sample(seq(num_total),num_test) # random sample of test
indicies
train_inds <- setdiff(seq(num_total),test_inds) # training data
indicies
hawks_train <- hawks_total %>% filter(row_number() %in% train_inds) #
train data
hawks_test <- hawks_total %>% filter(row_number() %in% test_inds) #
test data
Now, train a linear discriminant analysis model to carry out the classification task described above. Compute and report the train error and the test error.
(Q3) (*Optional) Implement your own linear discriminant analysis model. You should write a function called “generate_lda_model” which has two arguments, namely “training_data” (a data frame containing the training data) and “y_col_name” (a string specifying the name of the column that presents the labels). The remaining columns of “training_data” represents the feature vectors. The function “generate_lda_model” returns a function that is your linear discriminant analysis model trained on the training data.
2. Logistic regression
(Q1) Describe the probabilistic model which underpins logistic regression.
(Q2)
Recall that the sigmoid function S: ℝ → (0,1) is defined by s(z) = 1/(1 + e一Z). Generate the following plot which displays the sigmoid function:
(Q3) Now train a logistic regression model to predict whether a hawk belongs to either the “Sharp-shinned” or the “Cooper’s” species of hawks, based on a four-dimensional feature vector containing the weight, and the lengths of the wing, the tail and the hallux (similar to what you did in the Linear discriminant analysis question above. You can reuse the training test split in that question). Compute and report both the training error and the test error.
(Q4) (*optional)
Consider the following formula for the log-likelihood of the weights w ∈ ℝd and bias wo ∈ ℝ, given data D = ((x1, Y1), ⋯ , (xn, Yn)):
Show that and use it to demonstrate the following formulas for the derivatives:
Explain the role the above formula has in training a logistic regression model.
You can learn more about the glmnet approach to logistic regression here:
https://glmnet.stanford.edu/articles/glmnet.html#logistic-regression
3. Basic concepts in regularisation
Regularisation refers to the general technique within supervised learning of modifying an objective in some way so as to reduce the gap between test error and training error. Typically, this will increase the error on the training data. However, by reducing the gap between test and training errors, we can often improve performance (on unseen data).
Examples include
1) l2 regularisation in the context of ridge regression for regression
2) l1 or l2 regularised logistic regression for linear classification
(Q1)
Let’s review some key concepts relevant to regularisation. Write down your explanation of each of the following concepts.
1. Hyper-parameter (and give an example of a hyper-parameter)
2. Validation data
3. The train-validation-test split
(Q2) What is the Euclidean (l2) norm and what is the l1 norm of a vector?
(Q3) The Ridge regression method and the Lasso method are both for learning a linear regression model from the data, by minimising an objective function. Describe what kinds of terms are included in their objective functions. What is the difference between the two objective functions?
4. An investigation into ridge regression for high-dimensional regression
In this question we consider a high-dimensional regression problem. We consider a problem of predicting the melting point of a chemical compound from a relatively high-dimensional feature vector of chemical descriptors.
To do this we shall use data from the ““QSARdata”” data library. Begin by checking if the ““QSARdata”” library has been installed.
Next load the ““QSARdata”” library and load the ““MeltingPoint”” data set.
library(QSARdata)
data(MeltingPoint)
You will find a data frame called ““MP_Descriptors”” . The rows of the data frame correspond to different examples of chemical compounds and the columns correspond to various chemical descriptors. In addition you will find a vector called ““MP_Outcome”” which contains the corresponding melting point for each of the examples.
(Q1)
Begin by combining the data-frame of feature vectors ““MP_Descriptors”” together with the column vector of melting points ““MP_Outcome”” . Combine these together into a single data frame entitled ““mp_data_total”” as follows.
mp_data_total<-MP_Descriptors %>%
mutate(melting_pt=MP_Outcome)
How many variables are in your data frame? How many examples?
(Q2)
Next carry out a train-validate-test split of the ““mp_data_total”” data frame. You should use about 50% of the data to train the algorithm, about 25% to validate and about 25% to test
(Q3)
Our goal is to find a linear regression model φw,wo : ℝd → ℝ by φw,wo (X) = WXT + Wo which estimates the melting point based upon the chemical descriptors. Create a function which takes as input training data (a matrix of features for training data and a vector of labels for training data), validation data (a matrix of features for validation data and a vector of labels for validation data) and a hyper-parameter λ . The function should train a ridge regression model using the specified value of λ, then compute the validation error and output a single number corresponding to the validation error. Your function should not be specific to this particular regression problem, but should apply to ridge regression problems in general. You can use the “glmnet” library within your function.
(Q4)
Next generate a sequence of candidate hyper-parameters called ““lambdas”” . The sequence should begin with 10一5 and increase geometrically in multiples of 1.25 and should be of length 70. That is, ““lambdas”” should contain the numbers
(Q5)
Now use your function to estimate the mean squared error on the validation data for a ridge regression model for the problem of predicting the melting point based on the chemical descriptors. Consider all of the hyperparameter values within your vector ““lambdas”” . Store the results of this procedure in a data frame.
(Q6)
Plot the validation error as a function of the hyper-parameter λ. Use the “scale x continuous()” function to plot the λ coordinate on a logarithmc scale.
Your plot may look like this:
(Q7)
Now use your results data frame. to determine the hyper-parameter λ with the lowest validation error. Retrain your ridge regression model with your selected value of the hyper-parameter λ and estimate the test error by computing the mean squared error on the test data.
(Q8) (*optional)
Does the test error computed above lead to a biased estimate of the mean squared error? Why can’t we use the mean squared error on validation data for the ridge regression model with the selected hyper-parameter as an estimate of the mean squared error on test data? Observe that the ridge regression model with the selected choice of λ has actually been trained twice: It was trained in order to compute the validation error, and then again to compute the test error. Comment on the computational efficiency of this procedure. Could it be easily improved? What memory implications would this have?