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.
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.
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).
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
library(afex)
# library(emmeans) # for ANCOVA
# #optional packages:
# library(psych)
library(car) # for ANCOVA
# library(compute.es) # for ANCOVA
# library(msmeans) # outdated package replaced by emmeans
# Just Experiment 1
d <-read_csv("data/S1_voice_level_Final.csv")
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>, …
# DT::datatable(d)
#Acoustic variable names
#‘form’ = formant disperson/measure of resonance
#‘intense’ = intensity/loudness
#‘pitch’= pitch
#‘_s’ = script, acoustics recorded while saying negotiation script
#‘_r’ = rainbow, acoustics recorded during rainbow passage
#‘mean’ = average value of acoustic during a clip
#‘var’ = variance for acoustic during a clip
#‘MD’ = mean-deviated or centered
#‘Z’ = standardized version of acoustic
filtered_d <- d |>
select(c("voice","sex","plev","pitch_rmean","pitch_smean","pitch_rvar","pitch_svar","intense_rmean","intense_smean","intense_rvar","intense_svar","form_rmean","form_smean","form_rvar","form_svar"))
# Males recoded to 1, females recoded to 2
# filtered_d$sex <- ifelse(filtered_d$sex == 'M', "1", "2")
# filtered_d$sex <- as.numeric(filtered_d$sex)
filtered_d
## # A tibble: 161 × 15
## voice sex plev pitch_rmean pitch_smean pitch_rvar pitch_svar intense_rmean
## <dbl> <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 1 M 1 109. 116. 2482. 858. 58.7
## 2 2 M 1 133. 146. 1231. 1404. 63.9
## 3 3 M 1 108. 106. 229. 706. 53.0
## 4 4 M 1 108. 102. 1744. 450. 59.8
## 5 5 M 1 115. 122. 946. 478. 57.5
## 6 6 M 1 127. 137. 420. 623. 58.4
## 7 7 M 1 115. 106. 2267. 313. 57.9
## 8 8 M -1 101. 105. 2898. 1120. 57.3
## 9 9 M 1 162. 161. 697. 444. 58.1
## 10 10 M -1 118. 132. 185. 640. 58.4
## # ℹ 151 more rows
## # ℹ 7 more variables: intense_smean <dbl>, intense_rvar <dbl>,
## # intense_svar <dbl>, form_rmean <dbl>, form_smean <dbl>, form_rvar <dbl>,
## # form_svar <dbl>
any(is.na(filtered_d))
## [1] FALSE
# Prepare data for ANCOVA
# Calculate the means for each condition by subsetting the data
# condition_values <- unique(filtered_d$plev) # Get unique condition values
# Initialize an empty vector to store the means
# condition_means <- numeric(length(condition_values))
# filtered_d['diff_pitchmean'] = filtered_d['pitch_smean'] - filtered_d['pitch_rmean']
# filtered_d['diff_pitchvar'] = filtered_d['pitch_svar'] - filtered_d['pitch_rvar']
# filtered_d['diff_intensemean'] = filtered_d['intense_smean'] - filtered_d['intense_rmean']
# filtered_d['diff_intensevar'] = filtered_d['intense_svar'] - filtered_d['intense_rvar']
# filtered_d['diff_formmean'] = filtered_d['form_smean'] - filtered_d['form_rmean']
# filtered_d['diff_formvar'] = filtered_d['form_svar'] - filtered_d['form_rvar']
#
# filtered_d |>
# group_by(plev) |>
# summarize(
# mean_pitch = mean(diff_pitchmean),
# mean_pitchvar = mean(diff_pitchvar),
# mean_intense = mean(diff_intensemean),
# mean_intensevar = mean(diff_intensevar),
# mean_pitchs = mean(diff_formmean),
# mean_pitchvars = mean(diff_formvar)
# )
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:
# 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.
# Run the ANCOVA for each of the six hierarchy-based acoustic variables
pitchmean_model <- lm(pitch_smean ~ plev * sex + pitch_rmean, data=filtered_d)
ancova_pitch_result <- Anova(pitchmean_model, type="III")
# Fit an ANCOVA model using aov_ez
# pitchmean_model <- aov_ez(id = "voice", dv = "pitch_smean", data = filtered_d, between = c("plev", "sex"), covariate = c("pitch_rmean"))
pitchvar_model <- lm(pitch_svar ~ plev * sex + pitch_rvar, data=filtered_d)
intensemean_model <- lm(intense_smean ~ plev * sex + intense_rmean, data=filtered_d)
intensevar_model <- lm(intense_svar ~ plev * sex + intense_rvar, data=filtered_d)
formmean_model <- lm(form_smean ~ plev * sex + form_rmean, data=filtered_d)
formvar_model <- lm(form_svar ~ plev * sex + form_rvar, data=filtered_d)
# Display the summary of the ANCOVA for each variable
print(ancova_pitch_result) # summary(pitchmean_model) #anova(pitchmean_model)
## Anova Table (Type III tests)
##
## Response: pitch_smean
## Sum Sq Df F value Pr(>F)
## (Intercept) 379 1 4.4744 0.03600 *
## plev 590 1 6.9640 0.00916 **
## sex 396 1 4.6687 0.03224 *
## pitch_rmean 43334 1 511.2208 < 2e-16 ***
## plev:sex 222 1 2.6156 0.10783
## Residuals 13223 156
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(pitchvar_model)
##
## Call:
## lm(formula = pitch_svar ~ plev * sex + pitch_rvar, data = filtered_d)
##
## Residuals:
## Min 1Q Median 3Q Max
## -1133.0 -428.0 -94.1 292.2 3651.7
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1344.79878 110.60029 12.159 < 2e-16 ***
## plev -127.07528 72.56901 -1.751 0.0819 .
## sexM -820.49821 104.40916 -7.858 5.94e-13 ***
## pitch_rvar 0.34360 0.04246 8.092 1.55e-13 ***
## plev:sexM 30.80017 103.06722 0.299 0.7655
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 651.1 on 156 degrees of freedom
## Multiple R-squared: 0.5111, Adjusted R-squared: 0.4985
## F-statistic: 40.77 on 4 and 156 DF, p-value: < 2.2e-16
summary(intensemean_model)
##
## Call:
## lm(formula = intense_smean ~ plev * sex + intense_rmean, data = filtered_d)
##
## Residuals:
## Min 1Q Median 3Q Max
## -7.3607 -1.7424 -0.1007 1.5687 6.8850
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 18.25180 3.37665 5.405 2.38e-07 ***
## plev 0.06325 0.31684 0.200 0.842041
## sexM -1.66687 0.46567 -3.580 0.000459 ***
## intense_rmean 0.72368 0.05958 12.146 < 2e-16 ***
## plev:sexM 0.54133 0.45074 1.201 0.231577
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.847 on 156 degrees of freedom
## Multiple R-squared: 0.4987, Adjusted R-squared: 0.4858
## F-statistic: 38.8 on 4 and 156 DF, p-value: < 2.2e-16
summary(intensevar_model)
##
## Call:
## lm(formula = intense_svar ~ plev * sex + intense_rvar, data = filtered_d)
##
## Residuals:
## Min 1Q Median 3Q Max
## -126.73 -22.63 -5.05 20.21 130.57
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 57.56283 11.96221 4.812 3.50e-06 ***
## plev 7.12182 4.31875 1.649 0.101
## sexM 27.53904 6.19014 4.449 1.63e-05 ***
## intense_rvar 0.64991 0.05874 11.065 < 2e-16 ***
## plev:sexM -0.99304 6.13391 -0.162 0.872
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 38.86 on 156 degrees of freedom
## Multiple R-squared: 0.4709, Adjusted R-squared: 0.4573
## F-statistic: 34.71 on 4 and 156 DF, p-value: < 2.2e-16
summary(formmean_model)
##
## Call:
## lm(formula = form_smean ~ plev * sex + form_rmean, data = filtered_d)
##
## Residuals:
## Min 1Q Median 3Q Max
## -112.665 -52.619 -25.835 -0.345 229.149
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 458.6088 194.9988 2.352 0.019929 *
## plev -0.7387 9.2891 -0.080 0.936722
## sexM -45.2304 13.3128 -3.398 0.000863 ***
## form_rmean 0.5362 0.1500 3.576 0.000466 ***
## plev:sexM 2.0558 13.1922 0.156 0.876362
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 83.59 on 156 degrees of freedom
## Multiple R-squared: 0.1537, Adjusted R-squared: 0.132
## F-statistic: 7.082 on 4 and 156 DF, p-value: 2.892e-05
summary(formvar_model)
##
## Call:
## lm(formula = form_svar ~ plev * sex + form_rvar, data = filtered_d)
##
## Residuals:
## Min 1Q Median 3Q Max
## -17761 -8627 -4289 2159 59762
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.362e+04 5.726e+03 5.870 2.53e-08 ***
## plev -1.051e+03 1.517e+03 -0.693 0.489
## sexM -1.186e+03 2.403e+03 -0.493 0.622
## form_rvar 1.542e-01 9.431e-02 1.635 0.104
## plev:sexM 6.183e+02 2.159e+03 0.286 0.775
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 13650 on 156 degrees of freedom
## Multiple R-squared: 0.02184, Adjusted R-squared: -0.00324
## F-statistic: 0.8708 on 4 and 156 DF, p-value: 0.4829
# Calculate the pitch means for each condition
# for (i in 1:length(condition_values)) {
# condition_means[i] <- mean(filtered_d$pitch_smean[filtered_d$plev == condition_values[i]])
# }
#
# Create a 1x2 output table
# output_table <- data.frame(
# "Condition" = condition_values,
# "Pitch Mean" = condition_means
# )
#
# Print the output table
# print(output_table)
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
Were you able to reproduce the results you attempted to reproduce? If not, what part(s) were you unable to reproduce?
No. I couldn’t even make Table 4 even after 3 hours.
How difficult was it to reproduce your results?
I tried different ways of computing means, including taking the difference between postmanipulation and baseline, as well as running various ANCOVAs with lm, aov, and aov_ez. Even when reading the supplementary materials on OSF, I could not figure out if the ANCOVAs were performed using the raw values, MD values, or Z-transformed values without resorting to trial and error.
What aspects made it difficult? What aspects made it easy?
Difficulties were detailed above, I didn’t find it easy at all.