TaggerTestResults.Rmd

---
title: "Tagger test"
author: "Elen Le Foll"
date: "02/10/2021"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)

library(DescTools)
library(caret)
library(here)
library(paletteer)
library(readxl)
library(svglite)
library(tidyverse)

# Set the random number generator seed for reproducibility.
set.seed(13)

```

# Data import 

These chunks import the data directly from the Excel files in which I did the manual tag check and corrections. All warning messages can safely be ignored.

```{r import-evaluation-files}

importEval3 <- function(file, fileID, register, corpus) {
  Tag1 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag1, Tag1Gold) %>% 
  rename(Tag = Tag1, TagGold = Tag1Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)
  
  Tag2 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag2, Tag2Gold) %>% 
  rename(Tag = Tag2, TagGold = Tag2Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

Tag3 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag3, Tag3Gold) %>% 
  rename(Tag = Tag3, TagGold = Tag3Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

output <- rbind(Tag1, Tag2, Tag3) %>% 
  mutate(across(where(is.factor), str_remove_all, pattern = fixed(" "))) %>% # Removes all white spaces which are found in the excel files
  filter(!is.na(Output)) %>% 
  mutate_if(is.character, as.factor)

}

importEval4 <- function(file, fileID, register, corpus) {
  Tag1 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag1, Tag1Gold) %>% 
  rename(Tag = Tag1, TagGold = Tag1Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)
  
  Tag2 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag2, Tag2Gold) %>% 
  rename(Tag = Tag2, TagGold = Tag2Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

Tag3 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag3, Tag3Gold) %>% 
  rename(Tag = Tag3, TagGold = Tag3Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

Tag4 <- file %>% 
  add_column(FileID = fileID, Register = register, Corpus = corpus) %>%
  select(FileID, Corpus, Register, Output, Tokens, Tag4, Tag4Gold) %>% 
  rename(Tag = Tag4, TagGold = Tag4Gold, Token = Tokens) %>% 
  mutate(Evaluation = ifelse(is.na(TagGold), TRUE, FALSE)) %>% 
  mutate(TagGold = ifelse(is.na(TagGold), as.character(Tag), as.character(TagGold))) %>%
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

output <- rbind(Tag1, Tag2, Tag3, Tag4) %>% 
  mutate(across(where(is.factor), str_remove_all, pattern = fixed(" "))) %>% # Removes all white spaces which are found in the excel files
  filter(!is.na(Tag)) %>% 
  mutate_if(is.character, as.factor)

}

importEval <- function(file, fileID, register, corpus) { 
  if(sum(!is.na(file$Tag4)) > 0) {
    output = importEval4(file = file, fileID = fileID, register = register, corpus = corpus)
  }
  else{
    output = importEval3(file = file, fileID = fileID, register = register, corpus = corpus)
  }
}

BNCBERe39 <- importEval(file = read_excel(here("evaluation", "BNCBERe39.xlsx")), fileID = "BNCBERe39", register = "internet", corpus = "BNC2014")

BNCBAcbH_m1 <- importEval(file = read_excel(here("evaluation", "BNCBAcbH_m1.xlsx")), fileID = "BNCBAcbH_m1", register = "academic", corpus = "BNC2014")

BNCBAcjS6 <- importEval(file = read_excel(here("evaluation", "BNCBAcjS6.xlsx")), fileID = "BNCBAcjS6", register = "academic", corpus = "BNC2014")

BNCBAcjM105 <- importEval(file = read_excel(here("evaluation", "BNCBAcjM105.xlsx")), fileID = "BNCBAcjM105", register = "academic", corpus = "BNC2014")

BNCBAcjM102 <- importEval(file = read_excel(here("evaluation", "BNCBAcjM102.xlsx")), fileID = "BNCBAcjM102", register = "academic", corpus = "BNC2014")

BNCBEBl8 <- importEval(file = read_excel(here("evaluation", "BNCBEBl8.xlsx")), fileID = "BNCBEBl8", register = "internet", corpus = "BNC2014")

BNCBEEm10 <- importEval(file = read_excel(here("evaluation", "BNCBEEm10.xlsx")), fileID = "BNCBEEm10", register = "internet", corpus = "BNC2014")

BNCBFict_b2 <- importEval(file = read_excel(here("evaluation", "BNCBFict_b2.xlsx")), fileID = "BNCBFict_b2", register = "fiction", corpus = "BNC2014")

