<jupyter_start><jupyter_text>Section F. Logistic Regression vs. Bayesian Classifier  Question 7 [Discriminative vs Generative Models, 25 Marks]  Library<jupyter_code># library used
library(mvtnorm) # generates multivariate Gaussian sampels and calculate the densities
library(ggplot2) # draw ggplots
library(reshape2) # reshape dataframes<jupyter_output><empty_output><jupyter_text>Import data<jupyter_code># read csv
train_1e <- read.csv("./Task1E_train.csv")
test_1e <- read.csv("./Task1E_test.csv")
# subset data into x and y
train.data_all <- train_1e[,-3]
train.label_all <- train_1e[,3]
test.data_all <- test_1e[,-3]
test.label_all <- test_1e[,3]<jupyter_output><empty_output><jupyter_text>1. Create two models  Bayesian classifier<jupyter_code># a dataframe to record error of bc model
bc_error  <- data.frame("size" = seq(5,500,5), "train" = rep(0,100), "test" = rep(0,100))
# number of data points use in each run, increase by 5 each time
N <- seq(5,500,5)
# initialize count
c <- 1 
# run bc model
for (n in N){
    # subset data according to the number of data points use in this run
    train.data <- train.data_all[1:n,]
    train.label <- train.label_all[1:n]

    # 1. calculating the class probabilities
    p0.hat <- sum(train.label==1)/nrow(train.data) # total number of samples in class 0 divided by the total nmber of training data
    p1.hat <- sum(train.label==-1)/nrow(train.data) # or simply 1 - p1.hat

    # 2. estimate the class means
    mu0.hat <- colMeans(train.data[train.label==1,])
    mu1.hat <- colMeans(train.data[train.label==-1,])

    # 3. calculate class variance matrices
    # class covariance matrices:
    sigma0.hat <- var(train.data[train.label==1,])
    sigma1.hat <- var(train.data[train.label==-1,])

    # shared covariance matrix:
    sigma.hat <- p0.hat * sigma0.hat + p1.hat * sigma1.hat 

    # 4. calculate the postoriors based on the estimated parameters
    posterior0 <- p0.hat*dmvnorm(x=train.data_all, mean=mu0.hat, sigma=sigma.hat)
    posterior1 <- p1.hat*dmvnorm(x=train.data_all, mean=mu1.hat, sigma=sigma.hat)

    # 5. calculate predictions and assign to class
    train.predict <- ifelse(posterior0 > posterior1, 1, -1)
    test.predict <- ifelse(p0.hat*dmvnorm(x=test.data_all, mean=mu0.hat, sigma=sigma.hat) > p1.hat*dmvnorm(x=test.data_all, mean=mu1.hat, sigma=sigma.hat), 1, -1)
    
    # calculate accuracy:
    bc_error[c,2] <- p0.hat * sum(train.label_all!=train.predict)  # record train error
    bc_error[c,3] <- p0.hat * sum(test.label_all!=test.predict) # record test error
    c <- c + 1 # update count
}<jupyter_output><empty_output><jupyter_text>Logistic Regression<jupyter_code># auxiliary function that predicts class labels
predict <- function(w, X, c0, c1){
    sig <- sigmoid(w, X)
    return(ifelse(sig>0.5, c1,c0))
}
    
# auxiliary function that calculate a cost function
cost <- function (w, X, T, c0){
    sig <- sigmoid(w, X)
    return(sum(ifelse(T==c0, 1-sig, sig)))
}
# 1. define sigmoid function (=p(C1|X))
sigmoid <- function(w, x){
    return(1.0/(1.0+exp(-w%*%t(cbind(1,x)))))    
}

# 2. Initialize the weight vector, learning rate, stopping threshold, etc.
tau.max <- 1000 # maximum number of iterations
eta <- 0.01 # learning rate
epsilon <- 0.01 # a threshold on the cost (to terminate the process)
terminate <- FALSE

# a dataframe to record error of bc model
Lr_error  <- data.frame("size" = seq(5,500,5), "train" = rep(0,100), "test" = rep(0,100))
# number of data points use in each run, increase by 5 each time
N <- seq(5,500,5)
# initialize count
c <- 1 ## Just a few name/type conversion to make the rest of the code easy to follow



# run SGD to learn the weight vector, + 5 data points each time in trainning
for (n in N){
    # subset train data according to number of data points going to use
    train.data <- train.data_all[1:n,]
    train.label <- train.label_all[1:n]

    X <- as.matrix(train.data) # rename just for conviniance
    T <- train.label # rename just for conviniance
    W <- matrix(,nrow=tau.max, ncol=(ncol(X)+1)) # to be used to store the estimated coefficients

    set.seed(1343)
    W[1,] <- runif(ncol(W)) # initial weight (any better idea?)

    # project data using the sigmoid function (just for convenient)
    Y <- sigmoid(W[1,],X)
    
    costs <- data.frame('tau'=1:tau.max)  # to be used to trace the cost in each iteration
    costs[1, 'train'] <- cost(W[1,],as.matrix(train.data_all),train.label_all, 1) # initial costs
    costs[1,'test'] <- cost(W[1,],as.matrix(test.data_all),test.label_all, 1)


    # shuffle data:
    train.len  <- nrow(train.data)
    train.index <- sample(1:train.len, train.len, replace = FALSE)
    X <- X[train.index,]
    T <- T[train.index]
    tau <- 1 # iteration counter
    # for each datapoint:
    for (i in 1:train.len){

        # calcualte prediction
        Y <- sigmoid(W[tau,],X)

        # Update the weights
        W[(tau+1),] <- W[tau,] - eta * (Y[i]-T[i]) * cbind(1, t(X[i,]))

        # record the cost of train and test data:
        costs[(tau+1), 'train'] <- cost(W[tau,],as.matrix(train.data_all),train.label_all, 1)
        costs[(tau+1), 'test'] <- cost(W[tau,],as.matrix(test.data_all),test.label_all, 1)

        # update the counter:
        tau <- tau + 1
        u <- u + 1
        # decrease learning rate:
        eta = eta * 0.999
    }

# record errors 
    Lr_error[c,'train'] <- costs[tau,"train"]
    Lr_error[c,'test'] <- costs[tau,"test"]
    c  <-  c + 1
}

# merge two error dataframe together for plotting
bc_error[,"type"]  <- "BC"
Lr_error[,"type"] <- "LR"
model_error <- rbind(bc_error,Lr_error)<jupyter_output><empty_output><jupyter_text>2. Plot the errors of the models<jupyter_code># plot training error vs trainning size for the two model
ggplot(data = model_error, aes(x = size, y = train, color = type)) + geom_line() + 
scale_color_discrete(guide = guide_legend(title = 'Model'))+
ggtitle('Training error vs training size')+
theme_minimal()

# plot testing error vs testing size for the two model
ggplot(data = model_error, aes(x = size, y = test, color = type)) + geom_line() + 
scale_color_discrete(guide = guide_legend(title = 'Model'))+
ggtitle('Training error vs training size')+
theme_minimal()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.14431673052362706, 'fraction_numerical': 0.026500638569604086, 'mean_word_length': 3.6788648244958924, 'pattern_counts': {'":': 0, '<': 85, '<?xml version=': 0, '>': 24, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '5', 'dir': 'jupyter-structured-clean-dedup', 'id': '577126', 'license': 'MIT', 'n_tokens_mistral': 2398, 'n_tokens_neox': 2201, 'n_words': 608, 'path': 'Machine Learning study notes/Linear Models for Classification and regression/src/31436285_assessment_1_q7.ipynb', 'repo_name': 'alanwzy/my-projects', 'size': '6264'}