R version 4.4.1 (2024-06-14 ucrt) -- "Race for Your Life"
Copyright (C) 2024 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

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.

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
> #
> # 9   Logistic Models for Binary and Binomial Responses . . . . . . . . .117
> #
> #     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-2024 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.11/
To change this use: options(MLwiN_path="<path to MLwiN>")
> # 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)
> 
> # Change contrasts if wish to avoid warning indicating that, by default,
> # specified contrasts for ordered predictors will be ignored by runMLwiN
> # (they will be fitted as "contr.treatment" regardless of this setting). To
> # enable specified contrasts, set allowcontrast to TRUE (this will be the
> # default in future package releases). NB at the end of this script, the
> # specification for contrasts is changed back.
> my_contrasts <- options("contrasts")$contrasts
> options(contrasts = c(unordered = "contr.treatment",
+                       ordered = "contr.treatment"))
> 
> # As an alternative to changing contrasts, can instead use C() to specify
> # contrasts for ordered predictors in formula object, e.g.:
> 
> # (mymodel1 <- runMLwiN(logit(use) ~ 1 + C(lc, "contr.treatment"),
> #                      D = "Binomial",
> #                      data = bang,
> #                      allowcontrast = TRUE))
> 
> # 9.1 Introduction and description of the example data . . . . . . . . . 117
> 
> data(bang, package = "R2MLwiN")
> summary(bang)
     woman           district            use      
 Min.   :   1.0   Min.   : 1.00   Not_using:1728  
 1st Qu.: 717.5   1st Qu.:14.00   Using    :1139  
 Median :1434.0   Median :29.00                   
 Mean   :1434.0   Mean   :29.25                   
 3rd Qu.:2150.5   3rd Qu.:45.00                   
 Max.   :2867.0   Max.   :61.00                   
                       use4                 lc            age          
 Sterilization           : 302   None        : 774   Min.   :-14.0000  
 Modern_reversible_method: 555   One_child   : 517   1st Qu.: -8.0000  
 Traditional_method      : 282   Two_children: 461   Median : -2.0000  
 Not_using_contraception :1728   Three_plus  :1115   Mean   : -0.3279  
                                                     3rd Qu.:  6.0000  
                                                     Max.   : 19.0000  
   urban                       educ         hindu          d_lit       
 Rural:2063   None               :1806   Muslim:2480   Min.   :0.0000  
 Urban: 804   Lower_primary      : 357   Hindu : 387   1st Qu.:0.0850  
              Upper_primary      : 265                 Median :0.1100  
              Secondary_and_above: 439                 Mean   :0.1115  
                                                       3rd Qu.:0.1400  
                                                       Max.   :0.3000  
     d_pray            cons  
 Min.   :0.1000   Min.   :1  
 1st Qu.:0.2900   1st Qu.:1  
 Median :0.4100   Median :1  
 Mean   :0.4253   Mean   :1  
 3rd Qu.:0.5500   3rd Qu.:1  
 Max.   :0.7800   Max.   :1  
> 
> # 9.2 Single-level logistic regression . . . . . . . . . . . . . . . . . 119
> 
> # Link functions . . . . . . . . . . . . . . . . . . . . . . . . . . . . 119
> 
> # Interpretation of coeficients . . . . . . . . . . . . . . . . . . . . .120
> 
> # Fitting a single-level logit model in MLwiN . . . . . . . . . . . . . .121
> 
> addmargins(with(bang, table(lc, use)))
              use
lc             Not_using Using  Sum
  None               584   190  774
  One_child          283   234  517
  Two_children       234   227  461
  Three_plus         627   488 1115
  Sum               1728  1139 2867
> 
> (mymodel1 <- runMLwiN(logit(use) ~ 1 + lc, D = "Binomial", data = bang))
/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.11)  multilevel model (Binomial) 
Estimation algorithm:  IGLS MQL1        Elapsed time : 0.09s 
Number of obs:  2867 (from total 2867)        The model converged after 4 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc
Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                    Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept        -1.12288     0.08348   -13.45    3.05e-41   ***     -1.28650    -0.95926 
lcOne_child       0.93275     0.12156     7.67   1.676e-14   ***      0.69450     1.17100 
lcTwo_children    1.09250     0.12509     8.73   2.466e-18   ***      0.84733     1.33768 
lcThree_plus      0.87224     0.10302     8.47   2.523e-17   ***      0.67033     1.07416 
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 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> if (!require(car)) {
+   warning("car package required to use linearHypothesis() function")
+ } else {
+   linearHypothesis(mymodel1, "FP_lcOne_child = FP_lcTwo_children")
+ }
Loading required package: car
Loading required package: carData
Linear hypothesis test

