COP Recurrence and Emotions

Parkinson, and Elder groups

cat("\014")     # clean terminal



rm(list = ls()) # clean workspace
try(dev.off(), silent = TRUE) # close all plots
library(afex)
library(emmeans)
library(ggplot2)
library(ggridges)
library(ggdist)
library(dplyr)
library(reshape2)
library(GGally)
library(forcats)

theme_set(
  theme_minimal()
  )
a_posteriori <- function(afex_aov, sig_level = .05) {
  factors  <- as.list(rownames(afex_aov$anova_table))
  for (j in 1:length(factors)) {
    if (grepl(":", factors[[j]])) {
      factors[[j]] <- unlist(strsplit(factors[[j]], ":"))
    }
  }
  p_values <- afex_aov$anova_table$`Pr(>F)`
  for (i in 1:length(p_values)) {
    if (p_values[i] <= sig_level) {
      print(emmeans(afex_aov, factors[[i]], contr = "pairwise"))
    }
  }
}

master_dir     <- '~/Insync/OneDrive_shared/Fondecyt_Emociones_Estabilometria'
data_dir       <- paste(master_dir, 'data',  sep = '/')
emo_data_name  <- paste(data_dir, 'emo_data_clean_fix.csv', sep='/')
emo_data_clean_fix <- read.csv(emo_data_name, header = TRUE)
emo_data_clean_fix <- emo_data_clean_fix[(emo_data_clean_fix$Group == 'Elder' | emo_data_clean_fix$Group == 'Parkinson'), ]
emo_data_clean_fix$Group   <- factor(emo_data_clean_fix$Group, levels = c("Parkinson", "Elder"))
emo_data_clean_fix$Task    <- factor(emo_data_clean_fix$Task, levels = c("Unpleasant", "Neutral", "Pleasant"))
emo_data_clean_fix$Emotion <- factor(emo_data_clean_fix$Emotion)
emo_data_clean_fix$ID      <- factor(emo_data_clean_fix$ID)
emo_data_clean_fix$num_ID  <- factor(emo_data_clean_fix$num_ID)
outliers_x <- c()
outliers_y <- c()
rqa_x <- c('recurrence_rate_x', 'determinism_x', 'ave_diag_len_x', 'longest_diag_x', 'diag_entropy_x', 'laminarity_x',
           'trapping_time_x', 'longest_vertical_x', 'rec_time1_x', 'rec_time2_x', 'rec_per_dens_entr_x', 'clustering_x', 'transitivity_x')
rqa_y <- c('recurrence_rate_y', 'determinism_y', 'ave_diag_len_y', 'longest_diag_y', 'diag_entropy_y', 'laminarity_y',
           'trapping_time_y', 'longest_vertical_y', 'rec_time1_y', 'rec_time2_y', 'rec_per_dens_entr_y', 'clustering_y', 'transitivity_y')

1 RQA Measures

60 seconds

1.1 COP X

Delay: 23, by Mutual Information (MI) elbow

Embedding dimensions: 3, by False-Nearest Neighbors (FNN) elbow

options(width = 100)
rqa_x_pairs <- ggpairs(emo_data_clean_fix,
                       columns = rqa_x,
                       aes(colour = Group, alpha = .50),
                       progress = FALSE,
                       lower = list(continuous = wrap("points")))
suppressWarnings(print(rqa_x_pairs))

