get pain relationship, controlling for cue, cuetype, expect
model.stim <- lmer(event04_actual_angle ~
STIM_linear +
CUE_high_gt_low + STIM_quadratic+ EXPECT_demean +
EXPECT_cmc +
ses +
(1|sub), data = df
)
# CUE_high_gt_low+STIM+EXPECT_demean
sjPlot::tab_model(model.stim,
title = "Multilevel-modeling: \nlmer(NPSpos ~ CUE + STIM + EXPECT_demean + SES + (1| sub), data = pvc)",
CSS = list(css.table = '+font-size: 12;'))
Multilevel-modeling:
lmer(NPSpos ~ CUE + STIM + EXPECT_demean + SES + (1| sub), data = pvc)
|
event04_actual_angle
|
Predictors
|
Estimates
|
CI
|
p
|
(Intercept)
|
71.33
|
68.56 – 74.10
|
<0.001
|
STIM linear
|
29.89
|
28.43 – 31.34
|
<0.001
|
CUE high gt low
|
-2.04
|
-3.55 – -0.53
|
0.008
|
STIM quadratic
|
1.34
|
0.07 – 2.62
|
0.039
|
EXPECT demean
|
0.30
|
0.28 – 0.33
|
<0.001
|
EXPECT cmc
|
0.93
|
0.84 – 1.02
|
<0.001
|
sesses-03
|
-7.22
|
-8.86 – -5.59
|
<0.001
|
sesses-04
|
-7.45
|
-9.07 – -5.83
|
<0.001
|
Random Effects
|
σ2
|
366.84
|
τ00sub
|
157.35
|
ICC
|
0.30
|
N sub
|
96
|
Observations
|
4004
|
Marginal R2 / Conditional R2
|
0.645 / 0.751
|
ntrials = 12
lowintens = 48;
medintens = 49;
highintens = 50;
stim <- 48:50
painmean = 30 # average pain; arbitrary, on a 0 - 100 scale
painslope = fixEffect_expect['STIM_linear',1] # rise in pain per unit change in stim (per degree)
painslope_stan = 0.33621048
stdCoef.merMod <- function(object) {
sdy <- sd(getME(object,"y"))
sdx <- apply(getME(object,"X"), 2, sd)
sc <- fixef(object)*sdx/sdy
se.fixef <- coef(summary(object))[,"Std. Error"]
se <- se.fixef*sdx/sdy
return(data.frame(stdcoef=sc, stdse=se))
}
stdCoef.merMod(model.stim)
## stdcoef stdse
## (Intercept) 0.00000000 0.000000000
## STIM_linear 0.31648413 0.007848709
## CUE_high_gt_low -0.02640086 0.009953683
## STIM_quadratic 0.01621649 0.007849931
## EXPECT_demean 0.22157672 0.009947333
## EXPECT_cmc 0.69471095 0.034434185
## sesses-03 -0.08777975 0.010153298
## sesses-04 -0.09126342 0.010120143
# # library(limma)
# S <- rep(stim,times=ntrials) # stim
# C <- rep(rep(c(1,-1), each = 3), times = 6) #cue
# E <- painslope * (C + rnorm(length(C))) + painmean # pseudo nociception
# Szscore <- (S - mean(S)) / sd(S)
df$S <- as.numeric(mapvalues(df$stimintensity,
from = c("low", "med", "high"), c(48, 49, 50)))
df$C <- as.numeric(mapvalues(df$cuetype,
from = c("cuetype-low", "cuetype-high"), c(-1, 1)))
df$E <- painslope * (df$C + rnorm(length(df$C))) + painmean
df$Szscore <- (df$S - mean(df$S, na.rm = TRUE)) / sd(df$S)
df$Pcalib = df$Szscore * painslope + painmean + rnorm(length(df$C))
model.stim2pain <- lmer(Pcalib ~ S + (1|sub), df)
b_stim2pain = fixef(model.stim2pain)[2] #0.4126089 #36.5757
df$Sprime = df$Szscore * b_stim2pain + painmean # subjective pain experience, converted to a scale of 0-180, in order to match expectation ratings
# df$Sprime = df$S * b_stim2pain
df <- df %>%
group_by(sub) %>%
mutate(E = as.numeric(E)) %>%
mutate(avg_E = mean(E, na.rm = TRUE)) %>%
mutate(E_demean = E - avg_E) %>%
mutate(E_cmc = avg_E - mean(avg_E))
Lineplots Original
## Automatically converting the following non-factors to factors: cue_name
Lineplots P.assim
## Automatically converting the following non-factors to factors: cue_name
P.assim ~ demeaned_expect _ cue _ stim
Lineplots P.pe
## Automatically converting the following non-factors to factors: cue_name
P.pe ~ demeaned_expect _ cue _ stim
Lineplots P.adapt
## Automatically converting the following non-factors to factors: cue_name
P.adapt ~ demeaned_expect _ cue _ stim