```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
library(reticulate)
use_virtualenv("~/.virtualenvs/r-reticulate")
#py_install("keras")
if (!require('tidyverse')) install.packages("tidyverse")            # data wrangling
if (!require('textstem')) install.packages("textstem")            # data wrangling
if (!require('text2vec')) install.packages("text2vec")              # word embeddings
if (!require('uwot')) install.packages("uwot") # dimensionality reduction
if (!require('umap')) install.packages("umap") # visualization
if (!require('conText')) install.packages("conText") # context-specific description and embedding regression
if (!require('proxy')) install.packages("proxy")                    # for sentence similarity
if (!require('caret')) install.packages("caret")                    # for evaluation metrics (machine learning)
if (!require('keras3')) install.packages("keras3")                  # for machine learning
if (!require('quanteda')) install.packages("quanteda")              # functions for text processing
if (!require("smotefamily")) install.packages("smotefamily")        # to handle imbalanced data -> generates synthetic data
#if (!require("torch")) install.packages("torch")
#install_torch("lantern")
```


# Word embeddings and party manifestos
## 1. Tidying text
This is an essential task for word embeddings as well since reducing the number of tokens increase computational efficiency and gives us more meaningful relationships between words. 

```{r}
df <- readRDS("./data_prep.RDS")
df <- subset(df, issue!="NA")
df$id <- 1:nrow(df)

# pre-process text: remove punctuation, numbering, superfluous whitespace, linebreaks 
df$text_prep <- tolower(df$sentence_context) %>% 
  paste0(" ", .) %>%
  str_replace_all(., "\\d+", " ") %>% 
  str_replace_all(., "[[:punct:]]", " ") %>% 
  str_replace_all(., "[[:symbol:]]", " ") %>%  
  str_replace_all(., "\n", " ") 

# tokenize corpus into unigrams
toks <- quanteda::tokens(df$text_prep) %>% tokens_remove(quanteda::stopwords(language="en"))
toks[[1]][1:5]

# further pruning of vocabulary -> only terms that appear at least 5 times 
# first dfm transformation to extract those features before removing others from tokenized object
features <- dfm(toks, verbose = FALSE) %>% dfm_trim(min_termfreq = 5) %>% featnames()

toks_sub <- tokens_select(toks, features, padding=T)


# add meta information to dfm 
docvars(toks, "text_id") <- df$id
docvars(toks, "country") <- df$country
docvars(toks, "issue") <- df$issue
docvars(toks, "parfam") <- df$parfam_n
docvars(toks, "year") <- df$year

# for embedding regression, let's add a variable that distinguishes between pre- and post 2015
df$prepost <- ifelse(df$year<2015, "pre", "post")
table(df$prepost)

docvars(toks, "prepost") <- df$prepost
```

## 2. Load pre-trained embeddings 

