Project - Memory of pain (2018) – Data Analysis
Luis Anunciacao
June 26, 2021
#Global options
knitr::opts_chunk$set(echo = TRUE,
message = FALSE,
warning = FALSE,
include = TRUE, #<- show codes
cache = FALSE)
knitr::opts_knit$set(progress = TRUE, verbose = TRUE)
# auto format (kable)
options(kableExtra.auto_format = FALSE)Note These chunks present all analytical routine used in the manuscript “The Lighter Side of Pain: Do Positive Affective States Predict Memory of Pain Induced by Running a Marathon?”. You can reach out to me at luisfca@puc-rio.br
Thank you.
Last update: November 19, 2021
1 Data set
#get file
load(url("https://osf.io/z2xu8/download"))2 Data analysis
3 Descriptive results
3.1 Participants
view(summarytools::dfSummary(ds))3.2 Table 1
#sex 1 female 2 male
#run alone 1 yes 2 no
ds %>%
#filter(!is.na(pain_t2)) %>%
select(sex, age,run_alone, athlete, music, first_marathon, sum_positive, sum_negative, pain_t1, pain_t2) %>%
tableby(fct_rev(sex) ~ .,., control = tableby.control(cat.stats=c("Nmiss", "countpct"))) %>% #reverse factor order because of the table 1
summary() #fucking table with no percentage for missing... ##
##
## | | male (N=63) | female (N=45) | Total (N=108) | p value|
## |:---------------------------|:---------------:|:---------------:|:---------------:|-------:|
## |**age** | | | | 0.333|
## | N-Miss | 16 | 16 | 32 | |
## | Mean (SD) | 41.213 (9.851) | 39.172 (6.934) | 40.434 (8.858) | |
## | Range | 24.000 - 69.000 | 27.000 - 53.000 | 24.000 - 69.000 | |
## |**run_alone** | | | | 0.846|
## | N-Miss | 12 | 6 | 18 | |
## | no | 18 (35.3%) | 13 (33.3%) | 31 (34.4%) | |
## | yes | 33 (64.7%) | 26 (66.7%) | 59 (65.6%) | |
## |**athlete** | | | | 0.679|
## | N-Miss | 9 | 4 | 13 | |
## | no | 10 (18.5%) | 9 (22.0%) | 19 (20.0%) | |
## | yes | 44 (81.5%) | 32 (78.0%) | 76 (80.0%) | |
## |**music** | | | | 0.427|
## | N-Miss | 24 | 12 | 36 | |
## | no | 26 (66.7%) | 19 (57.6%) | 45 (62.5%) | |
## | yes | 13 (33.3%) | 14 (42.4%) | 27 (37.5%) | |
## |**first_marathon** | | | | 0.029|
## | N-Miss | 3 | 1 | 4 | |
## | no | 45 (75.0%) | 24 (54.5%) | 69 (66.3%) | |
## | yes | 15 (25.0%) | 20 (45.5%) | 35 (33.7%) | |
## |**sum_positive** | | | | 0.222|
## | Mean (SD) | 37.079 (10.258) | 39.511 (10.001) | 38.093 (10.177) | |
## | Range | 0.000 - 50.000 | 0.000 - 50.000 | 0.000 - 50.000 | |
## |**sum_negative** | | | | 0.463|
## | Mean (SD) | 11.968 (3.873) | 11.400 (4.064) | 11.731 (3.945) | |
## | Range | 0.000 - 28.000 | 0.000 - 31.000 | 0.000 - 31.000 | |
## |**pain_t1** | | | | 0.894|
## | Mean (SD) | 6.675 (2.773) | 6.600 (2.973) | 6.644 (2.844) | |
## | Range | 0.000 - 10.000 | 0.000 - 10.000 | 0.000 - 10.000 | |
## |**pain_t2** | | | | 0.205|
## | N-Miss | 31 | 21 | 52 | |
## | Mean (SD) | 5.406 (2.500) | 6.292 (2.629) | 5.786 (2.571) | |
## | Range | 0.000 - 10.000 | 0.000 - 10.000 | 0.000 - 10.000 | |ds %>%
#filter(!is.na(pain_t2)) %>%
mutate(sex = if_else(sex=="male","_male","female")) %>% #just for ordering
select(sex, age,athlete, run_alone, first_marathon, sum_positive, sum_negative, pain_t1, pain_t2) %>%
compareGroups::compareGroups(sex ~ ., data = ., include.miss = TRUE) %>% compareGroups::createTable()##
## --------Summary descriptives table by 'sex'---------
##
## _________________________________________________
## _male female p.overall
## N=63 N=45
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## age 41.2 (9.85) 39.2 (6.93) 0.294
## athlete: 0.639
## no 10 (15.9%) 9 (20.0%)
## yes 44 (69.8%) 32 (71.1%)
## 'Missing' 9 (14.3%) 4 (8.89%)
## run_alone: 0.721
## no 18 (28.6%) 13 (28.9%)
## yes 33 (52.4%) 26 (57.8%)
## 'Missing' 12 (19.0%) 6 (13.3%)
## first_marathon: 0.078
## no 45 (71.4%) 24 (53.3%)
## yes 15 (23.8%) 20 (44.4%)
## 'Missing' 3 (4.76%) 1 (2.22%)
## sum_positive 37.1 (10.3) 39.5 (10.0) 0.221
## sum_negative 12.0 (3.87) 11.4 (4.06) 0.467
## pain_t1 6.67 (2.77) 6.60 (2.97) 0.895
## pain_t2 5.41 (2.50) 6.29 (2.63) 0.209
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯ds %>%
filter(!is.na(pain_t2)) %>%
mutate(sex = if_else(sex=="male","_male","female")) %>% #just for ordering
select(sex, age, athlete, run_alone, first_marathon, sum_positive, sum_negative, pain_t1, pain_t2) %>%
compareGroups::compareGroups(sex ~ ., data = ., include.miss = TRUE) %>% compareGroups::createTable()##
## --------Summary descriptives table by 'sex'---------
##
## _________________________________________________
## _male female p.overall
## N=32 N=24
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯
## age 40.9 (10.1) 40.6 (8.11) 0.924
## athlete: 0.107
## no 5 (15.6%) 6 (25.0%)
## yes 22 (68.8%) 18 (75.0%)
## 'Missing' 5 (15.6%) 0 (0.00%)
## run_alone: 0.195
## no 11 (34.4%) 9 (37.5%)
## yes 14 (43.8%) 14 (58.3%)
## 'Missing' 7 (21.9%) 1 (4.17%)
## first_marathon: 0.151
## no 20 (62.5%) 12 (50.0%)
## yes 9 (28.1%) 12 (50.0%)
## 'Missing' 3 (9.38%) 0 (0.00%)
## sum_positive 37.7 (11.2) 41.4 (7.88) 0.154
## sum_negative 11.0 (3.73) 11.9 (2.69) 0.330
## pain_t1 6.70 (2.37) 7.17 (2.50) 0.485
## pain_t2 5.41 (2.50) 6.29 (2.63) 0.209
## ¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯Check difference in proportions
CrossTable(ds$sex)##
##
## Cell Contents
## |-------------------------|
## | N |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 108
##
##
## | female | male |
## |-----------|-----------|
## | 45 | 63 |
## | 0.417 | 0.583 |
## |-----------|-----------|
##
##
##
## chisq.test(table(ds$sex))##
## Chi-squared test for given probabilities
##
## data: table(ds$sex)
## X-squared = 3, df = 1, p-value = 0.08326Check difference in age
t.test(age ~ sex, var.equal = T, data = ds)##
## Two Sample t-test
##
## data: age by sex
## t = -0.97516, df = 74, p-value = 0.3327
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -6.209423 2.128719
## sample estimates:
## mean in group female mean in group male
## 39.17241 41.21277Check differences in run alone
CrossTable(ds$sex, ds$run_alone, chisq = T)##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 90
##
##
## | ds$run_alone
## ds$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 13 | 26 | 39 |
## | 0.014 | 0.007 | |
## | 0.333 | 0.667 | 0.433 |
## | 0.419 | 0.441 | |
## | 0.144 | 0.289 | |
## -------------|-----------|-----------|-----------|
## male | 18 | 33 | 51 |
## | 0.011 | 0.006 | |
## | 0.353 | 0.647 | 0.567 |
## | 0.581 | 0.559 | |
## | 0.200 | 0.367 | |
## -------------|-----------|-----------|-----------|
## Column Total | 31 | 59 | 90 |
## | 0.344 | 0.656 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.03762905 d.f. = 1 p = 0.8461899
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 2.151774e-30 d.f. = 1 p = 1
##
## Check differences in run alone
CrossTable(ds$sex, ds$athlete, chisq = T)##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 95
##
##
## | ds$athlete
## ds$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 9 | 32 | 41 |
## | 0.078 | 0.020 | |
## | 0.220 | 0.780 | 0.432 |
## | 0.474 | 0.421 | |
## | 0.095 | 0.337 | |
## -------------|-----------|-----------|-----------|
## male | 10 | 44 | 54 |
## | 0.059 | 0.015 | |
## | 0.185 | 0.815 | 0.568 |
## | 0.526 | 0.579 | |
## | 0.105 | 0.463 | |
## -------------|-----------|-----------|-----------|
## Column Total | 19 | 76 | 95 |
## | 0.200 | 0.800 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.171635 d.f. = 1 p = 0.6786628
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 0.02413618 d.f. = 1 p = 0.8765389
##
## Check differences in first marathon
CrossTable(ds$sex, ds$first_marathon, chisq = T)##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 104
##
##
## | ds$first_marathon
## ds$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 24 | 20 | 44 |
## | 0.924 | 1.821 | |
## | 0.545 | 0.455 | 0.423 |
## | 0.348 | 0.571 | |
## | 0.231 | 0.192 | |
## -------------|-----------|-----------|-----------|
## male | 45 | 15 | 60 |
## | 0.677 | 1.335 | |
## | 0.750 | 0.250 | 0.577 |
## | 0.652 | 0.429 | |
## | 0.433 | 0.144 | |
## -------------|-----------|-----------|-----------|
## Column Total | 69 | 35 | 104 |
## | 0.663 | 0.337 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 4.756635 d.f. = 1 p = 0.02918556
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 3.88465 d.f. = 1 p = 0.04872941
##
## Check differences in run listen to music
CrossTable(ds$sex, ds$music, chisq = T)##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 72
##
##
## | ds$music
## ds$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 19 | 14 | 33 |
## | 0.128 | 0.213 | |
## | 0.576 | 0.424 | 0.458 |
## | 0.422 | 0.519 | |
## | 0.264 | 0.194 | |
## -------------|-----------|-----------|-----------|
## male | 26 | 13 | 39 |
## | 0.108 | 0.181 | |
## | 0.667 | 0.333 | 0.542 |
## | 0.578 | 0.481 | |
## | 0.361 | 0.181 | |
## -------------|-----------|-----------|-----------|
## Column Total | 45 | 27 | 72 |
## | 0.625 | 0.375 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.630303 d.f. = 1 p = 0.4272442
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 0.3020979 d.f. = 1 p = 0.5825702
##
## 3.3 Participants at t2
ds_t2 %>% count(sex)## # A tibble: 2 x 2
## sex n
## <chr> <int>
## 1 female 24
## 2 male 33ds %>%
filter(!is.na(pain_t2)) %>%
select(sex, age,run_alone, athlete, music, first_marathon, sum_positive, sum_negative, pain_t1, pain_t2) %>%
arsenal::tableby(fct_rev(sex) ~ .,., control = arsenal::tableby.control(cat.stats=c("count"))) %>%
summary()##
##
## | | male (N=32) | female (N=24) | Total (N=56) | p value|
## |:---------------------------|:---------------:|:---------------:|:---------------:|-------:|
## |**age** | | | | 0.927|
## | N-Miss | 7 | 7 | 14 | |
## | Mean (SD) | 40.920 (10.140) | 40.647 (8.108) | 40.810 (9.266) | |
## | Range | 24.000 - 63.000 | 27.000 - 53.000 | 24.000 - 63.000 | |
## |**run_alone** | | | | 0.732|
## | no | 11 | 9 | 20 | |
## | yes | 14 | 14 | 28 | |
## |**athlete** | | | | 0.574|
## | no | 5 | 6 | 11 | |
## | yes | 22 | 18 | 40 | |
## |**music** | | | | 1.000|
## | no | 13 | 13 | 26 | |
## | yes | 7 | 7 | 14 | |
## |**first_marathon** | | | | 0.160|
## | no | 20 | 12 | 32 | |
## | yes | 9 | 12 | 21 | |
## |**sum_positive** | | | | 0.174|
## | Mean (SD) | 37.719 (11.229) | 41.417 (7.885) | 39.304 (10.023) | |
## | Range | 0.000 - 50.000 | 22.000 - 50.000 | 0.000 - 50.000 | |
## |**sum_negative** | | | | 0.352|
## | Mean (SD) | 11.031 (3.729) | 11.875 (2.692) | 11.393 (3.323) | |
## | Range | 0.000 - 17.000 | 10.000 - 21.000 | 0.000 - 21.000 | |
## |**pain_t1** | | | | 0.482|
## | Mean (SD) | 6.703 (2.365) | 7.167 (2.496) | 6.902 (2.411) | |
## | Range | 1.000 - 10.000 | 0.000 - 10.000 | 0.000 - 10.000 | |
## |**pain_t2** | | | | 0.205|
## | Mean (SD) | 5.406 (2.500) | 6.292 (2.629) | 5.786 (2.571) | |
## | Range | 0.000 - 10.000 | 0.000 - 10.000 | 0.000 - 10.000 | |Check difference in proportions
ds %>%
filter(!is.na(pain_t2)) %>%
{CrossTable(.$sex)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 56
##
##
## | female | male |
## |-----------|-----------|
## | 24 | 32 |
## | 0.429 | 0.571 |
## |-----------|-----------|
##
##
##
## ds %>%
filter(!is.na(pain_t2)) %>%
{chisq.test(table(.$sex))}##
## Chi-squared test for given probabilities
##
## data: table(.$sex)
## X-squared = 1.1429, df = 1, p-value = 0.285Check difference in age
ds %>%
filter(!is.na(pain_t2)) %>%
{t.test(age ~ sex, var.equal = T, data = .)}##
## Two Sample t-test
##
## data: age by sex
## t = -0.092558, df = 40, p-value = 0.9267
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -6.232813 5.686930
## sample estimates:
## mean in group female mean in group male
## 40.64706 40.92000Check differences in run alone
ds %>%
filter(!is.na(pain_t2)) %>%
{CrossTable(.$sex, .$run_alone, chisq = T)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 48
##
##
## | .$run_alone
## .$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 9 | 14 | 23 |
## | 0.036 | 0.025 | |
## | 0.391 | 0.609 | 0.479 |
## | 0.450 | 0.500 | |
## | 0.188 | 0.292 | |
## -------------|-----------|-----------|-----------|
## male | 11 | 14 | 25 |
## | 0.033 | 0.023 | |
## | 0.440 | 0.560 | 0.521 |
## | 0.550 | 0.500 | |
## | 0.229 | 0.292 | |
## -------------|-----------|-----------|-----------|
## Column Total | 20 | 28 | 48 |
## | 0.417 | 0.583 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.1168696 d.f. = 1 p = 0.7324548
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 0.002385093 d.f. = 1 p = 0.9610489
##
## Check differences in run alone
ds %>%
filter(!is.na(pain_t2)) %>%
{CrossTable(.$sex, .$athlete, chisq = T)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 51
##
##
## | .$athlete
## .$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 6 | 18 | 24 |
## | 0.131 | 0.036 | |
## | 0.250 | 0.750 | 0.471 |
## | 0.545 | 0.450 | |
## | 0.118 | 0.353 | |
## -------------|-----------|-----------|-----------|
## male | 5 | 22 | 27 |
## | 0.116 | 0.032 | |
## | 0.185 | 0.815 | 0.529 |
## | 0.455 | 0.550 | |
## | 0.098 | 0.431 | |
## -------------|-----------|-----------|-----------|
## Column Total | 11 | 40 | 51 |
## | 0.216 | 0.784 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.3155303 d.f. = 1 p = 0.5743062
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 0.04869792 d.f. = 1 p = 0.8253447
##
## check differences in listen to music
ds %>%
filter(!is.na(pain_t2)) %>%
{CrossTable(.$sex, .$music, chisq = T)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 40
##
##
## | .$music
## .$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 13 | 7 | 20 |
## | 0.000 | 0.000 | |
## | 0.650 | 0.350 | 0.500 |
## | 0.500 | 0.500 | |
## | 0.325 | 0.175 | |
## -------------|-----------|-----------|-----------|
## male | 13 | 7 | 20 |
## | 0.000 | 0.000 | |
## | 0.650 | 0.350 | 0.500 |
## | 0.500 | 0.500 | |
## | 0.325 | 0.175 | |
## -------------|-----------|-----------|-----------|
## Column Total | 26 | 14 | 40 |
## | 0.650 | 0.350 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0 d.f. = 1 p = 1
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0 d.f. = 1 p = 1
##
## Check differences in first marathon
ds %>%
filter(!is.na(pain_t2)) %>%
{CrossTable(.$sex, .$first_marathon, chisq = T)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 53
##
##
## | .$first_marathon
## .$sex | no | yes | Row Total |
## -------------|-----------|-----------|-----------|
## female | 12 | 12 | 24 |
## | 0.428 | 0.652 | |
## | 0.500 | 0.500 | 0.453 |
## | 0.375 | 0.571 | |
## | 0.226 | 0.226 | |
## -------------|-----------|-----------|-----------|
## male | 20 | 9 | 29 |
## | 0.354 | 0.540 | |
## | 0.690 | 0.310 | 0.547 |
## | 0.625 | 0.429 | |
## | 0.377 | 0.170 | |
## -------------|-----------|-----------|-----------|
## Column Total | 32 | 21 | 53 |
## | 0.604 | 0.396 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 1.974446 d.f. = 1 p = 0.1599768
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 1.261253 d.f. = 1 p = 0.261414
##
## 3.4 Reliability of NRS
Test retest reliability
cor.test(ds$pain_t1, ds$pain_t2, use = "complete.obs")##
## Pearson's product-moment correlation
##
## data: ds$pain_t1 and ds$pain_t2
## t = 4.4587, df = 54, p-value = 4.208e-05
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.2962362 0.6878335
## sample estimates:
## cor
## 0.5187372ds %>%
filter(!is.na(pain_t2)) %>%
select(pain_t1, pain_t2) %>%
psych::ICC(.)## Call: psych::ICC(x = .)
##
## Intraclass correlation coefficients
## type ICC F df1 df2 p lower bound upper bound
## Single_raters_absolute ICC1 0.45 2.6 55 56 2.0e-04 0.26 0.61
## Single_random_raters ICC2 0.47 3.1 55 55 1.9e-05 0.27 0.64
## Single_fixed_raters ICC3 0.52 3.1 55 55 1.9e-05 0.34 0.66
## Average_raters_absolute ICC1k 0.62 2.6 55 56 2.0e-04 0.41 0.76
## Average_random_raters ICC2k 0.64 3.1 55 55 1.9e-05 0.42 0.78
## Average_fixed_raters ICC3k 0.68 3.1 55 55 1.9e-05 0.50 0.80
##
## Number of subjects = 56 Number of Judges = 23.5 Reliability of the PANAS
#positive
ds %>%
select(productive, strong, enthusiastic, inspired, vigorous, stimulated, determined, powerfull, dynamic, interested) %>%
psych::alpha(.)##
## Reliability analysis
## Call: psych::alpha(x = .)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.88 0.88 0.9 0.43 7.6 0.017 4 0.77 0.41
##
## lower alpha upper 95% confidence boundaries
## 0.84 0.88 0.91
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## productive 0.87 0.88 0.89 0.44 7.2 0.019 0.013 0.42
## strong 0.86 0.87 0.88 0.42 6.5 0.020 0.011 0.40
## enthusiastic 0.86 0.87 0.87 0.42 6.6 0.020 0.008 0.40
## inspired 0.86 0.87 0.88 0.42 6.4 0.020 0.010 0.40
## vigorous 0.87 0.87 0.89 0.43 6.9 0.019 0.014 0.40
## stimulated 0.86 0.87 0.88 0.42 6.5 0.020 0.012 0.40
## determined 0.87 0.88 0.89 0.44 7.1 0.019 0.013 0.42
## powerfull 0.87 0.88 0.89 0.44 7.1 0.019 0.013 0.43
## dynamic 0.87 0.88 0.89 0.45 7.3 0.018 0.012 0.44
## interested 0.87 0.88 0.89 0.45 7.3 0.019 0.013 0.44
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## productive 105 0.65 0.65 0.59 0.55 3.5 1.18
## strong 104 0.77 0.76 0.75 0.68 3.9 1.24
## enthusiastic 102 0.74 0.76 0.76 0.68 4.2 0.97
## inspired 101 0.78 0.79 0.78 0.72 4.2 0.96
## vigorous 105 0.71 0.69 0.65 0.60 3.5 1.45
## stimulated 102 0.77 0.78 0.76 0.72 4.2 0.86
## determined 105 0.65 0.65 0.59 0.56 4.4 0.97
## powerfull 104 0.65 0.66 0.61 0.57 4.2 1.04
## dynamic 101 0.63 0.62 0.57 0.52 3.8 1.21
## interested 104 0.65 0.63 0.57 0.54 3.9 1.09
##
## Non missing response frequency for each item
## 1 2 3 4 5 miss
## productive 0.07 0.15 0.24 0.33 0.21 0.03
## strong 0.07 0.09 0.13 0.27 0.44 0.04
## enthusiastic 0.03 0.04 0.08 0.37 0.48 0.06
## inspired 0.01 0.07 0.10 0.32 0.50 0.06
## vigorous 0.15 0.11 0.12 0.27 0.34 0.03
## stimulated 0.01 0.05 0.07 0.43 0.44 0.06
## determined 0.05 0.01 0.04 0.33 0.57 0.03
## powerfull 0.02 0.08 0.12 0.28 0.51 0.04
## dynamic 0.06 0.12 0.16 0.33 0.34 0.06
## interested 0.05 0.07 0.16 0.40 0.32 0.04#negative
ds %>%
select(distressed, apreensivo, amedrontado, incomodado, angustiado, afraid, humiliated, irritable, scared, nervous) %>%
psych::alpha()##
## Reliability analysis
## Call: psych::alpha(x = .)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.76 0.78 0.82 0.26 3.5 0.033 1.2 0.36 0.22
##
## lower alpha upper 95% confidence boundaries
## 0.7 0.76 0.83
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## distressed 0.76 0.77 0.82 0.27 3.3 0.033 0.036 0.20
## apreensivo 0.74 0.75 0.80 0.25 3.0 0.037 0.030 0.22
## amedrontado 0.72 0.73 0.77 0.24 2.8 0.040 0.024 0.20
## incomodado 0.78 0.78 0.83 0.28 3.6 0.031 0.034 0.25
## angustiado 0.75 0.76 0.80 0.26 3.1 0.035 0.037 0.20
## afraid 0.74 0.75 0.79 0.25 3.1 0.037 0.026 0.22
## humiliated 0.77 0.78 0.82 0.28 3.5 0.034 0.029 0.22
## irritable 0.75 0.76 0.79 0.26 3.2 0.034 0.032 0.22
## scared 0.72 0.74 0.78 0.24 2.8 0.039 0.024 0.20
## nervous 0.72 0.73 0.77 0.23 2.7 0.039 0.031 0.19
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## distressed 105 0.50 0.49 0.39 0.32 1.3 0.70
## apreensivo 105 0.65 0.60 0.55 0.49 1.4 0.85
## amedrontado 104 0.76 0.71 0.72 0.63 1.2 0.62
## incomodado 104 0.46 0.40 0.27 0.24 1.4 0.77
## angustiado 105 0.55 0.58 0.51 0.42 1.1 0.53
## afraid 103 0.59 0.60 0.56 0.50 1.2 0.56
## humiliated 104 0.36 0.42 0.32 0.25 1.0 0.24
## irritable 104 0.52 0.55 0.51 0.37 1.2 0.55
## scared 105 0.69 0.68 0.68 0.60 1.2 0.62
## nervous 105 0.71 0.72 0.71 0.61 1.2 0.56
##
## Non missing response frequency for each item
## 1 2 3 4 5 miss
## distressed 0.80 0.12 0.05 0.03 0.00 0.03
## apreensivo 0.72 0.16 0.07 0.04 0.01 0.03
## amedrontado 0.91 0.06 0.01 0.00 0.02 0.04
## incomodado 0.73 0.20 0.03 0.03 0.01 0.04
## angustiado 0.91 0.05 0.02 0.02 0.00 0.03
## afraid 0.89 0.09 0.00 0.01 0.01 0.05
## humiliated 0.97 0.02 0.01 0.00 0.00 0.04
## irritable 0.88 0.09 0.02 0.02 0.00 0.04
## scared 0.92 0.04 0.01 0.02 0.01 0.03
## nervous 0.89 0.08 0.03 0.00 0.01 0.03ds %>%
select(distressed:interested) %>% #get items
psych::alpha(., check.keys = T)##
## Reliability analysis
## Call: psych::alpha(x = ., check.keys = T)
##
## raw_alpha std.alpha G6(smc) average_r S/N ase mean sd median_r
## 0.87 0.87 0.92 0.25 6.8 0.017 4.4 0.48 0.23
##
## lower alpha upper 95% confidence boundaries
## 0.83 0.87 0.9
##
## Reliability if an item is dropped:
## raw_alpha std.alpha G6(smc) average_r S/N alpha se var.r med.r
## distressed- 0.86 0.87 0.92 0.26 6.6 0.017 0.030 0.23
## productive 0.86 0.86 0.91 0.25 6.4 0.018 0.029 0.23
## apreensivo- 0.87 0.87 0.92 0.26 6.8 0.016 0.027 0.25
## amedrontado- 0.86 0.87 0.91 0.26 6.6 0.017 0.026 0.24
## incomodado- 0.87 0.88 0.92 0.27 7.0 0.017 0.028 0.26
## angustiado- 0.87 0.87 0.91 0.26 6.6 0.017 0.030 0.23
## strong 0.85 0.86 0.91 0.25 6.2 0.019 0.026 0.22
## enthusiastic 0.85 0.86 0.91 0.24 6.1 0.019 0.026 0.21
## inspired 0.85 0.86 0.91 0.24 6.2 0.019 0.025 0.23
## afraid- 0.87 0.87 0.92 0.26 6.7 0.017 0.027 0.23
## vigorous 0.86 0.87 0.91 0.25 6.5 0.018 0.027 0.23
## stimulated 0.85 0.86 0.91 0.24 6.0 0.019 0.027 0.20
## humiliated- 0.87 0.87 0.92 0.26 6.8 0.017 0.029 0.24
## determined 0.86 0.87 0.91 0.26 6.5 0.018 0.027 0.22
## irritable- 0.86 0.87 0.91 0.25 6.5 0.017 0.029 0.23
## powerfull 0.86 0.87 0.91 0.26 6.5 0.018 0.027 0.23
## dynamic 0.86 0.87 0.91 0.25 6.4 0.018 0.029 0.23
## scared- 0.87 0.87 0.91 0.26 6.6 0.017 0.026 0.23
## nervous- 0.86 0.86 0.91 0.25 6.3 0.018 0.029 0.23
## interested 0.86 0.87 0.92 0.25 6.5 0.018 0.029 0.23
##
## Item statistics
## n raw.r std.r r.cor r.drop mean sd
## distressed- 105 0.46 0.49 0.44 0.39 4.7 0.70
## productive 105 0.63 0.58 0.56 0.55 3.5 1.18
## apreensivo- 105 0.31 0.38 0.35 0.23 4.6 0.85
## amedrontado- 104 0.45 0.50 0.49 0.37 4.8 0.62
## incomodado- 104 0.28 0.30 0.24 0.18 4.6 0.77
## angustiado- 105 0.42 0.49 0.46 0.36 4.9 0.53
## strong 104 0.72 0.67 0.66 0.65 3.9 1.24
## enthusiastic 102 0.75 0.73 0.73 0.70 4.2 0.97
## inspired 101 0.74 0.69 0.70 0.69 4.2 0.96
## afraid- 103 0.37 0.45 0.42 0.33 4.8 0.56
## vigorous 105 0.62 0.54 0.51 0.51 3.5 1.45
## stimulated 102 0.77 0.76 0.76 0.73 4.2 0.86
## humiliated- 104 0.36 0.40 0.36 0.31 5.0 0.24
## determined 105 0.58 0.53 0.50 0.49 4.4 0.97
## irritable- 104 0.48 0.53 0.51 0.42 4.8 0.55
## powerfull 104 0.59 0.53 0.51 0.51 4.2 1.04
## dynamic 101 0.61 0.57 0.55 0.53 3.8 1.21
## scared- 105 0.40 0.49 0.48 0.36 4.8 0.62
## nervous- 105 0.56 0.63 0.62 0.51 4.8 0.56
## interested 104 0.60 0.54 0.51 0.52 3.9 1.09
##
## Non missing response frequency for each item
## 1 2 3 4 5 miss
## distressed 0.80 0.12 0.05 0.03 0.00 0.03
## productive 0.07 0.15 0.24 0.33 0.21 0.03
## apreensivo 0.72 0.16 0.07 0.04 0.01 0.03
## amedrontado 0.91 0.06 0.01 0.00 0.02 0.04
## incomodado 0.73 0.20 0.03 0.03 0.01 0.04
## angustiado 0.91 0.05 0.02 0.02 0.00 0.03
## strong 0.07 0.09 0.13 0.27 0.44 0.04
## enthusiastic 0.03 0.04 0.08 0.37 0.48 0.06
## inspired 0.01 0.07 0.10 0.32 0.50 0.06
## afraid 0.89 0.09 0.00 0.01 0.01 0.05
## vigorous 0.15 0.11 0.12 0.27 0.34 0.03
## stimulated 0.01 0.05 0.07 0.43 0.44 0.06
## humiliated 0.97 0.02 0.01 0.00 0.00 0.04
## determined 0.05 0.01 0.04 0.33 0.57 0.03
## irritable 0.88 0.09 0.02 0.02 0.00 0.04
## powerfull 0.02 0.08 0.12 0.28 0.51 0.04
## dynamic 0.06 0.12 0.16 0.33 0.34 0.06
## scared 0.92 0.04 0.01 0.02 0.01 0.03
## nervous 0.89 0.08 0.03 0.00 0.01 0.03
## interested 0.05 0.07 0.16 0.40 0.32 0.044 H1
4.1 Descriptive results
ds %>%
select(pain_t1, pain_t2) %>% descr()## Descriptive Statistics
## ds
## N: 108
##
## pain_t1 pain_t2
## ----------------- --------- ---------
## Mean 6.64 5.79
## Std.Dev 2.84 2.57
## Min 0.00 0.00
## Q1 5.00 4.00
## Median 7.00 6.00
## Q3 9.00 8.00
## Max 10.00 10.00
## MAD 2.97 2.97
## IQR 4.00 4.00
## CV 0.43 0.44
## Skewness -0.79 -0.55
## SE.Skewness 0.23 0.32
## Kurtosis -0.23 -0.26
## N.Valid 108.00 56.00
## Pct.Valid 100.00 51.854.2 Graph results
ds %>%
filter(!is.na(pain_t2)) %>%
select(id, email,pain_t1, pain_t2)%>%
pivot_longer(cols = -c(id,email),
names_to = "Time",
values_to = "Result")%>%
mutate(Time = if_else(Time == "pain_t1", "First measurement","Second measurement")) %>%
ggboxplot(., "Time", "Result") +
stat_compare_means(method = "t.test", label = "p.signif", paired = TRUE, var.equal = TRUE)
4.3 Table 2
In the manuscript, this is the table 2.
t.test(ds$pain_t1, ds$pain_t2,
alternative = "greater",
paired = T)##
## Paired t-test
##
## data: ds$pain_t1 and ds$pain_t2
## t = 3.4123, df = 55, p-value = 0.0006075
## alternative hypothesis: true difference in means is greater than 0
## 95 percent confidence interval:
## 0.5688612 Inf
## sample estimates:
## mean of the differences
## 1.116071ds %>%
filter(!is.na(pain_t2)) %>%
{ effsize::cohen.d(.$pain_t1, .$pain_t2, paired = TRUE)}##
## Cohen's d
##
## d estimate: 0.4473575 (small)
## 95 percent confidence interval:
## lower upper
## 0.1748532 0.7198618Graph
ds %>%
filter(!is.na(pain_t2)) %>%
select(id,pain_t1, pain_t2) %>%
pivot_longer(-id) %>%
mutate(name = as.factor(if_else(name == "pain_t1", "T1","T2"))) %>%
ggplot(., aes(x=name, y=value, group=1)) +
stat_summary(fun = mean, geom = "line") +
stat_summary(fun.data = mean_se, geom = "errorbar", width=0.2) +
theme_bw() +
ylim(6,9)
#https://stackoverflow.com/questions/10357768/plotting-lines-and-the-group-aesthetic-in-ggplot2
#https://kohske.wordpress.com/2010/12/27/faq-geom_line-doesnt-draw-lines/Negative affect. M = 11.73, SD = 3.94 Positive affect. M = 38.09, SD = 10.18
ds %>%
select(sum_positive, sum_negative) %>% descr()## Descriptive Statistics
## ds
## N: 108
##
## sum_negative sum_positive
## ----------------- -------------- --------------
## Mean 11.73 38.09
## Std.Dev 3.94 10.18
## Min 0.00 0.00
## Q1 10.00 34.00
## Median 11.00 40.00
## Q3 13.00 45.00
## Max 31.00 50.00
## MAD 1.48 8.90
## IQR 3.00 11.00
## CV 0.34 0.27
## Skewness 1.38 -1.61
## SE.Skewness 0.23 0.23
## Kurtosis 7.91 3.45
## N.Valid 108.00 108.00
## Pct.Valid 100.00 100.00aov(pain_t2 ~ pain_t1 + sex, ds) %>% summary()## Df Sum Sq Mean Sq F value Pr(>F)
## pain_t1 1 97.79 97.79 19.923 4.25e-05 ***
## sex 1 5.48 5.48 1.116 0.296
## Residuals 53 260.16 4.91
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 52 observations deleted due to missingnessTest the difference between positive and negative affect
t.test(ds$sum_positive, ds$sum_negative, paired = T)##
## Paired t-test
##
## data: ds$sum_positive and ds$sum_negative
## t = 25.868, df = 107, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 24.34093 28.38130
## sample estimates:
## mean of the differences
## 26.36111ds %>%
{ effsize::cohen.d(.$sum_positive, .$sum_negative, paired = TRUE)}##
## Cohen's d
##
## d estimate: 3.364 (large)
## 95 percent confidence interval:
## lower upper
## 2.702566 4.025435To check if this second group (T2) is virtually the same as the respondents in T1
ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
arrange(desc(dropout)) %>%
{t.test(sum_positive ~ dropout, var.equal = T, data = .)}##
## Two Sample t-test
##
## data: sum_positive by dropout
## t = 1.2873, df = 106, p-value = 0.2008
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -1.35854 6.38876
## sample estimates:
## mean in group n mean in group y
## 39.30357 36.78846ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
{t.test(sum_negative ~ dropout, var.equal = T, data = .)}##
## Two Sample t-test
##
## data: sum_negative by dropout
## t = -0.9251, df = 106, p-value = 0.357
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.2105510 0.8039576
## sample estimates:
## mean in group n mean in group y
## 11.39286 12.09615ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
{t.test(sum_negative ~ dropout, var.equal = T, data = .)}##
## Two Sample t-test
##
## data: sum_negative by dropout
## t = -0.9251, df = 106, p-value = 0.357
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.2105510 0.8039576
## sample estimates:
## mean in group n mean in group y
## 11.39286 12.09615ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
{t.test(age ~ dropout, var.equal = T, data = .)}##
## Two Sample t-test
##
## data: age by dropout
## t = 0.40825, df = 74, p-value = 0.6843
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.255624 4.933495
## sample estimates:
## mean in group n mean in group y
## 40.80952 39.97059ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
{gmodels::CrossTable(.$sex, .$dropout,chisq = T)}##
##
## Cell Contents
## |-------------------------|
## | N |
## | Chi-square contribution |
## | N / Row Total |
## | N / Col Total |
## | N / Table Total |
## |-------------------------|
##
##
## Total Observations in Table: 108
##
##
## | .$dropout
## .$sex | n | y | Row Total |
## -------------|-----------|-----------|-----------|
## female | 24 | 21 | 45 |
## | 0.019 | 0.021 | |
## | 0.533 | 0.467 | 0.417 |
## | 0.429 | 0.404 | |
## | 0.222 | 0.194 | |
## -------------|-----------|-----------|-----------|
## male | 32 | 31 | 63 |
## | 0.014 | 0.015 | |
## | 0.508 | 0.492 | 0.583 |
## | 0.571 | 0.596 | |
## | 0.296 | 0.287 | |
## -------------|-----------|-----------|-----------|
## Column Total | 56 | 52 | 108 |
## | 0.519 | 0.481 | |
## -------------|-----------|-----------|-----------|
##
##
## Statistics for All Table Factors
##
##
## Pearson's Chi-squared test
## ------------------------------------------------------------
## Chi^2 = 0.0678179 d.f. = 1 p = 0.7945408
##
## Pearson's Chi-squared test with Yates' continuity correction
## ------------------------------------------------------------
## Chi^2 = 0.004238619 d.f. = 1 p = 0.9480907
##
## ds %>%
mutate(dropout = if_else(is.na(pain_t2),"y","n")) %>%
{chisq.test(table(.$sex, .$dropout))}##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: table(.$sex, .$dropout)
## X-squared = 0.0042386, df = 1, p-value = 0.94815 H2
5.1 Distribution Plot
ds %>%
select(sum_positive, sum_negative) %>%
pivot_longer(everything()) %>%
ggplot(., aes(x = value, fill = name)) +
geom_density(alpha = .5)
5.2 Count plot
ds %>%
filter(!is.na(pain_t2)) %>%
ggplot(data=., aes(x=panas_group)) +
geom_bar() +
geom_text(stat='count', aes(label=..count..), vjust=-1)
ds %>% filter(!is.na(pain_t2)) %>% select(pain_t2,panas_group) %>%
arrange(panas_group) %>%
count(panas_group)## # A tibble: 2 x 2
## panas_group n
## <chr> <int>
## 1 negative or neutral 2
## 2 positive 545.3 Scatter plot
ds %>%
filter(!is.na(pain_t2)) %>%
select(id,pain_t1,pain_t2,panas_group) %>% #select target variables
pivot_longer(cols = -c(id,panas_group),
names_to = "Time",
values_to = "Result") %>% #to long format
mutate(panas_group = fct_recode(factor(panas_group), "Negative" = "0", "Positive" = "1"),
Time = fct_recode(factor(Time), "First measurement" = "pain_t1", "Second measurement" = "pain_t2")) %>% #transform panas group and time into factors and then re-arrange its names
ggplot(., aes(x = Time, y = Result, color = panas_group, group = panas_group)) +
stat_summary(geom = "line", fun = "mean") +
geom_smooth(method="lm", level=0.25) +
#ylim(0,10)
labs(color = "Affective group") +
theme_bw() + theme(legend.position = "bottom")
5.4 Individual plot
ds %>%
filter(!is.na(pain_t2)) %>%
select(id,pain_t1,pain_t2,panas_group) %>% #select target variables
pivot_longer(cols = -c(id,panas_group),
names_to = "Time",
values_to = "Result") %>% #to long format
mutate(panas_group = fct_recode(factor(panas_group), "Negative" = "0", "Positive" = "1"),
Time = fct_recode(factor(Time), "First measurement" = "pain_t1", "Second measurement" = "pain_t2")) %>% #transform panas group and time into factors and then re-arrange its names
ggplot(., aes(x = Time, y = Result, color = id, group = id, linetype = panas_group)) +
geom_line() + geom_point() + theme_bw() + scale_linetype(name = "Sex")
5.5 Table 3
In the manuscript, this is the table 3.
5.5.1 Pre-registered hypothesis
5.5.2 continuous predictor via linear mixed models (used in post review)
#https://web.stanford.edu/class/psych252/section/Mixed_models_tutorial.html
afex::set_sum_contrasts() #https://stats.stackexchange.com/questions/249884/conflicting-results-of-summary-and-anova-for-a-mixed-model-with-interactions
mod_lmm_ancova <- ds %>%
filter(!is.na(pain_t2)) %>%
mutate(id = as.factor(id)) %>%
select(id,pain_t1,pain_t2,sum_positive) %>% #select target variables
pivot_longer(cols = -c(id,sum_positive),
names_to = "Time",
values_to = "Result") %>%
lmer(Result ~ Time * sum_positive + (1 | id),.) #probabilistic modelanova(mod_lmm_ancova)## Type III Analysis of Variance Table with Satterthwaite's method
## Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
## Time 2.25656 2.25656 1 54 0.7397 0.3936
## sum_positive 2.93533 2.93533 1 54 0.9621 0.3310
## Time:sum_positive 0.00323 0.00323 1 54 0.0011 0.9742summary(mod_lmm_ancova)## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Result ~ Time * sum_positive + (1 | id)
## Data: .
##
## REML criterion at convergence: 516
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.1625 -0.5782 0.1266 0.5469 1.6345
##
## Random effects:
## Groups Name Variance Std.Dev.
## id (Intercept) 3.190 1.786
## Residual 3.051 1.747
## Number of obs: 112, groups: id, 56
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 7.4699997 1.1842963 53.9999999 6.308 5.46e-08 ***
## Time1 0.5792785 0.6735538 54.0000001 0.860 0.394
## sum_positive -0.0286551 0.0292134 53.9999999 -0.981 0.331
## Time1:sum_positive -0.0005405 0.0166148 54.0000001 -0.033 0.974
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Time1 sm_pst
## Time1 0.000
## sum_positiv -0.970 0.000
## Tm1:sm_pstv 0.000 -0.970 0.000MuMIn::r.squaredGLMM(mod_lmm_ancova)## R2m R2c
## [1,] 0.05966178 0.5403515.5.3 continuous predictor via repeated measures
afex_ancova <- ds %>%
filter(!is.na(pain_t2)) %>% #REPEATED ANOVA doest not work with empty cells
select(id,pain_t1,pain_t2,sum_positive) %>% #select target variables
pivot_longer(cols = -c(id,sum_positive), #use these variables as index
names_to = "Time",
values_to = "Result") %>%
ez.glm(id="id",
dv="Result",
data = .,
within=c("Time"),
covariate = "sum_positive",
factorize = FALSE)
afex_ancova## Anova Table (Type 3 tests)
##
## Response: Result
## Effect df MSE F ges p.value
## 1 sum_positive 1, 54 9.43 0.96 .013 .331
## 2 Time 1, 54 3.05 0.74 .003 .394
## 3 sum_positive:Time 1, 54 3.05 0.00 <.001 .974
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1ds %>%
filter(!is.na(pain_t2)) %>% #REPEATED ANOVA doest not work with empty cells
select(id,pain_t1,pain_t2,sum_positive) %>% #select target variables
pivot_longer(cols = -c(id,sum_positive), #use these variables as index
names_to = "Time",
values_to = "Result") %>%
write.csv(.,"ds_pain_long.csv", row.names = F)5.5.4 Categorical predictor via lmm (not used)
#https://web.stanford.edu/class/psych252/section/Mixed_models_tutorial.html
ds %>% #filter(!is.na(pain_t2)) %>%
mutate(id = as.factor(id)) %>%
select(id,pain_t1,pain_t2,panas_group) %>% #select target variables
pivot_longer(cols = -c(id,panas_group),
names_to = "Time",
values_to = "Result") %>%
lmer(Result ~ Time * panas_group + (1 | id),.) %>%
summary()## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: Result ~ Time * panas_group + (1 | id)
## Data: .
##
## REML criterion at convergence: 776.7
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.0549 -0.5027 0.1038 0.5516 1.5132
##
## Random effects:
## Groups Name Variance Std.Dev.
## id (Intercept) 4.589 2.142
## Residual 3.228 1.797
## Number of obs: 164, groups: id, 108
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) 6.7738 0.6512 122.1757 10.402 <2e-16 ***
## Time1 0.8005 0.4284 63.8755 1.869 0.0662 .
## panas_group1 0.7023 0.6512 122.1757 1.079 0.2829
## Time1:panas_group1 0.3234 0.4284 63.8755 0.755 0.4531
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Time1 pns_g1
## Time1 -0.354
## panas_grop1 0.916 -0.334
## Tm1:pns_gr1 -0.334 0.927 -0.354(In case someone wants to check a repeated measures anova)
5.5.5 categorical predictor via repeated measures
afex_rm_btw <- ds %>%
filter(!is.na(pain_t2)) %>% #REPEATED ANOVA doest not work with empty cells
select(id,pain_t1,pain_t2,panas_group) %>% #select target variables
pivot_longer(cols = -c(id,panas_group), #use these variables as index
names_to = "Time",
values_to = "Result") %>%
afex::aov_ez(
id = "id",
dv = "Result",
within = "Time",
between = "panas_group")
afex_rm_btw$Anova##
## Type III Repeated Measures MANOVA Tests: Pillai test statistic
## Df test stat approx F num Df den Df Pr(>F)
## (Intercept) 1 0.56982 71.528 1 54 1.815e-11 ***
## panas_group 1 0.00345 0.187 1 54 0.66737
## Time 1 0.05296 3.020 1 54 0.08796 .
## panas_group:Time 1 0.00492 0.267 1 54 0.60753
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 16 Exploratory analyzes
6.1 Babel (1)
pain intensity experienced, positive and negative affect as predictors, and recalled pain as a dependent variable.
lm(pain_t2 ~ pain_t1 + sum_positive + sum_negative, data = ds) %>%
olsrr::ols_regress()## Model Summary
## --------------------------------------------------------------
## R 0.531 RMSE 2.239
## R-Squared 0.282 Coef. Var 38.705
## Adj. R-Squared 0.241 MSE 5.015
## Pred R-Squared 0.136 MAE 1.804
## --------------------------------------------------------------
## RMSE: Root Mean Square Error
## MSE: Mean Square Error
## MAE: Mean Absolute Error
##
## ANOVA
## ------------------------------------------------------------------
## Sum of
## Squares DF Mean Square F Sig.
## ------------------------------------------------------------------
## Regression 102.662 3 34.221 6.824 6e-04
## Residual 260.767 52 5.015
## Total 363.429 55
## ------------------------------------------------------------------
##
## Parameter Estimates
## ----------------------------------------------------------------------------------------
## model Beta Std. Error Std. Beta t Sig lower upper
## ----------------------------------------------------------------------------------------
## (Intercept) 1.939 1.701 1.140 0.260 -1.475 5.353
## pain_t1 0.543 0.126 0.509 4.298 0.000 0.289 0.796
## sum_positive -0.023 0.033 -0.089 -0.703 0.485 -0.088 0.043
## sum_negative 0.088 0.098 0.114 0.900 0.372 -0.108 0.284
## ----------------------------------------------------------------------------------------6.2 Reviewer (2)
I miss results regarding if the participant were running alone or not, and if the participan regarded it self as an athlete
lm(pain_t2 ~ pain_t1 + run_alone + athlete, data = ds) %>%
olsrr::ols_regress()## Model Summary
## --------------------------------------------------------------
## R 0.557 RMSE 2.314
## R-Squared 0.310 Coef. Var 39.565
## Adj. R-Squared 0.261 MSE 5.353
## Pred R-Squared 0.163 MAE 1.781
## --------------------------------------------------------------
## RMSE: Root Mean Square Error
## MSE: Mean Square Error
## MAE: Mean Absolute Error
##
## ANOVA
## ------------------------------------------------------------------
## Sum of
## Squares DF Mean Square F Sig.
## ------------------------------------------------------------------
## Regression 101.100 3 33.700 6.295 0.0013
## Residual 224.835 42 5.353
## Total 325.935 45
## ------------------------------------------------------------------
##
## Parameter Estimates
## ---------------------------------------------------------------------------------------
## model Beta Std. Error Std. Beta t Sig lower upper
## ---------------------------------------------------------------------------------------
## (Intercept) 2.202 1.067 2.064 0.045 0.049 4.356
## pain_t1 0.570 0.140 0.539 4.073 0.000 0.288 0.853
## run_alone1 -0.171 0.356 -0.063 -0.482 0.633 -0.889 0.546
## athlete1 0.379 0.419 0.117 0.905 0.371 -0.466 1.224
## ---------------------------------------------------------------------------------------6.3 All predictors
ds %>%
select(pain_t2 , pain_t1 , sum_positive , sum_negative , sex , age , run_alone , athlete , first_marathon) %>%
str()## tibble [108 x 9] (S3: tbl_df/tbl/data.frame)
## $ pain_t2 : num [1:108] NA 8 7 NA 3 NA 4 NA 3 3 ...
## $ pain_t1 : num [1:108] 10 6 7 9 6 8 5 8 7 7 ...
## $ sum_positive : num [1:108] 36 34 40 16 37 32 42 31 32 0 ...
## $ sum_negative : num [1:108] 11 14 11 20 10 13 11 14 10 0 ...
## $ sex : Factor w/ 2 levels "female","male": 2 2 2 2 2 2 2 2 2 2 ...
## $ age : num [1:108] 69 63 58 56 55 54 51 51 50 50 ...
## $ run_alone : Factor w/ 2 levels "no","yes": 2 2 2 2 2 2 2 NA NA NA ...
## $ athlete : Factor w/ 2 levels "no","yes": 2 2 2 2 2 1 2 NA 2 NA ...
## $ first_marathon: Factor w/ 2 levels "no","yes": 1 1 1 1 1 2 2 1 1 NA ...6.4 Data mining
library(glmulti)
fit <- glmulti(pain_t2 ~ pain_t1 + sum_positive + sum_negative + sex + age + run_alone + athlete + first_marathon,
data = ds,
#xr = c("sex", "run_alone", "athlete", "first_marathon"),
level = 1, # No interaction considered
method = "g", # Exhaustive approach
crit = "aicc", # AIC as criteria
confsetsize = 5, # Keep 5 best models
plotty = F,
report = F, # No plot or interim reports
fitfunction = "lm", # lm function
marginality = T) #<-- Don't leave out the main effect## TASK: Genetic algorithm in the candidate set.
## Initialization...
## Algorithm started...
## Improvements in best and average IC have bebingo en below the specified goals.
## Algorithm is declared to have converged.
## Completed.summary(fit)## $name
## [1] "glmulti.analysis"
##
## $method
## [1] "g"
##
## $fitting
## [1] "lm"
##
## $crit
## [1] "aicc"
##
## $level
## [1] 1
##
## $marginality
## [1] TRUE
##
## $confsetsize
## [1] 5
##
## $bestic
## [1] 161.9397
##
## $icvalues
## [1] 161.9397 162.2118 163.8248 163.9958 164.0557
##
## $bestmodel
## [1] "pain_t2 ~ 1 + run_alone + first_marathon + pain_t1 + sum_positive + "
## [2] " sum_negative + age"
##
## $modelweights
## [1] 0.3370044 0.2941417 0.1313106 0.1205499 0.1169934
##
## $generations
## [1] 80
##
## $elapsed
## [1] 0.0102844
##
## $includeobjects
## [1] TRUEsummary(fit)$bestmodel## [1] "pain_t2 ~ 1 + run_alone + first_marathon + pain_t1 + sum_positive + "
## [2] " sum_negative + age"#Check models
top <- weightable(fit)
top <- top[top$aicc <= min(top$aicc) + 2,]
top## model
## 1 pain_t2 ~ 1 + run_alone + first_marathon + pain_t1 + sum_positive + sum_negative + age
## 2 pain_t2 ~ 1 + run_alone + athlete + first_marathon + pain_t1 + sum_positive + sum_negative + age
## 3 pain_t2 ~ 1 + sex + run_alone + athlete + first_marathon + pain_t1 + sum_positive + sum_negative + age
## aicc weights
## 1 161.9397 0.3370044
## 2 162.2118 0.2941417
## 3 163.8248 0.1313106summary(fit@objects[[1]])##
## Call:
## fitfunc(formula = as.formula(x), data = data)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2376 -1.1728 -0.1367 1.1670 3.9465
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.52671 4.27117 0.592 0.5589
## run_alone1 -0.73988 0.38017 -1.946 0.0617 .
## first_marathon1 -0.90557 0.41877 -2.162 0.0393 *
## pain_t1 0.71263 0.15161 4.700 6.29e-05 ***
## sum_positive -0.13932 0.06417 -2.171 0.0385 *
## sum_negative 0.41147 0.17854 2.305 0.0288 *
## age -0.01700 0.03600 -0.472 0.6404
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.01 on 28 degrees of freedom
## (73 observations deleted due to missingness)
## Multiple R-squared: 0.588, Adjusted R-squared: 0.4997
## F-statistic: 6.66 on 6 and 28 DF, p-value: 0.000185#http://www.metafor-project.org/doku.php/tips:model_selection_with_glmulti_and_mumin`6.5 Report
options('contrasts')options(contrasts=c('contr.treatment','contr.poly')) #attention!!!!!!!!!
mod_data_mining <- lm(pain_t2 ~ factor(run_alone) + factor(first_marathon) + pain_t1 + sum_positive + sum_negative + age , ds)olsrr::ols_regress(mod_data_mining)## Model Summary
## --------------------------------------------------------------
## R 0.767 RMSE 2.010
## R-Squared 0.588 Coef. Var 36.271
## Adj. R-Squared 0.500 MSE 4.042
## Pred R-Squared 0.366 MAE 1.498
## --------------------------------------------------------------
## RMSE: Root Mean Square Error
## MSE: Mean Square Error
## MAE: Mean Absolute Error
##
## ANOVA
## ------------------------------------------------------------------
## Sum of
## Squares DF Mean Square F Sig.
## ------------------------------------------------------------------
## Regression 161.511 6 26.918 6.66 2e-04
## Residual 113.175 28 4.042
## Total 274.686 34
## ------------------------------------------------------------------
##
## Parameter Estimates
## ------------------------------------------------------------------------------------------------------
## model Beta Std. Error Std. Beta t Sig lower upper
## ------------------------------------------------------------------------------------------------------
## (Intercept) 0.881 4.317 0.204 0.840 -7.963 9.725
## factor(run_alone)yes 1.480 0.760 0.255 1.946 0.062 -0.078 3.037
## factor(first_marathon)yes 1.811 0.838 0.307 2.162 0.039 0.096 3.527
## pain_t1 0.713 0.152 0.620 4.700 0.000 0.402 1.023
## sum_positive -0.139 0.064 -0.289 -2.171 0.039 -0.271 -0.008
## sum_negative 0.411 0.179 0.316 2.305 0.029 0.046 0.777
## age -0.017 0.036 -0.058 -0.472 0.640 -0.091 0.057
## ------------------------------------------------------------------------------------------------------#jtools::summ(mod_data_mining)
jtools::j_summ(mod_data_mining)| Observations | 35 (73 missing obs. deleted) |
| Dependent variable | pain_t2 |
| Type | OLS linear regression |
| F(6,28) | 6.66 |
| R² | 0.59 |
| Adj. R² | 0.50 |
| Est. | S.E. | t val. | p | |
|---|---|---|---|---|
| (Intercept) | 0.88 | 4.32 | 0.20 | 0.84 |
| factor(run_alone)yes | 1.48 | 0.76 | 1.95 | 0.06 |
| factor(first_marathon)yes | 1.81 | 0.84 | 2.16 | 0.04 |
| pain_t1 | 0.71 | 0.15 | 4.70 | 0.00 |
| sum_positive | -0.14 | 0.06 | -2.17 | 0.04 |
| sum_negative | 0.41 | 0.18 | 2.30 | 0.03 |
| age | -0.02 | 0.04 | -0.47 | 0.64 |
| Standard errors: OLS |
car::Anova(mod_data_mining, type = 3)## Anova Table (Type III tests)
##
## Response: pain_t2
## Sum Sq Df F value Pr(>F)
## (Intercept) 0.168 1 0.0417 0.83974
## factor(run_alone) 15.309 1 3.7875 0.06174 .
## factor(first_marathon) 18.901 1 4.6761 0.03928 *
## pain_t1 89.299 1 22.0931 6.291e-05 ***
## sum_positive 19.055 1 4.7144 0.03854 *
## sum_negative 21.468 1 5.3114 0.02881 *
## age 0.901 1 0.2230 0.64041
## Residuals 113.175 28
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 16.6 Individual plots
p1 <- ds %>%
filter(!is.na(pain_t2), !is.na(first_marathon)) %>%
ggplot(., aes(x = first_marathon, y = pain_t2)) +
geom_boxplot() +
ggsignif::geom_signif(
comparisons = list(c("no", "yes")),
annotation = c("p < 0.05")) +
labs(x = "First Marathon", y = "Recalled Pain") +
theme_minimal() +
scale_x_discrete(label = stringr::str_to_sentence(levels(ds$first_marathon)))p2 <- ds %>%
filter(!is.na(pain_t2)) %>%
ggplot(., aes(x = pain_t1, y = pain_t2)) +
#geom_jitter() +
geom_smooth(method = "lm", color = "black", fill = "navy") +
labs(x = "Experienced Pain", y = "") +
theme_minimal() + geom_rug()p3 <- ds %>%
filter(!is.na(pain_t2)) %>%
ggplot(., aes(x = sum_positive, y = pain_t2)) +
#geom_jitter() +
geom_smooth(method = "lm", color = "black", fill = "navy") +
labs(x = "Positive Affects", y = "Recalled Pain") +
theme_minimal() + geom_rug()p4 <- ds %>%
filter(!is.na(pain_t2)) %>%
ggplot(., aes(x = sum_negative, y = pain_t2)) +
#geom_jitter(aes(color = pain_t2)) +
geom_smooth(method = "lm", color = "black", fill = "navy") +
labs(x = "Negative Affects", y = "") +
theme_minimal() + geom_rug()library(patchwork)
patchwork <- p1 + p2 + p3 + p4 +
plot_layout(ncol = 2)
patchwork
#ggsave(filename = "dove_figure_1.tiff",
# plot = print(patchwork), device='tiff', dpi=300)6.7 All efects plot
plot(effects::allEffects( lm(pain_t2 ~ run_alone + first_marathon + pain_t1 + sum_positive + sum_negative + age, ds)))
Done (March 22, 2021) [Submission] Done (June 26, 2021) [After first review] Done (November 19, 2021 [After acceptation]) done!