BNCBFict_m54 <- importEval(file = read_excel(here("evaluation", "BNCBFict_m54.xlsx")), fileID = "BNCBFict_m54", register = "fiction", corpus = "BNC2014")

BNCBFict_e27 <- importEval(file = read_excel(here("evaluation", "BNCBFict_e27.xlsx")), fileID = "BNCBFict_e27", register = "fiction", corpus = "BNC2014")

BNCBEFor32 <- importEval(file = read_excel(here("evaluation", "BNCBEFor32.xlsx")), fileID = "BNCBEFor32", register = "internet", corpus = "BNC2014")

BNCBESm3 <- importEval(file = read_excel(here("evaluation", "BNCBESm3.xlsx")), fileID = "BNCBESm3", register = "internet", corpus = "BNC2014")

BNCBMass16 <- importEval(file = read_excel(here("evaluation", "BNCBMass16.xlsx")), fileID = "BNCBMass16", register = "news", corpus = "BNC2014")

BNCBMass23 <- importEval(file = read_excel(here("evaluation", "BNCBMass23.xlsx")), fileID = "BNCBMass23", register = "news", corpus = "BNC2014")

BNCBReg111 <- importEval(file = read_excel(here("evaluation", "BNCBReg111.xlsx")), fileID = "BNCBReg111", register = "news", corpus = "BNC2014")

BNCBReg750 <- importEval(file = read_excel(here("evaluation", "BNCBReg750.xlsx")), fileID = "BNCBReg750", register = "news", corpus = "BNC2014")

BNCBSer486 <- importEval(file = read_excel(here("evaluation", "BNCBSer486.xlsx")), fileID = "BNCBSer486", register = "news", corpus = "BNC2014")

BNCBSer562 <- importEval(file = read_excel(here("evaluation", "BNCBSer562.xlsx")), fileID = "BNCBSer562", register = "news", corpus = "BNC2014")

BNCBEsocFb <- importEval(file = read_excel(here("evaluation", "BNCBEsocFb.xlsx")), fileID = "BNCBEsocFb", register = "internet", corpus = "BNC2014")

S2DD <- importEval(file = read_excel(here("evaluation", "S2DD.xlsx")), fileID = "S2DD", register = "spoken", corpus = "BNC2014")

S3AV <- importEval(file = read_excel(here("evaluation", "S3AV.xlsx")), fileID = "S3AV", register = "spoken", corpus = "BNC2014")

SEL5 <- importEval(file = read_excel(here("evaluation", "SEL5.xlsx")), fileID = "SEL5", register = "spoken", corpus = "BNC2014")

SVLK <- importEval(file = read_excel(here("evaluation", "SVLK.xlsx")), fileID = "SVLK", register = "spoken", corpus = "BNC2014")

SZXQ <- importEval(file = read_excel(here("evaluation", "SZXQ.xlsx")), fileID = "SZXQ", register = "spoken", corpus = "BNC2014")

BNC2014Eval <- rbind(BNCBAcbH_m1, BNCBAcjS6, BNCBAcjM105, BNCBAcjM102, BNCBEBl8, BNCBEEm10, BNCBFict_b2, BNCBFict_m54, BNCBFict_e27, BNCBEFor32, BNCBERe39, BNCBESm3, BNCBMass16, BNCBMass23, BNCBReg111, BNCBReg750, BNCBSer486, BNCBSer562, BNCBEsocFb, S2DD, S3AV, SEL5, SVLK, SZXQ) 

summary(BNC2014Eval)

#saveRDS(BNC2014Eval, here("data", "MFTE_Evaluation_BNC2014_Results.rds")) # Last saved 30 Oct 2021

#write.csv(BNC2014Eval, here("data", "MFTE_Evaluation_BNC2014_Results.csv")) # Last saved 30 Oct 2021

```

## Quick import

```{r quick-import}

BNC2014Eval <- readRDS(here("data", "MFTE_Evaluation_BNC2014_Results.rds")) # 
summary(BNC2014Eval)

# Total number of tags manually checked
nrow(BNC2014Eval)

# Number of tags evaluated per file
BNC2014Eval %>% group_by(FileID) %>% count(.) %>% arrange(desc(n))

# Tagger evaluation
summary(BNC2014Eval$Evaluation)

# Number of UNCLEAR tokens
BNC2014Eval %>% filter(TagGold == "UNCLEAR") %>% count()

BinomCI(1335, 31311,
        conf.level = 0.95,
        sides = "two.sided",
        method = "wilsoncc") * 100

# Number of tags per feature
BNC2014Eval %>% group_by(TagGold) %>% count() %>% arrange(-n) %>% as.data.frame()


