2.3 Alternative explanations for association

On page 46, Hayes produced a couple correlations. Here’s how to get them from base R.

cor(glbwarm$sex, glbwarm$negemot)
#> [1] -0.117
cor(glbwarm$sex, glbwarm$govact)
#> [1] -0.0986

Again, if we wanted to get full Bayesian estimates, we’d fit an intercept-only multivariate model.

brm(data = glbwarm, family = gaussian,
    cbind(negemot, govact, sex) ~ 1,
    chains = 4, cores = 4) %>% 
  
  print(digits = 3)
#>  Family: MV(gaussian, gaussian, gaussian) 
#>   Links: mu = identity; sigma = identity
#>          mu = identity; sigma = identity
#>          mu = identity; sigma = identity 
#> Formula: negemot ~ 1 
#>          govact ~ 1 
#>          sex ~ 1 
#>    Data: glbwarm (Number of observations: 815) 
#> Samples: 4 chains, each with iter = 2000; warmup = 1000; thin = 1;
#>          total post-warmup samples = 4000
#> 
#> Population-Level Effects: 
#>                   Estimate Est.Error l-95% CI u-95% CI Eff.Sample  Rhat
#> negemot_Intercept    3.560     0.056    3.451    3.670       4000 1.000
#> govact_Intercept     4.588     0.050    4.489    4.686       4000 1.000
#> sex_Intercept        0.488     0.018    0.453    0.522       4000 0.999
#> 
#> Family Specific Parameters: 
#>                        Estimate Est.Error l-95% CI u-95% CI Eff.Sample  Rhat
#> sigma_negemot             1.531     0.038    1.461    1.607       4000 0.999
#> sigma_govact              1.363     0.034    1.299    1.431       4000 0.999
#> sigma_sex                 0.501     0.012    0.477    0.526       4000 1.000
#> rescor(negemot,govact)    0.575     0.023    0.528    0.618       4000 0.999
#> rescor(negemot,sex)      -0.117     0.034   -0.185   -0.049       4000 0.999
#> rescor(govact,sex)       -0.098     0.034   -0.164   -0.032       4000 1.000
#> 
#> Samples were drawn using sampling(NUTS). For each parameter, Eff.Sample 
#> is a crude measure of effective sample size, and Rhat is the potential 
#> scale reduction factor on split chains (at convergence, Rhat = 1).

For our purposes, the action is in the ‘rescor(\(i\), \(j\))’ portions of the ‘Family Specific Parameters’ section.

Anyway, if you wanted to get all the Pearson’s correlations among the glbwarm variables, rather than piecewise cor() approach, you could use the lowerCor() function from the psych package.

psych::lowerCor(glbwarm[, 1:7], digits = 3)
#>          govact posemt negemt idelgy age    sex    partyd
#> govact    1.000                                          
#> posemot   0.043  1.000                                   
#> negemot   0.578  0.128  1.000                            
#> ideology -0.418 -0.029 -0.349  1.000                     
#> age      -0.097  0.042 -0.057  0.212  1.000              
#> sex      -0.099  0.074 -0.117  0.133  0.166  1.000       
#> partyid  -0.360 -0.036 -0.324  0.619  0.154  0.109  1.000