Hypothesis:
FP_lcOne_child - FP_lcTwo_children = 0

Model 1: restricted model
Model 2: logit(use) ~ 1 + lc

  Res.Df Df  Chisq Pr(>Chisq)
1   2863                     
2   2862  1 1.5481     0.2134
> 
> # A probit model . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126
> 
> (mymodel2 <- runMLwiN(probit(use) ~ 1 + lc, D = "Binomial", data = bang))
/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.11)  multilevel model (Binomial) 
Estimation algorithm:  IGLS MQL1        Elapsed time : 0.1s 
Number of obs:  2867 (from total 2867)        The model converged after 4 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
probit(use) ~ 1 + lc
Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                    Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept        -0.68879     0.04915   -14.01   1.299e-44   ***     -0.78513    -0.59245 
lcOne_child       0.56972     0.07396     7.70   1.328e-14   ***      0.42476     0.71468 
lcTwo_children    0.66976     0.07631     8.78    1.69e-18   ***      0.52018     0.81933 
lcThree_plus      0.53191     0.06195     8.59   8.977e-18   ***      0.41049     0.65332 
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 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> (mymodel3 <- runMLwiN(logit(use) ~ 1 + lc + age, D = "Binomial", data = bang))
/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.11)  multilevel model (Binomial) 
Estimation algorithm:  IGLS MQL1        Elapsed time : 0.09s 
Number of obs:  2867 (from total 2867)        The model converged after 4 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age
Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                    Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept        -1.25598     0.09776   -12.85   8.821e-38   ***     -1.44758    -1.06438 
lcOne_child       0.99131     0.12376     8.01    1.15e-15   ***      0.74874     1.23388 
lcTwo_children    1.22356     0.13480     9.08    1.12e-19   ***      0.95935     1.48777 
lcThree_plus      1.11655     0.13824     8.08   6.649e-16   ***      0.84560     1.38750 
age              -0.01629     0.00609    -2.67    0.007514   **      -0.02823    -0.00435 
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 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> # 9.3 A two-level random intercept model . . . . . . . . . . . . . . . . 128
> 
> # Model specification . . . . . . . . . . . . . . . . . . . . . . . . . .128
> 
> # Estimation procedures . . . . . . . . . . . . . . . . . . . . . . . . .128
> 
> # Fitting a two-level random intercept model in MLwiN . . . . . . . . . .129
> 
> (mymodel4 <- runMLwiN(logit(use) ~ 1 + lc + age + (1 | district), D = "Binomial", data = bang))
/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

Convergence achieved
ECHO 0
Execution completed


-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
MLwiN (version: 3.11)  multilevel model (Binomial) 
          N min     mean max N_complete min_complete mean_complete max_complete
