diff --git a/code/09_user-navigation.R b/code/09_user-navigation.R
index 453cd5b..f90410d 100644
--- a/code/09_user-navigation.R
+++ b/code/09_user-navigation.R
@@ -1,17 +1,16 @@
-# 09_user_navigation.R
+# 09_user-navigation.R
 #
 # content: (1) Read data
 #           (1.1) Read log event data
 #           (1.2) Extract additional infos for clustering
 #          (2) Clustering
 #          (3) Fit tree
-#          (3) Investigate variants
+#          (4) Investigate variants
 #
 # input:  results/haum/event_logfiles_2024-02-21_16-07-33.csv
-# output: results/haum/event_logfiles_pre-corona_with-clusters_cases.csv
-#         results/haum/dattree.csv
+# output: results/haum/eventlogs_pre-corona_case-clusters.csv
 #
-# last mod: 2024-02-27
+# last mod: 2024-03-06
 
 
 # setwd("C:/Users/nwickelmaier/Nextcloud/Documents/MDS/2023ss/60100_master_thesis/analysis/code")
@@ -23,33 +22,11 @@ library(factoextra)
 
 #--------------- (1.1) Read log event data ---------------
 
-dat0 <- read.table("results/haum/event_logfiles_2024-02-21_16-07-33.csv",
-                   colClasses = c("character", "character", "POSIXct",
-                                  "POSIXct", "character", "integer",
-                                  "numeric", "character", "character",
-                                  rep("numeric", 3), "character",
-                                  "character", rep("numeric", 11),
-                                  "character", "character"),
-                   sep = ";", header = TRUE)
+load("results/haum/eventlogs_pre-corona_cleaned.RData")
 
-dat0$event <- factor(dat0$event, levels = c("move", "flipCard", "openTopic",
-                                          "openPopup"))
-dat0$topic <- factor(dat0$topic)
-
-dat0$weekdays <- factor(weekdays(dat0$date.start),
-                        levels = c("Montag", "Dienstag", "Mittwoch",
-                                   "Donnerstag", "Freitag", "Samstag",
-                                   "Sonntag"),
-                        labels = c("Monday", "Tuesday", "Wednesday",
-                                   "Thursday", "Friday", "Saturday",
-                                   "Sunday"))
-
-
-# Select data pre Corona
-dat <- dat0[as.Date(dat0$date.start) < "2020-03-13", ]
-dat <- dat[dat$path != 106098, ]
-
-rm(dat0)
+# Select one year to handle number of cases
+dat <- dat[as.Date(dat$date.start) > "2018-12-31" &
+           as.Date(dat$date.start) < "2020-01-01", ]
 
 #--------------- (1.2) Extract additional infos for clustering ---------------
 
