1 How are pro-environmental behaviours related to working hours? (secondary data)

Data selection:

file <- list.files("../data", pattern = ".dta", full.names = T)
data <- read_dta(file) %>% 
  select(age, sex, ethn, health, edu, income1, cohab, child, sclass, incsat, contains("sat"), contains("comm"), leissat, hhwkh, penvbav = index, j_jbsat, wkhn, wkh) %>% 
  drop_na() %>% 
  mutate(across(c(wkhn, wkh), log),
         across(c(sex, ethn, edu, sclass), factor))

TODO (Alex):

couldn’t find variables for commute
check which variables are categorical (atm I’m assuming: sex, ethn, edu, sclass)

The working hours data (wkhn, wkh) are assumed to come from a multivariate lognormal distribution \(\text{log}\mathcal{N}(\mu, \Sigma^2_e)\) with a variance-covariance matrix of the standard deviations \(\Sigma^2_e\) of both outcome variables and a mean \(\mu\) that is decomposed as

\[ \mu = \alpha + \beta_\text{age} \cdot \text{age} + \beta_\text{age2} \cdot \text{age}^2 + \beta_\text{gender} \cdot \text{gender} + \beta_\text{health} \cdot \text{health} + \beta_\text{education} \cdot \text{education}\\ + \beta_\text{income} \cdot \text{income} + β_\text{martial} \cdot \text{martial status} + \beta_\text{children} \cdot \text{children}\\ + \beta_\text{envbeh} \cdot \text{environmental behaviour} + \beta_\text{jobsat} \cdot \text{job satisfaction} + \beta_\text{leisat} \cdot \text{leisure satisfaction}\\ + \beta_\text{care} \cdot \text{care} + \beta_\text{hh work} \cdot \text{hh work} + \beta_\text{stress} \cdot \text{stress} + \beta_\text{family} \cdot \text{family value} \]

Predictors were added to the model incrementally as shown in the code here:

Models were compared using F-tests (base on the log-likelihood). The results of the model comparisons can be found in Table 1.1. Adding hhwkh as predictor to the model with the demographic predictor variables significantly improved the model fit (model 2). Adding leissat further improved the model fit (model 3), as well as adding incsat and j_jbsat (model 4). Finally adding penvbav improved the model fit (model 5).

Table 1.1: Model comparisons. Where two values are shown the left value corresponds to `wkh` and the right value corresponds to `wkhn`.
Model	\(R^2\)	\(R^2_\text{adj}\)	Log-likelihood	AIC	F-statistic	\(P(\geq F)\)
M1	.25, .23	.25, .23	-733.4, -721.7	1466.7, 1443.4	–	–
M2	.25, .24	.25, .24	-732.6, -720.2	1465.1, 1440.4	18.7	< .05
M3	.27, .24	.27, .24	-723, -715.9	1446, 1431.8	70	< .05
M4	.27, .24	.27, .24	-726.7, -719.8	1453.5, 1439.6	3.1	< .05
M5	.27, .25	.27, .24	-728.3, -720.5	1456.6, 1441	11.1	< .05

The regression coefficients of the full multivariate model M5 are shown in Table 1.2.

Table 1.2: Model coefffcients. Shown are the maximum likelihood estimate with 95% confidence intervals (CI) along with the null-hypothesis test of the t-statistic.
	Outcome variable: `wkh`			Outcome variable: `wkhn`
