Initialize R with the following commands:

options(contrasts=c("contr.sum","contr.poly"))
load(url("http://pnb.mcmaster.ca/bennett/psy710/datasets/p2-2023.rda") )

1 Spot detection thresholds

The data frame spot contains data from an experiment that used a 2(age) x 2 (luminance) x 3 (duration) factorial design to investigate how visual sensitivity (for a spot presented on a uniform, grey background) varies with age (young vs. old), background luminance (low vs high), and stimulus duration (d1 (short), d2 (medium), d3 (long)). The dependent variable was contrast sensitivity. Twelve participants were assigned to each condition.

Confirm the data are balanced. Then use a factorial ANOVA to evaluate the main effects of age, luminance, duration, and all of the interactions. List the ANOVA table.

summary(spot)
##       sID       sensitivity        age     luminance duration
##  s1     :  1   Min.   : 60.21   young:72   low :72   d1:48   
##  s2     :  1   1st Qu.:241.63   old  :72   high:72   d2:48   
##  s3     :  1   Median :298.53                        d3:48   
##  s4     :  1   Mean   :300.00                                
##  s5     :  1   3rd Qu.:367.84                                
##  s6     :  1   Max.   :540.83                                
##  (Other):138
with(spot,tapply(sensitivity,list(age,luminance,duration),length)) # balanced
## , , d1
## 
##       low high
## young  12   12
## old    12   12
## 
## , , d2
## 
##       low high
## young  12   12
## old    12   12
## 
## , , d3
## 
##       low high
## young  12   12
## old    12   12
xtabs(~duration+luminance+age,data=spot)
## , , age = young
## 
##         luminance
## duration low high
##       d1  12   12
##       d2  12   12
##       d3  12   12
## 
## , , age = old
## 
##         luminance
## duration low high
##       d1  12   12
##       d2  12   12
##       d3  12   12

spot.aov.01 <- aov(sensitivity~age*luminance*duration,data=spot)
summary(spot.aov.01)
##                         Df  Sum Sq Mean Sq F value   Pr(>F)    
## age                      1    4642    4642   0.568 0.452385    
## luminance                1   34728   34728   4.250 0.041213 *  
## duration                 2  128289   64144   7.850 0.000601 ***
## age:luminance            1    5208    5208   0.637 0.426093    
## age:duration             2    1322     661   0.081 0.922330    
## luminance:duration       2   21601   10800   1.322 0.270174    
## age:luminance:duration   2  104588   52294   6.400 0.002225 ** 
## Residuals              132 1078615    8171                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Evaluate all pairwise differences among the marginal means of duration.

TukeyHSD(spot.aov.01,which="duration")
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = sensitivity ~ age * luminance * duration, data = spot)
## 
## $duration
##             diff        lwr       upr     p adj
## d2-d1 -57.864363 -101.60354 -14.12519 0.0059485
## d3-d1   9.769321  -33.96986  53.50850 0.8570377
## d3-d2  67.633683   23.89451 111.37286 0.0010335
library(emmeans)
spot.emm <- emmeans(spot.aov.01,specs=~duration)
contrast(spot.emm,method="pairwise",adjust="tukey")
##  contrast estimate   SE  df t.ratio p.value
##  d1 - d2     57.86 18.5 132   3.136  0.0059
##  d1 - d3     -9.77 18.5 132  -0.529  0.8570
##  d2 - d3    -67.63 18.5 132  -3.665  0.0010
## 
## Results are averaged over the levels of: age, luminance 
## P value adjustment: tukey method for comparing a family of 3 estimates

Evaluate the simple luminance x duration interaction at each level of age. (N.B. You do not have to recalculate the \(F\) and \(p\) values).

# with aov
MS.error <- 8171
df.error <- 132
aov.old.01 <- aov(sensitivity~luminance*duration,data=subset(spot,age=="old"))
aov.young.01 <- aov(sensitivity~luminance*duration,data=subset(spot,age=="young"))
summary(aov.old.01)
##                    Df Sum Sq Mean Sq F value  Pr(>F)   
## luminance           1  33417   33417   4.053 0.04818 * 
## duration            2  67882   33941   4.116 0.02067 * 
## luminance:duration  2 105206   52603   6.379 0.00293 **
## Residuals          66 544211    8246                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(aov.young.01)
##                    Df Sum Sq Mean Sq F value Pr(>F)  
## luminance           1   6519    6519   0.805 0.3728  
## duration            2  61729   30865   3.812 0.0271 *
## luminance:duration  2  20983   10491   1.296 0.2806  
## Residuals          66 534403    8097                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(F.adj <- c(52603,10491)/MS.error)
## [1] 6.437768 1.283931
(p.adj <- 1-pf(F.adj,df1=2,df2=df.error))
## [1] 0.002148849 0.280382052
(p.adj < 0.05) # simple interaction significant only in old group
## [1]  TRUE FALSE

