Pre-processing Data
## [1] 9
##
## Pearson's product-moment correlation
##
## data: mf.stop$Stop.male and mf.stop$Stop.female
## t = 4.7099, df = 97, p-value = 8.261e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.2558004 0.5794841
## sample estimates:
## cor
## 0.4314261
## Warning in ICC::ICCest(ID.relation, value, data = mf.stop.l): 'x' has been
## coerced to a factor
## $ICC
## [1] 0.4183118
##
## $LowerCI
## [1] 0.2422021
##
## $UpperCI
## [1] 0.5678736
##
## $N
## [1] 99
##
## $k
## [1] 2
##
## $varw
## [1] 8478.709
##
## $vara
## [1] 6097.329
Obtain ICCs for tolerance time and other survey data
ICC::ICCest(ID.relation,value, data=mf.stop.l) #calculate the ICC Stop time
## Warning in ICC::ICCest(ID.relation, value, data = mf.stop.l): 'x' has been
## coerced to a factor
## $ICC
## [1] 0.4183118
##
## $LowerCI
## [1] 0.2422021
##
## $UpperCI
## [1] 0.5678736
##
## $N
## [1] 99
##
## $k
## [1] 2
##
## $varw
## [1] 8478.709
##
## $vara
## [1] 6097.329
ICC::ICCest(ID.relation,SelfEsteem, data=mf.relate.l.2) ##self-esteem
## Warning in ICC::ICCest(ID.relation, SelfEsteem, data = mf.relate.l.2): 'x'
## has been coerced to a factor
## $ICC
## [1] 0.2738551
##
## $LowerCI
## [1] 0.0597044
##
## $UpperCI
## [1] 0.4641515
##
## $N
## [1] 80
##
## $k
## [1] 2
##
## $varw
## [1] 16.3
##
## $vara
## [1] 6.14731
ICC::ICCest(ID.relation,CTS, data=mf.relate.l.2) ##current thoughts
## Warning in ICC::ICCest(ID.relation, CTS, data = mf.relate.l.2): 'x' has
## been coerced to a factor
## $ICC
## [1] 0.1036855
##
## $LowerCI
## [1] -0.1166582
##
## $UpperCI
## [1] 0.314582
##
## $N
## [1] 80
##
## $k
## [1] 2
##
## $varw
## [1] 133.0625
##
## $vara
## [1] 15.39264
ICC::ICCest(ID.relation,TT, data=mf.relate.l.2) ##total Trust
## Warning in ICC::ICCest(ID.relation, TT, data = mf.relate.l.2): 'x' has been
## coerced to a factor
## $ICC
## [1] 0.3330225
##
## $LowerCI
## [1] 0.1243261
##
## $UpperCI
## [1] 0.5137323
##
## $N
## [1] 80
##
## $k
## [1] 2
##
## $varw
## [1] 33.275
##
## $vara
## [1] 16.61424
Plots for each relationship (Key below)
- Fo = Female Observer Rating
- Fe = Female Experience Rating
- Mo = Male Observer Rating
- Me = Male Experience Rating
- Orange = First in the order
- Blue = Second in the order
Model Fit
#Basic Multilevel model
library(nlme)
CP.long$Pain.obsv <- as.numeric(CP.long$Pain.obsv)
## Warning: NAs introduced by coercion
CP.long <- na.omit(CP.long)
#CP.long <- CP.long[!is.na(CP.long$Pain.exp),] #remove NAs from only the predictor variable
####------this model is if I include time 1 and zero as controls
#model1 <- lme(fixed=Pain.exp~Pain.obsv + t.0+t.1, random=~1 + Pain.obsv | as.factor(ID.relation), data=CP.long)
#summary(model1)
#subset the data so that you remove the 0 and 30 seconds
CP.long.sub <- CP.long[!CP.long$time.e %in% c("EP.0.x","EP.30.x"),]
CP.long.sub$Pain.exp <- CP.long.sub$Pain.exp-mean(CP.long.sub$Pain.exp)
CP.long.sub$Pain.obsv <- CP.long.sub$Pain.obsv-mean(CP.long.sub$Pain.obsv)
model1 <- lme(fixed=Pain.exp~Pain.obsv, random=~1 + Pain.obsv | as.factor(ID.relation), data=CP.long.sub)
summary(model1)
## Linear mixed-effects model fit by REML
## Data: CP.long.sub
## AIC BIC logLik
## 1050.727 1072.91 -519.3636
##
## Random effects:
## Formula: ~1 + Pain.obsv | as.factor(ID.relation)
## Structure: General positive-definite, Log-Cholesky parametrization
## StdDev Corr
## (Intercept) 2.0103449 (Intr)
## Pain.obsv 0.3593478 0.365
## Residual 0.9549050
##
## Fixed effects: Pain.exp ~ Pain.obsv
## Value Std.Error DF t-value p-value
## (Intercept) 0.456664 0.27564356 237 1.656719 0.0989
## Pain.obsv 0.300801 0.07202603 237 4.176282 0.0000
## Correlation:
## (Intr)
## Pain.obsv 0.221
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.40763837 -0.56361266 0.05026319 0.54866684 2.63257238
##
## Number of Observations: 300
## Number of Groups: 62
model2 <- lme(fixed=Pain.exp~1, random=~1 | as.factor(ID.relation), data=CP.long.sub)
summary(model2)
## Linear mixed-effects model fit by REML
## Data: CP.long.sub
## AIC BIC logLik
## 1089.396 1100.497 -541.6979
##
## Random effects:
## Formula: ~1 | as.factor(ID.relation)
## (Intercept) Residual
## StdDev: 1.980183 1.13424
##
## Fixed effects: Pain.exp ~ 1
## Value Std.Error DF t-value p-value
## (Intercept) 0.5406304 0.2661333 238 2.031427 0.0433
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.2758970 -0.6278236 0.1335636 0.5718746 2.4381145
##
## Number of Observations: 300
## Number of Groups: 62
ICC.randomintercepts <- 1.98^2/(1.98^2+1.13^2)
ICC.randomintercepts
## [1] 0.7543147
CP.long.dif <- CP.long.sub
CP.long.dif$dif <- CP.long.dif$Pain.exp-CP.long.dif$Pain.obsv
model3 <- lme(fixed=dif~1, random=~1 | as.factor(ID.relation), data=CP.long.dif)
summary(model3)
## Linear mixed-effects model fit by REML
## Data: CP.long.dif
## AIC BIC logLik
## 1245.473 1256.574 -619.7363
##
## Random effects:
## Formula: ~1 | as.factor(ID.relation)
## (Intercept) Residual
## StdDev: 2.356141 1.500968
##
## Fixed effects: dif ~ 1
## Value Std.Error DF t-value p-value
## (Intercept) 0.4641934 0.3203022 238 1.449236 0.1486
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.49096362 -0.60367782 0.06255887 0.53811237 3.50516662
##
## Number of Observations: 300
## Number of Groups: 62
ICC.randomintercepts.dif <- 2.35^2/(2.35^2+1.50^2)
ICC.randomintercepts.dif
## [1] 0.7105179
CP.long.dif$dif.abs <- abs(CP.long.dif$dif)
model4 <- lme(fixed=dif.abs~1, random=~1 | as.factor(ID.relation), data=CP.long.dif)
summary(model4)
## Linear mixed-effects model fit by REML
## Data: CP.long.dif
## AIC BIC logLik
## 1050.412 1061.513 -522.2058
##
## Random effects:
## Formula: ~1 | as.factor(ID.relation)
## (Intercept) Residual
## StdDev: 1.379574 1.131352
##
## Fixed effects: dif.abs ~ 1
## Value Std.Error DF t-value p-value
## (Intercept) 2.156446 0.194289 238 11.09917 0
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -2.159040543 -0.631213897 0.009164038 0.447891089 3.261622687
##
## Number of Observations: 300
## Number of Groups: 62
ICC.randomintercepts.dif.abs <- 1.38^2/(1.38^2+1.14^2)
ICC.randomintercepts.dif.abs
## [1] 0.594382