Predictor	Est. with 95% CI	t-value	\(P(\geq \left\| t \right\|)\)	Est. with 95% CI	t-value	\(P(\geq \left\| t \right\|)\)
Intercept	3.8 [3.7, 3.9]	66.7	< .05	3.7 [3.6, 3.8]	65.5	< .05
penvbav	-0.0021 [-0.0041, -0.00019]	-2.1	< .05	-0.0037 [-0.0056, -0.0017]	-3.7	< .05
incsat	0.0015 [-0.0068, 0.0097]	0.3	0.73	-0.0031 [-0.011, 0.0051]	-0.7	0.46
jbsat	-0.0071 [-0.015, 0.00096]	-1.7	0.08	-0.0032 [-0.011, 0.0048]	-0.8	0.44
leissat	-0.039 [-0.047, -0.031]	-9.4	< .05	-0.027 [-0.035, -0.019]	-6.6	< .05
hhwkh	-0.0044 [-0.006, -0.0028]	-5.5	< .05	-0.0048 [-0.0063, -0.0032]	-6.0	< .05
age 1	-0.32 [-1.4, 0.71]	-0.6	0.54	1.2 [0.21, 2.3]	2.4	< .05
age 2	-12 [-13, -11]	-23.4	< .05	-12 [-13, -11]	-24.5	< .05
sex2	-0.092 [-0.12, -0.068]	-7.4	< .05	-0.097 [-0.12, -0.073]	-7.9	< .05
ethn2	-0.044 [-0.071, -0.017]	-3.2	< .05	-0.037 [-0.064, -0.01]	-2.7	< .05
health	-0.0016 [-0.014, 0.011]	-0.2	0.8	0.00018 [-0.012, 0.012]	0.0	0.98
edu2	-0.021 [-0.044, 0.0025]	-1.8	0.08	-0.0096 [-0.033, 0.014]	-0.8	0.42
income1	0.00007 [0.000063, 0.000077]	18.8	< .05	0.000056 [0.000049, 0.000063]	15.1	< .05
cohab	-0.06 [-0.083, -0.036]	-4.9	< .05	-0.067 [-0.09, -0.043]	-5.5	< .05
child	-0.21 [-0.24, -0.18]	-12.9	< .05	-0.2 [-0.23, -0.17]	-12.5	< .05
sclass2	0.016 [-0.036, 0.068]	0.6	0.55	0.05 [-0.0016, 0.1]	1.9	0.06
sclass3	0.068 [0.022, 0.11]	2.9	< .05	0.071 [0.026, 0.12]	3.1	< .05
sclass4	0.0081 [-0.043, 0.059]	0.3	0.76	0.045 [-0.0055, 0.096]	1.7	0.08
sclass6	0.11 [0.051, 0.17]	3.6	< .05	0.13 [0.07, 0.19]	4.3	< .05
sclass7	-0.07 [-0.12, -0.017]	-2.6	< .05	-0.076 [-0.13, -0.023]	-2.8	< .05
sclass8	-0.063 [-0.12, -0.0041]	-2.1	< .05	-0.035 [-0.094, 0.024]	-1.2	0.25

2 What factors influence individuals working hours? (primary data)

This code does the initial data wrangling of the input file.

file <- list.files("../data", pattern = ".xlsx", full.names = T)
data <- read_xlsx(file, sheet = 3) %>% 
  clean_names() %>% 
  select(-starts_with("org"), -ends_with("index"), -starts_with("q21"), -starts_with("q22")) %>% 
  rename_with(.fn = ~str_replace(., "q3_", "oep_")) %>% # organisations environment policy
  rename_with(.fn = ~str_replace(., "q7_", "pwh_")) %>% # prefered working hours 
  rename_with(.fn = ~str_replace(., "q8_", "envat_")) %>% # envattindex 
  rename_with(.fn = ~str_replace(., "q9_", "envbeh_")) %>% # environmental behaviour -> 1 is NA!
  rename_with(.fn = ~str_replace(., "q12_", "hhdiv_")) %>% # household division
  rename(pwh_2 = longerhoursifincome,
         pwh_3 = shorterhoursincome) %>% 
  mutate(across(c(industry, employed, workornot, 
                  maininc, family, housdv, age, 
                  sex, ethn, marstat, edu, income, 
                  own, value, savings), factor),
         # Reverse code items 
         across(starts_with("envbeh_"), ~ifelse(. == 6, NA, .)),
         across(c(envat_6, envat_7, envbeh_1, envbeh_3), ~re_code(., 1:5, 5:1)),
         across(c(pwh_3, pwh_5, pwh_6), ~re_code(., 1:5, 5:1)),
         across(sex, ~recode(., 
                             `1` = "female", 
                             `2` = "male", 
                             `3` = "non-binary"))) %>% 
  filter(sex != "non-binary") %>% # only one participant
  select(-value)

