Sex differences
This section demonstrates several different ways to estimate sex differences in personality. We begin with typical univariate differences in the broad Big Five domains, followed by multivariate differences in both domains and using personality items as a multivariate set. We will pay special attention to regularized methods because this approach enables us to examine both sex differences and the underlying distributions that give rise to these mean differences.
Univariate differences in trait domain sum scores using the d_pooled_sd -function
The d_pooled_sd function is used to calculate Cohen’s d with a pooled standard deviation across the male and female groups for each Big Five trait separately. The function takes several arguments, including the analyzed Big Five trait (var), sex as the binary grouping variable (group.var), the values that represent these groups (group.values), and a logical value indicating whether statistical inference is desired (infer). All arguments except for statistical inference are also included in the multivariate function D_regularized, which is presented further below.
d.N<-d_pooled_sd(data = US.dat,var = "N",
group.var = "gender.cat",group.values = c("Male","Female"),
infer=T)
d.E<-d_pooled_sd(data = US.dat,var = "E",
group.var = "gender.cat",group.values = c("Male","Female"),
infer=T)
d.O<-d_pooled_sd(data = US.dat,var = "O",
group.var = "gender.cat",group.values = c("Male","Female"),
infer=T)
d.A<-d_pooled_sd(data = US.dat,var = "A",
group.var = "gender.cat",group.values = c("Male","Female"),
infer=T)
d.C<-d_pooled_sd(data = US.dat,var = "C",
group.var = "gender.cat",group.values = c("Male","Female"),
infer=T)
# compile to a table
uni.d<-rbind(d.N,d.E,d.O,d.A,d.C)
rownames(uni.d)<-c("N","E","O","A","C")
round(uni.d,2)#> n.Male n.Female m.Male m.Female sd.Male sd.Female pooled.sd diff D
#> N 2945 5739 2.79 3.16 0.88 0.86 0.87 -0.37 -0.42
#> E 2945 5739 3.00 3.06 0.96 0.96 0.96 -0.07 -0.07
#> O 2945 5739 4.07 3.88 0.58 0.63 0.62 0.18 0.30
#> A 2945 5739 3.68 4.03 0.75 0.68 0.71 -0.35 -0.50
#> C 2945 5739 3.38 3.44 0.72 0.75 0.74 -0.07 -0.09
#> t_stat df p
#> N -18.45 5845.38 0
#> E -3.09 5919.79 0
#> O 13.44 6458.37 0
#> A -21.33 5455.81 0
#> C -3.95 6119.69 0The column labeled “D” indicates Cohen’s d. Females score higher in Neuroticism (N; d = -0.42) and Agreeableness (A; d = -0.50), while males score higher in Openness (O; d = 0.30). Differences in Extraversion (E) and Conscientiousness (C) are statistically significant (p < .001), but their magnitudes are negligible (d < |.10|).
The average absolute difference across the Big Five traits is therefore d = 0.27. However, this average absolute difference does not provide information about the overall distance between males and females. For examining questions regarding the overall distance in a certain domain, such as personality, one should use multivariate methods
Other columns in the output are:
n.Male and n.Female = sample sizes
m.Male and m.Female = mean levels
sd.Male and sd.Female = standard deviations
pooled.sd = pooled standard deviation
diff = the raw (non-standardized) difference
t_stat = t statistic
df = Welch modified degrees of freedom
p = p-value
Multivariate difference using trait sum scores as variables
Although the univariate differences presented above estimate sex differences in each broad Big Five trait separately, these single Cohen’s d estimates do not provide a clear picture of the overall distance between males and females. To understand overall distances, multivariate methods that use multiple variables and account for their covariability (intercorrelations) should be used. These methods provide D values, which can be interpreted in the same way as Cohen’s d, that is, as standardized distances between males and females.
Mahalanobis’ D
Mahalanobis’ D is the standardized distance between the male and female group centroids in multivariate space. To calculate Mahalanobis’ D, you need the correlation matrix of the multivariate set and a vector of standardized mean differences. The ordering of the variables should be the same in both inputs.
# let's use the previously obtained univariate d-values
d_vector<-uni.d[,"D"]
# calculate intercorrelations
cor_mat<-cor(US.dat[,c("N","E","O","A","C")])
# test that the ordering of these variables matches
names(d_vector)==colnames(cor_mat)#> [1] TRUE TRUE TRUE TRUE TRUE#> N E O A C
#> N 1.00 -0.25 -0.10 -0.10 -0.26
#> E -0.25 1.00 0.15 0.31 0.10
#> O -0.10 0.15 1.00 0.10 0.05
#> A -0.10 0.31 0.10 1.00 0.17
#> C -0.26 0.10 0.05 0.17 1.00#> [,1]
#> [1,] 0.7681502Mahalanobis’ D is, of course, larger than any of the univariate d’s for single variables, here the standardized distance between the group centroids of males and females in the multivariate space is D = 0.77
Regularized D with D_regularized -function
In the D_regularized function, a predictive approach is used to calculate the regularized variant of D. This approach involves predicting the binary sex variable using a multivariate set of variables with regularized logistic regression. The predicted values for each individual represent the logits of probabilities of being male or female, and are hereafter referred to as femininity-masculinity scores (FM).
FM scores can be used in the same way as any other variable to calculate standardized distances between the mean FM scores of males and females with pooled standard deviations. In addition to estimating mean differences, there are other informative features of the multivariate/FM distributions.
The D_regularized function uses separate partitions of the original data for regulation (i.e., training) and prediction (i.e., testing). During regulation, the beta-coefficients for variables in the multivariate set are penalized through a specified function, such as a elastic net penalty with 10-fold cross-validation (default). The resulting coefficient weights are then used in the prediction part of the data as measures of FM for each individual. The arguments out (a logical value indicating whether predictions are made out-of-bag) and size (the size of the training partition of the data, equal for males and females) are used to specify that a partition is required. The multivariate set is supplied with the mv.vars argument.
It is also possible to request pcc (probability of correctly classified), auc (area under the receiver operating characteristics), and frequencies indicating the proportions of the sample that fall within certain cutoffs of predicted probabilities from 0 to 1 for each group (pred.prob) by using these arguments within D_regularized call.
# set seed number for reproducibility
set.seed(37127)
# decide the size of the regularization (training data fold)
# here, half of the smaller of the male and female groups is used
size=round(min(table(US.dat$gender.cat))/2)
# run D_regularized
D_out<-
D_regularized(data=US.dat,
mv.vars=c("N","E","O","A","C"),
group.var="gender.cat",
group.values = c("Male","Female"),
out = T,
size=size,
pcc = T,
auc = T,
pred.prob = T)
# print the output
round(D_out$D,2)#> n.Male n.Female m.Male m.Female sd.Male sd.Female pooled.sd diff D
#> [1,] 1473 4267 0.34 -0.3 0.86 0.79 0.81 0.64 0.79
#> pcc.Male pcc.Female pcc.total auc
#> [1,] 0.63 0.67 0.66 0.71The output shows that the sample sizes are smaller compared to the univariate estimates because they only reflect the samples in the prediction/testing dataset. The standardized multivariate difference between males and females is estimated at D = 0.79, which is very similar in magnitude to the multivariate sex difference indexed by Mahalanobis’ D.
There are additional features in the output that may be of interest. Similar to univariate d estimates, descriptive statistics for both groups on the FM-score metric are provided, as well as an estimate of the pooled standard deviation used to calculate D. The probability of correctly classified individuals (pcc) and the area under the receiver operating characteristics (auc) are also included, which are commonly used indices of predictive performance. In this example, they indicate the performance of the broad Big Five domains in predicting sex. The predicted probabilities of being male across certain cutoffs (defaults were used in this example) for both groups can be obtained from the P.table part of the output.
#>
#> [0,0.2) [0.2,0.4) [0.4,0.6) [0.6,0.8) [0.8,1]
#> Female 7.2 39.0 37.0 14.1 2.7
#> Male 1.2 18.1 37.1 32.6 11.0The table shows that majority of males are nevertheless not very clearly classified as males based on the Big Five personality traits (11% in .80-1 range), and same is true for classifying females as females (7% in 0-0.20 range).
To examine the regularized beta weights obtained during the training phase, the coefficients function can be called on the cv.mod object of the output. The s argument specifies the value of the penalty parameter lambda used to calculate predictions. The default value of “lambda.min” is used to obtain the values where the mean cross-validated error is minimized, while “lambda.1se” uses the largest value of lambda such that the error is within one standard error of the minimum. It is also possible to plot the cross-validation curve to visualize how “lambda.min” and “lambda.1se” were obtained. For more information on regularized regression conducted with the glmnet package, please refer to https://glmnet.stanford.edu/articles/glmnet.html.
#> 6 x 1 sparse Matrix of class "dgCMatrix"
#> s1
#> (Intercept) 3.37297381
#> N -0.54865910
#> E -0.02070834
#> O 0.64983681
#> A -0.88854735
#> C -0.23751055#> 6 x 1 sparse Matrix of class "dgCMatrix"
#> s1
#> (Intercept) 2.3670214
#> N -0.3900418
#> E .
#> O 0.4626928
#> A -0.6930196
#> C -0.1053781
Overlap
Differences and similarities are two sides of the same coin. In the case of sex differences, visual inspection of the distributions of FM-scores can reveal both. Specifically, by assessing the overlap between the distributions of FM-scores of males and females, we can observe differences between the two groups. However, we can also identify similarities by examining the degree of overlap between the distributions.
# Obtain the D-estimate
D<-unname(D_out$D[,"D"])
# Overlap relative to a single distribution (OVL)
OVL<-2*pnorm((-D/2))
# Proportion of overlap relative to the joint distribution
OVL2<-OVL/(2-OVL)
# non-parametric overlap with overlapping package
np.overlap<-
overlap(x = list(D_out$pred.dat[
D_out$pred.dat$group=="Male","pred"],
D_out$pred.dat[
D_out$pred.dat$group=="Female","pred"]),
plot=T)
# this corresponds to Proportion of overlap relative to the joint distribution (OVL2)
np.OVL2<-unname(np.overlap$OV)
# from which Proportion of overlap relative to a single distribution (OVL) is approximated at
np.OVL<-(2*np.OVL2)/(1+np.OVL2)
# comparison of different overlap-estimates
round(cbind(OVL,np.OVL,OVL2,np.OVL2),2)#> OVL np.OVL OVL2 np.OVL2
#> [1,] 0.69 0.82 0.53 0.7In this example, the non-parametric overlap (np) estimates to these logistic distributions are almost identical to the parametric overlap estimates.
Bootstrap confidence intervals for Regularized D
To obtain confidence intervals for D-estimates using regularized techniques, a bootstrap approach can be employed. This involves estimating D with multiple bootstrap samples drawn with replacement from the original sample and using the distribution of these D values to construct a percentile confidence interval around the estimate. In this example, only 10 redraws with replacement are used, but in practice, at least 1000 or 5000 redraws should be used for more accurate results.
# number of redraws
reps=10
# placeholder for D-estimates
D_out_boot<-rep(NA,reps)
# repeat the procedure to obtain a distribution
t1<-Sys.time()
for (i in 1:reps){
# draw replaced samples with the same size as the original data
boot.dat<-sample_n(US.dat,size=nrow(US.dat),replace=T)
# run D_regularized for each sample
D_out_boot[i]<-
D_regularized(data=boot.dat,
mv.vars=c("N","E","O","A","C"),
group.var="gender.cat",
group.values = c("Male","Female"),
out = T,
size=size)$D[,"D"]
}
t2<-Sys.time()
t2-t1#> Time difference of 6.333763 secs#> 2.5% 97.5%
#> 0.77 0.83After obtaining the bootstrap confidence interval, it can be reported that the multivariate sex difference D = 0.79, 95% CI [0.77, 0.83].
Multivariate difference using item scores as variables
Recently, there have been numerous studies in the field of personality psychology emphasizing the significance of narrow personality traits in predicting criteria, beyond the commonly used broad trait domains such as the Big Five (e.g., Seeboth & Mõttus, 2018). These narrow personality traits are usually either facets or nuances that are commonly modeled as being part of a certain Big Five trait, although their unique predictive utility suggests that they have trait-like properties (stability, cross-informant agreement etc.) that go beyond the broad Big Five domains (e.g., Mõttus et al., 2019). Consequently, narrower traits can be regarded as distinct traits. Thus, it is interesting to investigate whether using item scores as the variable set in multivariate estimation would provide better predictions of sex and different estimates of sex differences (see Ilmarinen et al., 2022).
However, due to the increased number of variables (from five to fifty in this dataset), regularization methods and different data partitions for training and testing are highly recommended to prevent overfitting of the data. For illustrative purposes, Mahalanobis’ D is also calculated below.
Mahalanobis’ D
# place-holder vector for univariate Cohen d's for items
d_vector_items<-rep(NA,length(per.items))
# loop through the items and obtain d for each item
for (i in 1:length(per.items)){
d_vector_items[i]<-
d_pooled_sd(data=US.dat,
var=per.items[i],
group.var="gender.cat",
group.values=c("Male","Female"))[,"D"]
}
cor_mat_items<-cor(US.dat[,per.items])
# Calculate D
D_maha_items<-
sqrt(t(d_vector_items) %*% solve(cor_mat_items) %*% d_vector_items)
D_maha_items#> [,1]
#> [1,] 0.9212846Item-based Mahalanobis’ D is somewhat larger than D (Mahalanobis’ or regularized) based on Big Five domains. The standardized distance between the group centroids of males and females in a space comprised of personality items is D = 0.92
Regularized D with D_regularized -function
It is straightforward to use D_regularized -function even with large variable sets (especially when the variable names are save in a vector first; here “per.items”). Separate partitions of data are again used for regularization and prediction. Domain-based results are printed for comparison.
D_items_out<-
D_regularized(data=US.dat,
mv.vars=per.items,
group.var="gender.cat",
group.values = c("Male","Female"),
out = T,
size=size,pcc = T,auc = T,pred.prob = T)
# Item-based D
round(D_items_out$D,2)#> n.Male n.Female m.Male m.Female sd.Male sd.Female pooled.sd diff D
#> [1,] 1473 4267 0.46 -0.4 0.9 0.86 0.87 0.86 0.99
#> pcc.Male pcc.Female pcc.total auc
#> [1,] 0.7 0.69 0.69 0.76#> n.Male n.Female m.Male m.Female sd.Male sd.Female pooled.sd diff D
#> [1,] 1473 4267 0.34 -0.3 0.86 0.79 0.81 0.64 0.79
#> pcc.Male pcc.Female pcc.total auc
#> [1,] 0.63 0.67 0.66 0.71With items as the multivariate set, D = 0.99 is slightly larger than with Big Five trait domains, D = 0.79.
Differences in predictive performance can also be illustrated from comparing the predicted probabilities.
#>
#> [0,0.2) [0.2,0.4) [0.4,0.6) [0.6,0.8) [0.8,1]
#> Female 11.6 38.5 33.4 14.3 2.2
#> Male 1.6 15.1 32.9 35.4 15.0#>
#> [0,0.2) [0.2,0.4) [0.4,0.6) [0.6,0.8) [0.8,1]
#> Female 7.2 39.0 37.0 14.1 2.7
#> Male 1.2 18.1 37.1 32.6 11.0The analysis reveals that the distribution of the FM-scores shows more individuals in the extremes. However, despite this, the overall predictive performance based on the pcc and auc values is not substantially improved. It is noteworthy that over 30% of the data falls in the 40%-60% range, indicating that their personality profiles are rather indifferent in regards to sex.
Overlap
Overlaps give similar information but in another format:
# Obtain the D-estimate
D_items<-unname(D_items_out$D[,"D"])
# Overlap relative to a single distribution (OVL)
OVL_items<-2*pnorm((-D_items/2))
# Proportion of overlap relative to the joint distribution
OVL2_items<-OVL_items/(2-OVL_items)
# non-parametric overlap
np.overlap_items<-
overlap(x = list(D_items_out$pred.dat[
D_items_out$pred.dat$group=="Male","pred"],
D_items_out$pred.dat[
D_items_out$pred.dat$group=="Female","pred"]),
plot=T)
# this corresponds to Proportion of overlap relative to the joint distribution (OVL2)
np.OVL2_items<-unname(np.overlap_items$OV)
# from which Proportion of overlap relative to a single distribution (OVL) is approximated at
np.OVL_items<-(2*np.OVL2_items)/(1+np.OVL2_items)
# compare the different overlap-estimates
round(cbind(OVL_items,np.OVL_items,
OVL2_items,np.OVL2_items),2)#> OVL_items np.OVL_items OVL2_items np.OVL2_items
#> [1,] 0.62 0.77 0.45 0.63#> OVL np.OVL OVL2 np.OVL2
#> [1,] 0.69 0.82 0.53 0.7The overlap between FM-distributions of males and females is somewhat smaller when individual items are used as the multivariate set rather than broad Big Five domains.