```

# Analysis

In this chunk, I calculate the recall and precision rates of each feature, ignoring unclear tokens.

```{r recall-precision-f1}

data <- BNC2014Eval %>% 
  filter(TagGold != "UNCLEAR") %>% 
  mutate(Tag = factor(Tag, levels = union(levels(Tag), levels(TagGold)))) %>% # Ensure that the factor levels are the same for the next caret operation
  mutate(TagGold = factor(TagGold, levels = union(levels(Tag), levels(TagGold))))

# Spot gold tag corrections that are not actually errors
data[data$Tag==data$TagGold & data$Evaluation == FALSE,] %>% as.data.frame()

nrow(data) # Number of tags checked
head(data) # Check sanity of data
summary(data) # Check sanity of data

cm <- caret::confusionMatrix(data$Tag, data$TagGold) # Create confusion matrix
cm$overall # Note that is not very representative because it includes tags which are not intended for use in MDA studies, e.g., LS and FW, or which are part of the evaluation process, e.g., NULL and UNCLEAR.

# Quick summary of results: recall, precision and f1
cm$byClass[,5:7]

# Generate a better formatted results table: recall, precision and f1
confusion_matrix <- cm$table
total <- sum(confusion_matrix)
number_of_classes <- nrow(confusion_matrix)
correct <- diag(confusion_matrix)
# sum all columns
total_actual_class <- apply(confusion_matrix, 2, sum)
# sum all rows
total_pred_class <- apply(confusion_matrix, 1, sum)
# Precision = TP / all that were predicted as positive
precision <- correct / total_pred_class
# Recall = TP / all that were actually positive
recall <- correct / total_actual_class
# F1
f1 <- (2 * precision * recall) / (precision + recall)
# create data frame to output results
results <- data.frame(precision, recall, f1)
results


```

Removing features that have too low recall/precision

```{r feature removal}

# Problematic features?
results %>% filter(recall < 0.8 | precision < 0.8 | f1 < 0.8)

BNC2014Eval2 <- BNC2014Eval %>% 
  mutate(Tag = ifelse(Tag == "PHC", "CC", as.character(Tag))) %>% 
  mutate(TagGold = ifelse(TagGold == "PHC", "CC", as.character(TagGold))) %>% 
  mutate(Tag = ifelse(Tag == "QLIKE", "LIKE", as.character(Tag))) %>% 
  mutate(TagGold = ifelse(TagGold == "QLIKE", "LIKE", as.character(TagGold))) %>% 
  mutate_if(is.character, as.factor) %>% 
  mutate(Evaluation = ifelse(as.character(Tag) == as.character(TagGold), TRUE, FALSE))

data <- BNC2014Eval2 %>% 
  filter(TagGold != "UNCLEAR") %>% 
  mutate(Tag = factor(Tag, levels = union(levels(Tag), levels(TagGold)))) %>% # Ensure that the factor levels are the same for the next caret operation
  mutate(TagGold = factor(TagGold, levels = union(levels(Tag), levels(TagGold))))

# Spot gold tag corrections that are not actually errors (should return nothing if all is well)
data[data$Tag==data$TagGold & data$Evaluation == FALSE,] %>% as.data.frame()

head(data) # Check sanity of data
summary(data) # Check sanity of data

#saveRDS(data, here("data", "MFTE_Evaluation_BNC2014_Results_merged.rds")) # Last saved 31 Oct 2021

#write.csv(data, here("data", "MFTE_Evaluation_BNC2014_Results_merged.csv")) # Last saved 31 Oct 2021

```

## Comparing tagger accuracy across different registers

````{r register-based-accuracy}

registerEval <- function(data, register) {
  d <- data %>% filter(Register==register)
  cm <- caret::confusionMatrix(d$Tag, d$TagGold) 
  return(cm$overall)
  #return(cm$byClass[,5:7])
}

registerEval(data, "internet")
registerEval(data, "news")
registerEval(data, "academic")
registerEval(data, "spoken")
registerEval(data, "fiction")

```

# Visualising tagger errors per register