A descriptive overview of the data can be found in Table 2.1 with values separated by sex. Also shown is Cronbach’s \(\alpha\) for all scale items suggesting moderate to high internal consistency. Higher means for males were observed for wkh, and wkhchoose, with a higher oep mean for females. Positive correlations were observed for wkhchoose by wkh, as well as oep by envat, and envat by envbeh. Negative correlations were found between hhdiv and each envat and envbeh

Table 2.1: Descriptive summary by sex (females = 49; males = 29) with means and standard deviations (SD) and Cronbach’s \(\alpha\) with 95% Confidence interval, and all bivariate correlations.
	Means (SD)			Correlations
Variable	Females	Males	Cronbach’s \(\alpha\)	wkh	wkhchoose	pwh	oep	envat	envbeh
wkh	34.5 (12.7)	39.8 (11.9)
wkhchoose	32.6 (10.4)	37.3 (11.7)		.89*
pwh	2.9 (0.6)	2.8 (0.7)	.72 [.63, .82]	.42*	.20
oep	2.3 (0.9)	1.8 (0.6)	.84 [.79, .90]	.01	.05	.03
envat	1.9 (0.6)	1.8 (0.6)	.69 [.59, .80]	-.12	.03	-.13	.25*
envbeh	2.4 (0.6)	2.2 (0.5)	.70 [.60, .80]	-.12	-.04	-.05	.14	.46*
hhdiv	4.5 (0.8)	4.3 (0.8)	.90 [.86, .93]	.08	-.01	-.08	.06	-.45*	-.34*
Note:
* indicates significance on a level of \(\alpha\) = .05.

The following series of linear regression models was fitted. Continuous predictor variables were centred to reduce multicollinearity.

m_wkh1 <- lm(wkh_val ~ wkh_vars, data = data_wkh)
m_wkh2 <- lm(wkh_val ~ wkh_vars + oep + envat + envbeh + hhdiv, data = data_wkh)
m_wkh3 <- lm(wkh_val ~ wkh_vars * (oep + envat + envbeh + hhdiv), data = data_wkh)
m_wkh4 <- lm(wkh_val ~ sex + wkh_vars * (oep + envat + envbeh + hhdiv), data = data_wkh)
m_wkh5 <- lm(wkh_val ~ sex:wkh_vars + sex + wkh_vars * (oep + envat + envbeh + hhdiv), data = data_wkh)
m_wkh6 <- lm(wkh_val ~ sex:wkh_vars + (sex + wkh_vars) * (oep + envat + envbeh + hhdiv), data = data_wkh)
m_wkh7 <- lm(wkh_val ~ sex * wkh_vars * (oep + envat + envbeh + hhdiv), data = data_wkh)

Models were compared using F-test (log-likelihood based). The results of the model comparisons can be found in Table 2.2. The best fitting model was found to be model 6. In other words model fit does not just improve by adding sex to the model (model 4) but also the two way interactions of sex and the set of variables that represents the respondents’ environmental attitude. The interaction was sex and working hours variable per se does not improve the model fit (model 5). However adding interactions for each sex and working hours variable with the environmental attitude variables (model 6) increases the model fit significantly. Finally adding the three-way interaction with both sex and environmental attitude variables (model 7) does not show evidence for a significantly better model fit.

Table 2.2: Model comparisons.
Model	\(R^2\)	\(R^2_\text{adj}\)	Log-likelihood	AIC	F-statistic	\(P(\geq F)\)
M1	.72	0.72	-861.8	1731.6	–	–
M2	.73	0.72	-861.0	1738.1	0.4	.82
M3	.73	0.71	-860.0	1752.1	0.3	.98
M4	.74	0.72	-856.4	1746.8	7.3	< .05
M5	.74	0.72	-854.4	1746.9	1.9	.15
M6	.76	0.73	-846.1	1738.1	3.9	< .05
M7	.77	0.73	-841.7	1745.4	1	.46

