DEM 7283 - Homewk10 - Multi Level Regression Analysis

#Load data for analysis#sub-setting and Re-codding variables for analysis purposes

brfss <- readRDS("brfss_177.rds")
# Cleaning the  variable names for space, underscore & Uppercase Characters
renam<-names(brfss)
newnames<-tolower(gsub(pattern = "_",replacement =  "",x =  renam))
names(brfss)<-newnames
#get state names and abb.
sta <- read_csv("Stateinfo.csv")

brfss2 <- merge(brfss, sta, by.x = "state")

homewk3 <- brfss2 %>% 
  dplyr::select(state,Stateabbr,statenam,llcpwt,dispcode,seqno,psu,bmi5cat,racegr3,educa,ststr,smoke100,exerany2,sex,ageg,menthlth,smokday2,sleptim1) %>%  
  mutate( smokers=ifelse(smoke100==1 & smokday2%in%c(1,2),"Current",
                 ifelse(smoke100==1 & smokday2==3, "Former", 
                ifelse(smoke100==2, "Never", NA)))) %>% 
  
mutate(bmi = car::Recode(bmi5cat,recodes="1='Underweight';2='Normal';3='Overweight'; 4='Obese';else=NA",as.factor=T ),
       educa = car::Recode(educa,recodes="1:2='<HS';3:4='HS';5:6='>HS';else=NA", as.factor=T ),
       racegr3 = car::Recode(racegr3,recodes="1='NH-White';2='NH-Black';3:4='Other';5='Hispanic';else=NA",as.factor=T ), 
       ageg = car::Recode(ageg,recodes="1='18-24';2='25-30';3='35-44';4='45-54';5='55-64';else='65+'" ), 
       Physicala = car::Recode(exerany2,recodes="='Yes';2='No';else=NA",as.factor=T ), sex = car::Recode(sex,recodes="1='Male';2='Female';else=NA",as.factor=T ),
       sleepdurat=car::Recode(sleptim1, recodes="77=NA ;99= NA"),
        mentald=ifelse(menthlth<31, "Yes",
                       ifelse(menthlth==88, "No",NA))) 

#table(is.na(homewk3$sleepdurat), homewk3$state)
 
#Data containing only variables used in the analysis
 
Final_data <- homewk3 %>% 
 dplyr:: select(state,Stateabbr,statenam,llcpwt,dispcode,seqno,psu,smokers,bmi,educa,ageg,Physicala, sleepdurat,sex,racegr3,mentald,ststr) %>% 
 filter(complete.cases(.))

#  Survey Design
options(survey.lonely.psu = "adjust")

des<-svydesign(ids=~1, strata=~ststr, weights=~llcpwt, data = Final_data )

Data Checking

head(data.frame(name=table(Final_data$statenam), 
                id=unique(Final_data$state)))

#How many total states are in the data?
length(table(Final_data$state))

## [1] 7

Testing Group Differences.

fit.an<-lm( sleepdurat~as.factor(state),
           Final_data)
anova(fit.an)

Scaling the Outcome variable

Final_data$sleepduratz<-as.numeric(scale(Final_data$sleepdurat,
                              center=T,
                              scale=T))

Higher Level Predictor model

library(lme4)
library(lmerTest)
library(arm)
library(sjstats)
library(sjPlot)
fit<-lmer(sleepduratz~(1|state),
           data=Final_data)
arm::display(fit,
             detail=T)

## lmer(formula = sleepduratz ~ (1 | state), data = Final_data)
##             coef.est coef.se t value
## (Intercept) -0.01     0.02   -0.37  
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  state    (Intercept) 0.03    
##  Residual             1.00    
## ---
## number of obs: 11432, groups: state, 7
## AIC = 32451.4, DIC = 32432.4
## deviance = 32438.9

performance::icc(fit)

## # Intraclass Correlation Coefficient
## 
##      Adjusted ICC: 0.001
##   Conditional ICC: 0.001

The fixed effect term is less than zero (-0.01) and since the T-statistics is very small; less than 2, we could say that it is not significantly different from zero (average sleep duration between states is not significantly different). The variation at the state level is 0.03 while the individual level variation is fixed at 1. The intraclass correlation Coefficient (ICC) is very close to zero (0.001); individuals in each states are very much different considering their average sleep duration.

Individual Level Predictor

fit2<-lmer(sleepduratz~bmi+ racegr3+ ageg+ educa+sex+mentald+(1|state),
           data=Final_data,
           na.action=na.omit)
arm::display(fit2, detail=T)

## lmer(formula = sleepduratz ~ bmi + racegr3 + ageg + educa + sex + 
##     mentald + (1 | state), data = Final_data, na.action = na.omit)
##                 coef.est coef.se t value
## (Intercept)      0.17     0.08    2.06  
## bmiObese        -0.06     0.02   -2.45  
## bmiOverweight   -0.04     0.02   -1.48  
## bmiUnderweight   0.02     0.07    0.34  
## racegr3NH-Black -0.17     0.05   -3.16  
## racegr3NH-White -0.07     0.04   -1.92  
## racegr3Other    -0.10     0.05   -1.91  
## ageg25-30       -0.14     0.06   -2.35  
## ageg35-44       -0.14     0.06   -2.43  
## ageg45-54       -0.15     0.05   -2.82  
## ageg55-64       -0.03     0.05   -0.54  
## ageg65+          0.18     0.05    3.44  
## educa>HS        -0.02     0.06   -0.35  
## educaHS         -0.01     0.06   -0.14  
## sexMale         -0.05     0.02   -2.54  
## mentaldYes      -0.10     0.02   -5.22  
## 
## Error terms:
##  Groups   Name        Std.Dev.
##  state    (Intercept) 0.05    
##  Residual             0.99    
## ---
## number of obs: 11432, groups: state, 7
## AIC = 32250.3, DIC = 32049.3
## deviance = 32131.8

performance::icc(fit2)

## # Intraclass Correlation Coefficient
## 
##      Adjusted ICC: 0.003
##   Conditional ICC: 0.003

Compared to the Normal BMI categories,the Overweight and Obese categories on average have a lower sleep duration,holding all else constant. However, compared to the normal BMI category the underweight category has a higher average sleep duration, holding all else constatnt.
The age group 65+ has a higher average sleep duration compared to age group 18-24, holding all else constant.
The fixed effect term is greater than zero (0.17) and since the T-statistics is larger than 2, we could say that it is significantly different from zero (average sleep duration between states is significantly different hold all else constant). The variation at the state level is 0.05 while the individual level variation is 0.99. The intraclass correlation Coefficient (ICC) is very close to zero (0.003); individuals in each states are very much different considering their average sleep duration holding all other predictors constant.

Likelihood Ratio Test

anova(fit, fit2)

Yes, the inclusion of the individual level variables increased the overall model fit