1 Introduction

The data set contains four variables:

Cohesion (COHES),
Leadership Climate (LEAD),
Well-Being (WBEING) and
Work Hours (HRS).

Each of these variables has two variants - a group mean version that replicates each group mean for every individual in the group, and a within-group version where the group mean is subtracted from each individual response (i.e., a group-mean centered or demeaned variable). The group mean version is designated with a G. (e.g., G.HRS), and the within-group version is designated with a W. (e.g., W.HRS).

2 Data

Questionnaire data were analyzed from 7,382 enlisted personnel from 99 U.S. Army companies

data(bh1996, package="multilevel")
dta <- bh1996

str(dta)

'data.frame':   7382 obs. of  13 variables:
 $ GRP     : num  1 1 1 1 1 1 1 1 1 1 ...
 $ COHES   : num  3.75 3.25 3.38 3.75 4 ...
 $ G.COHES : num  2.93 2.93 2.93 2.93 2.93 ...
 $ W.COHES : num  0.821 0.321 0.446 0.821 1.071 ...
 $ LEAD    : num  3.18 3 3.64 3.36 3.55 ...
 $ G.LEAD  : num  2.93 2.93 2.93 2.93 2.93 ...
 $ W.LEAD  : num  0.2527 0.0709 0.7072 0.4345 0.6163 ...
 $ HRS     : num  12 11 12 9 7 8 9 8 12 11 ...
 $ G.HRS   : num  11 11 11 11 11 ...
 $ W.HRS   : num  1.027 0.027 1.027 -1.973 -3.973 ...
 $ WBEING  : num  2.11 3.33 2.11 4.39 1.72 ...
 $ G.WBEING: num  2.79 2.79 2.79 2.79 2.79 ...
 $ W.WBEING: num  -0.682 0.54 -0.682 1.596 -1.071 ...

3 Tables

# Number of companies and data in each company
n_distinct(dta$GRP)

[1] 99

with(dta, table(GRP)) |> quantile()

   0%   25%   50%   75%  100% 
 15.0  43.5  64.0  94.0 226.0

dta <- dta %>% group_by(GRP) %>% mutate(nc=n()) %>% ungroup()

with(dta, table(nc)) |> quantile()

    0%    25%    50%    75%   100% 
 15.00  60.50  94.00 166.25 324.00

此筆資料從99間公司收集共7,382筆資料，每家公司收集的資料數範圍為15-226。

4 Visualization

# 以Well-Being (WBEING)中位數排序,對每家公司的WBEING值分別作圖
ggplot(dta, aes(x=reorder(factor(GRP), WBEING, median), WBEING)) +
  geom_boxplot() +
  coord_flip()+
  labs(x="Group ID",
       y="Well-being score")

5 Null model - random intercepts only

\[{\rm WBEING}_{ij}= b_{0i} + \epsilon_{ij},~ i = 1, \ldots, 99,\] \[b_{0i} \sim N(\mu_{b0}, \sigma^2_{b0}),\epsilon_{it} \sim N(0, \sigma^2).\]