```{r, fig.width = 8, fig.height = 8}
exclude_tags <- c("NULL", "UNCLEAR")
min_n <- 250
jitter_dist <- 0.2
opacity <- 0.3

data %>%
  filter(
    !(Tag %in% exclude_tags),
    !(TagGold %in% exclude_tags)
  ) %>%
  add_count(Tag, name = "n_tagged") %>%
  add_count(TagGold, name = "n_tagged_gold") %>%
  filter(
    n_tagged >= min_n,
    n_tagged_gold >= min_n
  ) ->
  data_filtered

tags_remaining <- union(
  unique(data_filtered$Tag),
  unique(data_filtered$TagGold)
)

data_filtered %>%
  mutate(
    Tag = factor(Tag, levels = tags_remaining),
    TagGold = factor(TagGold, levels = tags_remaining)
  ) ->
  data_filtered

data_filtered %>%
  ggplot(aes(x = TagGold, y = Tag)) +
  theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) +
  coord_fixed() +
  scale_x_discrete(drop = FALSE) +
  scale_y_discrete(drop = FALSE) +
  geom_abline(slope = 1, intercept = 0) +
  geom_jitter(
    width = jitter_dist,
    height = jitter_dist,
    alpha = opacity
  ) ->
  fig

for(i in unique(data_filtered$Register)){
  print((
    fig %+% filter(data_filtered, Register == i)) +
      ggtitle(i)
  )
}
```

## Comparing tagger accuracy across individual files

Though this is not very informative because the individual test files really are quite short.

````{r file-based-accuracy}

fileEval <- function(data, file) {
  d <- data %>% filter(FileID==file) %>% 
    # Ensure that the factor levels are the same for the next caret operation
    mutate(Tag = factor(Tag, levels = union(levels(Tag), levels(TagGold)))) %>% 
  mutate(TagGold = factor(TagGold, levels = union(levels(Tag), levels(TagGold))))
  cm <- caret::confusionMatrix(d$Tag, d$TagGold) 
  return(cm$overall)
  #return(cm$byClass[,5:7])
}

levels(data$FileID)

fileEval(data, "BNCBAcbH_m1")
fileEval(data, "BNCBESm3")
fileEval(data, "BNCBMass16")
fileEval(data, "BNCBEsocFb")
fileEval(data, "S2DD")
fileEval(data, "SZXQ")


```

## Compute accuracy metrics per feature

The three accuracy metrics are recall, precision and F1 score.

```{r accuracy-per-feature}

cm <- caret::confusionMatrix(data$Tag, data$TagGold) 
cm$overall
cm$byClass[,5:7]

confusion_matrix <- cm$table
total <- sum(confusion_matrix)
number_of_classes <- nrow(confusion_matrix)
correct <- diag(confusion_matrix)
total_actual_class <- apply(confusion_matrix, 2, sum)
total_pred_class <- apply(confusion_matrix, 1, sum)
# Precision = TP / all that were predicted as positive
precision <- correct / total_pred_class
# Recall = TP / all that were actually positive
recall <- correct / total_actual_class
# F1
f1 <- (2 * precision * recall) / (precision + recall)
# create data frame to output results
results <- data.frame(precision, recall, f1)
results

#write.csv(results, here("data", "MFTEAccuracyResults.csv"))

```

## Compute accuracy metrics with bootstrapping

This next chunk is based on the method and code presented in Picoral et al. (2021).

```{r bootstrapped-recall-precision-f1}
library(boot)
library(caret)

# Save a bit of computation time by filtering the data before we start.
exclude_tags <- c("NULL", "UNCLEAR", "PRP")
data %>%
  select(Tag, TagGold) %>%
  filter(
    !(Tag %in% exclude_tags),
    !(TagGold %in% exclude_tags)
  ) ->
  data_filtered

# Now unify the factor levels of the two relevant columns.
tags_remaining <- union(
  unique(data_filtered$Tag),
  unique(data_filtered$TagGold)
)
data_filtered %>%
  mutate(
    Tag = factor(Tag, levels = tags_remaining),
    TagGold = factor(TagGold, levels = tags_remaining)
  ) ->
  data_filtered

# Function for calculating the statistics.
# We can simplify this a bit, and make it return all three statistics.
get_measure_for_feature <- function(data, indices, measure, feature){
  data <- data[indices, ]
  confusion <- confusionMatrix(data$Tag, data$TagGold)
  statistics <- confusion$byClass
  return(statistics[paste("Class:", feature), measure])
}

# An example
example_results <- boot(
  data = data_filtered, 
  statistic = get_measure_for_feature,
  measure = "Recall",
  feature = "ABLE", 
  R = 10
) # This would have to be increased obviously but sticking to a low number for now to reduce waiting time.

