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 > # > # 13 Fitting Models to Repeated Measures Data . . . . . . . . . . . . . 191 > # > # 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) > > > # 13.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . .191 > > # 13.2 A basic model . . . . . . . . . . . . . . . . . . . . . . . . . . 194 > > data(reading1, package = "R2MLwiN") > summary(reading1) id age1 read1 age2 Min. : 1.0 Min. :-10.000 Min. :-10.00000 Min. :-10.0000 1st Qu.:149.0 1st Qu.:-10.000 1st Qu.:-10.00000 1st Qu.: -1.5804 Median :363.0 Median : -2.500 Median : 4.26330 Median : -1.4704 Mean :342.6 Mean : -4.834 Mean : 0.01828 Mean : -2.1946 3rd Qu.:520.5 3rd Qu.: -2.380 3rd Qu.: 4.72640 3rd Qu.: -1.3604 Max. :751.0 Max. : -1.700 Max. : 7.50500 Max. : -0.6804 read2 age3 read3 age4 Min. :-10.000 Min. :-10.0000 Min. :-10.000 Min. :-10.0000 1st Qu.: 4.894 1st Qu.: -0.6104 1st Qu.: 5.743 1st Qu.:-10.0000 Median : 5.490 Median : -0.5004 Median : 6.531 Median : 0.4196 Mean : 4.302 Mean : -1.6561 Mean : 4.581 Mean : -2.4247 3rd Qu.: 6.214 3rd Qu.: -0.4004 3rd Qu.: 6.981 3rd Qu.: 0.5396 Max. : 7.150 Max. : 0.2896 Max. : 9.907 Max. : 1.2496 read4 age5 read5 age6 Min. :-10.000 Min. :-10.000 Min. :-10.000 Min. :-10.000 1st Qu.:-10.000 1st Qu.:-10.000 1st Qu.:-10.000 1st Qu.:-10.000 Median : 7.031 Median : 1.390 Median : 7.603 Median :-10.000 Mean : 2.728 Mean : -2.607 Mean : 1.980 Mean : -3.089 3rd Qu.: 7.853 3rd Qu.: 1.520 3rd Qu.: 9.068 3rd Qu.: 4.210 Max. : 9.908 Max. : 2.250 Max. : 10.420 Max. : 4.510 read6 Min. :-10.0000 1st Qu.:-10.0000 Median :-10.0000 Mean : 0.3622 3rd Qu.: 11.4330 Max. : 13.8690 > > reading1[reading1 == -10] <- NA > > summary(reading1) id age1 read1 age2 Min. : 1.0 Min. :-2.71 Min. :3.893 Min. :-2.0104 1st Qu.:149.0 1st Qu.:-2.51 1st Qu.:4.263 1st Qu.:-1.5404 Median :363.0 Median :-2.43 Median :4.541 Median :-1.4504 Mean :342.6 Mean :-2.41 Mean :4.720 Mean :-1.4372 3rd Qu.:520.5 3rd Qu.:-2.34 3rd Qu.:5.004 3rd Qu.:-1.3504 Max. :751.0 Max. :-1.70 Max. :7.505 Max. :-0.6804 NA's :130 NA's :130 NA's :36 read2 age3 read3 age4 Min. :4.341 Min. :-1.0604 Min. :4.955 Min. :-0.0404 1st Qu.:5.107 1st Qu.:-0.5704 1st Qu.:6.080 1st Qu.: 0.3996 Median :5.575 Median :-0.4754 Median :6.643 Median : 0.4896 Mean :5.690 Mean :-0.4608 Mean :6.670 Mean : 0.4869 3rd Qu.:6.299 3rd Qu.:-0.3804 3rd Qu.:7.093 3rd Qu.: 0.5696 Max. :7.150 Max. : 0.2896 Max. :9.907 Max. : 1.2496 NA's :36 NA's :51 NA's :51 NA's :113 read4 age5 read5 age6 Min. :5.490 Min. :0.9596 Min. : 5.462 Min. :3.660 1st Qu.:6.928 1st Qu.:1.3996 1st Qu.: 7.744 1st Qu.:4.130 Median :7.545 Median :1.4896 Median : 8.617 Median :4.210 Mean :7.620 Mean :1.4838 Mean : 8.610 Mean :4.207 3rd Qu.:8.367 3rd Qu.:1.5696 3rd Qu.: 9.519 3rd Qu.:4.290 Max. :9.908 Max. :2.2496 Max. :10.420 Max. :4.510 NA's :113 NA's :145 NA's :145 NA's :209 read6 Min. : 6.451 1st Qu.:10.575 Median :11.544 Mean :11.300 3rd Qu.:12.319 Max. :13.869 NA's :209 > > reading <- reshape(reading1, idvar = "student", timevar = "id", varying = c("read1", "age1", "read2", "age2", "read3", + "age3", "read4", "age4", "read5", "age5", "read6", "age6"), sep = "", direction = "long") > > reading <- reading[c("student", "id", "age", "read")] > reading <- reading[order(reading$student, reading$id), ] > > colnames(reading) <- c("student", "occasion", "age", "reading") > rownames(reading) <- NULL > > > summary(reading) student occasion age reading Min. : 1 Min. :1.0 Min. :-2.7104 Min. : 3.893 1st Qu.:102 1st Qu.:2.0 1st Qu.:-1.4704 1st Qu.: 5.497 Median :204 Median :3.5 Median :-0.4304 Median : 6.756 Mean :204 Mean :3.5 Mean : 0.0000 Mean : 7.125 3rd Qu.:306 3rd Qu.:5.0 3rd Qu.: 1.3096 3rd Qu.: 8.279 Max. :407 Max. :6.0 Max. : 4.5096 Max. :13.869 NA's :684 NA's :684 > > head(reading, 5) student occasion age reading 1 1 1 -2.6504 6.3009 2 1 2 -1.6104 6.8948 3 1 3 -0.6204 7.7685 4 1 4 0.3396 8.8804 5 1 5 1.3396 10.0820 > > tab <- aggregate(reading ~ occasion, reading, function(x) c(N = length(x), mean = mean(x), sd = sd(x))) > tab <- rbind(tab, c(NA, NA)) > tab$reading[7, ] <- c(length(na.omit(reading$reading)), mean(na.omit(reading$reading)), sd(na.omit(reading$reading))) > rownames(tab)[7] <- "Total" > tab occasion reading.N reading.mean reading.sd 1 1 277.0000000 4.7200072 0.6099982 2 2 371.0000000 5.6899981 0.7399990 3 3 356.0000000 6.6699997 0.8699961 4 4 294.0000000 7.6200003 0.9900053 5 5 262.0000000 8.6099938 1.1199879 6 6 198.0000000 11.3000247 1.4700124 Total NA 1758.0000000 7.1254008 2.1544282 > > tab <- aggregate(age ~ occasion, reading, function(x) c(N = length(x), mean = mean(x), sd = sd(x))) > tab <- rbind(tab, c(NA, NA)) > tab$age[7, ] <- c(length(na.omit(reading$age)), mean(na.omit(reading$age)), sd(na.omit(reading$age))) > rownames(tab)[7] <- "Total" > tab occasion age.N age.mean age.sd 1 1 2.770000e+02 -2.410003e+00 1.449145e-01 2 2 3.710000e+02 -1.437192e+00 1.539733e-01 3 3 3.560000e+02 -4.607933e-01 1.571502e-01 4 4 2.940000e+02 4.868789e-01 1.269735e-01 5 5 2.620000e+02 1.483837e+00 1.307948e-01 6 6 1.980000e+02 4.206671e+00 1.150215e-01 Total NA 1.758000e+03 9.554091e-06 1.943460e+00 > > (mymodel1 <- runMLwiN(reading ~ 1 + (1 | student) + (1 | occasion), data = reading)) /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 (Normal) N min mean max N_complete min_complete mean_complete max_complete student 407 6 6 6 407 1 4.31941 6 Estimation algorithm: IGLS Elapsed time : 0.09s Number of obs: 1758 (from total 2442) The model converged after 3 iterations. Log likelihood: -3842.9 Deviance statistic: 7685.7 --------------------------------------------------------------------------------------------------- The model formula: reading ~ 1 + (1 | student) + (1 | occasion) Level 2: student Level 1: occasion --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept 7.11491 0.05303 134.17 0 *** 7.01097 7.21885 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 0.07761 0.08309 --------------------------------------------------------------------------------------------------- The random part estimates at the occasion level: Coef. Std. Err. var_Intercept 4.56149 0.17245 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > > # 13.3 A linear growth curve model . . . . . . . . . . . . . . . . . . . 201 > > (mymodel2 <- runMLwiN(reading ~ 1 + age + (1 | student) + (1 | occasion), estoptions = list(startval = list(FP.b = mymodel1@FP, + FP.v = mymodel1@FP.cov, RP.b = mymodel1@RP, RP.v = mymodel1@RP.cov)), data = reading)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved JOIN 7.11491148444759 0 '_FP_b' JOIN 0.0028122591557062 0 0 '_FP_v' JOIN 0.0776119399796296 4.56148954336666 '_RP_b' JOIN 0.00690456124179504 -0.0060668859565199 0.0297395862616497 '_RP_v' TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 iteration 4 iteration 5 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Normal) N min mean max N_complete min_complete mean_complete max_complete student 407 6 6 6 407 1 4.31941 6 Estimation algorithm: IGLS Elapsed time : 0.14s Number of obs: 1758 (from total 2442) The model converged after 6 iterations. Log likelihood: -1897.8 Deviance statistic: 3795.6 --------------------------------------------------------------------------------------------------- The model formula: reading ~ 1 + age + (1 | student) + (1 | occasion) Level 2: student Level 1: occasion --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept 7.11710 0.04128 172.41 0 *** 7.03619 7.19801 age 0.99676 0.00722 138.07 0 *** 0.98261 1.01091 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 0.60298 0.04844 --------------------------------------------------------------------------------------------------- The random part estimates at the occasion level: Coef. Std. Err. var_Intercept 0.30724 0.01181 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > (mymodel3 <- runMLwiN(reading ~ 1 + age + (1 + age | student) + (1 | occasion), estoptions = list(resi.store = TRUE, + startval = list(FP.b = mymodel2@FP, FP.v = mymodel2@FP.cov, RP.b = mymodel2@RP, RP.v = mymodel2@RP.cov)), data = reading)) /nogui option ignored ECHO 0 Echoing is ON BATC 1 Batch mode is ON MAXI 2 STAR iteration 0 iteration 1 Convergence not achieved JOIN 7.11710098696708 0.996758938786093 '_FP_b' JOIN 0.00170404929632015 1.11939939725668e-05 5.21184505439158e-05 '_FP_v' JOIN 0.602981523429859 0 0 0.307237233390386 '_RP_b' JOIN 0.00234657001927917 0 0 0 0 0 -3.83095220098178e-05 0 0 0.000139481811276294 '_RP_v' TOLE 2 MAXI 20 NEXT iteration 2 iteration 3 iteration 4 iteration 5 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Normal) N min mean max N_complete min_complete mean_complete max_complete student 407 6 6 6 407 1 4.31941 6 Estimation algorithm: IGLS Elapsed time : 0.22s Number of obs: 1758 (from total 2442) The model converged after 6 iterations. Log likelihood: -1604.7 Deviance statistic: 3209.4 --------------------------------------------------------------------------------------------------- The model formula: reading ~ 1 + age + (1 + age | student) + (1 | occasion) Level 2: student Level 1: occasion --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept 7.11726 0.04335 164.19 0 *** 7.03230 7.20222 age 0.99506 0.01219 81.63 0 *** 0.97117 1.01896 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 0.68263 0.05302 cov_Intercept_age 0.12324 0.01237 var_age 0.03735 0.00386 --------------------------------------------------------------------------------------------------- The random part estimates at the occasion level: Coef. Std. Err. var_Intercept 0.16101 0.00705 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > u0 <- mymodel3@residual$lev_2_resi_est_Intercept > u0std <- (u0 - mean(u0))/sd(u0) > > u1 <- mymodel3@residual$lev_2_resi_est_age > u1std <- (u1 - mean(u1))/sd(u1) > > plot(u0std, u1std, asp = 1) > > e0 <- na.omit(mymodel3@residual$lev_1_resi_est_Intercept) > e0std <- (e0 - mean(e0))/sd(e0) > e0rank <- rank(e0) > e0uniform <- (e0rank - 0.5)/length(e0rank) > e0nscore <- qnorm(e0uniform) > > plot(e0std, e0nscore, asp = 1) > > # 13.4 Complex level 1 variation . . . . . . . . . . . . . . . . . . . . 204 > > (mymodel4 <- runMLwiN(reading ~ 1 + age + (1 + age | student) + (1 + age | occasion), data = reading)) /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 iteration 8 Convergence achieved ECHO 0 Execution completed -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- MLwiN (version: 3.05) multilevel model (Normal) N min mean max N_complete min_complete mean_complete max_complete student 407 6 6 6 407 1 4.31941 6 Estimation algorithm: IGLS Elapsed time : 0.36s Number of obs: 1758 (from total 2442) The model converged after 9 iterations. Log likelihood: -1588.5 Deviance statistic: 3177 --------------------------------------------------------------------------------------------------- The model formula: reading ~ 1 + age + (1 + age | student) + (1 + age | occasion) Level 2: student Level 1: occasion --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept 7.11696 0.04390 162.12 0 *** 7.03092 7.20300 age 0.99505 0.01232 80.75 0 *** 0.97089 1.01920 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 0.70489 0.05437 cov_Intercept_age 0.12797 0.01278 var_age 0.03546 0.00404 --------------------------------------------------------------------------------------------------- The random part estimates at the occasion level: Coef. Std. Err. var_Intercept 0.12268 0.00825 cov_Intercept_age 0.00067 0.00273 var_age 0.01451 0.00305 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > # 13.5 Repeated measures modelling of non-linear polynomial growth . . . 205 > > (mymodel5 <- runMLwiN(reading ~ 1 + age + I(age^2) + (1 + age + I(age^2) | student) + (1 + age | occasion), estoptions = list(resi.store = TRUE), + data = reading)) /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 (Normal) N min mean max N_complete min_complete mean_complete max_complete student 407 6 6 6 407 1 4.31941 6 Estimation algorithm: IGLS Elapsed time : 0.44s Number of obs: 1758 (from total 2442) The model converged after 8 iterations. Log likelihood: -1566 Deviance statistic: 3132 --------------------------------------------------------------------------------------------------- The model formula: reading ~ 1 + age + I(age^2) + (1 + age + I(age^2) | student) + (1 + age | occasion) Level 2: student Level 1: occasion --------------------------------------------------------------------------------------------------- The fixed part estimates: Coef. Std. Err. z Pr(>|z|) [95% Conf. Interval] Intercept 7.11503 0.04629 153.72 0 *** 7.02431 7.20575 age 0.99454 0.01259 79.00 0 *** 0.96987 1.01921 I(age^2) 0.00084 0.00316 0.26 0.7916 -0.00537 0.00704 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 0.76685 0.06009 cov_Intercept_age 0.14083 0.01390 var_age 0.03939 0.00427 cov_Intercept_I(age^2) -0.01435 0.00306 cov_age_I(age^2) -0.00178 0.00077 var_I(age^2) 0.00133 0.00031 --------------------------------------------------------------------------------------------------- The random part estimates at the occasion level: Coef. Std. Err. var_Intercept 0.12908 0.00857 cov_Intercept_age -0.00608 0.00267 var_age 0.00039 0.00392 -*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- > > l2varfn <- mymodel5@RP["RP2_var_Intercept"] + 2 * mymodel5@RP["RP2_cov_Intercept_age"] * mymodel5@data$age + mymodel5@RP["RP2_var_age"] * + mymodel5@data$age^2 + 2 * mymodel5@RP["RP2_cov_Intercept_I(age^2)"] * mymodel5@data$age^2 + 2 * mymodel5@RP["RP2_cov_age_I(age^2)"] * + mymodel5@data$age * mymodel5@data[["I(age^2)"]] + mymodel5@RP["RP2_var_I(age^2)"] * mymodel5@data[["I(age^2)"]]^2 > > l1varfn <- mymodel5@RP["RP1_var_Intercept"] + 2 * mymodel5@RP["RP1_cov_Intercept_age"] * mymodel5@data$age + mymodel5@RP["RP1_var_age"] * + mymodel5@data$age^2 > > totvarfn <- l2varfn + l1varfn > > plot(mymodel5@data$age, totvarfn) > > xb <- predict(mymodel5) > > u0 <- mymodel5@residual$lev_2_resi_est_Intercept > u1 <- mymodel5@residual$lev_2_resi_est_age > u2 <- mymodel5@residual[["lev_2_resi_est_I(age^2)"]] > > yhat <- xb + u0[mymodel5@data$student] + u1[mymodel5@data$student] * mymodel5@data$age + u2[mymodel5@data$student] * + mymodel5@data[["I(age^2)"]] > > plot(mymodel5@data$age, yhat, type = "n") > lines(mymodel5@data$age[mymodel5@data$student == 1], yhat[mymodel5@data$student == 1], col = 1) > lines(mymodel5@data$age[mymodel5@data$student == 2], yhat[mymodel5@data$student == 2], col = 2) > lines(mymodel5@data$age[mymodel5@data$student == 3], yhat[mymodel5@data$student == 3], col = 3) > lines(mymodel5@data$age[mymodel5@data$student == 4], yhat[mymodel5@data$student == 4], col = 4) > > plot(mymodel5@data$age, yhat, type = "n") > lines(mymodel5@data$age[mymodel5@data$student == 10], yhat[mymodel5@data$student == 10], col = 1) > lines(mymodel5@data$age[mymodel5@data$student == 11], yhat[mymodel5@data$student == 11], col = 2) > lines(mymodel5@data$age[mymodel5@data$student == 12], yhat[mymodel5@data$student == 12], col = 3) > lines(mymodel5@data$age[mymodel5@data$student == 13], yhat[mymodel5@data$student == 13], col = 4) > lines(mymodel5@data$age[mymodel5@data$student == 14], yhat[mymodel5@data$student == 14], col = 4) > > # Chapter learning outcomes . . . . . . . . . . . . . . . . . . . . . . .209 > > ############################################################################ > > proc.time() user system elapsed 4.03 0.51 6.87