Homework1
## Loading required package: MASS
## Loading required package: gld
## Loading required package: mvtnorm
## Loading required package: lattice
## Loading required package: ggplot2
##
## Attaching package: 'PairedData'
## The following object is masked from 'package:base':
##
## summary
##
## Attaching package: 'dplyr'
## The following object is masked from 'package:MASS':
##
## select
## The following objects are masked from 'package:stats':
##
## filter, lag
## The following objects are masked from 'package:base':
##
## intersect, setdiff, setequal, union
## Prior Post
## 1 83.8 95.2
## 2 83.3 94.3
## 3 86.0 91.5
## 4 82.5 91.9
## 5 86.7 100.3
## 6 79.6 76.7
## 'data.frame': 17 obs. of 2 variables:
## $ Prior: num 83.8 83.3 86 82.5 86.7 79.6 76.9 94.2 73.4 80.5 ...
## $ Post : num 95.2 94.3 91.5 91.9 100.3 ...
## Prior Post
## 83.22941 90.49412
## Prior Post
## 5.02 8.48
## Prior Post
## Prior 25.16721 22.88268
## Post 22.88268 71.82684
## Prior Post
## Prior 1.000 0.538
## Post 0.538 1.000
##
## Welch Two Sample t-test
##
## data: Anorexia$Prior and Anorexia$Post
## t = -3.0414, df = 25.986, p-value = 0.005324
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -12.17472 -2.35469
## sample estimates:
## mean of x mean of y
## 83.22941 90.49412
##
## Paired t-test
##
## data: Anorexia$Prior and Anorexia$Post
## t = -4.1849, df = 16, p-value = 0.0007003
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -10.94471 -3.58470
## sample estimates:
## mean of the differences
## -7.264706
## Therapy Weight
## 1 Prior 83.8
## 2 Prior 83.3
## 3 Prior 86.0
## 4 Prior 82.5
## 5 Prior 86.7
## 6 Prior 79.6

##
## Call:
## lm(formula = Weight ~ Therapy, data = dtaL)
##
## Residuals:
## Min 1Q Median 3Q Max
## -15.294 -2.454 1.106 4.004 11.106
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 90.494 1.689 53.578 < 2e-16 ***
## TherapyPrior -7.265 2.389 -3.041 0.00467 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 6.964 on 32 degrees of freedom
## Multiple R-squared: 0.2242, Adjusted R-squared: 0.2
## F-statistic: 9.25 on 1 and 32 DF, p-value: 0.004673
## Df Sum Sq Mean Sq F value Pr(>F)
## Therapy 1 448.6 448.6 9.25 0.00467 **
## Residuals 32 1551.9 48.5
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Subject
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 16 1142 71.38
##
## Error: Subject:Therapy
## Df Sum Sq Mean Sq F value Pr(>F)
## Therapy 1 448.6 448.6 17.51 7e-04 ***
## Residuals 16 409.8 25.6
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Df Sum Sq Mean Sq F value Pr(>F)
## Therapy 1 448.6 448.6 17.513 0.0007 ***
## Subject 16 1142.1 71.4 2.787 0.0240 *
## Residuals 16 409.8 25.6
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Df Sum Sq Mean Sq
## Therapy 1 448.6 448.6
## Sbj 32 1551.9 48.5