district 60   3 47.78333 173         60            3      47.78333          173
Estimation algorithm:  IGLS MQL1        Elapsed time : 0.11s 
Number of obs:  2867 (from total 2867)        The model converged after 5 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age + (1 | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                    Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept        -1.36711     0.12338   -11.08   1.557e-28   ***     -1.60893    -1.12529 
lcOne_child       0.98998     0.12643     7.83   4.869e-15   ***      0.74218     1.23777 
lcTwo_children    1.27523     0.13816     9.23   2.711e-20   ***      1.00443     1.54603 
lcThree_plus      1.21568     0.14245     8.53   1.413e-17   ***      0.93648     1.49487 
age              -0.01878     0.00625    -3.00    0.002659   **      -0.03102    -0.00653 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                  Coef.   Std. Err. 
var_Intercept   0.27409     0.07138 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> (mymodel5 <- runMLwiN(logit(use) ~ 1 + lc + age + (1 | district), D = "Binomial", estoptions = list(nonlinear = c(N = 1,
+   M = 2), startval = list(FP.b = mymodel4@FP, FP.v = mymodel4@FP.cov, RP.b = mymodel4@RP, RP.v = mymodel4@RP.cov)),
+   data = bang))
/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 -1.3671107908183 0.989975145198351 1.27523297062681 1.21567634144111 -0.0187756435513459 '_FP_b'
JOIN 0.0152221068761287 -0.00857801335032769 0.0159844374382819 -0.010151723167991 0.00868106697570791 0.0190895131628167 -0.0124154753668934 0.0096459226360998 0.012442006063399 0.0202921736477393 0.000326589121711455 -0.000146849092291606 -0.000325928858650572 -0.000596278835798142 3.9048334047441e-05 '_FP_v'
JOIN 0.274088563383711 1 '_RP_b'
JOIN 0.0050956166989184 0 0 '_RP_v'
TOLE 2
MAXI 20
NEXT
iteration 2
iteration 3
iteration 4

Convergence achieved
ECHO 0
Execution completed


-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
MLwiN (version: 3.11)  multilevel model (Binomial) 
          N min     mean max N_complete min_complete mean_complete max_complete
district 60   3 47.78333 173         60            3      47.78333          173
Estimation algorithm:  IGLS PQL2        Elapsed time : 0.14s 
Number of obs:  2867 (from total 2867)        The model converged after 5 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age + (1 | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                    Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept        -1.46602     0.12791   -11.46   2.059e-30   ***     -1.71671    -1.21532 
lcOne_child       1.06285     0.12882     8.25   1.575e-16   ***      0.81037     1.31533 
lcTwo_children    1.37010     0.14167     9.67   4.014e-22   ***      1.09242     1.64778 
lcThree_plus      1.30391     0.14595     8.93   4.104e-19   ***      1.01786     1.58997 
age              -0.02005     0.00640    -3.13    0.001747   **      -0.03260    -0.00749 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                  Coef.   Std. Err. 
var_Intercept   0.30776     0.07899 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
>  if (!require(car)) {
+   warning("car package required to use linearHypothesis() function")
+ } else {
+   linearHypothesis(mymodel5, "RP2_var_Intercept = 0")
+ }
Linear hypothesis test

Hypothesis:
RP2_var_Intercept = 0

Model 1: restricted model
Model 2: logit(use) ~ 1 + lc + age + (1 | district)

  Res.Df Df  Chisq Pr(>Chisq)    
1   2861                         
2   2860  1 15.181  9.769e-05 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> 
> # Variance partition coeficient . . . . . . . . . . . . . . . . . . . . .131
> 
> set.seed(1)
> 
> invlogit <- function(x) exp(x)/(1 + exp(x))
> 
> u <- sqrt(coef(mymodel5)["RP2_var_Intercept"]) * qnorm(runif(5000))
> 
> p1 <- invlogit(coef(mymodel5)["FP_Intercept"] + u)
> 
> p2 <- invlogit(coef(mymodel5)["FP_Intercept"] + coef(mymodel5)["FP_lc3plus"] + coef(mymodel5)["FP_age"] * -9.7 + u)
> 
> p3 <- invlogit(coef(mymodel5)["FP_Intercept"] + coef(mymodel5)["FP_age"] * 15.3 + u)
> 
> v1 <- p1 * (1 - p1)
> lev2var1 <- sd(p1)^2
> lev1var1 <- mean(v1)
> 
> v2 <- p2 * (1 - p2)
> lev2var2 <- sd(p2)^2
> lev1var2 <- mean(v2)
> 
> v3 <- p3 * (1 - p3)
> lev2var3 <- sd(p3)^2
> lev1var3 <- mean(v3)
> 
> cat(paste0("VPC = ", lev2var1/(lev2var1 + lev1var1)))
VPC = 0.0491754301807682> 
> cat(paste0("VPC for a young women with 3+ children (low probability use) = ", lev2var2/(lev2var2 + lev1var2)))
VPC for a young women with 3+ children (low probability use) = NA> 
> cat(paste0("VPC for an old woman with no children (high probability use) = ", lev2var3/(lev2var3 + lev1var3)))
VPC for an old woman with no children (high probability use) = 0.0419272279288916> 
> # Adding further explanatory variables . . . . . . . . . . . . . . . . . 134
> 
> table(bang$educ)

               None       Lower_primary       Upper_primary Secondary_and_above 
               1806                 357                 265                 439 
> 
> (mymodel6 <- runMLwiN(logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 | district), D = "Binomial", estoptions = list(nonlinear = c(N = 1,
+   M = 2), startval = list(FP.b = mymodel5@FP, FP.v = mymodel5@FP.cov, RP.b = mymodel5@RP, RP.v = mymodel5@RP.cov)),
+   data = bang))
/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 -1.46601726374528 1.06284833330899 1.37009902258045 1.30391318250182 -0.0200481959375968 0 0 0 0 0 '_FP_b'
JOIN 0.0163603695419892 -0.00890762958205187 0.0165943547340283 -0.0106419958708061 0.00897007438476784 0.0200716252826363 -0.0130375955794902 0.0100149271005426 0.0130197255583241 0.0213006854805154 0.000347319944589179 -0.000157647285641369 -0.000346447992241024 -0.00063025959287392 4.10233254113274e-05 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 '_FP_v'
JOIN 0.307757058623237 1 '_RP_b'
JOIN 0.00623907692791556 0 0 '_RP_v'
TOLE 2
MAXI 20
NEXT
iteration 2
iteration 3
iteration 4

Convergence achieved
ECHO 0
Execution completed


-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
MLwiN (version: 3.11)  multilevel model (Binomial) 
          N min     mean max N_complete min_complete mean_complete max_complete
district 60   3 47.78333 173         60            3      47.78333          173
Estimation algorithm:  IGLS PQL2        Elapsed time : 0.17s 
Number of obs:  2867 (from total 2867)        The model converged after 5 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                             Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept                 -2.05249     0.13819   -14.85   6.704e-50   ***     -2.32334    -1.78164 
lcOne_child                1.15141     0.13413     8.58    9.12e-18   ***      0.88853     1.41429 
lcTwo_children             1.51227     0.14734    10.26   1.026e-24   ***      1.22349     1.80105 
lcThree_plus               1.50191     0.15271     9.83   7.972e-23   ***      1.20260     1.80122 
age                       -0.01736     0.00665    -2.61    0.009033   **      -0.03039    -0.00433 
urbanUrban                 0.53306     0.10482     5.09   3.665e-07   ***      0.32762     0.73850 
educLower_primary          0.24654     0.12836     1.92     0.05478   .       -0.00505     0.49812 
educUpper_primary          0.72433     0.14380     5.04   4.731e-07   ***      0.44248     1.00618 
educSecondary_and_above    1.17020     0.12716     9.20   3.506e-20   ***      0.92096     1.41944 
hinduHindu                 0.43282     0.12765     3.39   0.0006973   ***      0.18263     0.68301 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                  Coef.   Std. Err. 
var_Intercept   0.23364     0.06534 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> # 9.4 A two-level random coeficient model . . . . . . . . . . . . . . . .135
> 
> (mymodel7 <- runMLwiN(logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 + urban | district), D = "Binomial",
+   estoptions = list(nonlinear = c(N = 1, M = 2), startval = list(FP.b = mymodel6@FP, FP.v = mymodel6@FP.cov, RP.b = mymodel6@RP,
+     RP.v = mymodel6@RP.cov)), data = bang))
/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 -2.0524878126349 1.1514091831992 1.51227328126239 1.50190931077976 -0.0173564094857298 0.533063022676533 0.246535751852934 0.724328418970923 1.17019898626514 0.432817014945708 '_FP_b'
JOIN 0.0190969057181895 -0.010249298461099 0.0179897690949509 -0.0123640566867248 0.00988274764015615 0.0217090961479245 -0.0151256474103302 0.0110367966927285 0.0143831177429424 0.0233214607718276 0.000362885872261872 -0.000160504805163484 -0.000363883885511234 -0.000670143247882038 4.41947424646325e-05 -0.00280269546137152 0.0005261144156233 0.000849328751980673 0.000702752807306554 -1.93563397557956e-05 0.0109869565822594 -0.00283420656434017 0.00036229338945785 0.000199158958934153 -0.000127364884113156 4.59999223868248e-05 -0.00059630716642084 0.0164773314339917 -0.00265716065127714 -0.000132304943914623 -1.94267258205956e-05 -0.00012744129452743 7.55732006752395e-05 -0.000986714124785436 0.00312301236770023 0.0206794937660382 -0.00433453673756542 0.00178163571915053 0.00218751328572626 0.00282048124886048 3.59899532831458e-05 -0.00290733241170694 0.00339389445448919 0.00374076332064882 0.0161709104496365 -0.00310427091118153 0.000291548581329297 0.000563339535587039 0.00123544577460758 -5.5691862556009e-05 6.03577341367192e-05 2.13428023572961e-05 -0.000294041097998737 0.000123287082941657 0.0162946578344111 '_FP_v'
JOIN 0.233641534263787 0 0 1 '_RP_b'
JOIN 0.00426981357741411 0 0 0 0 0 0 0 0 8.77180296005631e-36 '_RP_v'
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.11)  multilevel model (Binomial) 
          N min     mean max N_complete min_complete mean_complete max_complete
