Homework Exercise 3
Hank Lee
Sun Oct 04 11:30:56 2020
1 Description
Examine the channel effects of the EEG study by following the analysis in the left brow EMG example.
2 input data
pacman::p_load(car, tidyverse, lme4, GGally)## Installing package into 'C:/Users/HANK/Documents/R/win-library/4.0'
## (as 'lib' is unspecified)## Warning: unable to access index for repository http://cran.rstudio.com/src/contrib:
## 起始 '<html>…' 的列異常！## Warning: unable to access index for repository http://www.stats.ox.ac.uk/pub/RWin/src/contrib:
## 起始 '<html>…' 的列異常！## Warning: package 'GGally' is not available (for R version 4.0.2)## Warning: unable to access index for repository http://cran.rstudio.com/bin/windows/contrib/4.0:
## 起始 '<html>…' 的列異常！## Warning: unable to access index for repository http://www.stats.ox.ac.uk/pub/RWin/bin/windows/contrib/4.0:
## 起始 '<html>…' 的列異常！## Warning: 'BiocManager' not available. Could not check Bioconductor.
##
## Please use `install.packages('BiocManager')` and then retry.## Warning in p_install(package, character.only = TRUE, ...):## Warning in library(package, lib.loc = lib.loc, character.only = TRUE,
## logical.return = TRUE, : there is no package called 'GGally'## Warning in pacman::p_load(car, tidyverse, lme4, GGally): Failed to install/load:
## GGallydta5 <- read.table("C:/Users/HANK/Desktop/HOMEWORK/thetaEEG.txt", header=T)
str(dta5)## 'data.frame': 19 obs. of 10 variables:
## $ ID : int 1 2 3 4 5 6 7 8 9 10 ...
## $ Ch3 : num -3.54 5.72 0.52 0 2.07 1.67 9.13 -0.43 -0.56 1.28 ...
## $ Ch4 : num -3.11 5.07 -0.18 0.74 0.76 ...
## $ Ch5 : num -0.24 6.87 0.9 1.1 3.51 2.77 3.44 -0.31 -1.22 1.89 ...
## $ Ch6 : num 0.42 5.96 0.6 0.13 0.6 4.55 4.8 -0.61 0.67 1.77 ...
## $ Ch7 : num -0.49 8.2 1.27 0.19 3.71 1.8 0.48 -1.04 -0.97 1.83 ...
## $ Ch8 : num 2.13 4.87 1.28 0.07 1.86 4.79 1.63 -0.13 -0.98 0.91 ...
## $ Ch17: num -4.15 5.48 -0.95 0.8 1.49 4.51 9.94 -0.61 0 1.4 ...
## $ Ch18: num 2.87 5.57 1.74 0.25 3.11 3.24 1.34 -0.61 -1.22 1.1 ...
## $ Ch19: num 1.34 6.33 0.79 -0.66 1.8 3.99 1.53 -0.43 -0.91 -0.12 ...head(dta5)## ID Ch3 Ch4 Ch5 Ch6 Ch7 Ch8 Ch17 Ch18 Ch19
## 1 1 -3.54 -3.11 -0.24 0.42 -0.49 2.13 -4.15 2.87 1.34
## 2 2 5.72 5.07 6.87 5.96 8.20 4.87 5.48 5.57 6.33
## 3 3 0.52 -0.18 0.90 0.60 1.27 1.28 -0.95 1.74 0.79
## 4 4 0.00 0.74 1.10 0.13 0.19 0.07 0.80 0.25 -0.66
## 5 5 2.07 0.76 3.51 0.60 3.71 1.86 1.49 3.11 1.80
## 6 6 1.67 4.18 2.77 4.55 1.80 4.79 4.51 3.24 3.99cha <- as.factor(colnames(dta5[, -1]))
# linear model
m0 <- lm(as.matrix(dta5[, -1]) ~ 1)
m0_aov <- Anova(m0, idata = data.frame(cha), idesign = ~ cha)## Note: model has only an intercept; equivalent type-III tests substituted.summary(m0_aov)##
## Type III Repeated Measures MANOVA Tests:
##
## ------------------------------------------
##
## Term: (Intercept)
##
## Response transformation matrix:
## (Intercept)
## Ch3 1
## Ch4 1
## Ch5 1
## Ch6 1
## Ch7 1
## Ch8 1
## Ch17 1
## Ch18 1
## Ch19 1
##
## Sum of squares and products for the hypothesis:
## (Intercept)
## (Intercept) 1761.616
##
## Multivariate Tests: (Intercept)
## Df test stat approx F num Df den Df Pr(>F)
## Pillai 1 0.2035907 4.601445 1 18 0.045839 *
## Wilks 1 0.7964093 4.601445 1 18 0.045839 *
## Hotelling-Lawley 1 0.2556358 4.601445 1 18 0.045839 *
## Roy 1 0.2556358 4.601445 1 18 0.045839 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## ------------------------------------------
##
## Term: cha
##
## Response transformation matrix:
## cha1 cha2 cha3 cha4 cha5 cha6 cha7 cha8
## Ch3 0 0 0 1 0 0 0 0
## Ch4 0 0 0 0 1 0 0 0
## Ch5 0 0 0 0 0 1 0 0
## Ch6 0 0 0 0 0 0 1 0
## Ch7 0 0 0 0 0 0 0 1
## Ch8 -1 -1 -1 -1 -1 -1 -1 -1
## Ch17 1 0 0 0 0 0 0 0
## Ch18 0 1 0 0 0 0 0 0
## Ch19 0 0 1 0 0 0 0 0
##
## Sum of squares and products for the hypothesis:
## cha1 cha2 cha3 cha4 cha5 cha6
## cha1 5.7475 1.650000e-02 1.485000000 0.506000000 7.70000000 1.930500000
## cha2 0.0165 4.736842e-05 0.004263158 0.001452632 0.02210526 0.005542105
## cha3 1.4850 4.263158e-03 0.383684211 0.130736842 1.98947368 0.498789474
## cha4 0.5060 1.452632e-03 0.130736842 0.044547368 0.67789474 0.169957895
## cha5 7.7000 2.210526e-02 1.989473684 0.677894737 10.31578947 2.586315789
## cha6 1.9305 5.542105e-03 0.498789474 0.169957895 2.58631579 0.648426316
## cha7 3.0635 8.794737e-03 0.791526316 0.269705263 4.10421053 1.028984211
## cha8 -0.0165 -4.736842e-05 -0.004263158 -0.001452632 -0.02210526 -0.005542105
## cha7 cha8
## cha1 3.063500000 -1.650000e-02
## cha2 0.008794737 -4.736842e-05
## cha3 0.791526316 -4.263158e-03
## cha4 0.269705263 -1.452632e-03
## cha5 4.104210526 -2.210526e-02
## cha6 1.028984211 -5.542105e-03
## cha7 1.632889474 -8.794737e-03
## cha8 -0.008794737 4.736842e-05
##
## Multivariate Tests: cha
## Df test stat approx F num Df den Df Pr(>F)
## Pillai 1 0.5243665 1.515882 8 11 0.25588
## Wilks 1 0.4756335 1.515882 8 11 0.25588
## Hotelling-Lawley 1 1.1024593 1.515882 8 11 0.25588
## Roy 1 1.1024593 1.515882 8 11 0.25588
##
## Univariate Type III Repeated-Measures ANOVA Assuming Sphericity
##
## Sum Sq num Df Error SS den Df F value Pr(>F)
## (Intercept) 195.735 1 765.68 18 4.6014 0.04584 *
## cha 10.702 8 389.31 144 0.4948 0.85840
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Mauchly Tests for Sphericity
##
## Test statistic p-value
## cha 0.00037434 9.7367e-11
##
##
## Greenhouse-Geisser and Huynh-Feldt Corrections
## for Departure from Sphericity
##
## GG eps Pr(>F[GG])
## cha 0.31531 0.6559
##
## HF eps Pr(>F[HF])
## cha 0.3709213 0.6854028# wide format to long format
dta5L <- dta5 %>%
gather( key = channel, value = EMG, 2:10) %>%
mutate( channel = factor(channel)) %>%
mutate( ID = factor(ID))#
str(dta5L)## 'data.frame': 171 obs. of 3 variables:
## $ ID : Factor w/ 19 levels "1","2","3","4",..: 1 2 3 4 5 6 7 8 9 10 ...
## $ channel: Factor w/ 9 levels "Ch17","Ch18",..: 4 4 4 4 4 4 4 4 4 4 ...
## $ EMG : num -3.54 5.72 0.52 0 2.07 1.67 9.13 -0.43 -0.56 1.28 ...# set theme to black and white and save parameters
ot <- theme_set(theme_bw())3 Data Analysis
# boxplot
ggplot(dta5L, aes(reorder(channel, EMG, median), EMG)) +
geom_boxplot() +
geom_point(alpha = .3) +
geom_hline(yintercept = mean(dta5L$EMG)) +
labs(x ="Channel condition", y = "EMG amplitutde")
# univariate approach
summary(m1<- aov(EMG ~ channel + Error(ID/channel), data = dta5L))##
## Error: ID
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 18 765.7 42.54
##
## Error: ID:channel
## Df Sum Sq Mean Sq F value Pr(>F)
## channel 8 10.7 1.338 0.495 0.858
## Residuals 144 389.3 2.704# show tables of model output
model.tables(m1, se = T)## Tables of effects
##
## channel
## channel
## Ch17 Ch18 Ch19 Ch3 Ch4 Ch5 Ch6 Ch7 Ch8
## 0.3327 -0.2157 -0.0751 -0.1688 0.5196 -0.0325 0.0759 -0.2188 -0.2173
##
## Standard errors of effects
## channel
## 0.3772
## replic. 19model.tables(m1, "means")## Tables of means
## Grand mean
##
## 1.069883
##
## channel
## channel
## Ch17 Ch18 Ch19 Ch3 Ch4 Ch5 Ch6 Ch7 Ch8
## 1.4026 0.8542 0.9947 0.9011 1.5895 1.0374 1.1458 0.8511 0.8526# relevel factor channel by mean EMG #[c(3,2,1,4)]?
dta5L$channel <- factor(dta5L$channel, levels(dta5L$channel))# EMG over channel condition with CIs
ggplot(dta5L, aes(channel, EMG)) +
stat_summary(fun.data = "mean_cl_boot", geom = "pointrange", col = "firebrick") +
geom_hline(yintercept = mean(dta5L$EMG), linetype = "dashed") +
geom_point(alpha = .3, size = rel(.7)) +
labs(x = "channel condition", y = "EMG amplititude") ## Warning: Computation failed in `stat_summary()`:
## Hmisc package required for this function
# Patient by emg with cha colors
ggplot(dta5L, aes(reorder(ID, EMG, mean), EMG, color = channel)) +
geom_point() +
geom_hline(yintercept = mean(dta5L$EMG), linetype = "dashed") +
coord_flip() +
scale_color_manual(values = c("gray", "forestgreen", "orangered", "navy", "blue", "red", "coral", "green", "skyblue"))+
labs(y = "EMG amplitude", x = "Patient ID") +
theme(legend.position = c(.9, .2))
# mixed-effect approach
summary(m2 <- lmer(EMG ~ channel + (1 | ID), data = dta5L))## Linear mixed model fit by REML ['lmerMod']
## Formula: EMG ~ channel + (1 | ID)
## Data: dta5L
##
## REML criterion at convergence: 697
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.6502 -0.3882 -0.0002 0.4215 4.6394
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 4.426 2.104
## Residual 2.704 1.644
## Number of obs: 171, groups: ID, 19
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 1.4026 0.6126 2.290
## channelCh18 -0.5484 0.5335 -1.028
## channelCh19 -0.4079 0.5335 -0.765
## channelCh3 -0.5016 0.5335 -0.940
## channelCh4 0.1868 0.5335 0.350
## channelCh5 -0.3653 0.5335 -0.685
## channelCh6 -0.2568 0.5335 -0.481
## channelCh7 -0.5516 0.5335 -1.034
## channelCh8 -0.5500 0.5335 -1.031
##
## Correlation of Fixed Effects:
## (Intr) chnC18 chnC19 chnnC3 chnnC4 chnnC5 chnnC6 chnnC7
## channelCh18 -0.435
## channelCh19 -0.435 0.500
## channelCh3 -0.435 0.500 0.500
## channelCh4 -0.435 0.500 0.500 0.500
## channelCh5 -0.435 0.500 0.500 0.500 0.500
## channelCh6 -0.435 0.500 0.500 0.500 0.500 0.500
## channelCh7 -0.435 0.500 0.500 0.500 0.500 0.500 0.500
## channelCh8 -0.435 0.500 0.500 0.500 0.500 0.500 0.500 0.500#
confint(m2)## Computing profile confidence intervals ...## 2.5 % 97.5 %
## .sig01 1.4944091 2.9720041
## .sigma 1.4360421 1.7986384
## (Intercept) 0.2052226 2.6000405
## channelCh18 -1.5725629 0.4757208
## channelCh19 -1.4320366 0.6162471
## channelCh3 -1.5257208 0.5225629
## channelCh4 -0.8372997 1.2109839
## channelCh5 -1.3894050 0.6588787
## channelCh6 -1.2809839 0.7672997
## channelCh7 -1.5757208 0.4725629
## channelCh8 -1.5741418 0.4741418#
ranef(m2)## $ID
## (Intercept)
## 1 -1.49820126
## 2 4.62406308
## 3 -0.38071113
## 4 -0.72927649
## 5 0.96568759
## 6 2.27566903
## 7 3.70218575
## 8 -1.60121012
## 9 -1.44617658
## 10 0.14057616
## 11 1.24662086
## 12 0.11560432
## 13 -2.64586570
## 14 -1.24848279
## 15 -1.24744230
## 16 2.16017424
## 17 -0.09977786
## 18 -2.34100107
## 19 -1.99243572
##
## with conditional variances for "ID"# residual plot
plot(m2, form = resid(., type = "pearson") ~ fitted(.) | channel,
abline = 0, lty = 2, layout = c(4,1), pch = 20,
xlab = "Fitted values", ylab = "Residuals")
# variance emg by cha
v_cha <- dta5L %>%
group_by(channel) %>%
summarize( var_cha = var(EMG)) %>%
as.data.frame## `summarise()` ungrouping output (override with `.groups` argument)# augument data with variance of emg by channel
dta5L2 <- merge(dta5L, v_cha, by = "channel")
as.data.frame(dta5L2)## channel ID EMG var_cha
## 1 Ch17 1 -4.15 10.738754
## 2 Ch17 2 5.48 10.738754
## 3 Ch17 3 -0.95 10.738754
## 4 Ch17 4 0.80 10.738754
## 5 Ch17 5 1.49 10.738754
## 6 Ch17 6 4.51 10.738754
## 7 Ch17 7 9.94 10.738754
## 8 Ch17 8 -0.61 10.738754
## 9 Ch17 9 0.00 10.738754
## 10 Ch17 10 1.40 10.738754
## 11 Ch17 11 3.31 10.738754
## 12 Ch17 12 1.46 10.738754
## 13 Ch17 13 -1.40 10.738754
## 14 Ch17 14 -0.12 10.738754
## 15 Ch17 15 0.71 10.738754
## 16 Ch17 16 6.12 10.738754
## 17 Ch17 17 1.20 10.738754
## 18 Ch17 18 -1.77 10.738754
## 19 Ch17 19 -0.77 10.738754
## 20 Ch18 1 2.87 6.769626
## 21 Ch18 2 5.57 6.769626
## 22 Ch18 3 1.74 6.769626
## 23 Ch18 4 0.25 6.769626
## 24 Ch18 5 3.11 6.769626
## 25 Ch18 6 3.24 6.769626
## 26 Ch18 7 1.34 6.769626
## 27 Ch18 8 -0.61 6.769626
## 28 Ch18 9 -1.22 6.769626
## 29 Ch18 10 1.10 6.769626
## 30 Ch18 11 1.20 6.769626
## 31 Ch18 12 2.26 6.769626
## 32 Ch18 13 -1.28 6.769626
## 33 Ch18 14 0.98 6.769626
## 34 Ch18 15 0.35 6.769626
## 35 Ch18 16 2.83 6.769626
## 36 Ch18 17 -0.24 6.769626
## 37 Ch18 18 -0.12 6.769626
## 38 Ch18 19 -7.14 6.769626
## 39 Ch19 1 1.34 5.566882
## 40 Ch19 2 6.33 5.566882
## 41 Ch19 3 0.79 5.566882
## 42 Ch19 4 -0.66 5.566882
## 43 Ch19 5 1.80 5.566882
## 44 Ch19 6 3.99 5.566882
## 45 Ch19 7 1.53 5.566882
## 46 Ch19 8 -0.43 5.566882
## 47 Ch19 9 -0.91 5.566882
## 48 Ch19 10 -0.12 5.566882
## 49 Ch19 11 2.15 5.566882
## 50 Ch19 12 2.01 5.566882
## 51 Ch19 13 -1.59 5.566882
## 52 Ch19 14 0.06 5.566882
## 53 Ch19 15 -1.29 5.566882
## 54 Ch19 16 4.39 5.566882
## 55 Ch19 17 0.36 5.566882
## 56 Ch19 18 -3.66 5.566882
## 57 Ch19 19 2.81 5.566882
## 58 Ch3 1 -3.54 8.645877
## 59 Ch3 2 5.72 8.645877
## 60 Ch3 3 0.52 8.645877
## 61 Ch3 4 0.00 8.645877
## 62 Ch3 5 2.07 8.645877
## 63 Ch3 6 1.67 8.645877
## 64 Ch3 7 9.13 8.645877
## 65 Ch3 8 -0.43 8.645877
## 66 Ch3 9 -0.56 8.645877
## 67 Ch3 10 1.28 8.645877
## 68 Ch3 11 3.21 8.645877
## 69 Ch3 12 1.47 8.645877
## 70 Ch3 13 -2.14 8.645877
## 71 Ch3 14 -1.40 8.645877
## 72 Ch3 15 0.29 8.645877
## 73 Ch3 16 2.58 8.645877
## 74 Ch3 17 0.85 8.645877
## 75 Ch3 18 -1.71 8.645877
## 76 Ch3 19 -1.89 8.645877
## 77 Ch4 1 -3.11 12.348605
## 78 Ch4 2 5.07 12.348605
## 79 Ch4 3 -0.18 12.348605
## 80 Ch4 4 0.74 12.348605
## 81 Ch4 5 0.76 12.348605
## 82 Ch4 6 4.18 12.348605
## 83 Ch4 7 12.92 12.348605
## 84 Ch4 8 -1.59 12.348605
## 85 Ch4 9 0.92 12.348605
## 86 Ch4 10 0.92 12.348605
## 87 Ch4 11 3.41 12.348605
## 88 Ch4 12 2.38 12.348605
## 89 Ch4 13 -2.44 12.348605
## 90 Ch4 14 0.37 12.348605
## 91 Ch4 15 0.31 12.348605
## 92 Ch4 16 3.16 12.348605
## 93 Ch4 17 3.02 12.348605
## 94 Ch4 18 -1.40 12.348605
## 95 Ch4 19 0.76 12.348605
## 96 Ch5 1 -0.24 5.560076
## 97 Ch5 2 6.87 5.560076
## 98 Ch5 3 0.90 5.560076
## 99 Ch5 4 1.10 5.560076
## 100 Ch5 5 3.51 5.560076
## 101 Ch5 6 2.77 5.560076
## 102 Ch5 7 3.44 5.560076
## 103 Ch5 8 -0.31 5.560076
## 104 Ch5 9 -1.22 5.560076
## 105 Ch5 10 1.89 5.560076
## 106 Ch5 11 3.92 5.560076
## 107 Ch5 12 1.22 5.560076
## 108 Ch5 13 -2.01 5.560076
## 109 Ch5 14 -1.10 5.560076
## 110 Ch5 15 -0.13 5.560076
## 111 Ch5 16 2.67 5.560076
## 112 Ch5 17 -0.59 5.560076
## 113 Ch5 18 -1.83 5.560076
## 114 Ch5 19 -1.15 5.560076
## 115 Ch6 1 0.42 5.073637
## 116 Ch6 2 5.96 5.073637
## 117 Ch6 3 0.60 5.073637
## 118 Ch6 4 0.13 5.073637
## 119 Ch6 5 0.60 5.073637
## 120 Ch6 6 4.55 5.073637
## 121 Ch6 7 4.80 5.073637
## 122 Ch6 8 -0.61 5.073637
## 123 Ch6 9 0.67 5.073637
## 124 Ch6 10 1.77 5.073637
## 125 Ch6 11 0.85 5.073637
## 126 Ch6 12 1.71 5.073637
## 127 Ch6 13 -1.95 5.073637
## 128 Ch6 14 0.18 5.073637
## 129 Ch6 15 -0.89 5.073637
## 130 Ch6 16 3.77 5.073637
## 131 Ch6 17 2.17 5.073637
## 132 Ch6 18 -2.01 5.073637
## 133 Ch6 19 -0.95 5.073637
## 134 Ch7 1 -0.49 5.181343
## 135 Ch7 2 8.20 5.181343
## 136 Ch7 3 1.27 5.181343
## 137 Ch7 4 0.19 5.181343
## 138 Ch7 5 3.71 5.181343
## 139 Ch7 6 1.80 5.181343
## 140 Ch7 7 0.48 5.181343
## 141 Ch7 8 -1.04 5.181343
## 142 Ch7 9 -0.97 5.181343
## 143 Ch7 10 1.83 5.181343
## 144 Ch7 11 2.77 5.181343
## 145 Ch7 12 -0.25 5.181343
## 146 Ch7 13 0.31 5.181343
## 147 Ch7 14 -1.71 5.181343
## 148 Ch7 15 -0.65 5.181343
## 149 Ch7 16 1.26 5.181343
## 150 Ch7 17 0.65 5.181343
## 151 Ch7 18 0.24 5.181343
## 152 Ch7 19 -1.43 5.181343
## 153 Ch8 1 2.13 4.281098
## 154 Ch8 2 4.87 4.281098
## 155 Ch8 3 1.28 4.281098
## 156 Ch8 4 0.07 4.281098
## 157 Ch8 5 1.86 4.281098
## 158 Ch8 6 4.79 4.281098
## 159 Ch8 7 1.63 4.281098
## 160 Ch8 8 -0.13 4.281098
## 161 Ch8 9 -0.98 4.281098
## 162 Ch8 10 0.91 4.281098
## 163 Ch8 11 0.79 4.281098
## 164 Ch8 12 -1.52 4.281098
## 165 Ch8 13 -3.30 4.281098
## 166 Ch8 14 0.37 4.281098
## 167 Ch8 15 -1.06 4.281098
## 168 Ch8 16 3.61 4.281098
## 169 Ch8 17 1.25 4.281098
## 170 Ch8 18 -0.61 4.281098
## 171 Ch8 19 0.24 4.281098# weighted fit
summary(m3 <- lmer(EMG ~ channel + (1 | ID), data = dta5L2,
weights = 1/var_cha))## Linear mixed model fit by REML ['lmerMod']
## Formula: EMG ~ channel + (1 | ID)
## Data: dta5L2
## Weights: 1/var_cha
##
## REML criterion at convergence: 683.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.8607 -0.4157 -0.0211 0.4295 4.0287
##
## Random effects:
## Groups Name Variance Std.Dev.
## ID (Intercept) 3.9833 1.9958
## Residual 0.3712 0.6093
## Number of obs: 171, groups: ID, 19
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 1.4026 0.6477 2.166
## channelCh18 -0.5484 0.5849 -0.938
## channelCh19 -0.4079 0.5644 -0.723
## channelCh3 -0.5016 0.6154 -0.815
## channelCh4 0.1868 0.6716 0.278
## channelCh5 -0.3653 0.5643 -0.647
## channelCh6 -0.2568 0.5558 -0.462
## channelCh7 -0.5516 0.5577 -0.989
## channelCh8 -0.5500 0.5417 -1.015
##
## Correlation of Fixed Effects:
## (Intr) chnC18 chnC19 chnnC3 chnnC4 chnnC5 chnnC6 chnnC7
## channelCh18 -0.554
## channelCh19 -0.574 0.636
## channelCh3 -0.526 0.583 0.604
## channelCh4 -0.482 0.534 0.553 0.508
## channelCh5 -0.574 0.636 0.659 0.604 0.554
## channelCh6 -0.583 0.645 0.669 0.613 0.562 0.669
## channelCh7 -0.581 0.643 0.667 0.611 0.560 0.667 0.677
## channelCh8 -0.598 0.662 0.686 0.629 0.577 0.686 0.697 0.694# residual plot
plot(m3, form = resid(., type = "pearson") ~ fitted(.) | channel,
abline = 0, lty = 2, layout = c(4,1), pch = 20,
xlab = "Fitted values", ylab = "Residuals")
# residual plots - combined and compared
# dta5L_m3 <- fortify(m3, dta5L2) %>% select() %>% mutate(fit3 = .fitted, res3 = .scresid) %>% select(- contains("."))
dta5L_m3 <- dta5L2 %>%
mutate(fit3 = fitted(m3),res3 = residuals(m3))
# dta5L_m2 <- fortify(m2, dta3L) %>% select(1:4, 6) %>% mutate( fit2 = .fitted, res2 = .scresid) %>% select( - contains("."))
dta5L_m2 <- dta5L %>%
mutate(fit2 = fitted(m2),res2 = residuals(m2))
dta5_2m <- merge(dta5L_m2, dta5L_m3, by = c("ID","channel", "EMG"))
#
ggplot(dta5_2m, aes(fit2, res2)) +
geom_point(alpha = .3) +
geom_point(aes(fit3, res3), alpha = .6) +
geom_hline(yintercept = 0, linetype="dashed") +
facet_grid(. ~ channel) +
labs(x = "Fitted values", y = "Standardized Residuals") 