Homework2
##
## Attaching package: 'nlme'
## The following object is masked from 'package:dplyr':
##
## collapse
## Loading required package: Matrix
##
## Attaching package: 'Matrix'
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
##
## Attaching package: 'lme4'
## The following object is masked from 'package:nlme':
##
## lmList
## Grouped Data: effort ~ Type | Subject
## effort Type Subject
## 1 12 T1 1
## 2 15 T2 1
## 3 12 T3 1
## 4 10 T4 1
## 5 10 T1 2
## 6 14 T2 2
## Classes 'nffGroupedData', 'nfGroupedData', 'groupedData' and 'data.frame': 36 obs. of 3 variables:
## $ effort : num 12 15 12 10 10 14 13 12 7 14 ...
## $ Type : Factor w/ 4 levels "T1","T2","T3",..: 1 2 3 4 1 2 3 4 1 2 ...
## $ Subject: Ord.factor w/ 9 levels "8"<"5"<"4"<"9"<..: 8 8 8 8 9 9 9 9 6 6 ...
## - attr(*, "formula")=Class 'formula' language effort ~ Type | Subject
## .. ..- attr(*, ".Environment")=<environment: R_GlobalEnv>
## - attr(*, "labels")=List of 2
## ..$ x: chr "Type of stool"
## ..$ y: chr "Effort required to arise"
## - attr(*, "units")=List of 1
## ..$ y: chr "(Borg scale)"
## - attr(*, "FUN")=function (x)
## ..- attr(*, "source")= chr "function (x) mean(x, na.rm = TRUE)"
## - attr(*, "order.groups")= logi TRUE
## starting httpd help server ...
## done
| T1 |
8.555556 |
| T2 |
12.444444 |
| T3 |
10.777778 |
| T4 |
9.222222 |
| 8 |
8.25 |
| 5 |
8.50 |
| 4 |
9.25 |
| 9 |
10.00 |
| 6 |
10.25 |
| 3 |
10.75 |
| 7 |
10.75 |
| 1 |
12.25 |
| 2 |
12.25 |
## Fixed Effects:
## coef.est coef.se
## (Intercept) 8.56 0.58
## TypeT2 3.89 0.52
## TypeT3 2.22 0.52
## TypeT4 0.67 0.52
##
## Random Effects:
## Groups Name Std.Dev.
## Subject (Intercept) 1.33
## Residual 1.10
## ---
## number of obs: 36, groups: Subject, 9
## AIC = 133.1, DIC = 123.2
## deviance = 122.1
## Fixed Effects:
## coef.est coef.se
## TypeT1 8.56 0.58
## TypeT2 12.44 0.58
## TypeT3 10.78 0.58
## TypeT4 9.22 0.58
##
## Random Effects:
## Groups Name Std.Dev.
## Subject (Intercept) 1.33
## Residual 1.10
## ---
## number of obs: 36, groups: Subject, 9
## AIC = 133.1, DIC = 123.2
## deviance = 122.1
## Linear mixed-effects model fit by REML
## Data: ergoStool
## AIC BIC logLik
## 133.1308 141.9252 -60.56539
##
## Random effects:
## Formula: ~1 | Subject
## (Intercept) Residual
## StdDev: 1.332465 1.100295
##
## Fixed effects: effort ~ Type
## Value Std.Error DF t-value p-value
## (Intercept) 8.555556 0.5760123 24 14.853079 0.0000
## TypeT2 3.888889 0.5186838 24 7.497610 0.0000
## TypeT3 2.222222 0.5186838 24 4.284348 0.0003
## TypeT4 0.666667 0.5186838 24 1.285304 0.2110
## Correlation:
## (Intr) TypeT2 TypeT3
## TypeT2 -0.45
## TypeT3 -0.45 0.50
## TypeT4 -0.45 0.50 0.50
##
## Standardized Within-Group Residuals:
## Min Q1 Med Q3 Max
## -1.80200345 -0.64316591 0.05783115 0.70099706 1.63142054
##
## Number of Observations: 36
## Number of Groups: 9
## Approximate 95% confidence intervals
##
## Fixed effects:
## lower est. upper
## (Intercept) 7.3667247 8.5555556 9.744386
## TypeT2 2.8183781 3.8888889 4.959400
## TypeT3 1.1517114 2.2222222 3.292733
## TypeT4 -0.4038442 0.6666667 1.737177
## attr(,"label")
## [1] "Fixed effects:"
##
## Random Effects:
## Level: Subject
## lower est. upper
## sd((Intercept)) 0.749509 1.332465 2.368835
##
## Within-group standard error:
## lower est. upper
## 0.8292494 1.1002946 1.4599324
Homework3
## 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.99
## '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 ...
## Note: model has only an intercept; equivalent type-III tests substituted.
##
## 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: newhw3
##
## Response transformation matrix:
## newhw31 newhw32 newhw33 newhw34 newhw35 newhw36 newhw37 newhw38
## 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:
## newhw31 newhw32 newhw33 newhw34 newhw35
## newhw31 5.7475 1.650000e-02 1.485000000 0.506000000 7.70000000
## newhw32 0.0165 4.736842e-05 0.004263158 0.001452632 0.02210526
## newhw33 1.4850 4.263158e-03 0.383684211 0.130736842 1.98947368
## newhw34 0.5060 1.452632e-03 0.130736842 0.044547368 0.67789474
## newhw35 7.7000 2.210526e-02 1.989473684 0.677894737 10.31578947
## newhw36 1.9305 5.542105e-03 0.498789474 0.169957895 2.58631579
## newhw37 3.0635 8.794737e-03 0.791526316 0.269705263 4.10421053
## newhw38 -0.0165 -4.736842e-05 -0.004263158 -0.001452632 -0.02210526
## newhw36 newhw37 newhw38
## newhw31 1.930500000 3.063500000 -1.650000e-02
## newhw32 0.005542105 0.008794737 -4.736842e-05
## newhw33 0.498789474 0.791526316 -4.263158e-03
## newhw34 0.169957895 0.269705263 -1.452632e-03
## newhw35 2.586315789 4.104210526 -2.210526e-02
## newhw36 0.648426316 1.028984211 -5.542105e-03
## newhw37 1.028984211 1.632889474 -8.794737e-03
## newhw38 -0.005542105 -0.008794737 4.736842e-05
##
## Multivariate Tests: newhw3
## 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 *
## newhw3 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
## newhw3 0.00037434 9.7367e-11
##
##
## Greenhouse-Geisser and Huynh-Feldt Corrections
## for Departure from Sphericity
##
## GG eps Pr(>F[GG])
## newhw3 0.31531 0.6559
##
## HF eps Pr(>F[HF])
## newhw3 0.3709213 0.6854028
## '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 ...

