def cohens_d(X, Y):
return (
- (mean(X) + mean(Y)) /
+ (mean(X) - mean(Y)) /
sqrt(
(len(X)*variance(X) + len(Y)*variance(Y)) /
(len(X) + len(Y))
)
)
-def population_with_error(μ, σ, n):
+def population_with_error(μ, ε, n):
def trait():
return normal(μ, 1)
def measurement_error():
- return normal(0, σ)
+ return normal(0, ε)
return [trait() + measurement_error() for _ in range(n)]
# trait differs by 1 standard deviation
-adjusted_f = population_with_error(1, 0, 10000)
-adjusted_m = population_with_error(0, 0, 10000)
+true_f = population_with_error(1, 0, 10000)
+true_m = population_with_error(0, 0, 10000)
# as above, but with 0.5 standard units measurment error
measured_f = population_with_error(1, 0.5, 10000)
measured_m = population_with_error(0, 0.5, 10000)
-smart_d = cohens_d(adjusted_f, adjusted_m)
-print(smart_d) # 1.0193773432617055 — d≈1.0, as expected!
+true_d = cohens_d(true_f, true_m)
+print(true_d) # 1.0193773432617055 — d≈1.0, as expected!
naïve_d = cohens_d(measured_f, measured_m)
-print(naïve_d) # 0.8953395386313235
+print(naïve_d) # 0.8953395386313235 — deflated!
-def performance(g, σ_g, s, n):
+def performance(μ_g, σ_g, s, n):
def general_ability():
- return normal(g, σ_g)
+ return normal(μ_g, σ_g)
def special_ability():
return normal(s, 1)
return [general_ability() + special_ability() for _ in range(n)]
matched_m = performance(0, 0, 0, 10000)
population_d = cohens_d(population_f, population_m)
-print(population_d) # 0.7287587808164793
+print(population_d) # 0.7287587808164793 — deflated!
matched_d = cohens_d(matched_f, matched_m)
-print(matched_d) # 1.018362581243161
+print(matched_d) # 1.018362581243161 — as you would expect