Data: Each year, beginning at age 14, 82 teenagers completed a 4-item questionnaire assessing their alcohol consumption during the previous year. Using a 8-point scale (0 = “not al all, 8 =”every day") teenagers described the frequency with which they
drank beer or wine, drank hard liquor, had five or more drinks in a row, and got drunk. Two potential predictors of alcohol use are whether the teenager is a child of an alcoholic parent; and alcohol use among the teenager’s peers. The teenager used a 6-point scale to estimate the proportion of their friends who drank alcohol occasionally (item 1) or regularly (item 2). This was obtained during the first wave of data collection.
Column 1: Teenager ID
Column 2: Whether the teenager is a child of a alcohlic parent
Column 3: Sex (male = 1, female = 0)
Column 4: Number of year since age 14
Column 5: Alcohol use of the teenager (sqrt-root of mean of 6 items) Column 6: Alcohol use among the teenager’s peers (sqrt-root of mean of 2 items) Column 7: Alcoholic parenet variable centered
Column 8: Peer variable centered
# Loading package
pacman::p_load(tidyverse, afex, segmented, foreign, haven, lme4, lmerTest, ggplot2, Rmisc, car)
# data input
dta3 <- read.table("C:/Users/HANK/Desktop/HOMEWORK/alcohol_use.txt", header = T)
head(dta3)
## sid coa sex age14 alcuse peer cpeer ccoa
## 1 1 1 0 0 1.73205 1.26491 0.24691 0.549
## 2 1 1 0 1 2.00000 1.26491 0.24691 0.549
## 3 1 1 0 2 2.00000 1.26491 0.24691 0.549
## 4 2 1 1 0 0.00000 0.89443 -0.12357 0.549
## 5 2 1 1 1 0.00000 0.89443 -0.12357 0.549
## 6 2 1 1 2 1.00000 0.89443 -0.12357 0.549
dta3 <- dta3 %>%
mutate(sid = factor(sid),
coa = factor(coa),
sex = factor(sex),
age14 = factor(age14))
str(dta3)
## 'data.frame': 246 obs. of 8 variables:
## $ sid : Factor w/ 82 levels "1","2","3","4",..: 1 1 1 2 2 2 3 3 3 4 ...
## $ coa : Factor w/ 2 levels "0","1": 2 2 2 2 2 2 2 2 2 2 ...
## $ sex : Factor w/ 2 levels "0","1": 1 1 1 2 2 2 2 2 2 2 ...
## $ age14 : Factor w/ 3 levels "0","1","2": 1 2 3 1 2 3 1 2 3 1 ...
## $ alcuse: num 1.73 2 2 0 0 ...
## $ peer : num 1.265 1.265 1.265 0.894 0.894 ...
## $ cpeer : num 0.247 0.247 0.247 -0.124 -0.124 ...
## $ ccoa : num 0.549 0.549 0.549 0.549 0.549 0.549 0.549 0.549 0.549 0.549 ...
ggplot(dta3, aes(age14, alcuse, color = sex))+
geom_point()+
facet_wrap( ~ coa)+
labs(x = "Age (after 14)", y = "Alcohol use") +
theme(legend.position = "bottom")
summary(m1 <- lmer(alcuse ~ coa + peer*age14 + (1 | sid), data = dta3))
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: alcuse ~ coa + peer * age14 + (1 | sid)
## Data: dta3
##
## REML criterion at convergence: 622.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.46909 -0.65888 0.02567 0.51863 2.56829
##
## Random effects:
## Groups Name Variance Std.Dev.
## sid (Intercept) 0.3373 0.5808
## Residual 0.4834 0.6952
## Number of obs: 246, groups: sid, 82
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.284160 0.180572 166.151585 -1.574 0.117468
## coa1 0.565134 0.158782 79.000000 3.559 0.000633 ***
## peer 0.648244 0.139127 175.437485 4.659 6.25e-06 ***
## age141 0.341846 0.187127 160.000000 1.827 0.069592 .
## age142 0.849373 0.187127 160.000000 4.539 1.11e-05 ***
## peer:age141 -0.008533 0.149775 160.000000 -0.057 0.954637
## peer:age142 -0.302754 0.149775 160.000000 -2.021 0.044906 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) coa1 peer age141 age142 pr:141
## coa1 -0.297
## peer -0.734 -0.127
## age141 -0.518 0.000 0.438
## age142 -0.518 0.000 0.438 0.500
## peer:age141 0.422 0.000 -0.538 -0.814 -0.407
## peer:age142 0.422 0.000 -0.538 -0.407 -0.814 0.500
summary(m2 <- lmer(alcuse ~ coa + peer*age14 + (1 | age14) + (1 | sid), data = dta3))
## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, : Model is nearly unidentifiable: large eigenvalue ratio
## - Rescale variables?
## Warning in as_lmerModLT(model, devfun): Model may not have converged with 1
## eigenvalue close to zero: 1.9e-09
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: alcuse ~ coa + peer * age14 + (1 | age14) + (1 | sid)
## Data: dta3
##
## REML criterion at convergence: 622.3
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.46909 -0.65888 0.02567 0.51863 2.56829
##
## Random effects:
## Groups Name Variance Std.Dev.
## sid (Intercept) 0.33733 0.5808
## age14 (Intercept) 0.09478 0.3079
## Residual 0.48336 0.6952
## Number of obs: 246, groups: sid, 82; age14, 3
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.284160 0.356914 212.474160 -0.796 0.426830
## coa1 0.565134 0.158782 79.000062 3.559 0.000633 ***
## peer 0.648244 0.139127 175.437652 4.659 6.25e-06 ***
## age141 0.341846 0.473898 159.999963 0.721 0.471747
## age142 0.849373 0.473898 159.999963 1.792 0.074972 .
## peer:age141 -0.008533 0.149775 159.999964 -0.057 0.954637
## peer:age142 -0.302754 0.149775 159.999964 -2.021 0.044906 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) coa1 peer age141 age142 pr:141
## coa1 -0.150
## peer -0.371 -0.127
## age141 -0.664 0.000 0.173
## age142 -0.664 0.000 0.173 0.500
## peer:age141 0.214 0.000 -0.538 -0.322 -0.161
## peer:age142 0.214 0.000 -0.538 -0.161 -0.322 0.500
## convergence code: 0
## Model is nearly unidentifiable: large eigenvalue ratio
## - Rescale variables?