# with emmeans:
spot.aov.01 <- aov(sensitivity~age*luminance*duration,data=spot)
spot.emm.02 <- emmeans(spot.aov.01,specs=~luminance*duration|age)
joint_tests(spot.emm.02,by="age")
## age = young:
##  model term         df1 df2 F.ratio p.value
##  luminance            1 132   0.798  0.3734
##  duration             2 132   3.777  0.0254
##  luminance:duration   2 132   1.284  0.2804
## 
## age = old:
##  model term         df1 df2 F.ratio p.value
##  luminance            1 132   4.090  0.0452
##  duration             2 132   4.154  0.0178
##  luminance:duration   2 132   6.437  0.0021

Evaluate the simple-simple main effects of age in the low-luminance, long-duration condition, and in the high-luminance, long-duration (d3) condition.

MS.error <- 8171
df.error <- 132

# low luminance:
aov.age.low <- aov(sensitivity~age,data=subset(spot,duration=="d3"&luminance=="low"))
summary(aov.age.low)
##             Df Sum Sq Mean Sq F value Pr(>F)  
## age          1  35057   35057   3.262 0.0846 .
## Residuals   22 236426   10747                 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
(F.adj <- 35057/MS.error)
## [1] 4.290417
(p.adj <- 1-pf(F.adj,df1=2,df2=df.error)) # significant
## [1] 0.0156587

# high luminance:
aov.age.high <- aov(sensitivity~age,data=subset(spot,duration=="d3"&luminance=="high"))
summary(aov.age.high)
##             Df Sum Sq Mean Sq F value Pr(>F)
## age          1   8697    8697   1.412  0.247
## Residuals   22 135491    6159
(F.adj <- 8697/MS.error)
## [1] 1.064374
(p.adj <- 1-pf(F.adj,df1=2,df2=df.error)) # not significant
## [1] 0.3478852

2 Story Presentation and Organization

An experiment examined how comprehension of a statistics lecture is affected by the mode of presentation and its organization. Undergraduate students were presented with material on the analysis of variance, and their comprehension of the material was tested one day later. The lecture was presented in three modes: audio only (i.e., a podcast), visual only (i.e., a written document), and audiovisual (i.e., a combination of podcast and written document). In addition, the statistics material was either organized or disorganized. The data are stored in the data frame learning, which contains the factors sID (subject ID), mode, and material, and the dependent variable comprehension.

Use a 2 x 3 factorial ANOVA to analyze the effects of material and mode (and their interaction) on memory.

summary(learning)
##       sID     comprehension             mode            material 
##  s1     : 1   Min.   : 64.77   audio      :20   organized   :30  
##  s2     : 1   1st Qu.: 89.32   visual     :20   disorganized:30  
##  s3     : 1   Median :101.38   audiovisual:20                    
##  s4     : 1   Mean   :100.00                                     
##  s5     : 1   3rd Qu.:110.44                                     
##  s6     : 1   Max.   :132.19                                     
##  (Other):54                                                      
##                       group   
##  audio.organized         :10  
##  visual.organized        :10  
##  audiovisual.organized   :10  
##  audio.disorganized      :10  
##  visual.disorganized     :10  
##  audiovisual.disorganized:10  
## 
learn.aov.01 <- aov(comprehension~material*mode,data=learning)
summary(learn.aov.01)
##               Df Sum Sq Mean Sq F value   Pr(>F)    
## material       1   3532    3532  26.425 3.89e-06 ***
## mode           2   1249     625   4.674   0.0134 *  
## material:mode  2   1277     639   4.778   0.0123 *  
## Residuals     54   7217     134                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Analyze the simple main effect of material for each mode.

aov.material.audio <- aov(comprehension~material,data=subset(learning,mode=="audio"))
aov.material.visual <- aov(comprehension~material,data=subset(learning,mode=="visual"))
aov.material.audiovis <- aov(comprehension~material,data=subset(learning,mode=="audiovisual"))
summary(aov.material.audio)
##             Df Sum Sq Mean Sq F value   Pr(>F)    
## material     1   2950  2949.7   17.06 0.000629 ***
## Residuals   18   3113   172.9                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(aov.material.visual)
##             Df Sum Sq Mean Sq F value  Pr(>F)   
## material     1   1824  1824.1   14.61 0.00125 **
## Residuals   18   2247   124.9                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(aov.material.audiovis)
##             Df Sum Sq Mean Sq F value Pr(>F)
## material     1   34.9   34.93   0.339  0.568
## Residuals   18 1856.4  103.13

# if you recalculate F & p values...
MS.residuals <- 134
df.residuals <- 54
(F.adj <- c(2949.7, 1824.1, 34.93)/MS.residuals)
## [1] 22.0126866 13.6126866  0.2606716
(p.adj <- 1-pf(F.adj,df1=1,df2=df.residuals))
## [1] 1.888876e-05 5.234969e-04 6.117400e-01
p.adj < .05
## [1]  TRUE  TRUE FALSE