The coefficients of the best fitting model (model 6) can be found in Table 2.3. Importantly, we found evidence for interactions of sex and both envbeh and hhdiv but not for oep and envat.

Table 2.3: Model coefffcients. Shown are the maximum likelihood estimate with 95% confidence intervals (CI) along with the null-hypothesis test of the t-statistic
Predictor	Est. with 95% CI	t-value	\(P(\geq \left\| t \right\|)\)
(Intercept)	3.1 [0.4, 5.8]	2.2	< .05
sexmale	0.5 [-4.2, 5.2]	0.2	.83
wkh	-5 [-8.8, -1.2]	-2.6	< .05
wkhchoose	29.5 [25.7, 33.3]	15.2	< .05
oep	1.3 [-1.7, 4.3]	0.8	.40
envat	-3.4 [-8.6, 1.7]	-1.3	.19
envbeh	-2.5 [-7.2, 2.1]	-1.1	.28
hhdiv	-3.5 [-7.3, 0.3]	-1.8	.07
sexmale:wkh	5.8 [-0.7, 12.2]	1.8	.08
sexmale:wkhchoose	5.3 [-1.2, 11.8]	1.6	.11
sexmale:oep	-1.2 [-5.2, 2.8]	-0.6	.55
sexmale:envat	3.5 [-2.2, 9.2]	1.2	.23
sexmale:envbeh	7.5 [1.9, 13.1]	2.6	< .05
sexmale:hhdiv	6 [2, 10]	3.0	< .05
wkh:oep	1.4 [-2.6, 5.3]	0.7	.50
wkhchoose:oep	1.3 [-2.7, 5.3]	0.6	.52
wkh:envat	-1.8 [-8.1, 4.6]	-0.5	.58
wkhchoose:envat	0.7 [-5.6, 7.1]	0.2	.82
wkh:envbeh	-1.1 [-7.2, 5]	-0.4	.73
wkhchoose:envbeh	-0.6 [-6.7, 5.5]	-0.2	.85
wkh:hhdiv	0.7 [-3.8, 5.3]	0.3	.75
wkhchoose:hhdiv	0.3 [-4.3, 4.8]	0.1	.90

The predictions of model 6 are visualised in Figure 2.1.

Figure 2.1: Model predictions. Values on the x-axis are centred.

Figure 2.2 shows the estimated marginal means with 95% confidence intervals for the working hours variables presented by sex. The estimates were obtained from model 6, the model that included the interaction of the working hour variables and sex. We used pairwise comparisons with Bonferroni adjustment for multiple comparisons to test the differences between males and females for each working hour variable. Consistent with Figure 2.2 we found higher actual working hours (est. = -6.3, SE = 2.4, t(212) = -2.7, p < .05) and chosen working hours (est. = -5.8, SE = 2.4, t(212) = -2.5, p < .05) for males compared to females but no difference in the attitude towards working hour results (est. = -0.5, SE = 2.4, t(212) = -0.2, p = .83).

Figure 2.2: Estimated marginal mean and 95% confidence intervals (CI) for working hour values by sex estimated on the basis of model 6.

3 Session info

R version 4.0.3 (2020-10-10)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 18.04.6 LTS

Matrix products: default
BLAS:   /usr/lib/x86_64-linux-gnu/openblas/libblas.so.3
LAPACK: /usr/lib/x86_64-linux-gnu/libopenblasp-r0.2.20.so

locale:
 [1] LC_CTYPE=en_GB.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_GB.UTF-8        LC_COLLATE=en_GB.UTF-8    
 [5] LC_MONETARY=en_GB.UTF-8    LC_MESSAGES=en_GB.UTF-8   
 [7] LC_PAPER=en_GB.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] grid      stats     graphics  grDevices utils     datasets  methods  
[8] base     