district 60   3 47.78333 173         60            3      47.78333          173
Estimation algorithm:  IGLS PQL2        Elapsed time : 0.25s 
Number of obs:  2867 (from total 2867)        The model converged after 8 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 + urban | 
    district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                             Coef.   Std. Err.        z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept                 -2.09365     0.14823   -14.12   2.704e-45   ***     -2.38418    -1.80312 
lcOne_child                1.16597     0.13489     8.64    5.44e-18   ***      0.90159     1.43036 
lcTwo_children             1.52627     0.14841    10.28    8.31e-25   ***      1.23539     1.81715 
lcThree_plus               1.52278     0.15408     9.88   4.924e-23   ***      1.22079     1.82477 
age                       -0.01818     0.00670    -2.71     0.00667   **      -0.03131    -0.00504 
urbanUrban                 0.57420     0.13647     4.21    2.58e-05   ***      0.30674     0.84167 
educLower_primary          0.24518     0.12951     1.89     0.05834   .       -0.00866     0.49901 
educUpper_primary          0.73327     0.14533     5.05   4.523e-07   ***      0.44843     1.01811 
educSecondary_and_above    1.17969     0.12839     9.19   4.003e-20   ***      0.92804     1.43134 
hinduHindu                 0.50956     0.13292     3.83   0.0001263   ***      0.24904     0.77007 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                              Coef.   Std. Err. 
var_Intercept               0.35974     0.09867 
cov_Intercept_urbanUrban   -0.25801     0.11151 
var_urbanUrban              0.34896     0.17340 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> if (!require(car)) {
+   warning("car package required to use linearHypothesis() function")
+ } else {
+   linearHypothesis(mymodel7, "RP2_cov_Intercept_urbanUrban = 0")
+ }
Linear hypothesis test