# with emmeans (uses recalculated F & p values)
require(emmeans)
learn.aov.01 <- aov(comprehension~material*mode,data=learning)
learn.emm <- emmeans(learn.aov.01,specs=~material*mode)
joint_tests(learn.emm,by="mode")
## mode = audio:
##  model term df1 df2 F.ratio p.value
##  material     1  54  22.071  <.0001
## 
## mode = visual:
##  model term df1 df2 F.ratio p.value
##  material     1  54  13.649  0.0005
## 
## mode = audiovisual:
##  model term df1 df2 F.ratio p.value
##  material     1  54   0.261  0.6113

When evaluating simple main effects, we sometimes choose to recalculate the \(F\) and \(p\) values using values of MS-residuals and df-residuals from the main ANOVA. Why would we want to do that?

# we will dicuss this answer in class

Evaluate the assumptions of normality and homogeneity of variance.

shapiro.test(residuals(learn.aov.01))
## 
##  Shapiro-Wilk normality test
## 
## data:  residuals(learn.aov.01)
## W = 0.96808, p-value = 0.1174

learning$group <- interaction(learning$mode,learning$material)
bartlett.test(residuals(learn.aov.01)~learning$group)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  residuals(learn.aov.01) by learning$group
## Bartlett's K-squared = 2.172, df = 5, p-value = 0.8249

For a factorial design, it is possible to calculate an association strength for each main effect and interaction using eta-squared and partial eta-squared. For any given effect, will eta-squared or partial eta-squared be larger? Why?

# for an AxB design, eta-squared for A is the ratio of SS_A divided by SS_Total (i.e., SS_A+SS_B+SS_AxB+SS_residuals)
# however, partial eta-squared for A is SS_A divided by (SS_A + SS_residuals)
# so the denominator is smaller for partial eta-squared, and therefore partial eta-squared > eta-squared be

Calculate partial omega-squared and partial Cohen’s \(f\) for each main effect and interaction.

library(effectsize)
omega_squared(learn.aov.01)
## # Effect Size for ANOVA (Type I)
## 
## Parameter     | Omega2 (partial) |       95% CI
## -----------------------------------------------
## material      |             0.30 | [0.14, 1.00]
## mode          |             0.11 | [0.00, 1.00]
## material:mode |             0.11 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].
cohens_f(learn.aov.01)
## # Effect Size for ANOVA (Type I)
## 
## Parameter     | Cohen's f (partial) |      95% CI
## -------------------------------------------------
## material      |                0.70 | [0.45, Inf]
## mode          |                0.42 | [0.14, Inf]
## material:mode |                0.42 | [0.15, Inf]
## 
## - One-sided CIs: upper bound fixed at [Inf].

3 Anorexia study

The data frame anorexia contains the results from an experiment that examined the effectiveness of three types of therapy (behavioural, cognitive, and psychodynamic) for inducing weight loss (wt) in adolescent girls who were high or low in neuroticism. The factors therapy and neuroticism were crossed factorially, but the design was unbalanced.

Verify that the design is unbalanced.

xtabs(~therapy+neuroticism,data=anorexia)
##                neuroticism
## therapy         high low
##   behavioural      8   7
##   cognitve        10  10
##   psychodynamic    8   9

Use factorial ANOVA to evaluate the effects of neuroticism and therapy on weight loss (wt). Set up your ANOVA so that the main effects are evaluated with null hypothesis that the unweighted marginal means are equal.

wt.aov.01 <- aov(wt~neuroticism*therapy,data=anorexia)
drop1(wt.aov.01,.~.,test="F")
## Single term deletions
## 
## Model:
## wt ~ neuroticism * therapy
##                     Df Sum of Sq    RSS    AIC F value    Pr(>F)    
## <none>                           300.40 103.20                      
## neuroticism          1   108.033 408.43 117.18 16.5433 0.0001845 ***
## therapy              2     7.546 307.94 100.49  0.5778 0.5651615    
## neuroticism:therapy  2    57.184 357.58 108.26  4.3783 0.0181709 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# alternative method:
require(car)
Anova(wt.aov.01,type="3")
## Anova Table (Type III tests)
## 
## Response: wt
##                      Sum Sq Df  F value    Pr(>F)    
## (Intercept)         1269.26  1 194.3633 < 2.2e-16 ***
## neuroticism          108.03  1  16.5433 0.0001845 ***
## therapy                7.55  2   0.5778 0.5651615    
## neuroticism:therapy   57.18  2   4.3783 0.0181709 *  
## Residuals            300.40 46                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Within the framework of nested linear models, describe the difference between Type I and Type II sums-of-squares for the main effect of therapy.

# we will discuss this answer in class

Are Type II or Type III sums-of-squares more appropriate for evaluating the main effects of therapy and neuroticism? Why?

# we will discuss this answer in class

PSYCH 710 Practice Questions

1 Spot detection thresholds

2 Story Presentation and Organization

3 Anorexia study