Results

Demographics

# demographics 
min(tdata$Age)
## [1] 18
max(tdata$Age)
## [1] 64
mean(tdata$Age)
## [1] 32.1
sd(tdata$Age)
## [1] 12.51379
# 1 = male, 2 = female, 3 = other
table(tdata$Sex)
## 
##  1  2  3 
## 47 69  4

Graphs

myTheme <- theme(plot.title = element_text(face="bold", size = 22),
        axis.title.x = element_text(face = "bold", size = 20),
        axis.title.y = element_text(face = "bold", size = 20),
        axis.text.x = element_text(size = 18, angle = 0), 
        axis.text.y = element_text(size = 14, angle = 0),
        legend.text = element_text(size = 18),
        legend.title = element_text(face = "bold", size = 18),
        strip.text.x = element_text(size = 18),
        #panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(), 
        panel.background = element_blank(), 
        axis.line.x = element_line(colour = "black"), 
        axis.line.y = element_line(colour = "black"),
        axis.text = element_text(colour ="black"), 
        axis.ticks = element_line(colour ="black"))


library(see)
## Warning: package 'see' was built under R version 4.0.4
## first, turn sID into a factor
tdata_sub$sID <- factor(tdata_sub$sID)

pd <- position_dodge(width = 0.3)

tdata_sub$valueJitter <- jitter(tdata_sub$value, factor = 1, amount = 0.04)

theme_set(theme_light(base_size = 20, base_family = "Poppins"))

# new labes for the facets 
process.labs <- c("Process: Generative", "Process: Preventive")
names(process.labs) <- c("generative", "preventive")

g <- ggplot(tdata_sub, aes(x=variable, y=valueJitter, group = sID)) +
  guides(fill=FALSE)+
  facet_grid( ~ Process, labeller = labeller(Process =process.labs))+
  #ggtitle("Subjects' causal srength ratings") +
  scale_y_continuous(limits = c(-0.05, 1.05), breaks=seq(0, 1, 0.1), expand = c(0,0)) +
  scale_x_discrete(labels=c("Single-effect \n cause", "Multiple-effects \n cause")) +
  #stat_summary(fun.y = mean, geom = "bar", position = "dodge", colour = "black", alpha =0.5) +
  geom_violinhalf(aes(y = value, group = variable, fill = variable), color = NA, position=position_dodge(1), alpha = 0.2)+
  geom_line(position = pd, color = "black", size = 1, alpha=0.04) +
  geom_point(aes(color = variable), position = pd, alpha = 0.2) +
  stat_summary(aes(y = value,group=1), fun.data = mean_cl_boot, geom = "errorbar", width = 0, size = 1) +
  stat_summary(aes(y = value,group=1), fun.y=mean, colour="black", geom="line",group=1, size = 1.5, linetype = "solid", alpha = 1)+
  stat_summary(aes(y = value,group=1, fill = variable), fun.y=mean, geom="point", color = "black", shape = 22, size = 5, group=1, alpha = 1)+
  stat_summary(aes(y = value,group=1), fun.y=median, geom="point", color = "black", shape = 3, size = 4, group=1, alpha = 1, position = position_dodge(width = 0.5))+
  labs(x = "Entity", y = "Causal Strength Rating") +
  scale_color_manual(name = "Entity",values=c("#fc9272", "#3182bd"))+
  scale_fill_manual(name = "Entity",values=c("#fc9272", "#3182bd"))+
  theme(legend.position = "none")+
  myTheme
## Warning: `fun.y` is deprecated. Use `fun` instead.
## Warning: `fun.y` is deprecated. Use `fun` instead.

## Warning: `fun.y` is deprecated. Use `fun` instead.
g

#ggsave("results_lines.svg",width=11.5,height=6)
#ggsave("results_lines.pdf",width=11.5,height=6)

Descriptive Stats

## : single
## : generative
##       median         mean      SE.mean CI.mean.0.95          var      std.dev 
##   0.92000000   0.85025000   0.01599521   0.03167212   0.03070162   0.17521877 
##     coef.var 
##   0.20607912 
## ------------------------------------------------------------ 
## : multiple
## : generative
##       median         mean      SE.mean CI.mean.0.95          var      std.dev 
##   0.41000000   0.49525000   0.02156385   0.04269857   0.05579994   0.23622010 
##     coef.var 
##   0.47697143 
## ------------------------------------------------------------ 
## : single
## : preventive
## NULL
## ------------------------------------------------------------ 
## : multiple
## : preventive
## NULL
library(afex)
## ************
## Welcome to afex. For support visit: http://afex.singmann.science/
## - Functions for ANOVAs: aov_car(), aov_ez(), and aov_4()
## - Methods for calculating p-values with mixed(): 'KR', 'S', 'LRT', and 'PB'
## - 'afex_aov' and 'mixed' objects can be passed to emmeans() for follow-up tests
## - NEWS: library('emmeans') now needs to be called explicitly!
## - Get and set global package options with: afex_options()
## - Set orthogonal sum-to-zero contrasts globally: set_sum_contrasts()
## - For example analyses see: browseVignettes("afex")
## ************
## 
## Attaching package: 'afex'
## The following object is masked from 'package:lme4':
## 
##     lmer
library(emmeans)