@@ -107,7 +84,7 @@ datcase$morning <- aggregate(date.start ~ case, dat,
 
 dat_split <- split(dat, ~ case)
 
-time_minmax <- function(subdata) {
+time_minmax_ms <- function(subdata) {
   subdata$min_time <- min(subdata$timeMs.start)
     if (all(is.na(subdata$timeMs.stop))) {
       subdata$max_time <- NA
@@ -116,9 +93,9 @@ time_minmax <- function(subdata) {
     }
   subdata
 }
-# TODO: Export from package mtt
+# TODO: Move to helper file
 
-dat_list <- pbapply::pblapply(dat_split, time_minmax)
+dat_list <- pbapply::pblapply(dat_split, time_minmax_ms)
 dat_minmax <- dplyr::bind_rows(dat_list)
 
 datcase$min_time <- aggregate(min_time ~ case, dat_minmax, unique)$min_time
@@ -160,8 +137,6 @@ datcase$thema <- ifelse(is.na(datcase$thema), 0, datcase$thema)
 datcase$ntopics <- ifelse(is.na(datcase$ntopics), 0, datcase$ntopics)
 datcase$ntopiccards <- ifelse(is.na(datcase$ntopiccards), 0, datcase$ntopiccards)
 
-
-
 cor_mat <- cor(datcase[, -1], use = "pairwise")
 diag(cor_mat) <- NA
 heatmap(cor_mat)
@@ -216,26 +191,24 @@ get_centrality <- function(case, data) {
   net <- process_map(alog, render = FALSE)
   inet <- DiagrammeR::to_igraph(net)
 
-  c(igraph::centr_degree(inet, loops = FALSE)$centralization,
-    igraph::centr_degree(inet, loops = TRUE)$centralization,
-    igraph::centr_betw(inet)$centralization)
+  #c(igraph::centr_degree(inet, loops = FALSE)$centralization,
+  #  igraph::centr_degree(inet, loops = TRUE)$centralization,
+  #  igraph::centr_betw(inet)$centralization)
+  igraph::centr_betw(inet)$centralization
 }
 # TODO: Move to helper file
 
-# centrality <- lapply(dattree$case, get_centrality, data = dat)
-# centrality <- do.call(rbind, centrality)
-#
+centrality <- pbapply::pblapply(dattree$case, get_centrality, data = dat)
+centrality <- do.call(rbind, centrality)
+
 # save(centrality, file = "results/haum/tmp_centrality.RData")
-load("results/haum/tmp_centrality.RData")
+#load("results/haum/tmp_centrality.RData")
 
-#dattree$DegreeCentrality <- centrality[, 2]
-dattree$BetweenCentrality <- centrality[, 3]
-
-## Add average duration per item
+dattree$BetweenCentrality <- unlist(centrality)
 
+# Average duration per item
 dat_split <- split(dat[, c("item", "case", "path", "timeMs.start", "timeMs.stop")], ~ path)
-
-dat_list <- pbapply::pblapply(dat_split, time_minmax)
+dat_list <- pbapply::pblapply(dat_split, time_minmax_ms)
 dat_minmax <- dplyr::bind_rows(dat_list)
 
 tmp <- aggregate(min_time ~ path, dat_minmax, unique)
@@ -244,10 +217,23 @@ tmp$duration <- tmp$max_time - tmp$min_time
 tmp$case <- aggregate(case ~ path, dat_minmax, unique)$case
 
 dattree$AvDurItem <- aggregate(duration ~ case, tmp, mean)$duration
-#dattree$AvDurItem <- dattree$AvDurItem / datcase$duration
 
 rm(tmp)
 
+# Indicator variable if table was used as info terminal only
+dattree$InfocardOnly <- factor(datcase$infocardOnly, levels = 0:1,
+                               labels = c("no", "yes"))
+
+# Add pattern to datcase; loosely based on Bousbia et al. (2009)
+dattree$Pattern <- "Dispersion"
+dattree$Pattern <- ifelse(dattree$PathLinearity > 0.8 & dattree$Singularity > 0.8, "Scholar",
+                          dattree$Pattern)
+dattree$Pattern <- ifelse(dattree$PathLinearity <= 0.8 &
+                          dattree$BetweenCentrality > 0.5, "Star",
+                          dattree$Pattern)
+dattree$Pattern <- factor(dattree$Pattern)
+
+
 summary(dattree)
 
 plot(dattree[, -1], pch = ".")
@@ -262,19 +248,6 @@ hist(dattree$PathLinearity, breaks = 50, main = "")
 hist(dattree$Singularity, breaks = 50, main = "")
 hist(dattree$BetweenCentrality, breaks = 50, main = "")
 
-# Indicator variable if table was used as info terminal only
-dattree$InfocardOnly <- factor(datcase$infocardOnly, levels = 0:1, labels = c("no", "yes"))
-
-# Add pattern to datcase; loosely based on Bousbia et al. (2009)
-dattree$Pattern <- "Dispersion"
-dattree$Pattern <- ifelse(dattree$PathLinearity > 0.8 & dattree$Singularity > 0.8, "Scholar",
-                          dattree$Pattern)
-dattree$Pattern <- ifelse(dattree$PathLinearity <= 0.8 &
-                          dattree$BetweenCentrality > 0.5, "Star",
-                          dattree$Pattern)
-dattree$Pattern <- factor(dattree$Pattern)
-
-
 # Remove cases with extreme outliers
 # TODO: Do I want this???
 
@@ -298,9 +271,9 @@ hist(dattree$BetweenCentrality, breaks = 50, main = "")
 
 #--------------- (2) Clustering ---------------
 
-library(cluster)
+#library(cluster)
 
-df <- dattree[1:10000, -1] # remove case variable
+#df <- dattree[, -1] # remove case variable
 # TODO: Do I need to scale or does normalization also work?
 
 # Normalize Duration and Numbers
@@ -312,42 +285,51 @@ df <- dattree[1:10000, -1] # remove case variable
 # summary(df)
 
 # Look at collinearity
-cor_mat <- cor(df)
-diag(cor_mat) <- NA
-heatmap(cor_mat)
+# cor_mat <- cor(df)
+# diag(cor_mat) <- NA
+# heatmap(cor_mat)
 
 #df <- as.data.frame(scale(dattree[, -1]))
 
 
 #--------------- (2.2) Hierarchical clustering ---------------
 
-mat <- daisy(df, metric = "gower")
+dist_mat <- cluster::daisy(dattree[, -1], metric = "gower")
+
+# # "Flatten" with PCA
+# mm <- model.matrix( ~ ., df)[, -1]   # remove intercept
+# tmp <- as.data.frame(lapply(as.data.frame(mm), normalize))
+# pc <- prcomp(mm)
+# coor_2d <- as.data.frame(pc$x[, 1:2])
+# coor_3d <- as.data.frame(pc$x[, 1:3])
 
 # "Flatten" with MDS
-coor_2d <- as.data.frame(cmdscale(mat, k = 2))
-coor_3d <- as.data.frame(cmdscale(mat, k = 3))
+coor_2d <- as.data.frame(cmdscale(dist_mat, k = 2))
+coor_3d <- as.data.frame(cmdscale(dist_mat, k = 3))
+
+# Idea from
+# https://stats.stackexchange.com/questions/264912/mds-on-large-dataset-r-or-python
+# https://www.inf.uni-konstanz.de/exalgo/software/mdsj/
+write.table(as.matrix(dist_mat), file = "mds/dist_mat.txt", row.names = FALSE,
+            col.names = FALSE)
+# Run java script
+system("java -jar mdsj.jar -d2 mds/dist_mat.txt mds/mds_coor_2d.txt")
+system("java -jar mdsj.jar -d3 mds/dist_mat.txt mds/mds_coor_3d.txt")
+
+coor_2d_java <- read.table("mds/mds_coor_2d.txt", header = FALSE, sep = " ")
+plot(coor_2d_java)
+
 
 plot(coor_2d)
 rgl::plot3d(coor_3d)
 
-#mat <- dist(df)
-
-# https://uc-r.github.io/hc_clustering
 method <- c(average = "average", single = "single", complete = "complete",
-            ward = "ward.D2")
+            ward = "ward")
 
-hc_method <- function(x) {
-  hclust(mat, method = x)
-}
+hcs <- pbapply::pblapply(method, function(x) cluster::agnes(dist_mat, method = x))
+acs <- pbapply::sapply(hcs, function(x) x$ac)
 
-hcs <- lapply(method, hc_method)
-cds <- lapply(hcs, cophenetic)
-cors <- sapply(cds, cor, y = mat)
-
-# https://en.wikipedia.org/wiki/Cophenetic_correlation
-# https://stats.stackexchange.com/questions/195446/choosing-the-right-linkage-method-for-hierarchical-clustering
-
-hc <- hcs$average
+hc <- hcs$ward
 
 # Something like a scree plot (??)
 plot(rev(hc$height)[1:100], type = "b", pch = 16, cex = .5)
@@ -356,72 +338,31 @@ k <- 4
 
 mycols <- c("#78004B", "#FF6900", "#3CB4DC", "#91C86E")
 
-grp_hclust <- cutree(hc, k = k)
+cluster <- cutree(as.hclust(hc), k = k)
 
-table(grp_hclust)
+table(cluster)
 
-fviz_cluster(list(data = df, cluster = grp_hclust),
-             palette = mycols,
-             ellipse.type = "convex",
-             show.clust.cent = FALSE,
-             ggtheme = theme_bw())
-
-plot(coor_2d, col = mycols[grp_hclust])
+plot(coor_2d, col = mycols[cluster])
 legend("topleft", paste("Cl", 1:4), col = mycols, pch = 21)
-rgl::plot3d(coor_3d, col = mycols[grp_hclust])
+rgl::plot3d(coor_3d, col = mycols[cluster])
 
-table(datcase[grp_hclust == 1, "Pattern"])
-table(datcase[grp_hclust == 2, "Pattern"])
-table(datcase[grp_hclust == 3, "Pattern"])
-table(datcase[grp_hclust == 4, "Pattern"])
+table(dattree[cluster == 1, "Pattern"])
+table(dattree[cluster == 2, "Pattern"])
+table(dattree[cluster == 3, "Pattern"])
+table(dattree[cluster == 4, "Pattern"])
 
-aggregate(. ~ grp_hclust, df, mean)
+aggregate(. ~ cluster, df, mean)
 
 aggregate(cbind(duration, distance, scaleSize, rotationDegree, length,
-                nmove, nflipCard, nopenTopic, nopenPopup) ~ grp_hclust, datcase,
+                nmove, nflipCard, nopenTopic, nopenPopup) ~ cluster, datcase,
           mean)
 
-#--------------- (2.3) DBSCAN clustering ---------------
-
-library(dbscan)
-d1 <- dbscan(df, eps = 1, minPts = ncol(df) + 1)
-hullplot(df, d1)
-
-grp_db <- d1$cluster
-table(grp_db)
-
-kNNdistplot(df, k = ncol(df))
-abline(h = 0.2, col = "red")
-abline(h = 1, col = "red")
-
-fviz_cluster(list(data = df[grp_db != 0, ], cluster = grp_db[grp_db != 0]),
-             palette = mycols,
-             ellipse.type = "convex",
-             show.clust.cent = FALSE,
-             ggtheme = theme_bw())
-
-mycols <- c("black", mycols)
-
-plot(coor_2d, col = mycols[grp_db + 1])
-legend("topleft", paste("Cl", 0:4), col = mycols, pch = 21)
-rgl::plot3d(coor_3d, col = mycols[grp_db + 1])
-
-aggregate(. ~ grp_db, df, mean)
-
-table(datcase[grp_db == 0, "Pattern"])
-table(datcase[grp_db == 1, "Pattern"])
-table(datcase[grp_db == 2, "Pattern"])
-table(datcase[grp_db == 3, "Pattern"])
-table(datcase[grp_db == 4, "Pattern"])
-
 ### Look at selected cases ###########################################
-dattree[grp_db == 0, ]
-
 tmp <- dat
 tmp$start <- tmp$date.start
 tmp$complete <- tmp$date.stop
 
-alog <- activitylog(tmp[tmp$case == 15, ],
+alog <- activitylog(tmp[tmp$case == 24016, ],
                     case_id     = "case",
                     activity_id = "item",
                     resource_id = "path",
@@ -433,30 +374,30 @@ rm(tmp)
 
 ######################################################################
 
-res <- merge(dat, data.frame(case = dattree$case, grp_km, grp_hclust, grp_db),
+res <- merge(dat, data.frame(case = dattree$case, cluster),
              by = "case", all.x = TRUE)
 res <- res[order(res$fileId.start, res$date.start, res$timeMs.start), ]
 
-xtabs( ~ item + grp_db, res)
-aggregate(event ~ grp_db, res, table)
+xtabs( ~ item + cluster, res)
+aggregate(event ~ cluster, res, table)
 
 # Look at clusters
 par(mfrow = c(2, 2))
-vioplot::vioplot(duration ~ grp_db, res)
-vioplot::vioplot(distance ~ grp_db, res)
-vioplot::vioplot(scaleSize ~ grp_db, res)
-vioplot::vioplot(rotationDegree ~ grp_db, res)
+vioplot::vioplot(duration ~ cluster, res)
+vioplot::vioplot(distance ~ cluster, res)
+vioplot::vioplot(scaleSize ~ cluster, res)
+vioplot::vioplot(rotationDegree ~ cluster, res)
 
-aggregate(cbind(duration, distance, scaleSize, rotationDegree) ~ grp_db, res, mean)
-aggregate(cbind(duration, distance, scaleSize, rotationDegree) ~ grp_db, res, median)
+aggregate(cbind(duration, distance, scaleSize, rotationDegree) ~ cluster, res, mean)
+aggregate(cbind(duration, distance, scaleSize, rotationDegree) ~ cluster, res, median)
 
 write.table(res,
-            file = "results/haum/event_logfiles_pre-corona_with-clusters_cases.csv",
+            file = "results/haum/eventlogs_2019_case-clusters.csv",
             sep = ";",
             quote = FALSE,
             row.names = FALSE)
 
-save(res, mat, h1, h2, h3, h4, h5, c1, c2, c3, c4, c5, datcase, dattree, df,
+save(res, dist_mat, hcs, acs, datcase, dattree, 
      file = "results/haum/tmp_user-navigation.RData")
 
 #--------------- (3) Fit tree ---------------
@@ -464,28 +405,13 @@ save(res, mat, h1, h2, h3, h4, h5, c1, c2, c3, c4, c5, datcase, dattree, df,
 library(rpart)
 library(partykit)
 
-dattree_db <- dattree[grp_db != 0, -1]
-dattree_db$grp <- factor(grp_db[grp_db != 0])
-
-c1 <- rpart(grp ~ ., data = dattree_db, method = "class")
+c1 <- rpart(as.factor(cluster) ~ ., data = dattree[, -1], method = "class")
 plot(as.party(c1))
 
-c2 <- rpart(as.factor(grp_hclust) ~ ., data = dattree[, -1], method = "class")
-plot(as.party(c2))
-
 # with conditional tree
-c2 <- ctree(grp ~ ., data = dattree_db, alpha = 0.05)
+c2 <- ctree(as.factor(cluster) ~ ., data = dattree[, -1], alpha = 0)
 plot(c2)
 
-# with excluded points
-c5 <- ctree(factor(grp_db) ~ ., data = dattree[, -1], alpha = 0)
-plot(c5)
-
-# with excluded points
-c6 <- ctree(factor(grp_db) ~ ., data = df, alpha = 0)
-plot(c6)
-# --> just checking
-
 #--------------- (4) Investigate variants ---------------
 
 res$start <- res$date.start
@@ -501,10 +427,10 @@ trace_explorer(alog, n_traces = 25)
 # --> sequences of artworks are just too rare
 
 tr <- traces(alog)
-trace_length <- sapply(strsplit(tr$trace, ","), length)
+trace_length <- pbapply::pbsapply(strsplit(tr$trace, ","), length)
 tr[trace_length > 10, ]
 
-trace_varied <- sapply(strsplit(tr$trace, ","), function(x) length(unique(x)))
+trace_varied <- pbapply::pbsapply(strsplit(tr$trace, ","), function(x) length(unique(x)))
 tr[trace_varied > 1, ]
 table(tr[trace_varied > 2, "absolute_frequency"])
 table(tr[trace_varied > 3, "absolute_frequency"])
@@ -528,7 +454,7 @@ tr[trace_varied == 5 & trace_length > 50, ]
 # --> every variant exists only once, of course
 datcase[datcase$nitems == 5 & datcase$length > 50,]
 
-sapply(datcase[, -c(1, 9)], median)
+pbapply::pbsapply(datcase[, -c(1, 9)], median)
 
 #ex <- datcase[datcase$nitems == 4 & datcase$length == 15,]
 ex <- datcase[datcase$nitems == 5,]
@@ -569,32 +495,3 @@ for (case in cases) {
 
 }
 
-
-
-
-########################### TODO: Still need it?
-
-
-net <- process_map(alog, render = FALSE)
-#DiagrammeR::get_node_df(net)
-
-DiagrammeR::get_node_info(net)
-
-DiagrammeR::get_degree_distribution(net)
-
-DiagrammeR::get_degree_in(net)
-DiagrammeR::get_degree_out(net)
-DiagrammeR::get_degree_total(net)
-
-
-N <- DiagrammeR::count_nodes(net) - 2   # Do not count start and stop nodes
-
-dc <- DiagrammeR::get_degree_total(net)[1:N, "total_degree"] / (N - 1)
-
-inet <- DiagrammeR::to_igraph(net)
-igraph::centr_degree(inet, loops = FALSE)
-igraph::centr_betw(inet)
-igraph::centr_clo(inet)
-
-
-