For this exercise, please try to reproduce the results from Experiment 1 of the associated paper (Ko, Sadler & Galinsky, 2015). The PDF of the paper is included in the same folder as this Rmd file.

Methods summary:

A sense of power has often been tied to how we perceive each other’s voice. Social hierarchy is embedded into the structure of society and provides a metric by which others relate to one another. In 1956, the Brunswik Lens Model was introduced to examine how vocal cues might influence hierarchy. In “The Sound of Power: Conveying and Detecting Hierarchical Rank Through Voice,” Ko and colleagues investigated how manipulation of hierarchal rank within a situation might impact vocal acoustic cues. Using the Brunswik Model, six acoustic metrics were utilized (pitch mean & variability, loudness mean & variability, and resonance mean & variability) to isolate a potential contribution between individuals of different hierarchal rank. In the first experiment, Ko, Sadler & Galinsky examined the vocal acoustic cues of individuals before and after being assigned a hierarchal rank in a sample of 161 subjects (80 male). Each of the six hierarchy acoustic cues were analyzed with a 2 (high vs. low rank condition) x 2 (male vs. female) analysis of covariance, controlling for the baseline of the respective acoustic cue.


Target outcomes:

Below is the specific result you will attempt to reproduce (quoted directly from the results section of Experiment 1):

The impact of hierarchical rank on speakers’ acoustic cues. Each of the six hierarchy-based (i.e., postmanipulation) acoustic variables was submitted to a 2 (condition: high rank, low rank) × 2 (speaker’s sex: female, male) between-subjects analysis of covariance, controlling for the corresponding baseline acoustic variable. Table 4 presents the adjusted means by condition. Condition had a significant effect on pitch, pitch variability, and loudness variability. Speakers’ voices in the high-rank condition had higher pitch, F(1, 156) = 4.48, p < .05; were more variable in loudness, F(1, 156) = 4.66, p < .05; and were more monotone (i.e., less variable in pitch), F(1, 156) = 4.73, p < .05, compared with speakers’ voices in the low-rank condition (all other Fs < 1; see the Supplemental Material for additional analyses of covariance involving pitch and loudness). (from Ko et al., 2015, p. 6; emphasis added)

The adjusted means for these analyses are reported in Table 4 (Table4_AdjustedMeans.png, included in the same folder as this Rmd file).


Step 1: Load packages

library(tidyverse) # for data munging
library(knitr) # for kable table formating
library(haven) # import and export 'SPSS', 'Stata' and 'SAS' Files
library(readxl) # import excel files

# #optional packages:
library(psych)
library(car) # for ANCOVA
library(compute.es) # for ANCOVA
library(lsmeans) # for ANCOVA

Step 2: Load data

# Just Experiment 1
d <-read_csv("data/S1_voice_level_Final.csv")
# DT::datatable(d)
d
## # A tibble: 161 × 45
##    voice form_smean form_svar form_rmean form_rvar intense_smean intense_svar
##    <dbl>      <dbl>     <dbl>      <dbl>     <dbl>         <dbl>        <dbl>
##  1     1      1043.    38805.      1259.    68275.          65.5        157. 
##  2     2      1083.    35609.      1278.    54612.          65.5        301. 
##  3     3      1092.    25979.      1334.    55175.          53.6         71.2
##  4     4      1268.    71346.      1298.    74340.          60.7        201. 
##  5     5      1071.    38246.      1256.    67846.          60.1        204. 
##  6     6      1095.    40716.      1278.    76674.          60.1        150. 
##  7     7      1060.    38855.      1310.    88947.          60.5        127. 
##  8     8      1051.    28191.      1223.    52682.          57.3        207. 
##  9     9      1038.    34105.      1256.    68508.          57.9        199. 
## 10    10      1346.    69260.      1314.    55400.          60.3        149. 
## # ℹ 151 more rows
## # ℹ 38 more variables: intense_rmean <dbl>, intense_rvar <dbl>,
## #   pitch_smean <dbl>, pitch_svar <dbl>, pitch_rmean <dbl>, pitch_rvar <dbl>,
## #   pow <dbl>, age <dbl>, sex <chr>, race <chr>, native <chr>, feelpower <dbl>,
## #   plev <dbl>, vsex <dbl>, pitch_rmeanMD <dbl>, pitch_rvarMD <dbl>,
## #   intense_rmeanMD <dbl>, intense_rvarMD <dbl>, formant_rmeanMD <dbl>,
## #   formant_rvarMD <dbl>, pitch_smeanMD <dbl>, pitch_svarMD <dbl>, …