a1 <- aov_car(value ~ Order*Multiple_Effects*Target*variable + Error(sID/(variable)), tdata_sub)
## Contrasts set to contr.sum for the following variables: Order, Multiple_Effects, Target
a1
## Anova Table (Type 3 tests)
## 
## Response: value
##                                    Effect     df  MSE          F  ges p.value
## 1                                   Order 1, 108 0.04    9.12 ** .039    .003
## 2                        Multiple_Effects 1, 108 0.04       1.71 .008    .193
## 3                                  Target 2, 108 0.04       0.80 .007    .453
## 4                  Order:Multiple_Effects 1, 108 0.04     3.81 + .017    .054
## 5                            Order:Target 2, 108 0.04       0.67 .006    .516
## 6                 Multiple_Effects:Target 2, 108 0.04       0.39 .003    .680
## 7           Order:Multiple_Effects:Target 2, 108 0.04       0.68 .006    .508
## 8                                variable 1, 108 0.04 180.80 *** .462   <.001
## 9                          Order:variable 1, 108 0.04       1.84 .009    .178
## 10              Multiple_Effects:variable 1, 108 0.04       0.70 .003    .403
## 11                        Target:variable 2, 108 0.04       0.60 .006    .552
## 12        Order:Multiple_Effects:variable 1, 108 0.04     4.50 * .021    .036
## 13                  Order:Target:variable 2, 108 0.04       0.44 .004    .647
## 14       Multiple_Effects:Target:variable 2, 108 0.04       0.65 .006    .523
## 15 Order:Multiple_Effects:Target:variable 2, 108 0.04     3.27 * .030    .042
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
lmeModel <- lmer(value ~ variable + (1|sID), data=tdata_sub)
## boundary (singular) fit: see ?isSingular
# follow-up analysis 

ls1 <- lsmeans(a1, c("variable"))
## NOTE: Results may be misleading due to involvement in interactions
ls1
##  variable lsmean     SE  df lower.CL upper.CL
##  single    0.850 0.0184 216    0.814    0.887
##  multiple  0.495 0.0184 216    0.459    0.532
## 
## Results are averaged over the levels of: Order, Multiple_Effects, Target 
## Warning: EMMs are biased unless design is perfectly balanced 
## Confidence level used: 0.95
############### 
# a conditional analysis 

ls2 <- lsmeans(a1, c("variable")) # group means by between-condition
## NOTE: Results may be misleading due to involvement in interactions
ls2
##  variable lsmean     SE  df lower.CL upper.CL
##  single    0.850 0.0184 216    0.814    0.887
##  multiple  0.495 0.0184 216    0.459    0.532
## 
## Results are averaged over the levels of: Order, Multiple_Effects, Target 
## Warning: EMMs are biased unless design is perfectly balanced 
## Confidence level used: 0.95
# simple main effects 
pairs(ls2) # compares rep-measure differences separately for each between-factor level
##  contrast          estimate     SE  df t.ratio p.value
##  single - multiple    0.355 0.0264 108 13.446  <.0001 
## 
## Results are averaged over the levels of: Order, Multiple_Effects, Target
# interaction contrast 
pairs(pairs(ls2), by = NULL)
##  contrast  estimate SE df z.ratio p.value
##  (nothing)   nonEst NA NA NA      NA     
## 
## Results are averaged over the levels of: Order, Multiple_Effects, Target
#test(pairs(pairs(ls2), by = NULL), joint = TRUE) # This reproduces the F-Value of the ANOVA interaction


#lsmip(a1, High_Strength_Component ~ variable) # lsemans can also produce graphs
# compute the confidence interval for the singular causation differences in each between-subject condition

# sending

t <- qt(0.975, 116, lower.tail = TRUE, log.p = FALSE)
#t

effect <- "Mdiff"
Mdiff <- 0.355
SE <- 0.0266
CI <- SE*t
CI_low <- Mdiff - CI
CI_up <- Mdiff + CI

Mdiff
## [1] 0.355
CI_low
## [1] 0.3023153
CI_up
## [1] 0.4076847
# Plot 

myTheme <- theme(plot.title = element_text(face="bold", size = 22),
        axis.title.x = element_text(face = "bold", size = 20),
        axis.title.y = element_text(face = "bold", size = 20),
        axis.text.x = element_text(size = 18, angle = 0), 
        axis.text.y = element_text(size = 25, angle = 0),
        legend.text = element_text(size = 18),
        legend.title = element_text(face = "bold", size = 18),
        strip.text.x = element_text(size = 18),
        #panel.grid.major = element_blank(), 
        #panel.grid.minor = element_blank(), 
        #panel.background = element_blank(), 
        axis.line.x = element_line(colour = "black"), 
        axis.line.y = element_line(colour = "black"),
        axis.text = element_text(colour ="black"), 
        axis.ticks = element_line(colour ="black"))

theme_set(theme_light(base_size = 30, base_family = "Poppins"))