summary(emo_data_clean_fix[rqa_x])

 recurrence_rate_x determinism_x    ave_diag_len_x    longest_diag_x   diag_entropy_x 
 Min.   :0.02971   Min.   :0.9241   Min.   :  3.667   Min.   : 120.0   Min.   :1.728  
 1st Qu.:0.06680   1st Qu.:0.9994   1st Qu.: 23.338   1st Qu.: 288.0   1st Qu.:3.998  
 Median :0.09165   Median :0.9997   Median : 33.865   Median : 380.0   Median :4.393  
 Mean   :0.09689   Mean   :0.9991   Mean   : 39.207   Mean   : 475.3   Mean   :4.345  
 3rd Qu.:0.11710   3rd Qu.:0.9998   3rd Qu.: 47.400   3rd Qu.: 539.5   3rd Qu.:4.739  
 Max.   :0.77406   Max.   :1.0000   Max.   :630.694   Max.   :3406.0   Max.   :7.046  
  laminarity_x    trapping_time_x   longest_vertical_x  rec_time1_x      rec_time2_x    
 Min.   :0.9761   Min.   :  4.212   Min.   :  29.0     Min.   : 1.308   Min.   : 103.9  
 1st Qu.:0.9996   1st Qu.: 28.695   1st Qu.: 299.5     1st Qu.: 8.645   1st Qu.: 405.2  
 Median :0.9998   Median : 41.063   Median : 416.0     Median :11.045   Median : 473.8  
 Mean   :0.9995   Mean   : 48.032   Mean   : 492.8     Mean   :12.759   Mean   : 474.0  
 3rd Qu.:0.9999   3rd Qu.: 55.965   3rd Qu.: 592.5     3rd Qu.:15.039   3rd Qu.: 536.6  
 Max.   :1.0000   Max.   :804.363   Max.   :3410.0     Max.   :34.068   Max.   :1049.6  
 rec_per_dens_entr_x  clustering_x    transitivity_x  
 Min.   :0.4539      Min.   :0.4297   Min.   :0.4752  
 1st Qu.:0.6822      1st Qu.:0.5296   1st Qu.:0.5668  
 Median :0.7118      Median :0.5706   Median :0.6080  
 Mean   :0.7101      Mean   :0.5713   Mean   :0.6106  
 3rd Qu.:0.7455      3rd Qu.:0.6079   3rd Qu.:0.6457  
 Max.   :0.8038      Max.   :0.9388   Max.   :0.9502  

1.2 COP Y

Delay: 23, by Mutual Information (MI) elbow

Embedding dimensions: 3, by False-Nearest Neighbors (FNN) elbow

options(width = 100)
rqa_y_pairs <- ggpairs(emo_data_clean_fix,
                       columns = rqa_y,
                       aes(colour = Group, alpha = .50),
                       progress = FALSE,
                       lower = list(continuous = wrap("points")))
suppressWarnings(print(rqa_y_pairs))

summary(emo_data_clean_fix[rqa_y])

 recurrence_rate_y determinism_y    ave_diag_len_y    longest_diag_y   diag_entropy_y 
 Min.   :0.01939   Min.   :0.9822   Min.   :  4.874   Min.   :  72.0   Min.   :2.181  
 1st Qu.:0.05486   1st Qu.:0.9995   1st Qu.: 23.481   1st Qu.: 255.0   1st Qu.:4.009  
 Median :0.07329   Median :0.9997   Median : 31.659   Median : 329.0   Median :4.352  
 Mean   :0.07697   Mean   :0.9995   Mean   : 35.610   Mean   : 379.9   Mean   :4.335  
 3rd Qu.:0.09666   3rd Qu.:0.9999   3rd Qu.: 45.077   3rd Qu.: 442.5   3rd Qu.:4.740  
 Max.   :0.16422   Max.   :1.0000   Max.   :106.420   Max.   :1894.0   Max.   :5.607  
  laminarity_y    trapping_time_y  longest_vertical_y  rec_time1_y      rec_time2_y   
 Min.   :0.9930   Min.   :  6.59   Min.   :  95.0     Min.   : 6.165   Min.   :153.7  
 1st Qu.:0.9997   1st Qu.: 27.38   1st Qu.: 257.0     1st Qu.:10.460   1st Qu.:480.8  
 Median :0.9998   Median : 35.65   Median : 333.0     Median :13.813   Median :539.1  
 Mean   :0.9997   Mean   : 41.74   Mean   : 376.3     Mean   :15.478   Mean   :533.8  
 3rd Qu.:0.9999   3rd Qu.: 53.20   3rd Qu.: 462.0     3rd Qu.:18.444   3rd Qu.:596.3  
 Max.   :1.0000   Max.   :134.04   Max.   :1001.0     Max.   :52.205   Max.   :952.9  
 rec_per_dens_entr_y  clustering_y    transitivity_y  
 Min.   :0.5647      Min.   :0.4093   Min.   :0.4519  
 1st Qu.:0.7336      1st Qu.:0.5072   1st Qu.:0.5427  
 Median :0.7619      Median :0.5401   Median :0.5765  
 Mean   :0.7534      Mean   :0.5406   Mean   :0.5803  
 3rd Qu.:0.7813      3rd Qu.:0.5745   3rd Qu.:0.6157  
 Max.   :0.8300      Max.   :0.6890   Max.   :0.7352  