With embeddings, there are two different options. First, we can create our own embeddings locally. Here, we infer from the data we have how different words relate to each other. The main advantage of this approach is that our embeddings have in-domain knowledge. This may play a role, for instance, if some terms are used differently in other contexts. However, embeddings are most accurate if they are trained on large instances of data. This takes a lot of time and often requires much data. Therefore, we follow a second approach. Here, we download pre-trained embeddings from numberbatch (https://github.com/commonsense/conceptnet-numberbatch?tab=readme-ov-file) or GloVe (https://nlp.stanford.edu/projects/glove/) and retrieve them for our text.
```{r}
#path <- "./glove.6B/glove.6B.100d.txt" # use this on a less powerful machine
#path <- "./glove.42B.300d/glove.42B.300d.txt" # for GloVe embeddings

path <- "./numberbatch-en.txt" # for numberbatch, an ensemble embedding structure
embedding_vct <- data.table::fread(path, quote = "", skip=1, data.table = FALSE) %>% 
  as_tibble() # if using glove, remove skip=1 -> in numberbatch, the first line only indicate dimensionality

dim(embedding_vct)
head(embedding_vct)

embedding_vct[10000:10010,]
```

Every word in this embedding model is represented by a vector of 100/300 dimensions, other dimensionalities are available, too.

Do pre-trained embeddings contain all words of our MARPOR corpus?

Let us first tokenize the data. It's a bit different from before. We first tokenize but then create an iterator that goes through all tokens and stores them in a vocabulary object.
```{r}
# tokenize (use text2vec tokenizer here to create correct object type)
tokens <- word_tokenizer(df$text_prep)
# Create vocabulary. Terms will be unigrams (simple words).
it <- itoken(tokens)
vocab <- create_vocabulary(it, stopwords = c(stopwords(language="en"),"also", "s", "t", "d"))
vocab <- subset(vocab, term_count>5)
tail(vocab)
```

Lemmatization can help for later interpretability since similar words will be put into one vector.

```{r}
lemmatized_terms <- lemmatize_words(vocab$term)
vocab$term <- lemmatized_terms

vocab <- vocab %>%
  group_by(term) %>%
  summarise(term_count = sum(term_count, na.rm=T), 
            doc_count = sum(doc_count, na.rm=T))
```


Now, let's check which tokens are not part of pre-trained embeddings

```{r}
not_in_emb <- vocab %>%
  filter(!term %in% embedding_vct$V1)

print(paste0(nrow(not_in_emb), " words of our MARPOR corpus are not represented in the pre-trained embeddings."))
not_in_emb <- not_in_emb[order(not_in_emb$term_count, decreasing=T),]

print(paste0("The most frequently used words which are not captured by the pre-trained embeddings are ", paste0(head(not_in_emb$term), collapse=", ")))
```

As some of the terms are not represented in pre-trained corpus, we need to create our own matrix only with those terms that are represented.

```{r}
# fill matrix with embeddings -> embeddings that are not available, remain 0 in the matrix
# store V1 (terms) as row names
row_names <- embedding_vct %>%
  filter(V1 %in% vocab$term)  %>%
  select(V1)

# filter the embeddings to terms of our corpus
embeddings <- embedding_vct %>%
  filter(V1 %in% vocab$term) %>%
  select(-V1) %>%
  as.matrix()

# set terms as rownames
rownames(embeddings) <- row_names$V1

# add those words which do not exist in input embeddings with a 0
embeddings_na <- matrix(data = 0, nrow = nrow(not_in_emb), ncol = 300)

# set terms as rownames
rownames(embeddings_na) <- not_in_emb$term

# row bind available and not available embeddings
embeddings <- rbind(embeddings, embeddings_na)

saveRDS(embeddings, "./embeddings_mat.RDS")
```


Before doing analysis of our input data, let's do a very common illustration of the idea of word embeddings. 
Paris - France should equal which city - Germany? 
To execute this common example, we look for the cosine similarity of Germany given the aforementioned calculation.

```{r}
embeddings <- readRDS("./embeddings_mat.RDS")

which_capital = embeddings["paris", , drop = FALSE] -
  embeddings["france", , drop = FALSE] +
  embeddings["germany", , drop = FALSE]
capital_cos_sim = sim2(x = embeddings, y = which_capital, method = "cosine", norm = "l2")
head(sort(capital_cos_sim[,1], decreasing = TRUE), 10)
```
Yes, it's Berlin!

*Your Turn*
If you want, try out a similar logic yourself!

```{r}

```




```{r}
which_leaning = embeddings["democrat", , drop = FALSE] -
  embeddings["liberal", , drop = FALSE] +
  embeddings["conservative", , drop = FALSE]
leaning_cos_sim = sim2(x = embeddings, y = which_leaning, method = "cosine", norm = "l2")
head(sort(leaning_cos_sim[,1], decreasing = TRUE), 10)
```

We can explore their nearest neighbours. Which words have close embeddings to economy, which to migrants?

```{r}
embeddings <- readRDS("./embeddings_mat.RDS")


# which are close words to the term "migrant"?
find_nns(embeddings['migrant',], pre_trained = embeddings, N = 20)


# and which terms are close to "economy"?
find_nns(embeddings['economy',], pre_trained = embeddings, N = 10)
```

# Exploratory analysis 

We can also connect them with our dataset and look at which words are most often used in vicinity to a key term. For this purpose, we first create a so-called feature context matrix.

```{r}
# how are specific words used in these contexts?
# 1. create a feature context matrix
toks_fcm <- fcm(toks, context = "window", window = 10, count = "frequency", tri = FALSE) # important to set tri = FALSE
toks_fcm

# 2. create a transformation matrix (this gives important terms a high weight, terms like stopwords a low weight)
trans_mat <- compute_transform(x = toks_fcm, pre_trained = embeddings, weighting=500) # as weighting, set a high number

# 3. Let's further zero in on how the word "economy" is used usually
toks_context <- tokens_context(x = toks, pattern = "economic", window = 5L)
toks_context 

# 4. transform to document frequency matrix
context_dfm <- dfm(toks_context)
context_dfm

# 5. create document embedding matrix
context_dem <- dem(x = context_dfm, pre_trained = embeddings, transform = TRUE, transform_matrix = trans_mat, verbose = TRUE)
dim(context_dem)

# 6. now, we have several embeddings for our target word; to get a single estimate, aggregate column
context_wv <- colMeans(context_dem)
```


We can scrutinize how different actors are talking about a keyword of interest.
```{r}
context_wv_parfam <- dem_group(context_dem, groups = context_dem@docvars$parfam)
dim(context_wv_parfam) # 5 different party families


# nearest neighbours for different party families
context_nns <- nns(context_wv_parfam, pre_trained = embeddings, N = 10, candidates = context_wv_parfam@features, as_list = TRUE)
context_nns
```
That gives us a very nice descriptive idea of how different actors talk differently about a subject. The rad left is talking about "revolution", whereas the mainstream right emphasizes "growth". The Greens stress the environmental aspect. 


# visualization
Let's visualize relationships between words. As our embedding matrix has 300 dimensions, we first need to reduce dimensionality before we can visualize anything in a two-dimensional space. 

```{r}
# let's find the 50 nearest neighbors surrounding our target word
words <- find_nns(context_wv, pre_trained = embeddings, N = 50, candidates = context_dem@features)

# filter our vocabulary by these terms
vocab_sub <- vocab %>%
  filter(term %in% c(words)) 

# create data frame of embedding matrix
emb_df <- data.frame(embeddings)

# row names to column
emb_df <- emb_df %>% rownames_to_column()

# filter by 50 most frequent terms
emb_df_sub <- subset(emb_df, rowname %in% vocab_sub$term)

# convert back to matrix
emb_df_sub2 <- emb_df_sub %>% remove_rownames() %>% column_to_rownames("rowname") %>% as.matrix()


# now finally: visualizing
# reducing dimensionality 
visualization <- umap(emb_df_sub2, n_neighbors = 10, n_threads = 2) 


df_visualization  <- data.frame(word = rownames(emb_df_sub2),  
                  xpos = gsub(".+//", "", rownames(emb_df_sub2)),  
                  x = visualization$layout[, 1], y = visualization$layout[, 2],  
                  stringsAsFactors = FALSE) 

df_visualization %>%
  ggplot(aes(x,y, label=word, colour=word=="economic")) + 
  geom_point() + geom_text(hjust=0, vjust=0) + theme_light() +
  scale_color_manual(values=c("black", "red"), guide="none") 
```
We can learn quite something from this. Growth and innovation are quite close to each other, as social and environmental are. Both are on the opposite ends of the two-dimensional space, so there may be a trade-off between them. 

# Embedding regression
Let's first replicate the example of Rodriguez and Spirling to show how embedding regression (averaging between embeddings of a word by covariates) can help us to retrieve changing meaning of a word.
```{r}
set.seed(451)
model2 <- conText(formula = trump ~ prepost,
                 data = toks,
                 pre_trained = embeddings,
                 transform = TRUE, transform_matrix = trans_mat,
                 bootstrap = TRUE,
                 num_bootstraps = 100,
                 permute = TRUE, num_permutations = 10,
                 window = 10, 
                 verbose = T)

model2@normed_coefficients
```

We see: there is a significant difference in the word usage of "trump" pre and post the 2015. The substantive differences can be clarified by looking into the nearest neighbors before and after 2015. 

```{r}
post <- model2['(Intercept)',]
pre <- model2['(Intercept)',]  + model2['prepost_pre',]

# nearest neighbors
nns(rbind(post,pre), N = 40, pre_trained = embeddings, candidates = model2@features)
```

Before 2015, trump was mainly used as a verb, afterwards it's a description of the US president.


Sometimes, we want to make inference about the different usage of a word by specific groups. Following the example above, we know which terms different parties use when talking about the economy. But we do not know yet, whether there are substantial and significant differences in the words parties use. To get this single parameter, embedding regression breaks down the embeddings towards a focal word (here economy) to a single number and compares whether there are statistically significant difference. For instance, do the Greens use a significantly different vocabulary on the economy compared to the mainstream right? 

We've got all the ingredients we need already, let's feed them into conText, developed by Rodrigues et al. 2021. 

```{r}
set.seed(510)
model1 <- conText(formula = economic ~ parfam + year,
                 data = toks,
                 pre_trained = embeddings,
                 transform = TRUE, transform_matrix = trans_mat,
                 bootstrap = TRUE,
                 num_bootstraps = 100,
                 permute = TRUE, num_permutations = 10,
                 window = 10, 
                 verbose = T)


rownames(model1)
# (normed) coefficient table
model1@normed_coefficients

```

How have different party families talked about the economy?

```{r}
Greens_wv <- model1['(Intercept)',] + 2010*model1['year',]
MR_wv <- model1['(Intercept)',]  + model1['parfam_Mainstream Right',]  + 2010*model1['year',]
RL_wv <- model1['(Intercept)',]  + model1['parfam_Radical Left',]   + 2010*model1['year',]

# nearest neighbors
nns(rbind(MR_wv,RL_wv), N = 40, pre_trained = embeddings, candidates = model1@features)
```

# Document embeddings

```{r}
embeddings <- readRDS("./embeddings_mat.RDS")

text <- c("The President greets the press in Chicago", "Obama speaks to the media in Illinois", "The band gave a concert in Japan.", "Obama performed on stage in China.")
dfm_example <- quanteda::tokens(text) %>% dfm()

common_features <- intersect(colnames(dfm_example), rownames(embeddings))
dfm_example2 <- dfm_example[, common_features, drop = FALSE]
embeddings2 <- embeddings[common_features, , drop = FALSE]

text_mat <- as.matrix(dfm_example2)

# Matrix multiplication: dfm %*% embeddings
document_embeddings <- text_mat %*% embeddings2

# document similarity
simil(document_embeddings, method = "cosine")
```
This is, of course, only a toy example. We can also create an embedding matrix from our documents (similarly to the one we created above), and calculate the similarity of each sentence to each other.

```{r}
# feature co-occurrence matrix
toks_fcm <- fcm(toks, context = "window", window = 10, count = "frequency", tri = FALSE) # important to set tri = FALSE

# transformation matrix
trans_mat <- compute_transform(x = toks_fcm, pre_trained = embeddings, weighting=500) # as weighting, set a high number

# document feature matrix of our tokens
dfm_marpor <- dfm(toks)

# document-embedding matrix
dem_marpor <- dem(x = dfm_marpor, pre_trained = embeddings, transform = TRUE, transform_matrix = trans_mat, verbose = TRUE)


# consider, for example, the following statements
print(df$sentence_context[2])
print(df$sentence_context[25])
print(df$sentence_context[26])
print(df$sentence_context[27])
print(df$sentence_context[119])
print(df$sentence_context[120])

mat_marpor <- as.matrix(dem_marpor)

mat_marpor2 <- mat_marpor[c(2,25,26,27,119,120),]
simil(mat_marpor2, method="cosine")
```


# Supervised classification with embeddings

We are moving closer to a more accurate representation of text. In principle, we've got more information, so this should also help us in the predicting task we have started on Friday. Does it? 

Don't be afraid if you cannot follow every step here. This is meant to be a showcase session. At home, go through it line by line. Most of it will be easy to understand. 

The first part is just to ensure that we use exactly the same sample as yesterday. For your own application, this is not really important.
```{r}
embeddings <- readRDS("./embeddings_mat.RDS")

# Classification of so many categories is difficult, let's start simply. We want to distinguish statements on welfare from statements on the national way of life (immigration, patriotism, etc.)
df_sub <- subset(df, issue %in% c("504", "505", "601", "602"))

# for comparison, let's take the exact same ids as on Friday
df_sub2 <- readRDS("pred_df.RDS")  

df_sub <- subset(df_sub, id %in% df_sub2$id)
```


Now this is new: every input vector needs to be of the same size. There are two options to proceed: either we take the longest sequence as our default length or we take a length informed by summary measures (like the mean or median) as a benchmark. We decide for the latter since taking the longest sequence makes computation much slower. But be aware: this is a feature we could/should tweak if we want to improve performance.
```{r}
# count tokens
df_sub$count <- str_count(df_sub$text_prep, "\\w+")
summary(df_sub$count)

# on the basis of summary statistics, we need to decide how to prune sequences to n number of vectors -> here, we use the mean length of sequences but this is arbitrary and a potential hyperparameter
max_len <- round(median(df_sub$count,na.rm=T))
```


Now, we can create a new vocabulary (because we have a subset of our data frame by now) and prune our sequences to *max_len*.
```{r}
tokens <- word_tokenizer(df_sub$text_prep)
it <- itoken(tokens)
vocab <- create_vocabulary(it, stopwords = c(stopwords(language="en"),"also", "s", "t", "d"))
vocab <- subset(vocab, term_count>5 & term_count<20000)

# sequencing
vectorize_layer <- layer_text_vectorization(
  max_tokens = nrow(vocab),          
  output_sequence_length = max_len    
)
vectorize_layer %>% adapt(df_sub$text_prep)
features <- vectorize_layer(df_sub$text_prep)
features <- as.array(features)
head(features)
```


In machine learning, we want to train a model based on annotated data and evaluate its performance on test data. We always need to split our data into two subsets. 

**Note**: The best way is actually to split into three subsets: train, test and validation (left out) because train and test data can be mixed in the repeated training sequence.

```{r}
# create training data (let's use same IDs as yesterday)
index_train <- readRDS("./train_ids.RDS") %>% as.numeric
index_test <- readRDS("./test_ids.RDS")  %>% as.numeric

# extract features and labels (topics)
row_names_train <- as.numeric(row.names(df_sub)[with(df_sub, id %in% index_train)])
row_names_test <- as.numeric(row.names(df_sub)[with(df_sub, id %in% index_test)])


x_train <- features[row_names_train,]
x_test <- features[row_names_test,]

# define the labels of interest. Here, distinguish statements on welfare state from migration
df_sub$labels <- ifelse(df_sub$issue %in% c("504", "505"), 1, 0)
y_train <- subset(df_sub, id %in% index_train) 
y_train <- unlist(y_train[,"labels"])
y_test <- subset(df_sub, id %in% index_test) 
y_test <- unlist(y_test[,"labels"])
```


*Tuning your code*
In case you want to balance the data, there are different approaches. You could undersample majority classes, oversample minority classes or use SMOTE which creates synthetic minority classes.
Here, you'll see how SMOTE would work
```{r}
data_train <- data.frame(as.array(x_train), y_train)
smote_data <- SLS(data_train, target = data_train$y_train, K = 5, C = 5)
table(smote_data$data$y_train) # worked!

x_train_balanced <- as.matrix(smote_data$data[, 1:max_len])
y_train_balanced <- smote_data$data$y_train
```


Let's save the data for later usage (the transformer model).
```{r}
train <- data.frame(text_prepared = df_sub[df_sub$id %in% index_train,7], labels = y_train)
train$labels_text <- ifelse(train$labels==0, "migration", "welfare")
write.csv(train, "./data_transformers/training.csv")
test <- data.frame(text_prepared = df_sub[df_sub$id %in% index_test,7], labels = y_test)
test$labels_text <- ifelse(test$labels==0, "migration", "welfare")
write.csv(test, "./data_transformers/test.csv")
```


Now, we are at the point to define our model, to compile it, and to run it! 

A simple model (only using information from our embeddings)
```{r}
model <- keras_model_sequential() %>%
  layer_embedding(
    input_dim = dim(embeddings)[1],  # Assuming the number of words in your vocabulary matches the rows of t_mat
    input_length = max_len,
    output_dim = dim(embeddings)[2],  # Assuming the dimensionality of your embeddings
    weights = list(embeddings),
    trainable = F
  ) %>%
  layer_flatten() %>%
  layer_dense(units = 1, activation = "sigmoid")  

summary(model)

# Step 4: Compile your model
model %>% compile(
  loss = "binary_crossentropy",
  optimizer = optimizer_adam(learning_rate = 0.001),
  metrics = c('accuracy', metric_precision(), metric_recall())
)


stop_if_no_improvement <- callback_early_stopping(
  monitor = "val_loss",
  patience = 5,
  restore_best_weights = TRUE
)



# Step 5: Train your model
history <- model %>% fit(
  x = x_train,
  y = y_train,
  epochs = 20,
  batch_size = 128,
  validation_data = list(x_test, y_test),
  callbacks = list(stop_if_no_improvement)
)

print(paste0("Our best epoch is epoch ", (best_epoch <- which(history$metrics$val_loss == min(history$metrics$val_loss)))))
print(paste0("The smallest loss in our validation set is ", (loss <- history$metrics$val_loss[best_epoch] %>% round(3))))
print(paste0("The accuracy in our validation set is ", acc <- history$metrics$val_acc[best_epoch] %>% round(3)))
```
Okay, we see that there is quite a large error. The loss in our validation set is 0.23. Even the accuracy is not that great. 9% of our cases are wrongly predicted. That's more than it was on Friday with Naive Bayes!

If we look at the confusion matrix: 
```{r}
# predict
pred_results <- as.data.frame(predict(model, x_test))

# classify cases with probability > 0.5 as 1 
pred_results$pred <- ifelse(pred_results$V1>=0.5, 1, 0)

# combine predictions with true annotated data
pred_results$true <- y_test

confusionMatrix(as.factor(pred_results$pred), as.factor(pred_results$true))
```
That's not a good performance. Accuracy is ok but most sentences on migration were misclassified (, thus the low sensitivity).


The pre-trained embeddings treat each word separately and disregard structure sentence and inter-word relationships which might help to classify our data correctly. 

Therefore, we usually add hidden layers which can do that. Let's tune our deep learning model with these layers.
1. Dense layers provide a weight to each output.
2. Convolutional layers build n-grams of words which might actually help us to retrieve more contextual meaining. 
3. Our pre-trained embeddings can actually be fine-tuned by setting "trainable" to TRUE -> this makes them more domain specific. 
```{r}
# define model
model <- keras_model_sequential() %>%
  layer_embedding(
    input_dim = dim(embeddings)[1],  # Assuming the number of words in your vocabulary matches the rows of t_mat
    input_length = max_len,
    output_dim = dim(embeddings)[2],  # Assuming the dimensionality of your embeddings
    weights = list(embeddings),
    trainable = T
  ) %>%
  layer_conv_1d(filters = 128, kernel_size = 5, activation = 'relu') %>%
  layer_global_max_pooling_1d() %>%
  layer_dense(units = 128, activation = "relu") %>% # for dense layers between input and output, we often rely on "rectified linear units" that are non-linear and do not lead to vanishing gradient (weights become too small)
  layer_dropout(rate = 0.4) %>%  # dropout to reduce overfitting (nodes are randomly held out)
  layer_dense(units = 1, activation = "sigmoid")  

summary(model)

# Step 4: Compile your model
model %>% compile(
  loss = "binary_crossentropy",
  optimizer = optimizer_adam(learning_rate = 0.001),
  metrics = c('accuracy', metric_precision(), metric_recall())
)


stop_if_no_improvement <- callback_early_stopping(
  monitor = "val_loss",
  patience = 5,
  restore_best_weights = TRUE
)



# Step 5: Train your model
history <- model %>% fit(
  x = x_train,
  y = y_train,
  epochs = 10,
  batch_size = 128,
  validation_data = list(x_test, y_test),
  callbacks = list(stop_if_no_improvement)
)

print(paste0("Our best epoch is epoch ", (best_epoch <- which(history$metrics$val_loss == min(history$metrics$val_loss)))))
print(paste0("The smallest loss in our validation set is ", (loss <- history$metrics$val_loss[best_epoch] %>% round(3))))
print(paste0("The accuracy in our validation set is ", acc <- history$metrics$val_acc[best_epoch] %>% round(3)))
#saveRDS(model, "./model_e7b16d300.RDS")
```
Let's predict our test data with the more sophisticated deep learning model.
```{r}
# now predict probability for every speech
pred_results <- as.data.frame(predict(model, x_test))
```

...and show evaluation metrics.
```{r}
# transform probability into label: <0.5 -> 0, >0.5=1
pred_results$pred <- ifelse(pred_results$V1>=0.5, 1, 0)

# add true category, text and party
pred_results$true <- y_test
pred_results$issue_cmp <- df_sub$issue[index_test]
pred_results$sentence <- df_sub$sentence[index_test]
pred_results$party <- df_sub$partyname[index_test]

(confusion_matrix <- table(pred_results$pred, pred_results$true))

confusionMatrix(as.factor(pred_results$pred), as.factor(pred_results$true))
```
Now, the performance is way better than it was when using only quanteda's built-in functions. 

*Try it out*
Try building your own neural network. Can you beat the predictive quality?
```{r}

```




# Multi-class prediction 
```{r}
index_train <- readRDS("./train_ids.RDS") %>% as.numeric
index_test <- readRDS("./test_ids.RDS")  %>% as.numeric

# extract features and labels (topics)
row_names_train <- as.numeric(row.names(df_sub)[with(df_sub, id %in% index_train)])
row_names_test <- as.numeric(row.names(df_sub)[with(df_sub, id %in% index_test)])

x_train <- features[[row_names_train-1]]
x_test <- features[[row_names_test-1]]

# define the labels of interest. Here, distinguish statements on welfare state from migration
df_sub$labels[df_sub$issue=="504"] <- 0
df_sub$labels[df_sub$issue=="505"] <- 1
df_sub$labels[df_sub$issue=="601"] <- 2
df_sub$labels[df_sub$issue=="602"] <- 3

y_train <- subset(df_sub, id %in% index_train) 
y_train <- unlist(y_train[,"labels"]) 
y_test <- subset(df_sub, id %in% index_test) 
y_test <- unlist(y_test[,"labels"])

y_train <- to_categorical(y_train, num_classes = 4)
y_test <- to_categorical(y_test)

# define model
model2 <- keras_model_sequential() %>%
  layer_embedding(
    input_dim = dim(embeddings)[1],  # Assuming the number of words in your vocabulary matches the rows of t_mat
    input_length = max_len,
    output_dim = dim(embeddings)[2],  # Assuming the dimensionality of your embeddings
    weights = list(embeddings),
    trainable = T
  ) %>%
  layer_conv_1d(filters = 128, kernel_size = 5, activation = 'relu') %>%
  layer_global_max_pooling_1d() %>%
  layer_dense(units = 128, activation = "relu") %>% # for dense layers between input and output, we often rely on "rectified linear units" that are non-linear and do not lead to vanishing gradient (weights become too small)
  layer_dropout(rate = 0.4) %>%  # dropout to reduce overfitting (nodes are randomly held out)
  layer_dense(units = 4, activation = "softmax")  # this is the decisive parameter for multiclass prediction; often softmax is used as activation, the number of units must equal the number of categories

summary(model2)

# Step 4: Compile your model
model2 %>% compile(
  loss = "categorical_crossentropy", # for multiclass, we need to set a different loss function
  optimizer = optimizer_adam(learning_rate=0.0001),
  metrics = c("accuracy")
)


stop_if_no_improvement <- callback_early_stopping(
  monitor = "val_loss",
  patience = 5,
  restore_best_weights = TRUE
)


# Step 5: Train your model
history <- model2 %>% fit(
  x = x_train,
  y = y_train,
  epochs = 15,
  batch_size = 128,
  validation_data = list(x_test, y_test),
  callbacks = list(stop_if_no_improvement)
)

print(paste0("Our best epoch is epoch ", (best_epoch <- which(history$metrics$val_loss == min(history$metrics$val_loss)))))
print(paste0("The smallest loss in our validation set is ", (loss <- history$metrics$val_loss[best_epoch] %>% round(3))))
print(paste0("The accuracy in our validation set is ", acc <- history$metrics$val_acc[best_epoch] %>% round(3)))
```

```{r}
# now predict probability for every speech
pred_results2 <- as.data.frame(predict(model2, x_test))
colnames(pred_results2) <- c("504", "505", "601", "602")
pred_results2$pred <- colnames(pred_results2)[apply(pred_results2,1,which.max)]

# add true category, text and party
pred_results2$true <- df_sub$issue[df_sub$id %in% index_test]

confusionMatrix(as.factor(pred_results2$pred), as.factor(pred_results2$true))
```


## Addendum: Train your own local embeddings

In the following code, you'll see how you could train your own local embeddings. This is not always advisable. In some instances, the data is simply to small to retrieve (more) meaningful representations than pre-trained embeddings. Besides, you can also fine-tune the pre-trained embeddings by applying transfer-learning. 

Still, sometimes, it can be helpful (particularly if your text is very specific). 

First, we do a similar preparation as before when we transformed our words into a vocabulary. What we do add is a function to vectorize (transform our words into numbers) and create a tokens co-occurrence matrix. We did not do that before since we just loaded the numeric representation of pre-trained embeddings and connected it to our raw vocabulary. 
```{r}
tokens <- word_tokenizer(df$text_prep)
it <- itoken(tokens)
vocab <- create_vocabulary(it, stopwords = c(stopwords(language="en"),"also", "s", "t", "d"))

# now we can actually start to transform our tokens into vector representations
vectorizer <- vocab_vectorizer(vocab) # defines a function for vectorization

# define the term-co-occurrence matrix
## this step is crucial: we define the context window: the larger the window, the more context is captured (but the more computationally intense the process is)
tcm <- create_tcm(itoken(tokens), vectorizer, skip_grams_window = 10L)
```


Now, we can already estimate our own embeddings with glove. Here, we create a vocab*300 dimensional structure.
```{r}
glove = GlobalVectors$new(rank = 300, x_max = 10) 
wv_main <- glove$fit_transform(tcm, n_iter = 100, convergence_tol = 0.01)


# model learns two sets of word vectors - main and context -> main are the vectors for our target word whereas context are vectors for all words in context window (here n=10); convention: just average them
dim(wv_main)
wv_context <- glove$components
dim(wv_context)
word_vectors = wv_main + t(wv_context)
```


Let's do a short test whether it has worked.
```{r}
# similarity of words of interest
# which words appear often in relation to "economic"?
europe = word_vectors["economic", , drop = F]
cos_sim_europe = sim2(x = word_vectors, y = europe, method = "cosine", norm = "l2")
head(sort(cos_sim_europe[,1], decreasing = T), 10)
```


Does the equation between Paris, France and Germany still work?
```{r}
which_capital = word_vectors["paris", , drop = FALSE] -
  word_vectors["france", , drop = FALSE] +
  word_vectors["germany", , drop = FALSE]
capital_cos_sim = sim2(x = word_vectors, y = which_capital, method = "cosine", norm = "l2")
head(sort(capital_cos_sim[,1], decreasing = TRUE), 5)
```

Well, not really. And this could be due to the fact that our data has not enough examples of Berlin and Germany in one sentence. 