Week 9 Lab Activity / Homework

## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr     1.1.2     ✔ readr     2.1.4
## ✔ forcats   1.0.0     ✔ stringr   1.5.0
## ✔ ggplot2   3.4.2     ✔ tibble    3.2.1
## ✔ lubridate 1.9.2     ✔ tidyr     1.3.0
## ✔ purrr     1.0.1     
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
## Loading required package: Matrix
## 
## 
## Attaching package: 'Matrix'
## 
## 
## The following objects are masked from 'package:tidyr':
## 
##     expand, pack, unpack
## 
## 
## 
## Attaching package: 'psych'
## 
## 
## The following objects are masked from 'package:ggplot2':
## 
##     %+%, alpha
## 
## 
## 
## Attaching package: 'MASS'
## 
## 
## The following object is masked from 'package:dplyr':
## 
##     select
## 
## 
## 
## arm (Version 1.13-1, built: 2022-8-25)
## 
## 
## Working directory is C:/Users/brand/Desktop/UCR/PhD CCN/PSYC213_Spring2023/Lab/Week 9
## 
## 
## 
## Attaching package: 'arm'
## 
## 
## The following objects are masked from 'package:psych':
## 
##     logit, rescale, sim
Homework Question
Olympic Data Example

The Olympic data has seven judges’ ratings of seven figure skaters ( on two cirteria: “Technical Merit” and “Artistic Expression”) from the 1932 Winter Olympics in Lake Placid.

(a) Formulate the data into a 98 x 4 array where the first two columns are technical merit and artistic expression scores. The third column is a skater id and the fourth is a judge id.
olym <- read_excel("1932_Olympic_Fig_skt_data_v01.xlsx")

A =  olym %>% filter(Criterion == "Artistic Expression")
M =  olym %>% filter(Criterion == "Technical Merit")
A_l = A %>% pivot_longer(cols=c('judge 1', 'judge 2','judge 3','judge 4',
                                    'judge 5','judge 6','judge 7'),
                    names_to='judge',
                    values_to='Artistic Expression')

M_l = M %>% pivot_longer(cols=c('judge 1', 'judge 2','judge 3','judge 4',
                                    'judge 5','judge 6','judge 7'),
                    names_to='judge',
                    values_to='Technical Merit')
olympic <- data.frame(A_l , M_l)
olympic_1932 <- olympic %>% dplyr::select(Technical.Merit, Artistic.Expression,
                                  Pair, judge) %>% mutate(Merit = Technical.Merit, Art = Artistic.Expression, Skater = Pair,Technical.Merit = NULL, Artistic.Expression = NULL, Pair = NULL)
(b) Write the notation for a non-nested multilevel model(varying across skaters and judges) for the technical merit ratings and fit using lmer().
m1 <- lmer(Merit ~ 1 + (1 | Skater) + (1 |judge),data = olympic_1932)
display(m1)
## lmer(formula = Merit ~ 1 + (1 | Skater) + (1 | judge), data = olympic_1932)
## coef.est  coef.se 
##     5.15     0.20 
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  Skater   (Intercept) 0.41    
##  judge    (Intercept) 0.29    
##  Residual             0.33    
## ---
## number of obs: 49, groups: Skater, 7; judge, 7
## AIC = 68.5, DIC = 57.6
## deviance = 59.1
coef(m1)
## $Skater
##   (Intercept)
## 1    5.496361
## 2    5.509435
## 3    5.391773
## 4    5.130303
## 5    5.078009
## 6    5.078009
## 7    4.358967
## 
## $judge
##         (Intercept)
## judge 1    5.418946
## judge 2    5.036209
## judge 3    5.514630
## judge 4    4.892683
## judge 5    5.000327
## judge 6    4.868762
## judge 7    5.311301
## 
## attr(,"class")
## [1] "coef.mer"
(c) Fit the model in (b) using the artistic expression ratings.
m2 <- lmer(Art ~ 1 + (1 | Skater) + (1 | judge),data = olympic_1932)
display(m2)
## lmer(formula = Art ~ 1 + (1 | Skater) + (1 | judge), data = olympic_1932)
## coef.est  coef.se 
##     5.09     0.20 
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  Skater   (Intercept) 0.45    
##  judge    (Intercept) 0.28    
##  Residual             0.27    
## ---
## number of obs: 49, groups: Skater, 7; judge, 7
## AIC = 54.2, DIC = 43.4
## deviance = 44.8
coef(m2)
## $Skater
##   (Intercept)
## 1    5.411951
## 2    5.507015
## 3    5.479854
## 4    5.208242
## 5    4.882307
## 6    4.855146
## 7    4.298341
## 
## $judge
##         (Intercept)
## judge 1    5.212086
## judge 2    5.262339
## judge 3    5.375409
## judge 4    4.747242
## judge 5    4.810059
## judge 6    4.910566
## judge 7    5.325156
## 
## attr(,"class")
## [1] "coef.mer"
(d) Display your results for both outcomes graphically.
j_1 <- coef(m1)[[1]]$`(Intercept)`
s_1 <- coef(m1)[[2]]$`(Intercept)`