1.3 Diagonal Entropy X

options(width = 100)
diag_entropy_x_rep_anova <- aov_ez("ID", "diag_entropy_x", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data           <- diag_entropy_x_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

diag_entropy_x_rain <- ggplot(rep_anova_data, aes(y = Task, x = diag_entropy_x, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(diag_entropy_x_rain))

diag_entropy_x_afex_plot <-
  afex_plot(
    diag_entropy_x_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(diag_entropy_x_afex_plot))

nice(diag_entropy_x_rep_anova)

Anova Table (Type 3 tests)

Response: diag_entropy_x
      Effect           df  MSE      F   ges p.value
1      Group       1, 105 0.80   0.01 <.001    .920
2       Task 1.92, 201.73 0.21 2.76 +  .009    .068
3 Group:Task 1.92, 201.73 0.21   0.76  .002    .466
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(diag_entropy_x_rep_anova)

1.4 Recurrence time of 1st type X

options(width = 100)
rec_time1_x_rep_anova <- aov_ez("ID", "rec_time1_x", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data        <- rec_time1_x_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

rec_time1_x_rain <- ggplot(rep_anova_data, aes(y = Task, x = rec_time1_x, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(rec_time1_x_rain))

rec_time1_x_afex_plot <-
  afex_plot(
    rec_time1_x_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(rec_time1_x_afex_plot))

nice(rec_time1_x_rep_anova)

Anova Table (Type 3 tests)

Response: rec_time1_x
      Effect           df   MSE      F   ges p.value
1      Group       1, 105 61.04   0.02 <.001    .882
2       Task 2.00, 209.63 22.04 4.16 *  .016    .017
3 Group:Task 2.00, 209.63 22.04 2.78 +  .011    .064
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(rec_time1_x_rep_anova)

$emmeans
 Task       emmean    SE  df lower.CL upper.CL
 Unpleasant   12.4 0.541 105     11.3     13.5
 Neutral      12.1 0.536 105     11.0     13.1
 Pleasant     13.8 0.645 105     12.5     15.1

Results are averaged over the levels of: Group 
Confidence level used: 0.95 

$contrasts
 contrast              estimate    SE  df t.ratio p.value
 Unpleasant - Neutral     0.359 0.659 105   0.544  0.8495
 Unpleasant - Pleasant   -1.404 0.639 105  -2.197  0.0763
 Neutral - Pleasant      -1.763 0.639 105  -2.757  0.0187

Results are averaged over the levels of: Group 
P value adjustment: tukey method for comparing a family of 3 estimates 

1.5 Recurrence time of 2nd type X

options(width = 100)
rec_time2_x_rep_anova <- aov_ez("ID", "rec_time2_x", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data        <- rec_time2_x_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

rec_time2_x_rain <- ggplot(rep_anova_data, aes(y = Task, x = rec_time2_x, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(rec_time2_x_rain))

rec_time2_x_afex_plot <-
  afex_plot(
    rec_time2_x_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(rec_time2_x_afex_plot))

nice(rec_time2_x_rep_anova)

Anova Table (Type 3 tests)

Response: rec_time2_x
      Effect           df      MSE    F   ges p.value
1      Group       1, 105 31755.69 0.44  .003    .509
2       Task 1.95, 204.93  7130.88 0.55  .002    .571
3 Group:Task 1.95, 204.93  7130.88 0.30 <.001    .732
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(rec_time2_x_rep_anova)

1.6 Diagonal Entropy Y

options(width = 100)
diag_entropy_y_rep_anova <- aov_ez("ID", "diag_entropy_y", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data           <- diag_entropy_y_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

diag_entropy_y_rain <- ggplot(rep_anova_data, aes(y = Task, x = diag_entropy_y, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(diag_entropy_y_rain))

diag_entropy_y_afex_plot <-
  afex_plot(
    diag_entropy_y_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(diag_entropy_y_afex_plot))

nice(diag_entropy_y_rep_anova)

Anova Table (Type 3 tests)

Response: diag_entropy_y
      Effect           df  MSE      F   ges p.value
1      Group       1, 105 0.80   0.64  .005    .425
2       Task 1.93, 202.84 0.08 2.86 +  .004    .062
3 Group:Task 1.93, 202.84 0.08   0.32 <.001    .719
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(diag_entropy_y_rep_anova)

1.7 Recurrence time of 1st type Y

options(width = 100)
rec_time1_y_rep_anova <- aov_ez("ID", "rec_time1_y", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data        <- rec_time1_y_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

rec_time1_y_rain <- ggplot(rep_anova_data, aes(y = Task, x = rec_time1_y, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(rec_time1_y_rain))

rec_time1_y_afex_plot <-
  afex_plot(
    rec_time1_y_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(rec_time1_y_afex_plot))

nice(rec_time1_y_rep_anova)

Anova Table (Type 3 tests)

Response: rec_time1_y
      Effect           df    MSE       F   ges p.value
1      Group       1, 105 119.79    0.01 <.001    .924
2       Task 2.00, 209.57  14.13 6.93 **  .012    .001
3 Group:Task 2.00, 209.57  14.13    0.59  .001    .557
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(rec_time1_y_rep_anova)

$emmeans
 Task       emmean    SE  df lower.CL upper.CL
 Unpleasant   15.5 0.721 105     14.0     16.9
 Neutral      14.5 0.594 105     13.3     15.7
 Pleasant     16.4 0.729 105     15.0     17.9

Results are averaged over the levels of: Group 
Confidence level used: 0.95 

$contrasts
 contrast              estimate    SE  df t.ratio p.value
 Unpleasant - Neutral     0.984 0.529 105   1.860  0.1556
 Unpleasant - Pleasant   -0.942 0.512 105  -1.840  0.1616
 Neutral - Pleasant      -1.925 0.511 105  -3.767  0.0008

Results are averaged over the levels of: Group 
P value adjustment: tukey method for comparing a family of 3 estimates 

1.8 Recurrence time of 2nd type Y

options(width = 100)
rec_time2_y_rep_anova <- aov_ez("ID", "rec_time2_y", emo_data_clean_fix, within = c("Task"), between = c("Group"))

Warning: Missing values for 1 ID(s), which were removed before analysis:
FOSA_210
Below the first few rows (in wide format) of the removed cases with missing data.
  Contrasts set to contr.sum for the following variables: Group

rep_anova_data        <- rec_time2_y_rep_anova$data$long
xtabs(~ Task + Group, data = rep_anova_data)

            Group
Task         Parkinson Elder
  Unpleasant        60    47
  Neutral           60    47
  Pleasant          60    47

rec_time2_y_rain <- ggplot(rep_anova_data, aes(y = Task, x = rec_time2_y, color = Group, fill = Group)) +
  stat_halfeye(
    trim   = FALSE,
    .width = 0,
    justification = -.15,
    alpha  = .4,
    point_colour = NA) +
  geom_point(size = 2, alpha = .4, position = position_jitter(width = 0, height = .05))
suppressWarnings(print(rec_time2_y_rain))

rec_time2_y_afex_plot <-
  afex_plot(
    rec_time2_y_rep_anova,
    x = "Task",
    trace = "Group",
    error = "within",
    error_arg = list(width = .25),
    dodge = -.5,
    mapping = c("color"),
    point_arg = list(size = 4)
  )

Warning: Panel(s) show a mixed within-between-design.
Error bars do not allow comparisons across all means.
Suppress error bars with: error = "none"

suppressWarnings(print(rec_time2_y_afex_plot))

nice(rec_time2_y_rep_anova)

Anova Table (Type 3 tests)

Response: rec_time2_y
      Effect           df      MSE    F  ges p.value
1      Group       1, 105 25288.75 0.43 .003    .514
2       Task 1.87, 196.63  3943.24 0.94 .002    .388
3 Group:Task 1.87, 196.63  3943.24 0.60 .001    .541
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘+’ 0.1 ‘ ’ 1

Sphericity correction method: GG 

a_posteriori(rec_time2_y_rep_anova)