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

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 schools and students in each school
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組的總人數為7382人，中位數是64人，每組的人數落在15人至324人之間。

4 Visualization

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

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

-The estimated well-being score is 2.77 and the variances are 0.0358 and 0.7895 at the GRP and observation levels.Most of the variance comes from the observation, and the variance in the GRP level is smaller than observation level. These results indicate fairly significant variation in well-being score in group.

-ICC=0.0358/(0.0358+0.7895)=0.043

-At the group level, the intra-class correlation is 0.043. These values indicate the correlation in well-being score between two measurements in the same group. It means that 4.3% of the variation in well-being score can be attributed to attending the same group.

VarCorr(m0)

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

-在組的層次下，幸福感的標準差為0.1892，殘差的標準差為0.8885。

performance::icc(m0)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.043
  Unadjusted ICC: 0.043

-調整前後的ICC相同，表示每個測量值的的相似性為0.043，為較差的信度。

6 Work hours - individual and group

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.21367   22.19
HRS         -0.04646    0.00489   -9.51
G.HRS       -0.12693    0.01940   -6.54

-每增加1單位工作時間，幸福感下降0.04646；每增加1單位組間工作時間，幸福感下降0.12693。

-The estimated well-being score is 4.74 and the variances are 0.0135 and 0.7801 at the GRP and observation levels.Most of the variance comes from the observation, and the variance in the GRP level is smaller than observation level. These results indicate fairly significant variation in well-being score in group.

-The fixed-effect estimates yield an intercept of 4.74 and slope of -0.04646 and -0.12693. These values can be interpreted as the well-being score in groups, and the average rate of change per hours. Average well-being score clearly decreases significantly over time.

-ICC=0.0135/(0.0135+0.7801)=0.017

-At the group level, the intra-class correlation is 0.017. These values indicate the correlation in well-being score between two measurements in the same group. It means that 1.7% of the variation in well-being score can be attributed to attending the same group.

-About 1.2%(=(0.7895-0.7801)/0.7895) of variance in well-being scores can be attributed to differences in work hours among measurements attending different groups.

-About 3.8%(=(0.0358+0.7895)-(0.0135+0.7801)/0.0358+0.7895) of variance in well-being scores can be attributed to differences in work hours among pupils attending the same group.

VarCorr(m1)

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

-在組的層次下，幸福感標準差為0.1164，殘差的標準差為0.8832。

performance::icc(m1)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.017
  Unadjusted ICC: 0.016

-調整前後的ICC分別為0.017及0.016，相較於m0，ICC皆較小，表示增加固定變項降低觀測值之間的相似性，為較差的信度。

7 Leadership consideration - individual

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()

-在組別的層次下，不同領導氛圍下所獲得的幸福感分數。整體而言，領導氛圍越好，幸福感越高。

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.00924816 (tol = 0.002, component 1)

-每增加1單位工作時間，幸福感下降0.02848；每增加1單位領導氛圍，幸福感上升0.049454；每增加1單位不同(組間)工作時間，幸福感下降0.07057。

-The estimated well-being score is 2.46 and the variances are 0.1466, 0.0107 and 0.6413 at the GRP,LEAD and observation levels.Most of the variance comes from the observation, and the variance in the observation level is greater than GRP and LEAD level. These results indicate fairly no significant variation in well-being score in group and leadership climate.

-The fixed-effect estimates yield an intercept of 2.46 and slope of -0.02848,0.49454 and -0.07057. These values can be interpreted as the well-being score in groups, and the average rate of change per hours. Average well-being score clearly decreases significantly over time.However,the average rate of change per leadership climate. Average well-being score clearly increases significantly over time.

-ICC=0.1466/(0.1466+0.0107+0.6413)=0.18357

-At the group level, the intra-class correlation is 0.18357. These values indicate the correlation in well-being score between two measurements in the same group. It means that 18.4% of the variation in well-being score can be attributed to attending the same group.

-About 18.8%(=(0.7895-0.6413)/0.7895) of variance in well-being scores can be attributed to differences in work hours and leadership among measurements attending different groups.

-About-21.2%(=(0.0135+0.7801)-(0.1466+0.0107+0.6413)/0.0135+0.7801) of variance in well-being scores can be attributed to differences in leadership climate among pupils attending the same group.

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

-The intercepts and slopes have correlations of -0.97 at the GRP level. It can be interpreted as the correlation between well-being score around group and leadership climate.

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

-每增加1單位工作時間，幸福感下降0.02856；每增加1單位領導氛圍，幸福感上升0.12895；每增加1單位不同(組間)工作時間，幸福感下降0.17402；在不同工作時間的領導氛圍在每增加1單位時，會上升0.03217幸福感。

-The estimated well-being score is 3.64 and the variances are 0.13598, 0.00986 and 0.64129 at the GRP,LEAD and observation levels.Most of the variance comes from the observation, and the variance in the observation level is greater than GRP and LEAD level. These results indicate fairly significant variation in well-being score in group and leadership climate.

-The intercepts and slopes have correlations of -0.97 at the GRP level.It can be interpreted as the correlation between well-being score around group and leadership climate.

-The fixed-effect estimates yield an intercept of 3.64 and slope of -0.02856,0.12895,-0.17402 and 0.03217. These values can be interpreted as the well-being score in groups, and the average rate of change per hours. Average well-being score clearly decreases significantly over time.However,the average rate of change per leadership climate. Average well-being score clearly increases significantly over time.In addition,the average rate of change per leadership climate of each group.Average well-being score also increases significantly over time.

-The random effects indicate that the average well-being score varies by GRP, with a variance of 0.13598, while the rate of change has a variance of 0.00986. The correlation between intercept and slope is -0.97, suggesting that it can not be attributed to the leadership climate of each group.

-At the group level, the intra-class correlation is 0.185(=(0.13598+0.00986)/(0.13598+0.00986+0.64129)). These values indicate the correlation in well-being score between two measurements in the same leadership climate within the same G.HRS. It means that 18.5% of the variation in well-being score can be attributed to attending the same leadership climate within the same G.HRS.

-Gain well-being scores of any two measurements from the same GRP has a correlation of 0.173(=0.13598/(0.13598+0.00986+0.64129))

VarCorr(m3)

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

-在m3中領導氛圍與GRP的相關為-0.97

performance::icc(m3)

# Intraclass Correlation Coefficient

    Adjusted ICC: 0.030
  Unadjusted ICC: 0.024

-調整後ICC為0.03，未調整ICC為0.024，測量值之間的相似度不高，為較差的信度。

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

Huang En-Li

2022-11-27