# 僅估計公司層次的截距隨機效果
m0 <- lme4::lmer(WBEING ~ (1 | GRP), data=dta)
summary(m0, corr=FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: WBEING ~ (1 | GRP)
   Data: dta

REML criterion at convergence: 19347.3

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-3.322 -0.648  0.031  0.718  2.667 

Random effects:
 Groups   Name        Variance Std.Dev.
 GRP      (Intercept) 0.0358   0.189   
 Residual             0.7895   0.889   
Number of obs: 7382, groups:  GRP, 99

Fixed effects:
            Estimate Std. Error t value
(Intercept)   2.7743     0.0222     125

VarCorr(m0)

 Groups   Name        Std.Dev.
 GRP      (Intercept) 0.1892  
 Residual             0.8885

performance::icc(m0)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.043
  Unadjusted ICC: 0.043

- 所有公司的員工幸福感(well-being)估計平均值為2.7743，標準差為0.189
- 所有員工的幸福感(well-being)估計平均值為2.7743，標準差為0.906[sqrt(0.0358+0.785)]
- 同一間公司中，員工間幸福感的相關為0.043

6 Work hours - individual and group

\[{\rm WBEING}_{ij}= b_{0i} + \beta_{1} {\rm HRS}_{ij} + \beta_{2} {\rm G.HRS}_{i}+ \epsilon_{ij},~ i = 1, \ldots, 99,\] \[b_{0i} \sim N(\mu_{b0}, \sigma^2_{b0}),\epsilon_{it} \sim N(0, \sigma^2).\]

# 以個人及公司的HRS作為解釋變項，並估計公司層次的截距隨機效果
m1 <- lme4::lmer(WBEING ~ HRS + G.HRS + (1 | GRP) , data=dta)
summary(m1, corr=FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: WBEING ~ HRS + G.HRS + (1 | GRP)
   Data: dta

REML criterion at convergence: 19212.3

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-3.353 -0.650  0.038  0.713  2.709 

Random effects:
 Groups   Name        Variance Std.Dev.
 GRP      (Intercept) 0.0135   0.116   
 Residual             0.7801   0.883   
Number of obs: 7382, groups:  GRP, 99

Fixed effects:
            Estimate Std. Error t value
(Intercept)  4.74083    0.21368   22.19
HRS         -0.04646    0.00489   -9.51
G.HRS       -0.12693    0.01940   -6.54

VarCorr(m1)

 Groups   Name        Std.Dev.
 GRP      (Intercept) 0.1164  
 Residual             0.8832

performance::icc(m1)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.017
  Unadjusted ICC: 0.016

- 所有員工的幸福感(well-being)估計平均值為4.741，標準差為0.891[sqrt(0.0135+0.7801)]
- 公司平均增加一小時的工作時間(work hours)，員工的幸福感會下降0.127
- 員工每增加一小時的工作時間，員工的幸福感會下降0.046
- 同一間公司中，員工間幸福感的相關為0.016

7 Leadership consideration - individual

# 針對員工人數大於100的公司分別作圖(x=Leadership, y=Well-being)
ggplot(subset(dta, nc > 100), aes(LEAD, WBEING))+
  stat_smooth(method='lm', formula=y~x, se=FALSE,
              size=rel(.5), col=1)+
  geom_point(size=rel(.5), col=8, alpha=.5)+
  facet_wrap(. ~ GRP)+
  labs(x="Leadership",
       y="Well-being")+
  theme_minimal()

- 在所有公司中，隨著領導氣氛(leadership climate)的增加，員工的幸福感(well-being)也有上升的趨勢

\[{\rm WBEING}_{ij}= b_{0i} + \beta_{1} {\rm HRS}_{ij} + b_{2i} {\rm LEAD}_{ij} + \beta_{3} {\rm G.HRS}_{i}+ \epsilon_{ij},~ i = 1, \ldots, 99,\] \[(b_{0i},~ b_{2i})' \sim N({\pmb \mu_b},~~{\pmb \Sigma_b}),\epsilon_{it} \sim N(0, \sigma^2).\]

# 以個人HRS、公司HRS以及個人LEAD作為解釋變項，並估計公司層次的截距隨機效果及LEAD隨機斜率效果
m2 <- lme4::lmer(WBEING ~ HRS + LEAD + G.HRS + (LEAD | GRP), data=dta)
summary(m2, corr=FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: WBEING ~ HRS + LEAD + G.HRS + (LEAD | GRP)
   Data: dta

REML criterion at convergence: 17822.6

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-3.871 -0.656  0.041  0.697  3.958 

Random effects:
 Groups   Name        Variance Std.Dev. Corr 
 GRP      (Intercept) 0.1466   0.383         
          LEAD        0.0107   0.103    -0.97
 Residual             0.6413   0.801         
Number of obs: 7382, groups:  GRP, 99

Fixed effects:
            Estimate Std. Error t value
(Intercept)  2.46418    0.20754   11.87
HRS         -0.02848    0.00447   -6.37
LEAD         0.49454    0.01687   29.31
G.HRS       -0.07057    0.01782   -3.96
optimizer (nloptwrap) convergence code: 0 (OK)
Model failed to converge with max|grad| = 0.00923328 (tol = 0.002, component 1)

# 上述模型有警示模型無法收斂，因此對m2模型進行優化
library(optimx)
update(m2, control = lmerControl(optimizer= "optimx",
                                 optCtrl  = list(method="nlminb")))

Linear mixed model fit by REML ['lmerMod']
Formula: WBEING ~ HRS + LEAD + G.HRS + (LEAD | GRP)
   Data: dta
REML criterion at convergence: 17822.6
Random effects:
 Groups   Name        Std.Dev. Corr 
 GRP      (Intercept) 0.383         
          LEAD        0.103    -0.97
 Residual             0.801         
Number of obs: 7382, groups:  GRP, 99
Fixed Effects:
(Intercept)          HRS         LEAD        G.HRS  
     2.4641      -0.0285       0.4945      -0.0706

performance::icc(m2)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.030
  Unadjusted ICC: 0.024

- 所有員工的幸福感(well-being)估計平均值為2.464，標準差為0.888[sqrt(0.1466+0.6413)]
- 公司平均增加一小時的工作時間(work hours)，員工的幸福感會下降0.071
- 員工每增加一小時的工作時間，員工的幸福感會下降0.028
- 領導氣氛(leadership climate)每增加一單位，員工的幸福感會上升0.495，各公司間的標準差為0.103
- 截距隨機效果及領導氣氛隨機斜率效果間有很高的負相關(r=-.97)，顯示平均幸福感越低的公司，其員工越容易受到領導氣氛的影響而提升幸福感
- 同一間公司中，員工間幸福感的相關為0.024

\[{\rm WBEING}_{ij}= b_{0i} + \beta_{1} {\rm HRS}_{ij} + b_{2i} {\rm LEAD}_{ij} + \beta_{3} {\rm G.HRS}_{i} + \beta_{4} {\rm G.HRS}_{i}\times{\rm LEAD}_{ij} + \epsilon_{ij},~ i = 1, \ldots, 99,\] \[(b_{0i},~ b_{2i})' \sim N({\pmb \mu_b},~~{\pmb \Sigma_b}),\epsilon_{it} \sim N(0, \sigma^2).\]

# 以個人HRS、公司HRS、個人LEA及公司HRS與個人LEA交互作用作為解釋變項，並估計公司層次的截距隨機效果及LEAD隨機斜率效果
m3 <- lme4::lmer(WBEING ~ HRS + LEAD + G.HRS + LEAD:G.HRS + (LEAD | GRP), 
                 data=dta,
                 control = lmerControl(optimizer= "optimx",
                                       optCtrl  = list(method="nlminb")))
summary(m3, corr=FALSE)

Linear mixed model fit by REML ['lmerMod']
Formula: WBEING ~ HRS + LEAD + G.HRS + LEAD:G.HRS + (LEAD | GRP)
   Data: dta
Control: lmerControl(optimizer = "optimx", optCtrl = list(method = "nlminb"))

REML criterion at convergence: 17825.9

Scaled residuals: 
   Min     1Q Median     3Q    Max 
-3.837 -0.660  0.041  0.695  3.953 

Random effects:
 Groups   Name        Variance Std.Dev. Corr 
 GRP      (Intercept) 0.13598  0.3688        
          LEAD        0.00986  0.0993   -0.97
 Residual             0.64129  0.8008        
Number of obs: 7382, groups:  GRP, 99

Fixed effects:
            Estimate Std. Error t value
(Intercept)  3.64326    0.73255    4.97
HRS         -0.02856    0.00447   -6.39
LEAD         0.12895    0.21881    0.59
G.HRS       -0.17402    0.06415   -2.71
LEAD:G.HRS   0.03217    0.01919    1.68

VarCorr(m3)

 Groups   Name        Std.Dev. Corr 
 GRP      (Intercept) 0.36876       
          LEAD        0.09929  -0.97
 Residual             0.80081

performance::icc(m3)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.030
  Unadjusted ICC: 0.024

- 所有員工的幸福感(well-being)估計平均值為3.643，標準差為0.881[sqrt(0.13598+0.64129)]
- 公司平均增加一小時的工作時間(work hours)，員工的幸福感會下降0.174
- 員工每增加一小時的工作時間，員工的幸福感會下降0.029
- 領導氣氛(leadership climate)每增加一單位，員工的幸福感會上升0.129，各公司間的標準差為0.099
- 同一間公司中，員工間幸福感的相關為0.024
- 截距隨機效果及領導氣氛隨機斜率效果間有很高的負相關(r=-.97)，顯示平均幸福感越低的公司，其員工越容易受到領導氣氛的影響而提升幸福感
- 此處並未針對交互作用進行分析，若交互作用顯著，則要進一步檢驗交互作用效果

8 References

Bliese, P. D. & Halverson, R. R. (1996). Individual and nomothetic models of job stress: An examination of work hours, cohesion, and well-being. Journal of Applied Social Psychology, 26, 1171-1189.

Job Stress

Ching-Yuan Chang

2022-10-08