library(foreign)
library(lme4)
## Loading required package: Matrix
library(rockchalk)

setwd("/Users/yahyaalshehri/Desktop/PhD Courses/Multilevel Modeling/mlmus3 2")
dta <- read.dta("wheeze.dta")

rockchalk::summarize(dta)
## Numeric variables
##              id        var     smoking      y        age  
## min           1         2         0         0        -2   
## med         269         3.50      0         0        -0.50
## max         537         5         1         1         1   
## mean        269         3.50      0.35      0.15     -0.50
## sd          155.05      1.12      0.48      0.36      1.12
## skewness      0         0         0.64      1.94      0   
## kurtosis     -1.20     -1.36     -1.60      1.76     -1.36
## nobs       2148      2148      2148      2148      2148   
## nmissing      0         0         0         0         0
head(dta)
##   id var smoking y age
## 1  1   2       0 0  -2
## 2  1   3       0 0  -1
## 3  1   4       0 0   0
## 4  1   5       0 0   1
## 5  2   2       0 0  -2
## 6  2   3       0 0  -1
### the first model:

dicotomous.logistic.regression <- glm(y ~ age + smoking, family = binomial, data = dta)

summary(dicotomous.logistic.regression)
## 
## Call:
## glm(formula = y ~ age + smoking, family = binomial, data = dta)
## 
## Deviance Residuals: 
##     Min       1Q   Median       3Q      Max  
## -0.6685  -0.5909  -0.5608  -0.5045   2.0613  
## 
## Coefficients:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -1.88373    0.08384 -22.467   <2e-16 ***
## age         -0.11341    0.05408  -2.097   0.0360 *  
## smoking      0.27214    0.12347   2.204   0.0275 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 1829.1  on 2147  degrees of freedom
## Residual deviance: 1819.9  on 2145  degrees of freedom
## AIC: 1825.9
## 
## Number of Fisher Scoring iterations: 4
##### the second model
Variance.Componant.dicotomous.logistic.regression <- glmer(y ~ 1+ (1|id), family = binomial, data = dta)

summary(Variance.Componant.dicotomous.logistic.regression)
## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: binomial  ( logit )
## Formula: y ~ 1 + (1 | id)
##    Data: dta
## 
##      AIC      BIC   logLik deviance df.resid 
##   1602.8   1614.2   -799.4   1598.8     2146 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.1991 -0.1604 -0.1604 -0.1604  2.1547 
## 
## Random effects:
##  Groups Name        Variance Std.Dev.
##  id     (Intercept) 5.43     2.33    
## Number of obs: 2148, groups:  id, 537
## 
## Fixed effects:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept)  -3.1161     0.2412  -12.92   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##### the third model

Rando.intercept.dicotomous.logistic.regression <- glmer(y ~ age + smoking + (1|id), family = binomial, data = dta)

summary(Rando.intercept.dicotomous.logistic.regression)
## Generalized linear mixed model fit by maximum likelihood (Laplace
##   Approximation) [glmerMod]
##  Family: binomial  ( logit )
## Formula: y ~ age + smoking + (1 | id)
##    Data: dta
## 
##      AIC      BIC   logLik deviance df.resid 
##   1597.9   1620.6   -794.9   1589.9     2144 
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.4027 -0.1802 -0.1577 -0.1321  2.5176 
## 
## Random effects:
##  Groups Name        Variance Std.Dev.
##  id     (Intercept) 5.49     2.343   
## Number of obs: 2148, groups:  id, 537
## 
## Fixed effects:
##             Estimate Std. Error z value Pr(>|z|)    
## (Intercept) -3.37395    0.27498 -12.270   <2e-16 ***
## age         -0.17676    0.06797  -2.601   0.0093 ** 
## smoking      0.41478    0.28704   1.445   0.1485    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##         (Intr) age   
## age      0.227       
## smoking -0.419 -0.010