##
## 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
## 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
## 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. 19
dta3L$channel <- factor(dta3L$channel, levels(dta3L$channel))
ggplot(dta3L, aes(reorder(channel, EMG, mean),EMG)) +
stat_summary(fun.data = "mean_cl_boot", geom = "pointrange", col = "firebrick") +
geom_hline(yintercept = mean(dta3L$EMG), linetype = "dashed") +
geom_point(alpha = .3, size = rel(.7)) +
labs(x = "channel condition", y = "EMG amplititude")


## Linear mixed model fit by REML ['lmerMod']
## Formula: EMG ~ channel + (1 | ID)
## Data: dta3L
##
## 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
## 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
## $ID
## (Intercept)
## 1 -1.49820132
## 2 4.62406329
## 3 -0.38071115
## 4 -0.72927652
## 5 0.96568763
## 6 2.27566913
## 7 3.70218592
## 8 -1.60121019
## 9 -1.44617664
## 10 0.14057617
## 11 1.24662091
## 12 0.11560432
## 13 -2.64586581
## 14 -1.24848285
## 15 -1.24744236
## 16 2.16017433
## 17 -0.09977786
## 18 -2.34100118
## 19 -1.99243581
##
## with conditional variances for "ID"



## channel ID EMG var_chan
## 1 Ch17 1 -4.15 10.73875
## 2 Ch17 2 5.48 10.73875
## 3 Ch17 3 -0.95 10.73875
## 4 Ch17 4 0.80 10.73875
## 5 Ch17 5 1.49 10.73875
## 6 Ch17 6 4.51 10.73875
## Linear mixed model fit by REML ['lmerMod']
## Formula: EMG ~ channel + (1 | ID)
## Data: dta3L2
## Weights: 1/var_chan
##
## 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



Homework4
#load
dta4 <- read.table("D:/sheu/cognitive_task.txt", h=T)
#varPlot
VCA::varPlot(Score ~ Method/Class/ID,
Data=dta4,
YLabel=list(text="Score", side=2, cex=1),
MeanLine=list(var=c("Method", "Class"),
col=c("darkred", "salmon"), lwd=c(1, 2)))
## Warning in min(x): min 中沒有無漏失的引數; 回傳 Inf
## Warning in max(x): max 中沒有無漏失的引數;回傳 -Inf

##
## Error: Class
## Df Sum Sq Mean Sq F value Pr(>F)
## Method 1 112.5 112.50 6.459 0.044 *
## Residuals 6 104.5 17.42
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 24 18.5 0.7708
##
## Error: Method
## Df Sum Sq Mean Sq
## Method 1 112.5 112.5
##
## Error: Method:Klass
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 6 104.5 17.42
##
## Error: Within
## Df Sum Sq Mean Sq F value Pr(>F)
## Residuals 24 18.5 0.7708
## Computing bootstrap confidence intervals ...
##
## 3 warning(s): Model failed to converge with max|grad| = 0.00211488 (tol = 0.002, component 1) (and others)
## 2.5 % 97.5 %
## .sig01 0.8016683 3.299184
## .sigma 0.6328451 1.114697
## (Intercept) 1.5741244 5.804193
## MethodI2 0.9221775 6.606790