barchart <- ggplot()+
  myTheme+
  #guides(fill=FALSE)+
  #facet_wrap(~Latency + SampleSize, ncol=2)+
  #ggtitle("Mean difference (95% CI)") +
  #coord_cartesian(ylim=c(-0.1,2)) + 
  scale_y_continuous(limits = c(-0.1, 0.5), breaks=seq(-0.1, 0.5, 0.1), expand = c(0,0)) +
  scale_x_discrete(labels=c("r")) +
  #annotate("rect", xmin=1.7, xmax=2.3, ymin=0.95, ymax=1.05, color="#31a354", fill = "white", size = 1) +
  #stat_summary(fun.y=mean, colour="grey20", geom="point", shape = 21, size = 3)+
  #stat_summary(fun.y = mean, geom = "bar", position = "dodge", colour = "black")+
  #stat_summary(fun.data = mean_cl_boot, geom = "errorbar", position = position_dodge(width = 0.90), width = 0.2) +
  #geom_jitter(width = 0.3, height = 0.02, alpha = 0.6, colour = "red") +
  #ggtitle("Means (95% bootstr. CIs)") +
  #theme(axis.text.x = element_text(size = 10, angle = 0, hjust = 0.5))+
  labs(x= "", y = "Mean change") +
  #scale_color_manual(values=c("#005083",  "#f0b64d"))# +
  #scale_fill_manual(values=c("#969696",  "#969696"))
  #annotate("point", x = 1, y = 100, colour = "firebrick", size = 2)+
  #annotate(xmin = -Inf, xmax = Inf, ymin = 4.77-1.96*0.297, ymax = 4.77+1.96*0.297, geom = "rect", alpha = 0.2, fill = "firebrick")+
  #annotate(xmin = -Inf, xmax = Inf, ymin = 5.02-1.96*0.372, ymax = 5.02+1.96*0.372, geom = "rect", alpha = 0.2, fill = "blue")+
  #annotate(geom = "hline",yintercept = 100, y = 100, color = "red")+
  annotate("pointrange", x = 1, y = Mdiff, ymin = CI_low, ymax = CI_up, colour = "black", size = 1.5, shape = 24, fill = "darkgrey")+
  #annotate("pointrange", x = 2, y = 5.02, ymin = 5.02-1.96*0.372, ymax = 5.02+1.96*0.372, colour = "blue", size = 0.8, shape = 15)+
  #annotate("text", x = 0.5, y = 2.6, family = "Poppins", size = 6, color = "gray20", label = "Impfeffekt")+
  #geom_curve(aes(x = 0.5, y = 3, xend = 0.9, yend = 4),arrow = arrow(length = unit(0.03, "npc")),color = "gray20", curvature = +0.2)+
  #annotate("text", x = 1.8, y = 2.6, family = "Poppins", size = 6, color = "gray20", label = "Dosierungseffekt")+
  #geom_curve(aes(x = 1.8, y = 3, xend = 2, yend = 4),arrow = arrow(length = unit(0.03, "npc")),color = "gray20", curvature = +0.2)+
  annotate(geom = "hline",yintercept = 0, y = 0, color = "red", size = 1.2)+
  theme(plot.background = element_rect(
    fill = "white",
    colour = "black",
    size = 1
  ))
## Warning: Ignoring unknown aesthetics: y
barchart

#ggsave("delta_posGen.svg",width=2.5,height=4)
#ggsave("delta_posGen.pdf",width=2.5,height=4)

Compute Cohen’s ds and their SE, CI.

# since we have a repeated-meausres design, we now need the correlations of the ratings

cor.test(tdata$multiple_strength_rating, tdata$single_strength_rating)
## 
##  Pearson's product-moment correlation
## 
## data:  tdata$multiple_strength_rating and tdata$single_strength_rating
## t = -0.17482, df = 118, p-value = 0.8615
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1947710  0.1636218
## sample estimates:
##         cor 
## -0.01609148
library(dplyr) # for pipe operator
# using the esc package because it gives SE of the ES directly
library(esc)

# get means and sds
m1 <- tdata_sub %>%
          filter(variable == "single")%>%
          summarize(Mean1 = mean(value))

sd1 <- tdata_sub %>%
          filter(variable == "single")%>%
          summarize(SD1 = sd(value))


m2 <- tdata_sub %>%
          filter(variable == "multiple")%>%
          summarize(Mean2 = mean(value))

sd2 <- tdata_sub %>%
          filter(variable == "multiple")%>%
          summarize(SD2 = sd(value))



esc_mean_sd(
  grp1m = m1[,1], grp1sd = sd1[,1], grp1n = 360/2,
  grp2m = m2[,1], grp2sd = sd2[,1], grp2n = 360/2,
  r = -0.01609148,
  es.type = "d"
)
## 
## Effect Size Calculation for Meta Analysis
## 
##      Conversion: mean and sd (within-subject) to effect size d
##     Effect Size:   1.1978
##  Standard Error:   0.1145
##        Variance:   0.0131
##        Lower CI:   0.9735
##        Upper CI:   1.4222
##          Weight:  76.3131