Step 3: Tidy data

d_tidy <- d %>%
  pivot_longer(
    cols = -c(voice, pow, age, sex, race, native, feelpower, plev, vsex),
    names_to = "variable",
    values_to = "value"
  )
d_tidy
## # A tibble: 5,796 × 11
##    voice   pow   age sex   race  native feelpower  plev  vsex variable     value
##    <dbl> <dbl> <dbl> <chr> <chr> <chr>      <dbl> <dbl> <dbl> <chr>        <dbl>
##  1     1     1    19 M     X     Y              6     1    -1 form_smean  1.04e3
##  2     1     1    19 M     X     Y              6     1    -1 form_svar   3.88e4
##  3     1     1    19 M     X     Y              6     1    -1 form_rmean  1.26e3
##  4     1     1    19 M     X     Y              6     1    -1 form_rvar   6.83e4
##  5     1     1    19 M     X     Y              6     1    -1 intense_sm… 6.55e1
##  6     1     1    19 M     X     Y              6     1    -1 intense_sv… 1.57e2
##  7     1     1    19 M     X     Y              6     1    -1 intense_rm… 5.87e1
##  8     1     1    19 M     X     Y              6     1    -1 intense_rv… 1.74e2
##  9     1     1    19 M     X     Y              6     1    -1 pitch_smean 1.16e2
## 10     1     1    19 M     X     Y              6     1    -1 pitch_svar  8.58e2
## # ℹ 5,786 more rows

Step 4: Run analysis

Pre-processing

# Verify that there are 161 participants and check distribution of plev and sex
nrow(d)
## [1] 161
table(d$plev, d$vsex)
##     
##      -1  1
##   -1 38 41
##   1  42 40

Descriptive statistics

In the paper, the adjusted means by condition are reported (see Table 4, or Table4_AdjustedMeans.png, included in the same folder as this Rmd file). Reproduce these values below:

# Pitch Mean
model_pitch_mean <- lm(pitch_smean ~ pitch_rmean + plev * vsex, data = d)
lsm_pitch_mean <- lsmeans(model_pitch_mean, specs = "plev")

# Pitch Variability
model_pitch_var <- lm(pitch_svar ~ pitch_rvar + plev * vsex, data = d)
lsm_pitch_var <- lsmeans(model_pitch_var, specs = "plev")

# Loudness Mean
model_intensity_mean <- lm(intense_smean ~ intense_rmean + plev * vsex, data = d)
lsm_intensity_mean <- lsmeans(model_intensity_mean, specs = "plev")

# Loudness Variability 
model_intensity_var <- lm(intense_svar ~ intense_rvar + plev * vsex, data = d)
lsm_intensity_var <- lsmeans(model_intensity_var, specs = "plev")

# Resonance Mean 
model_formant_mean <- lm(form_smean ~ form_rmean + plev * vsex, data = d)
lsm_formant_mean <- lsmeans(model_formant_mean, specs = "plev")

# Resonance Variability 
model_formant_var <- lm(form_svar ~ form_rvar + plev * vsex, data = d)
lsm_formant_var <- lsmeans(model_formant_var, specs = "plev")


