R version 3.6.3 (2020-02-29) -- "Holding the Windsock" Copyright (C) 2020 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 (64-bit) R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. Natural language support but running in an English locale R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > ############################################################################ > # MLwiN User Manual > # > # 14 Multivariate Response Models . . . . . . . . . . . . . . . . . . . 211 > # > # Rasbash, J., Steele, F., Browne, W. J. and Goldstein, H. (2012). > # A User's Guide to MLwiN, v2.26. Centre for Multilevel Modelling, > # University of Bristol. > ############################################################################ > # R script to replicate all analyses using R2MLwiN > # > # Zhang, Z., Charlton, C., Parker, R, Leckie, G., and Browne, W.J. > # Centre for Multilevel Modelling, 2012 > # http://www.bristol.ac.uk/cmm/software/R2MLwiN/ > ############################################################################ > > library(R2MLwiN) R2MLwiN: A package to run models implemented in MLwiN from R Copyright 2013-2017 Zhengzheng Zhang, Christopher M. J. Charlton, Richard M. A. Parker, William J. Browne and George Leckie Support provided by the Economic and Social Research Council (ESRC) (Grants RES-149-25-1084, RES-576-25-0032 and ES/K007246/1) To cite R2MLwiN in publications use: Zhengzheng Zhang, Richard M. A. Parker, Christopher M. J. Charlton, George Leckie, William J. Browne (2016). R2MLwiN: A Package to Run MLwiN from within R. Journal of Statistical Software, 72(10), 1-43. doi:10.18637/jss.v072.i10 A BibTeX entry for LaTeX users is @Article{, title = {{R2MLwiN}: A Package to Run {MLwiN} from within {R}}, author = {Zhengzheng Zhang and Richard M. A. Parker and Christopher M. J. Charlton and George Leckie and William J. Browne}, journal = {Journal of Statistical Software}, year = {2016}, volume = {72}, number = {10}, pages = {1--43}, doi = {10.18637/jss.v072.i10}, } The MLwiN_path option is currently set to C:/Program Files/MLwiN v3.05/ To change this use: options(MLwiN_path="") > # MLwiN folder > mlwin <- getOption("MLwiN_path") > while (!file.access(mlwin, mode = 1) == 0) { + cat("Please specify the root MLwiN folder or the full path to the MLwiN executable:\n") + mlwin <- scan(what = character(0), sep = "\n") + mlwin <- gsub("\\", "/", mlwin, fixed = TRUE) + } > options(MLwiN_path = mlwin) > > > # 14.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . .211 > > data(gcsemv1, package = "R2MLwiN") > summary(gcsemv1) school student female agemths written Min. :20920 Min. : 1 Male : 777 Min. :186.0 Min. : 0.625 1st Qu.:60501 1st Qu.: 64 Female:1128 1st Qu.:205.0 1st Qu.:37.500 Median :68133 Median : 133 Median :510.0 Median :46.875 Mean :62128 Mean :1037 Mean :389.8 Mean :46.798 3rd Qu.:68411 3rd Qu.: 458 3rd Qu.:510.0 3rd Qu.:55.625 Max. :84772 Max. :5521 Max. :510.0 Max. :90.000 NA's :202 csework cons Min. : 9.259 Min. :1 1st Qu.: 62.963 1st Qu.:1 Median : 75.926 Median :1 Mean : 73.435 Mean :1 3rd Qu.: 86.111 3rd Qu.:1 Max. :100.000 Max. :1 NA's :180 > > # 14.2 Specifying a multivariate model . . . . . . . . . . . . . . . . . 212 > > # 14.3 Setting up the basic model . . . . . . . . . . . . . . . . . . . .214 > > (mymodel1 <- runMLwiN(c(written, csework) ~ 1 + (1 | student), D = "Multivariate Normal", estoptions = list(sort.ignore = TRUE), + data = gcsemv1)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved TOLE 2 MAXI 20 NEXT iteration 2 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Multivariate Normal) Estimation algorithm: IGLS Elapsed time : 0.29s Number of obs: 1905 (from total 1905) The model converged after 3 iterations. Log likelihood: -13903.9 Deviance statistic: 27807.9 --------------------------------------------------------------------------------------------------- The model formula: c(written, csework) ~ 1 + (1 | student) Level 1: student --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept_written 46.80278 0.32017 146.18 0 *** 46.17527 47.43030 Intercept_csework 73.36400 0.38821 188.98 0 *** 72.60313 74.12488 Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 --------------------------------------------------------------------------------------------------- The random part estimates at the student level: Coef. Std. Err. var_Intercept_written 178.71008 6.10784 cov_Intercept_written_Intercept_csework 102.31138 5.91751 var_Intercept_csework 265.44838 9.01669 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > (mymodel2 <- runMLwiN(c(written, csework) ~ 1 + female + (1 | school) + (1 | student), D = "Multivariate Normal", + data = gcsemv1)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Multivariate Normal) N min mean max N_complete min_complete mean_complete max_complete school 73 2 26.09589 104 73 2 26.09589 104 Estimation algorithm: IGLS Elapsed time : 0.31s Number of obs: 1905 (from total 1905) The model converged after 4 iterations. Log likelihood: -13400.2 Deviance statistic: 26800.5 --------------------------------------------------------------------------------------------------- The model formula: c(written, csework) ~ 1 + female + (1 | school) + (1 | student) Level 2: school Level 1: student --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept_written 49.45213 0.93384 52.96 0 *** 47.62183 51.28243 Intercept_csework 69.67166 1.17179 59.46 0 *** 67.37500 71.96831 femaleFemale_written -2.50295 0.56072 -4.46 8.052e-06 *** -3.60194 -1.40396 femaleFemale_csework 6.75139 0.67065 10.07 7.734e-24 *** 5.43694 8.06584 Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 --------------------------------------------------------------------------------------------------- The random part estimates at the school level: Coef. Std. Err. var_Intercept_written 46.81298 9.18733 cov_Intercept_written_Intercept_csework 24.87783 8.88036 var_Intercept_csework 75.16623 14.56485 --------------------------------------------------------------------------------------------------- The random part estimates at the student level: Coef. Std. Err. var_Intercept_written 124.63432 4.34983 cov_Intercept_written_Intercept_csework 73.00323 4.17829 var_Intercept_csework 180.09817 6.24580 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > mymodel2@RP["RP2_cov_Intercept_written_Intercept_csework"]/sqrt(mymodel2@RP["RP2_var_Intercept_written"] * mymodel2@RP["RP2_var_Intercept_csework"]) RP2_cov_Intercept_written_Intercept_csework 0.4193899 > > mymodel2@RP["RP1_cov_Intercept_written_Intercept_csework"]/sqrt(mymodel2@RP["RP1_var_Intercept_written"] * mymodel2@RP["RP1_var_Intercept_csework"]) RP1_cov_Intercept_written_Intercept_csework 0.4872688 > > # 14.4 A more elaborate model . . . . . . . . . . . . . . . . . . . . . .219 > > (mymodel3 <- runMLwiN(c(written, csework) ~ 1 + female + (1 + female | school) + (1 | student), D = "Multivariate Normal", + data = gcsemv1)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 iteration 4 iteration 5 iteration 6 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Multivariate Normal) N min mean max N_complete min_complete mean_complete max_complete school 73 2 26.09589 104 73 2 26.09589 104 Estimation algorithm: IGLS Elapsed time : 0.78s Number of obs: 1905 (from total 1905) The model converged after 7 iterations. Log likelihood: -13378.1 Deviance statistic: 26756.1 --------------------------------------------------------------------------------------------------- The model formula: c(written, csework) ~ 1 + female + (1 + female | school) + (1 | student) Level 2: school Level 1: student --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept_written 49.40112 0.99594 49.60 0 *** 47.44912 51.35312 Intercept_csework 69.30076 1.35678 51.08 0 *** 66.64152 71.96000 femaleFemale_written -2.47108 0.64394 -3.84 0.0001243 *** -3.73318 -1.20898 femaleFemale_csework 7.15671 1.13530 6.30 2.904e-10 *** 4.93157 9.38186 Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 --------------------------------------------------------------------------------------------------- The random part estimates at the school level: Coef. Std. Err. var_Intercept_written 54.54177 11.77079 cov_Intercept_written_Intercept_csework 38.89313 12.68564 var_Intercept_csework 104.69574 21.82904 cov_Intercept_written_femaleFemale_written -7.44933 5.68375 cov_Intercept_csework_femaleFemale_written -6.97482 7.29527 var_femaleFemale_written 5.28220 4.23883 cov_Intercept_written_femaleFemale_csework -21.10519 9.96561 cov_Intercept_csework_femaleFemale_csework -39.67547 14.69664 cov_femaleFemale_written_femaleFemale_csework 11.00389 6.29634 var_femaleFemale_csework 49.90979 14.60603 --------------------------------------------------------------------------------------------------- The random part estimates at the student level: Coef. Std. Err. var_Intercept_written 123.40991 4.37104 cov_Intercept_written_Intercept_csework 70.55722 4.10549 var_Intercept_csework 169.80446 5.99977 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > (mymodel4 <- runMLwiN(c(written, csework) ~ 1 + female + (1 + female[1] | school) + (1 | student), D = "Multivariate Normal", + estoptions = list(resi.store = TRUE), data = gcsemv1)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 iteration 4 iteration 5 iteration 6 iteration 7 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Multivariate Normal) N min mean max N_complete min_complete mean_complete max_complete school 73 2 26.09589 104 73 2 26.09589 104 Estimation algorithm: IGLS Elapsed time : 0.83s Number of obs: 1905 (from total 1905) The model converged after 8 iterations. Log likelihood: -13399.6 Deviance statistic: 26799.2 --------------------------------------------------------------------------------------------------- The model formula: c(written, csework) ~ 1 + female + (1 + female[1] | school) + (1 | student) Level 2: school Level 1: student --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept_written 49.51235 0.94927 52.16 0 *** 47.65181 51.37289 Intercept_csework 69.65891 1.17223 59.42 0 *** 67.36139 71.95644 femaleFemale_1 -2.58355 0.63117 -4.09 4.254e-05 *** -3.82063 -1.34648 femaleFemale_2 6.75772 0.67057 10.08 6.943e-24 *** 5.44342 8.07201 Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 --------------------------------------------------------------------------------------------------- The random part estimates at the school level: Coef. Std. Err. var_Intercept_written 48.19051 10.44075 cov_Intercept_written_Intercept_csework 23.72850 9.34337 var_Intercept_csework 75.23904 14.57476 cov_Intercept_written_femaleFemale_1 -2.54470 4.60831 cov_Intercept_csework_femaleFemale_1 2.00384 5.27357 var_femaleFemale_1 4.10047 3.41483 --------------------------------------------------------------------------------------------------- The random part estimates at the student level: Coef. Std. Err. var_Intercept_written 123.99295 4.37418 cov_Intercept_written_Intercept_csework 73.30694 4.18134 var_Intercept_csework 180.08374 6.24510 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > mymodel4@RP["RP2_cov_Intercept_written_Intercept_csework"]/sqrt(mymodel4@RP["RP2_var_Intercept_written"] * mymodel4@RP["RP2_var_Intercept_csework"]) RP2_cov_Intercept_written_Intercept_csework 0.394065 > > mymodel4@RP["RP2_cov_Intercept_written_femaleFemale_1"]/sqrt(mymodel4@RP["RP2_var_Intercept_written"] * mymodel4@RP["RP2_var_femaleFemale_1"]) RP2_cov_Intercept_written_femaleFemale_1 -0.1810253 > > mymodel4@RP["RP2_cov_Intercept_csework_femaleFemale_1"]/sqrt(mymodel4@RP["RP2_var_Intercept_csework"] * mymodel4@RP["RP2_var_femaleFemale_1"]) RP2_cov_Intercept_csework_femaleFemale_1 0.114084 > > u0 <- mymodel4@residual$lev_2_resi_est_Intercept.written > u1 <- mymodel4@residual$lev_2_resi_est_Intercept.csework > u2 <- mymodel4@residual$lev_2_resi_est_femaleFemale.1 > > plot(u0, u0, asp = 1) > plot(u0, u1, asp = 1) > plot(u0, u2, asp = 1) > plot(u1, u1, asp = 1) > plot(u1, u2, asp = 1) > plot(u2, u2, asp = 1) > > # 14.5 Multivariate models for discrete responses . . . . . . . . . . . .222 > > data(tutorial, package = "R2MLwiN") > > tutorial$binexam <- as.integer(tutorial$normexam > 0) > tutorial$binlrt <- as.integer(tutorial$standlrt > 0) > > (mymodel5 <- runMLwiN(c(logit(binexam), logit(binlrt)) ~ 1, D = c("Mixed", "Binomial", "Binomial"), estoptions = list(sort.ignore = TRUE), + data = tutorial)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Mixed) Estimation algorithm: IGLS MQL1 Elapsed time : 0.55s Number of obs: 4059 (from total 4059) The model converged after 4 iterations. Log likelihood: NA Deviance statistic: NA --------------------------------------------------------------------------------------------------- The model formula: c(logit(binexam), logit(binlrt)) ~ 1 Level 1: l1id --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept_binexam 0.04879 0.03140 1.55 0.1202 -0.01276 0.11034 Intercept_binlrt 0.06062 0.03141 1.93 0.05357 . -0.00093 0.12218 Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1 --------------------------------------------------------------------------------------------------- The random part estimates at the l1id level: Coef. Std. Err. var_bcons_1 1.00000 0.00000 cov_bcons_1_bcons_2 0.41914 0.01193 var_bcons_2 1.00000 0.00000 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > > # Chapter learning outcomes . . . . . . . . . . . . . . . . . . . . . . .224 > > ############################################################################ > > proc.time() user system elapsed 3.64 0.45 7.90