print(example_results)
```

In the following chunk, the bootstrapping is applied to every combination of feature and measure. This will take a long time, so printout indicates the update us on the progress. However, I do not recommend running this code because it is incredibly slow and the {boot} library appears to have some weird bugs that cause various errors.

```{r bootstrapped-CI, eval=FALSE}
n_samples <- 10 # This would obviously have to be increased to 1000+
statistics <- c("Precision", "Recall", "F1")

# Get a dataframe ready.
all_results <- expand_grid(
  tag = tags_remaining,
  statistic = statistics,
  lower = NA,
  upper = NA
)

# Crunch the numbers painfully slowly in a loop.
for(row in 1:nrow(all_results)){
  
  # Find out what feature and measure we are working with this time.
  current_feature <- as.character(all_results[row, "tag"])
  current_measure <- as.character(all_results[row, "statistic"])
  
  # Make a progress printout.
  cat(current_measure, "for", current_feature, ":", n_samples, "samples\n")
  flush.console()
  
  # Filter out the irrelevant data to save a bit of time,
  # then hand on to the boot and ci functions.
  data_filtered %>%
    filter((Tag == current_feature) | (TagGold == current_feature)) %>%
    boot(
      statistic = get_measure_for_feature,
      measure = current_measure,
      feature = current_feature,
      R = n_samples
    ) %>%
    boot.ci(type = "perc") ->
    result
  
  # If we got a valid result, put it into the data frame.
  if(!is.null(result$perc)){
    all_results[row, "lower"] <- result$perc[4]
    all_results[row, "upper"] <- result$perc[5]
  }

}

all_results

```

Given that the above chunk proved too slow to run in R, the code was "translated" to run in python. Many thanks to Luke Tudge who did this conversion for me. The script is included in this project's repository and is called `Bootstrapped_Accuracy.ipynb`.

The results of the python script are plotted in the following chunk.

```{r plot-accuracy-CI}

resultsCI <- read.csv(here("data", "MFTE_Evaluation_BNC2014_CIs.csv")) # As computed in Bootstrapped_Accuracy.ipynb.
head(resultsCI)

resultsCI <- resultsCI %>% 
  mutate(tag = as.factor(tag)) %>% 
  filter(tag %in% c(str_extract(tag, "[A-Z0-9]+"))) %>% # Remove all punctuation tags which are uninteresting here.
  droplevels(.) %>% 
  mutate(metric = factor(metric, levels = c("precision", "recall", "f1")))

ggplot(resultsCI, aes(y = reorder(tag, desc(tag)), x = value, group = metric, colour = n)) +
  geom_point() +
  geom_errorbar(aes(xmin=lower, xmax = upper)) +
  ylab("") +
  xlab("") +
  facet_wrap(~ metric) +
  scale_color_paletteer_c("harrypotter::harrypotter", trans = "log", breaks = c(50,5000), labels = c(50,5000), name = "Number of tokens manually evaluated\n") +
  theme_bw() +
  theme(legend.position = "bottom")

ggsave(here("plots", "TaggerAccuracyResults95CI.svg"), width = 8, height = 12)

```

## Obtaining full list of errors

```{r errors}

# Adding an error tag with the incorrectly assigned tag and underscore and then the correct "gold" label
errors <- BNC2014Eval2 %>% 
  filter(Evaluation=="FALSE") %>% 
  filter(TagGold != "UNCLEAR") %>% 
  mutate(Error = paste(Tag, TagGold, sep = " -> ")) 

# Total number of errors
nrow(errors) # 1199

FreqErrors <- errors %>% 
  count(Error) %>% 
  arrange(desc(n)) 

FreqErrors %>% 
  #group_by(Register) %>% 
  filter(n > 9) %>% 
  print.data.frame() 

errors %>% 
  filter(Error == "NN -> JJAT") %>% 
  select(-Output, -Corpus, -Tag, -TagGold) %>% 
  filter(grepl(x = Token, pattern = "[A-Z]+.")) %>% 
  print.data.frame() 

errors %>% 
  filter(Error %in% c("NN -> VB", "VB -> NN", "NN -> VPRT", "VPRT -> NN")) %>% 
  count(Token) %>% 
  arrange(desc(n)) %>% 
  print.data.frame() 

errors %>% 
  filter(Error == "NN -> JJPR") %>% 
  count(Token) %>% 
  filter(grepl(x = Token, pattern = "[A-Z]+.")) %>% 
  arrange(desc(n)) %>% 
  print.data.frame() 

errors %>% 
  filter(Error == "ACT -> NULL") %>% 
  count(Token) %>% 
  arrange(desc(n)) %>% 
  print.data.frame() 



```