Hypothesis:
RP2_cov_Intercept_urbanUrban = 0

Model 1: restricted model
Model 2: logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 + urban | 
    district)

  Res.Df Df  Chisq Pr(>Chisq)  
1   2854                       
2   2853  1 5.3534    0.02068 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> if (!require(car)) {
+   warning("car package required to use linearHypothesis() function")
+ } else {
+   linearHypothesis(mymodel7, "RP2_var_urbanUrban = 0")
+ }
Linear hypothesis test

Hypothesis:
RP2_var_urbanUrban = 0

Model 1: restricted model
Model 2: logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 + urban | 
    district)

  Res.Df Df  Chisq Pr(>Chisq)  
1   2854                       
2   2853  1 4.0499    0.04417 *
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> if (!require(car)) {
+   warning("car package required to use linearHypothesis() function")
+ } else {
+   linearHypothesis(mymodel7, c("RP2_cov_Intercept_urbanUrban = 0", "RP2_var_urbanUrban = 0"))
+ }
Linear hypothesis test

Hypothesis:
RP2_cov_Intercept_urbanUrban = 0
RP2_var_urbanUrban = 0

Model 1: restricted model
Model 2: logit(use) ~ 1 + lc + age + urban + educ + hindu + (1 + urban | 
    district)

  Res.Df Df  Chisq Pr(>Chisq)  