print("====High Rank Condition====")
## [1] "====High Rank Condition===="
cat("Pitch mean:", summary(lsm_pitch_mean)$lsmean[2], "\n")
## Pitch mean: 158.6098
cat("Pitch variability:", summary(lsm_pitch_var)$lsmean[2], "\n")
## Pitch variability: 1425.016
cat("Loudness mean:", summary(lsm_intensity_mean)$lsmean[2], "\n")
## Loudness mean: 59.33567
cat("Loudness variability:", summary(lsm_intensity_var)$lsmean[2], "\n")
## Loudness variability: 196.7301
cat("Resonance mean:", summary(lsm_formant_mean)$lsmean[2], "\n")
## Resonance mean: 1129.384
cat("Resonance variability:", summary(lsm_formant_var)$lsmean[2], "\n")
## Resonance variability: 42170.78
print("====Low Rank Condition====")
## [1] "====Low Rank Condition===="
cat("Pitch mean:", summary(lsm_pitch_mean)$lsmean[1], "\n")
## Pitch mean: 155.5227
cat("Pitch variability:", summary(lsm_pitch_var)$lsmean[1], "\n")
## Pitch variability: 1648.367
cat("Loudness mean:", summary(lsm_intensity_mean)$lsmean[1], "\n")
## Loudness mean: 58.66784
cat("Loudness variability:", summary(lsm_intensity_var)$lsmean[1], "\n")
## Loudness variability: 183.4795
cat("Resonance mean:", summary(lsm_formant_mean)$lsmean[1], "\n")
## Resonance mean: 1128.806
cat("Resonance variability:", summary(lsm_formant_var)$lsmean[1], "\n")
## Resonance variability: 43654.54

Inferential statistics

The impact of hierarchical rank on speakers’ acoustic cues. Each of the six hierarchy-based (i.e., postmanipulation) acoustic variables was submitted to a 2 (condition: high rank, low rank) × 2 (speaker’s sex: female, male) between-subjects analysis of covariance, controlling for the corresponding baseline acoustic variable. […] Condition had a significant effect on pitch, pitch variability, and loudness variability. Speakers’ voices in the high-rank condition had higher pitch, F(1, 156) = 4.48, p < .05; were more variable in loudness, F(1, 156) = 4.66, p < .05; and were more monotone (i.e., less variable in pitch), F(1, 156) = 4.73, p < .05, compared with speakers’ voices in the low-rank condition (all other Fs < 1; see the Supplemental Material for additional analyses of covariance involving pitch and loudness).

# reproduce the above results here

print("=== Pitch Mean ===")
## [1] "=== Pitch Mean ==="
Anova(model_pitch_mean, type = 2)
## Anova Table (Type II tests)
## 
## Response: pitch_smean
##             Sum Sq  Df  F value  Pr(>F)    
## pitch_rmean  43334   1 511.2208 < 2e-16 ***
## plev           379   1   4.4744 0.03599 *  
## vsex           460   1   5.4218 0.02117 *  
## plev:vsex      222   1   2.6156 0.10783    
## Residuals    13223 156                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
print("=== Pitch Variability ===")
## [1] "=== Pitch Variability ==="
Anova(model_pitch_var, type = 2)
## Anova Table (Type II tests)
## 
## Response: pitch_svar
##              Sum Sq  Df F value    Pr(>F)    
## pitch_rvar 27760551   1 65.4850 1.548e-13 ***
## plev        2006828   1  4.7340   0.03108 *  
## vsex       26142260   1 61.6676 6.133e-13 ***
## plev:vsex     37857   1  0.0893   0.76546    
## Residuals  66131879 156                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
print("=== Loudness Variability ===")
## [1] "=== Loudness Variability ==="
Anova(model_intensity_var, type = 2)
## Anova Table (Type II tests)
## 
## Response: intense_svar
##              Sum Sq  Df  F value    Pr(>F)    
## intense_rvar 184927   1 122.4323 < 2.2e-16 ***
## plev           7052   1   4.6687   0.03224 *  
## vsex          29860   1  19.7693 1.651e-05 ***
## plev:vsex        40   1   0.0262   0.87160    
## Residuals    235630 156                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Step 5: Reflection

Were you able to reproduce the results you attempted to reproduce? If not, what part(s) were you unable to reproduce?

Yes!

How difficult was it to reproduce your results?

Quite difficult.

What aspects made it difficult? What aspects made it easy?

The most difficult aspect was correctly identifying which variables were the baseline vs. post-manipulation measures. I looked through the codebook which helped me resolve the discrepancies and produce the correct results. I am also confused about tidying the data, since the way I did it ended up not being helpful for the computations needed.