other attached packages:
 [1] broom_1.0.5           ggthemes_4.2.4        psyntur_0.1.0        
 [4] readxl_1.4.3          janitor_2.0.1         haven_2.5.1          
 [7] modelr_0.1.10         gtable_0.3.4          gridExtra_2.3        
[10] rmdformats_1.0.1      patchwork_1.1.1       kableExtra_1.3.4.9000
[13] lubridate_1.9.2       forcats_1.0.0         stringr_1.5.1        
[16] dplyr_1.1.4           purrr_1.0.1           readr_2.1.4          
[19] tidyr_1.3.0           tibble_3.2.1          ggplot2_3.4.4        
[22] tidyverse_2.0.0      

loaded via a namespace (and not attached):
  [1] TH.data_1.0-10       minqa_1.2.4          colorspace_2.1-0    
  [4] ellipsis_0.3.2       rio_0.5.26           estimability_1.3    
  [7] snakecase_0.11.0     htmlTable_2.2.1      base64enc_0.1-3     
 [10] rstudioapi_0.15.0    farver_2.1.1         mvtnorm_1.2-4       
 [13] fansi_1.0.5          xml2_1.3.3           codetools_0.2-16    
 [16] splines_4.0.3        mnormt_2.1.1         cachem_1.0.8        
 [19] knitr_1.45           Formula_1.2-4        jsonlite_1.8.7      
 [22] nloptr_1.2.2.2       cluster_2.1.0        png_0.1-8           
 [25] compiler_4.0.3       httr_1.4.4           emmeans_1.7.2       
 [28] backports_1.4.1      Matrix_1.5-3         fastmap_1.1.1       
 [31] cli_3.6.1            htmltools_0.5.7      tools_4.0.3         
 [34] coda_0.19-4          glue_1.6.2           reshape2_1.4.4      
 [37] Rcpp_1.0.11          carData_3.0-4        cellranger_1.1.0    
 [40] jquerylib_0.1.4      vctrs_0.6.5          svglite_2.1.0       
 [43] nlme_3.1-157         psych_2.1.6          xfun_0.41           
 [46] openxlsx_4.2.5.2     lme4_1.1-27          rvest_1.0.3         
 [49] timechange_0.1.1     lifecycle_1.0.4      formula.tools_1.7.1 
 [52] zoo_1.8-12           MASS_7.3-57          scales_1.3.0        
 [55] hms_1.1.2            sandwich_3.0-2       parallel_4.0.3      
 [58] RColorBrewer_1.1-3   yaml_2.3.7           curl_5.0.0          
 [61] sass_0.4.7           rpart_4.1-15         reshape_0.8.9       
 [64] latticeExtra_0.6-29  stringi_1.8.2        highr_0.10          
 [67] checkmate_2.3.0      boot_1.3-25          zip_2.2.0           
 [70] operator.tools_1.6.3 rlang_1.1.2          pkgconfig_2.0.3     
 [73] systemfonts_1.0.4    evaluate_0.23        lattice_0.20-41     
 [76] labeling_0.4.3       htmlwidgets_1.6.2    tidyselect_1.2.0    
 [79] GGally_2.1.2         plyr_1.8.9           magrittr_2.0.3      
 [82] bookdown_0.24        R6_2.5.1             generics_0.1.3      
 [85] Hmisc_4.5-0          multcomp_1.4-17      pillar_1.9.0        
 [88] foreign_0.8-79       withr_2.5.2          mgcv_1.8-33         
 [91] survival_3.2-7       abind_1.4-5          nnet_7.3-19         
 [94] car_3.0-10           utf8_1.2.4           tzdb_0.4.0          
 [97] rmarkdown_2.25       jpeg_0.1-9           data.table_1.14.6   
[100] digest_0.6.33        webshot_0.5.4        xtable_1.8-4        
[103] munsell_0.5.0        viridisLite_0.4.2    ez_4.4-0            
[106] bslib_0.5.1

Pro-environmental behaviour and working time – Can work-time reduction support sustainable living?

2 What factors influence individuals working hours? (primary data)

3 Session info

4 References