1   2855                       
2   2853  2 5.4717    0.06484 .
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> 
> (mymodel8 <- runMLwiN(logit(use) ~ 1 + lc + age + urban + educ + hindu + d_lit + d_pray + (1 + urban | district),
+   D = "Binomial", estoptions = list(nonlinear = c(N = 1, M = 2), startval = list(FP.b = mymodel7@FP, FP.v = mymodel7@FP.cov,
+     RP.b = mymodel7@RP, RP.v = mymodel7@RP.cov)), data = bang))
/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 -2.09365114401487 1.16597447430765 1.52626872549308 1.52277975007742 -0.0181756030954108 0.574203474031062 0.245175315697699 0.733269026333489 1.17968787196353 0.509555790433903 0 0 '_FP_b'
JOIN 0.0219733196802088 -0.0104541888031348 0.0181959598975733 -0.0126078929973787 0.0100424237520441 0.0220256026701775 -0.0154396033053014 0.0112258857845742 0.0146092310342419 0.0237399309100986 0.000370979834478551 -0.000165444373814286 -0.000370526959534443 -0.000682628446709666 4.48870021916425e-05 -0.00767744730675773 0.000480161961342514 0.000870968622887474 0.000659233695422207 -2.36492406938521e-05 0.018622783474728 -0.00291720245689676 0.000407641708114746 0.00022877717133059 -7.54876972434869e-05 4.6508417010466e-05 -0.000705193071727363 0.0167726708681131 -0.0028148671892566 -2.73405942759159e-05 7.53379457881814e-05 -4.51502866793914e-05 7.59301974829547e-05 -0.00108779055559781 0.00324217544128519 0.0211210684327045 -0.00459151677866056 0.00186306059818434 0.00229861328230719 0.00297792844854554 3.349332574253e-05 -0.00294660060949399 0.00351102245403737 0.00388494224260141 0.0164852284249341 -0.00349614929466233 0.000407001163384805 0.000674570662564585 0.00135662192237635 -5.71220594618246e-05 -2.59130801425283e-05 2.33979892126268e-05 -0.000157509096792477 0.000341272982594661 0.0176668751018998 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 '_FP_v'
JOIN 0.359738458340248 -0.258012324937657 0.348955050974729 1 '_RP_b'
JOIN 0.00973512093334196 -0.00773262176575147 0.0124351943508086 0.00611415126486965 -0.0150050721588829 0.0300670987206531 0 0 0 0 '_RP_v'
TOLE 2
MAXI 20
NEXT
iteration 2
iteration 3
iteration 4
iteration 5

Convergence achieved
ECHO 0
Execution completed


-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
MLwiN (version: 3.11)  multilevel model (Binomial) 
          N min     mean max N_complete min_complete mean_complete max_complete
district 60   3 47.78333 173         60            3      47.78333          173
Estimation algorithm:  IGLS PQL2        Elapsed time : 0.21s 
Number of obs:  2867 (from total 2867)        The model converged after 6 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use) ~ 1 + lc + age + urban + educ + hindu + d_lit + d_pray + 
    (1 + urban | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
                             Coef.   Std. Err.       z    Pr(>|z|)         [95% Conf.   Interval] 
Intercept                 -1.72295     0.26327   -6.54   5.976e-11   ***     -2.23896    -1.20695 
lcOne_child                1.17020     0.13503    8.67   4.462e-18   ***      0.90555     1.43486 
lcTwo_children             1.53401     0.14861   10.32     5.6e-25   ***      1.24273     1.82529 
lcThree_plus               1.52828     0.15425    9.91   3.857e-23   ***      1.22595     1.83061 
age                       -0.01814     0.00670   -2.71    0.006761   **      -0.03127    -0.00501 
urbanUrban                 0.52822     0.13814    3.82   0.0001314   ***      0.25747     0.79897 
educLower_primary          0.23770     0.12995    1.83     0.06738   .       -0.01700     0.49239 
educUpper_primary          0.74232     0.14559    5.10   3.421e-07   ***      0.45697     1.02768 
educSecondary_and_above    1.19596     0.12895    9.27   1.779e-20   ***      0.94323     1.44869 
hinduHindu                 0.50955     0.13249    3.85   0.0001201   ***      0.24987     0.76922 
d_lit                      2.07495     1.70615    1.22      0.2239           -1.26905     5.41894 
d_pray                    -1.40837     0.53393   -2.64    0.008346   **      -2.45486    -0.36188 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                              Coef.   Std. Err. 
var_Intercept               0.30501     0.08818 
cov_Intercept_urbanUrban   -0.23343     0.10549 
var_urbanUrban              0.35162     0.17408 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> # 9.5 Modelling binomial data . . . . . . . . . . . . . . . . . . . . . .139
> 
> # Modelling district-level variation with district-level proportions . . 139
> 
> # Creating a district-level data set . . . . . . . . . . . . . . . . . . 140
> 
> if (!require(doBy)) {
+   warning("package doBy required to run this example")
+ } else {
+   bangshort <- summaryBy(use + cons ~ district + d_lit + d_pray, FUN = c(mean, sum), data = bang)
+   bangshort$use.sum <- NULL
+   colnames(bangshort) <- c("district", "d_lit", "d_pray", "use", "cons", "denom")
+   bangshort$use <- bangshort$use - 1
+ 
+   # Fitting the model . . . . . . . . . . . . . . . . . . . . . . . . . . .142
+ 
+   (mymodel9 <- runMLwiN(logit(use, denom) ~ 1 + d_lit + d_pray + (1 | district), D = "Binomial", data = bangshort))
+   print(mymodel9)
+ 
+   (mymodel10 <- runMLwiN(logit(use, denom) ~ 1 + d_lit + d_pray + (1 | district), D = "Binomial", estoptions = list(nonlinear = c(N = 1,
+     M = 2), startval = list(FP.b = mymodel9@FP, FP.v = mymodel9@FP.cov, RP.b = mymodel9@RP, RP.v = mymodel9@RP.cov)),
+     data = bangshort))
+ }
Loading required package: doBy
/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.11)  multilevel model (Binomial) 
          N min mean max N_complete min_complete mean_complete max_complete
