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(lsr)
# #optional packages:
library(psych)
library(car) # for ANCOVA
library(compute.es) # for ANCOVA
library(lsmeans) # for ANCOVA# Just Experiment 1
d <-read_csv("data/S1_voice_level_Final.csv")
# DT::datatable(d)full_data <- select(d, c("voice", "age", "plev", ends_with("mean"), ends_with("var"), -starts_with("Z")))# Pitch features
P_mean = as.factor(full_data$plev)
P_var = as.factor(full_data$plev)
PM <- data.frame(P_mean, select(full_data,c( "pitch_rmean","pitch_smean")))
PV <- data.frame(P_var,select(full_data,c("pitch_rvar","pitch_svar")))
# Loudness features
L_mean = as.factor(full_data$plev)
L_var = as.factor(full_data$plev)
LM <- data.frame(L_mean, select(full_data,c( "intense_rmean","intense_smean")))
LV <- data.frame(L_var,select(full_data,c("intense_rvar","intense_svar")))
# Resonance features
R_mean = as.factor(full_data$plev)
R_var = as.factor(full_data$plev)
RM <- data.frame(R_mean, select(full_data,c( "form_rmean","form_smean")))
RV <- data.frame(R_var,select(full_data,c("form_rvar","form_svar")))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 Adjusted Means
anovapm <- aov(PM$pitch_smean ~ P_mean + PM$pitch_rmean, data = PM)
adj_m = emmeans(anovapm, ~ P_mean)
pmean_means <- tapply(PM$pitch_smean, PM$P_mean, mean)
cov_rain_mean <- mean(PM$pitch_rmean)
sprintf("Mean of Pitch for High rank: %f ", pmean_means[2])## [1] "Mean of Pitch for High rank: 155.965473 "sprintf("Mean of Pitch for Low rank: %f ", pmean_means[1])## [1] "Mean of Pitch for Low rank: 158.239942 "sprintf("Adjusted Mean of Pitch for High rank: %f ", summary(adj_m)$emmean[2])## [1] "Adjusted Mean of Pitch for High rank: 158.706800 "sprintf("Adjusted Mean of Pitch for Low rank: %f ", summary(adj_m)$emmean[1])## [1] "Adjusted Mean of Pitch for Low rank: 155.394514 "sprintf("Effect condition for Pitch : %s ", summary(etaSquared(anovapm))[1,2])## [1] "Effect condition for Pitch : Min. :0.03064 "sprintf("------------------------------------------")## [1] "------------------------------------------"# Pitch Variance Adjusted Means
anovapv <- aov(PV$pitch_svar ~ P_var + PV$pitch_rvar, data = PV)
adj_pv = emmeans(anovapv, ~ P_var)
pvar_means <- tapply(PV$pitch_svar, PV$P_var, mean)
cov_rain_mean <- mean(PV$pitch_rvar)
sprintf("Mean of Pitch Variance for High rank: %f ", pvar_means[2])## [1] "Mean of Pitch Variance for High rank: 1402.644728 "sprintf("Mean of Pitch Variance for Low rank: %f ", pvar_means[1])## [1] "Mean of Pitch Variance for Low rank: 1677.755540 "sprintf("Adjusted Mean of Pitch Variance for High rank: %f ", summary(adj_pv)$emmean[2])## [1] "Adjusted Mean of Pitch Variance for High rank: 1417.621220 "sprintf("Adjusted Mean of Pitch Variance for Low rank: %f ", summary(adj_pv)$emmean[1])## [1] "Adjusted Mean of Pitch Variance for Low rank: 1662.210321 "sprintf("Effect condition for Pitch Variance : %s ", summary(etaSquared(anovapv))[1,2])## [1] "Effect condition for Pitch Variance : Min. :0.02539 "sprintf("------------------------------------------")## [1] "------------------------------------------"# Loundess Adjusted Means
anovalm <- aov(LM$intense_smean ~ L_mean + LM$intense_rmean, data = LM)
adj_lm = emmeans(anovalm, ~ L_mean)
lmean_means <- tapply(LM$intense_smean, LM$L_mean, mean)
cov_rain_mean <- mean(LM$intense_rmean)
sprintf("Mean of Loudness for High rank: %f ", lmean_means[2])## [1] "Mean of Loudness for High rank: 59.364195 "sprintf("Mean of Loudness for Low rank: %f ", lmean_means[1])## [1] "Mean of Loudness for Low rank: 58.665912 "sprintf("Adjusted Mean of Loudness for High rank: %f ", summary(adj_lm)$emmean[2])## [1] "Adjusted Mean of Loudness for High rank: 59.324878 "sprintf("Adjusted Mean of Loudness for Low rank: %f ", summary(adj_lm)$emmean[1])## [1] "Adjusted Mean of Loudness for Low rank: 58.706723 "sprintf("Effect condition for Loudness : %s ", summary(etaSquared(anovalm))[1,2])## [1] "Effect condition for Loudness : Min. :0.01101 "sprintf("------------------------------------------")## [1] "------------------------------------------"# Loudness Variance Adjusted Means
anovalv <- aov(LV$intense_svar ~ L_var + LV$intense_rvar, data = LV)
adj_lv <- emmeans(anovalv, ~ L_var)
lvar_means <- tapply(LV$intense_svar, LV$L_var, mean)
cov_rain_mean <- mean(LV$intense_rvar)
sprintf("Mean of Loudness Variance for High rank: %f ", lvar_means[2])## [1] "Mean of Loudness Variance for High rank: 198.464118 "sprintf("Mean of Loudness Variance for Low rank: %f ", lvar_means[1])## [1] "Mean of Loudness Variance for Low rank: 181.473818 "sprintf("Adjusted Mean of Loudness Variance for High rank: %f ", summary(adj_lv)$emmean[2])## [1] "Adjusted Mean of Loudness Variance for High rank: 197.132794 "sprintf("Adjusted Mean of Loudness Variance for Low rank: %f ", summary(adj_lv)$emmean[1])## [1] "Adjusted Mean of Loudness Variance for Low rank: 182.855699 "sprintf("Effect condition for Loudness Variance : %s ", summary(etaSquared(anovalv))[1,2])## [1] "Effect condition for Loudness Variance : Min. :0.02991 "sprintf("------------------------------------------")## [1] "------------------------------------------"# Resonance Adjusted Means
anovarm <- aov(RM$form_smean ~ R_mean + RM$form_rmean, data = RM)
adj_rm = emmeans(anovarm, ~ R_mean)
rmean_means <- tapply(RM$form_smean, RM$R_mean, mean)
cov_rain_mean <- mean(RM$form_rmean)
sprintf("Mean of Resonance for High rank: %f ", rmean_means[2])## [1] "Mean of Resonance for High rank: 1127.423151 "sprintf("Mean of Resonance for Low rank: %f ", rmean_means[1])## [1] "Mean of Resonance for Low rank: 1131.192650 "sprintf("Adjusted Mean of Resonance for High rank: %f ", summary(adj_rm)$emmean[2])## [1] "Adjusted Mean of Resonance for High rank: 1129.045306 "sprintf("Adjusted Mean of Resonance for Low rank: %f ", summary(adj_rm)$emmean[1])## [1] "Adjusted Mean of Resonance for Low rank: 1129.508894 "sprintf("Effect condition for Resonance : %s ", summary(etaSquared(anovarm))[1,2])## [1] "Effect condition for Resonance : Min. :7.360e-06 "sprintf("------------------------------------------")## [1] "------------------------------------------"# Resonance Variance Adjusted Means
anovarv <- aov(RV$form_svar ~ R_var + RV$form_rvar, data = RV)
adj_rv = emmeans(anovarv, ~ R_var)
rvar_means <- tapply(RV$form_svar, RV$R_var, mean)
cov_rain_mean <- mean(RV$form_rvar)
sprintf("Mean of Resonance Variance for High rank: %f ", rvar_means[2])## [1] "Mean of Resonance Variance for High rank: 42020.666111 "sprintf("Mean of Resonance Variance for Low rank: %f ", rvar_means[1])## [1] "Mean of Resonance Variance for Low rank: 43837.423735 "sprintf("Adjusted Mean of Resonance Variance for High rank: %f ", summary(adj_rv)$emmean[2])## [1] "Adjusted Mean of Resonance Variance for High rank: 42143.429922 "sprintf("Adjusted Mean of Resonance Variance for Low rank: %f ", summary(adj_rv)$emmean[1])## [1] "Adjusted Mean of Resonance Variance for Low rank: 43709.998008 "sprintf("Effect condition for Resonance Variance : %s ", summary(etaSquared(anovarv))[1,2])## [1] "Effect condition for Resonance Variance : Min. :0.003362 "sprintf("------------------------------------------")## [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. […] 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
# Pitch statistics
f_pm <- summary(anovapm)[[1]][1,4]
p_pm <- summary(anovapm)[[1]][1,5]
sprintf("--------------------------------")## [1] "--------------------------------"sprintf("Statistics for pitch")## [1] "Statistics for pitch"sprintf("F statistic: %f", f_pm)## [1] "F statistic: 2.365219"sprintf("p-value: %f", p_pm)## [1] "p-value: 0.126067"# Pitch variance statistics
f_pv <- summary(anovapv)[[1]][1,4]
p_pv <- summary(anovapv)[[1]][1,5]
sprintf("--------------------------------")## [1] "--------------------------------"sprintf("Statistics for pitch variance")## [1] "Statistics for pitch variance"sprintf("F statistic: %f", f_pv)## [1] "F statistic: 5.212299"sprintf("p-value: %f *", p_pv)## [1] "p-value: 0.023761 *"# Loudness statistics
f_lm <- summary(anovalm)[[1]][1,4]
p_lm <- summary(anovalm)[[1]][1,5]
sprintf("--------------------------------")## [1] "--------------------------------"sprintf("Statistics for loudness")## [1] "Statistics for loudness"sprintf("F statistic: %f", f_lm)## [1] "F statistic: 2.245975"sprintf("p-value: %f", p_lm)## [1] "p-value: 0.135958"# Loudness variance statistics
f_lv <- summary(anovalv)[[1]][1,4]
p_lv <- summary(anovalv)[[1]][1,5]
sprintf("--------------------------------")## [1] "--------------------------------"sprintf("Statistics for loudness variance")## [1] "Statistics for loudness variance"sprintf("F statistic: %f", f_lv)## [1] "F statistic: 6.911326"sprintf("p-value: %f **", p_lv)## [1] "p-value: 0.009411 **"Were you able to reproduce the results you attempted to reproduce? If not, what part(s) were you unable to reproduce?
I was not able to reproduce the exact F statistic and significance for p values. The adjusted means are very close to the results quoted in the experiment. The effect conditions, when rounded, match the eta squared given in the paper.
How difficult was it to reproduce your results?
This reproduction was much harder than the paper from Group A. The paper is not clear about the columns to be taken from the dataset for analysis and does not give too many numerical reference measures during the various steps. The different steps to reproducing the analysis could have been conveyed better in the paper. It took me 5 hours to reproduce the results
What aspects made it difficult? What aspects made it easy?
Difficult: I struggled with reproducing the exact F statistics and p-values. Initially, I manually implemented adjusted means and ran into errors. Later, I modified it to a function that gave me results closer to the original. Easy: The table was a concise method to compare my results with the paper’s results. The word document was also helfpul in interpretting the different columns in the data.