ggplot( mapping = aes( x = j_1, y = s_1))+geom_point()

j_2 <- coef(m2)[[1]]$`(Intercept)`
s_2 <- coef(m2)[[2]]$`(Intercept)`

ggplot( mapping = aes( x = j_2, y = s_2))+geom_point()

olym <- read_excel("1932_Olympic_Fig_skt_data_v01.xlsx")
olym_l <- pivot_longer(olym, starts_with("judge"), names_to = "judge", values_to = "score")
olym_sw <- pivot_wider(olym_l, names_from = Criterion, values_from=score)
## Technical Merit vs Artistic Expression 

tm_lmer <- lmer(`Technical Merit` ~ 1 + (1 | judge) + (1 | Pair), olym_sw)
display(tm_lmer)
## lmer(formula = `Technical Merit` ~ 1 + (1 | judge) + (1 | Pair), 
##     data = olym_sw)
## coef.est  coef.se 
##     5.15     0.20 
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  judge    (Intercept) 0.29    
##  Pair     (Intercept) 0.41    
##  Residual             0.33    
## ---
## number of obs: 49, groups: judge, 7; Pair, 7
## AIC = 68.5, DIC = 57.6
## deviance = 59.1
coef(tm_lmer)
## $judge
##         (Intercept)
## judge 1    5.418946
## judge 2    5.036209
## judge 3    5.514630
## judge 4    4.892683
## judge 5    5.000327
## judge 6    4.868762
## judge 7    5.311301
## 
## $Pair
##   (Intercept)
## 1    5.496361
## 2    5.509435
## 3    5.391773
## 4    5.130303
## 5    5.078009
## 6    5.078009
## 7    4.358967
## 
## attr(,"class")
## [1] "coef.mer"
Skaters<-coef(tm_lmer)$Pair
#y = reorder(rownames(coef(tm_lmer)$Pair, `(Intercept)`))
#x<-
olym_sw %>% ggplot(aes(x = `Technical Merit`, 
                               y = `Artistic Expression`)) +
  geom_point(aes(colour = as.factor(Pair),shape = as.factor(judge),size = 5))
## Warning: The shape palette can deal with a maximum of 6 discrete values because
## more than 6 becomes difficult to discriminate; you have 7. Consider
## specifying shapes manually if you must have them.
## Warning: Removed 7 rows containing missing values (`geom_point()`).

ggplot(coef(tm_lmer)$judge, aes(x = `(Intercept)`, 
                               y = reorder(rownames(coef(tm_lmer)$judge), `(Intercept)`))) +
  geom_point(size = 3) 

ae_lmer <- lmer(`Artistic Expression` ~ 1 + (1 | judge) + (1 | Pair), olym_sw)
display(ae_lmer)
## lmer(formula = `Artistic Expression` ~ 1 + (1 | judge) + (1 | 
##     Pair), data = olym_sw)
## coef.est  coef.se 
##     5.09     0.20 
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  judge    (Intercept) 0.28    
##  Pair     (Intercept) 0.45    
##  Residual             0.27    
## ---
## number of obs: 49, groups: judge, 7; Pair, 7
## AIC = 54.2, DIC = 43.4
## deviance = 44.8
plot(`Technical Merit` ~ `Artistic Expression`, olym_sw, type="p")