district 60   1    1   1         60            1             1            1
Estimation algorithm:  IGLS MQL1        Elapsed time : 0.08s 
Number of obs:  60 (from total 60)        The model converged after 4 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use, denom) ~ 1 + d_lit + d_pray + (1 | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
               Coef.   Std. Err.       z   Pr(>|z|)         [95% Conf.   Interval] 
Intercept   -0.39797     0.23395   -1.70    0.08892   .       -0.85650     0.06056 
d_lit        3.75968     1.61790    2.32    0.02014   *        0.58865     6.93071 
d_pray      -1.19571     0.50106   -2.39    0.01702   *       -2.17777    -0.21366 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                  Coef.   Std. Err. 
var_Intercept   0.21061     0.05890 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
/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 -0.397970286407069 3.75968303276531 -1.19571286717458 '_FP_b'
JOIN 0.054731962557466 -0.164897319586594 2.61760668722093 -0.0718793179017356 -0.308613346250095 0.251058756151548 '_FP_v'
JOIN 0.210612498304647 1 '_RP_b'
JOIN 0.00346930817161946 0 0 '_RP_v'
TOLE 2
MAXI 20
NEXT
iteration 2
iteration 3
iteration 4

Convergence achieved
ECHO 0
Execution completed


-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
MLwiN (version: 3.11)  multilevel model (Binomial) 
          N min mean max N_complete min_complete mean_complete max_complete
district 60   1    1   1         60            1             1            1
Estimation algorithm:  IGLS PQL2        Elapsed time : 0.11s 
Number of obs:  60 (from total 60)        The model converged after 5 iterations.
Log likelihood:      NA 
Deviance statistic:  NA 
--------------------------------------------------------------------------------------------------- 
The model formula:
logit(use, denom) ~ 1 + d_lit + d_pray + (1 | district)
Level 2: district     Level 1: l1id      
--------------------------------------------------------------------------------------------------- 
The fixed part estimates:  
               Coef.   Std. Err.       z   Pr(>|z|)         [95% Conf.   Interval] 
Intercept   -0.42665     0.24072   -1.77    0.07633   .       -0.89845     0.04515 
d_lit        3.99760     1.68806    2.37    0.01788   *        0.68906     7.30614 
d_pray      -1.25083     0.52205   -2.40    0.01658   *       -2.27402    -0.22763 
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1  
--------------------------------------------------------------------------------------------------- 
The random part estimates at the district level: 
                  Coef.   Std. Err. 
var_Intercept   0.22523     0.06234 
--------------------------------------------------------------------------------------------------- 
The random part estimates at the l1id level: 
                Coef.   Std. Err. 
var_bcons_1   1.00000     0.00000 
-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*-*- 
> 
> # Addendum: changing contrasts back to pre-existing . . . . . . . . . . . NA
> 
> # Following re-specification of contrast settings towards the start of this
> # script, change contrasts back to pre-existing:
> options(contrasts = my_contrasts)
> 
> ############################################################################
> 
> proc.time()
   user  system elapsed 
   4.45    0.65    8.67