1A: Loading data, setting up

library(prcr)
library(tidyverse)

ctk <- read_csv("~/downloads/ctk.csv")

esm_df <- haven::read_sav('/Users/joshuarosenberg/Dropbox/1_Research/SciMo/JRST\ MEP\ Paper\ Final\ Docs/esm_data.sav')

to_remove <- esm_df$stud_ID == "AXR050594" & esm_df$month == 10 & esm_df$day == 20 & esm_df$year == 8 & esm_df$signal == 2 # removing 10/20/8, ID = "AXR050594", second beeper
esm_df <- esm_df[!to_remove, ]

df <- esm_df %>%
    mutate(behavioral_scale = ifelse(!is.na(conc) & !is.na(hardwk), (conc + hardwk) / 2, 
                                     ifelse(is.na(conc), hardwk, conc)),
           cognitive_scale = ifelse(!is.na(imp_y) & !is.na(imp_fut), (imp_y + imp_fut) / 2, 
                                    ifelse(is.na(imp_y), imp_fut, imp_y)),
           affective_scale = ifelse(!is.na(enjoy) & !is.na(interest), (enjoy + interest) / 2, 
                                    ifelse(is.na(enjoy), interest, enjoy)))

df <- df[ctk$ctk, ]

p6 <- create_profiles(df,
                      behavioral_scale,
                      cognitive_scale,
                      affective_scale,
                      n_profiles = 6,
                      to_center = T,
                      to_scale = T)
knitr::opts_chunk$set(cache=TRUE)

cv_output <- cross_validate(df,
                            behavioral_scale,
                            cognitive_scale,
                            affective_scale,
                            n_profiles = 6,
                            to_center = T,
                            to_scale = T,
                            k = 30)

knitr::kable(cv_output)
k_iteration kappa percentage_agree
1 0.44 0.53
2 0.75 0.79
3 0.71 0.76
4 0.90 0.92
5 0.57 0.64
6 0.55 0.62
7 0.63 0.69
8 0.74 0.78
9 0.49 0.57
10 0.71 0.76
11 0.92 0.93
12 0.56 0.63
13 0.85 0.88
14 0.65 0.71
15 0.54 0.62
16 0.66 0.72
17 0.70 0.75
18 0.61 0.67
19 0.58 0.65
20 0.68 0.73
21 0.59 0.66
22 0.89 0.91
23 0.67 0.72
24 0.62 0.69
25 0.68 0.74
26 0.92 0.93
27 0.63 0.69
28 0.68 0.73
29 0.79 0.83
30 0.60 0.67

1B: Processing data

library(tidyverse)
library(lme4)

to_m <- p6$.data
dcd <- p6$data_with_dummy_codes
dff <- bind_cols(to_m, dcd)
dff$response_date <- paste0(dff$month, "-", dff$day)
dff$beep_ID <- paste0(dff$response_date, dff$signal, dff$teacher_ID)
df_ss <- select(dff, contains("cluster"), -cluster, teacher_ID, stud_ID, beep_ID)

2A: Null / Variance Components Models

knitr::opts_chunk$set(cache=TRUE)

m1 <- glmer(cluster_1 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"),
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

m2 <- glmer(cluster_2 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

m3 <- glmer(cluster_3 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

m4 <- glmer(cluster_4 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

m5 <- glmer(cluster_5 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

m6 <- glmer(cluster_6 ~ 1 +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_ss)

xm1 <- sjstats::icc(m1)
xm2 <- sjstats::icc(m2)
xm3 <- sjstats::icc(m3)
xm4 <- sjstats::icc(m4)
xm5 <- sjstats::icc(m5)
xm6 <- sjstats::icc(m6)

2B: ICCs for each null model with plot of profiles

plot(p6)

df_with_ICCs <- data.frame(profile_1 = xm1, 
                           profile_2 = xm2, 
                           profile_3 = xm3,
                           profile_4 = xm4,
                           profile_5 = xm5, 
                           profile_6 = xm6)

row.names(df_with_ICCs) <- c("Momentary RE", "Student RE", "Teacher RE")

knitr::kable(df_with_ICCs, digits = 4)
profile_1 profile_2 profile_3 profile_4 profile_5 profile_6
Momentary RE 0.0000000 0.12122591 0.08742863 0.006808985 0.05680992 0.04423836
Student RE 0.3951807 0.27105319 0.46877232 0.303127910 0.40472867 0.17004291
Teacher RE 0.0000000 0.04996376 0.00000000 0.005911475 0.01024088 0.01520584

3A: Preparation for analysis with MLMs

student_df <- haven::read_sav("~/documents/myscimo/scimo_student_survey.sav")
student_df <- rename(student_df, stud_ID = Stud_ID)
student_df <- select(student_df, -teacher_ID)

teacher_df <- haven::read_sav("~/documents/myscimo/scimo_teacher_survey.sav")

student_df$stud_ID <- as.vector(student_df$stud_ID)
df_ss$stud_ID <- as.vector(df_ss$stud_ID)

df_joined <- left_join(df_ss, student_df, by = "stud_ID")

# df_joined <- rename(df_joined, teacher_ID = teacher_ID1)

df_joined$teacher_ID <- as.character(df_joined$teacher_ID)

teacher_df$teacher_ID <- as.character(teacher_df$teacher_ID)

df_joined <- left_join(df_joined, teacher_df, by = "teacher_ID")

df$teacher_ID <- as.character(df$teacher_ID)

df_j <- bind_cols(df_joined, df)

3B: Constructing variables

# IPs

df_j$instructional_practice <- as.factor(df_j$instructional_practice)

df_j$instructional_practice <- car::recode(df_j$instructional_practice,
                                           "'1' = 'Lecture';
                                              c('2', '3') = 'Individual Work';
                                              c('4', '5') = 'Group Work';
                                              c('6', '8', '7') = 'Quiz and Test';
                                              c('9') = 'Discussion';
                                              c('11', '10') = 'Presentation';
                                              c('12', '13') = 'Video';
                                              c('15', '14', '16') = 'Laboratory';
                                              c('17') = 'Non-instructional';
                                              c('18') = NA")

df_j$instructional_practice <- as.character(df_j$instructional_practice)

df_j$instructional_practice <- ifelse((df_j$instructional_practice == "Discussion" | df_j$instructional_practice == "Non-instructional" |
                                           df_j$instructional_practice == "Presentation" | df_j$instructional_practice == "Video" | df_j$instructional_practice == "Group Work"), "Other", df_j$instructional_practice)

# For classroom activity: use df_j$instructional_practice

# Choices

df_j$ch_who[is.na(df_j$ch_who)] <- 0
df_j$ch_howdo[is.na(df_j$ch_howdo)] <- 0
df_j$ch_mat[is.na(df_j$ch_mat)] <- 0
df_j$ch_time[is.na(df_j$ch_time)] <- 0
df_j$ch_doing[is.na(df_j$ch_doing)] <- 0
df_j$ch_topic[is.na(df_j$ch_topic)] <- 0
df_j$ch_defin[is.na(df_j$ch_defin)] <- 0
df_j$ch_other[is.na(df_j$ch_other)] <- 0
df_j$ch_none[is.na(df_j$ch_none)] <- 0

df_j$ch_none <- ifelse(((df_j$ch_who == 0 & df_j$ch_howdo == 0 & df_j$ch_mat == 0 & df_j$ch_time == 0 & df_j$ch_doing == 0 &
                             df_j$ch_defin == 0 & df_j$ch_topic == 0 & df_j$ch_other == 0) & df_j$ch_none == 0), 1, df_j$ch_none)

df_j$ch_none <- ifelse(((df_j$ch_who == 1 | df_j$ch_howdo == 1 | df_j$ch_mat == 1 | df_j$ch_time == 1 | df_j$ch_doing |
                             df_j$ch_defin == 1 | df_j$ch_topic == 1 | df_j$ch_other == 1) & df_j$ch_none == 1), 0, df_j$ch_none)

df_j$ch_doing_defining_topic <- ifelse(df_j$ch_defin == 1 | df_j$ch_topic == 1 | df_j$ch_doing, 1, 0)

df_j$ch_defining_topic <- ifelse(df_j$ch_defin == 1 | df_j$ch_topic == 1, 1, 0)
df_j$framing <- ifelse(df_j$ch_defin == 1 | df_j$ch_topic == 1 | df_j$ch_doing, 1, 0)


df_j$any_choice <- ifelse(df_j$framing == 1 | df_j$ch_mat == 1 | df_j$ch_howdo == 1 | df_j$ch_time == 1 |
                              df_j$ch_who == 1 | df_j$ch_other == 1, 1, 0) # this is the same as ch_none above

# For any vs. none use: df_j$any_choice
# For Framing use: df_j$framing
# For who use: df_j$ch_who

# need race, gender, grade, interest at T1, competence at T1, prior grades

df_j$race <- as.factor(as.vector(df_j$race))

df_j$race <- forcats::fct_recode(df_j$race,
                                 Black = "3",
                                 White  = "4",
                                 Hispanic = "2",
                                 Other = "1",
                                 Other = "5",
                                 Other = "6")

df_j$race[df_j$race == "NaN"] <- NA

df_j$gender <- as.factor(df_j$gender)

df_j$gender <- forcats::fct_recode(df_j$gender,
                                   Male = "1",
                                   Female = "2")

df_j$grade <- as.factor(df_j$grade)

df_j$grade <- forcats::fct_recode(df_j$grade,
                                  ninth = "1",
                                  tenth = "2",
                                  eleventh = "3",
                                  twelth = "4")

df_j <- filter(df_j, gengrades != 9 & gengrades != 98)
df_j <- filter(df_j, instructional_practice != "NA")

df_j$gengrades <- as.integer(as.vector(df_j$gengrades))

df_j$interestsci <- as.vector(df_j$interestsci)
df_j$expectwell <- as.vector(df_j$expectwell)
# df_j$CumGPA_2010 <- as.numeric(as.vector(df_j$CumGPA_2010))
df_j$instructional_practice[is.na(df_j$instructional_practice)] <- NA
df_j$instructional_practice <- as.factor(df_j$instructional_practice)
df_j$instructional_practice <- forcats::fct_relevel(df_j$instructional_practice,
                                                    "Other")
df_j$race <- forcats::fct_relevel(df_j$race, "Black")

4A: Analysis with MLMs

knitr::opts_chunk$set(cache=TRUE)

m1 <- glmer(cluster_1 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"),
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m1_out <- broom::tidy(m1)
sjPlot::sjt.glmer(m1, show.re.var = F, show.ci = F, show.se = T)
    cluster_1
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.18 0.11 .004
gender (Female)   0.81 0.18 .346
race (Other)   2.52 1.23 .059
race (Hispanic)   1.55 0.57 .236
race (White)   0.59 0.23 .178
grade (tenth)   1.39 0.43 .281
grade (eleventh)   1.02 0.29 .948
grade (twelth)   7.23 4.95 .004
gengrades   0.94 0.07 .473
interestsci   1.14 0.16 .362
expectwell   0.65 0.11 .011
instructional_practice (Individual Work)   0.92 0.16 .634
instructional_practice (Laboratory)   0.84 0.14 .275
instructional_practice (Lecture)   1.22 0.20 .242
instructional_practice (NaN)   0.87 0.48 .793
instructional_practice (Quiz and Test)   1.14 0.19 .429
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.000
ICCstud_ID   0.331
ICCteacher_ID   0.001
Observations   3613
Deviance   2519.743
m2 <- glmer(cluster_2 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m2_out <- broom::tidy(m2)
sjPlot::sjt.glmer(m2, show.re.var = F, show.ci = F, show.se = T)
    cluster_2
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.07 0.04 <.001
gender (Female)   1.34 0.30 .189
race (Other)   0.53 0.29 .251
race (Hispanic)   1.07 0.40 .866
race (White)   0.85 0.33 .674
grade (tenth)   1.07 0.33 .833
grade (eleventh)   0.77 0.23 .388
grade (twelth)   0.57 0.44 .466
gengrades   1.13 0.09 .136
interestsci   0.76 0.11 .053
expectwell   0.94 0.16 .721
instructional_practice (Individual Work)   1.65 0.37 .025
instructional_practice (Laboratory)   0.57 0.14 .023
instructional_practice (Lecture)   1.44 0.33 .115
instructional_practice (NaN)   0.82 0.72 .825
instructional_practice (Quiz and Test)   5.27 1.14 <.001
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.071
ICCstud_ID   0.291
ICCteacher_ID   0.000
Observations   3613
Deviance   1842.132
m3 <- glmer(cluster_3 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m3_out <- broom::tidy(m3)
sjPlot::sjt.glmer(m3, show.re.var = F, show.ci = F, show.se = T)
    cluster_3
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.49 0.40 .378
gender (Female)   0.48 0.15 .016
race (Other)   1.60 1.07 .486
race (Hispanic)   0.67 0.34 .437
race (White)   0.47 0.25 .156
grade (tenth)   0.74 0.33 .499
grade (eleventh)   1.67 0.66 .196
grade (twelth)   0.25 0.29 .235
gengrades   0.67 0.08 <.001
interestsci   0.93 0.18 .709
expectwell   0.68 0.17 .115
instructional_practice (Individual Work)   0.56 0.16 .038
instructional_practice (Laboratory)   1.47 0.36 .119
instructional_practice (Lecture)   0.94 0.24 .818
instructional_practice (NaN)   0.43 0.42 .386
instructional_practice (Quiz and Test)   0.87 0.24 .619
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.082
ICCstud_ID   0.421
ICCteacher_ID   0.000
Observations   3613
Deviance   1433.118
m4 <- glmer(cluster_4 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m4_out <- broom::tidy(m4)
sjPlot::sjt.glmer(m4, show.re.var = F, show.ci = F, show.se = T)
    cluster_4
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.07 0.04 <.001
gender (Female)   1.05 0.21 .830
race (Other)   0.34 0.17 .030
race (Hispanic)   0.63 0.22 .175
race (White)   1.36 0.47 .374
grade (tenth)   1.22 0.34 .483
grade (eleventh)   1.15 0.29 .583
grade (twelth)   0.38 0.28 .192
gengrades   1.00 0.07 .995
interestsci   1.29 0.17 .055
expectwell   1.18 0.18 .252
instructional_practice (Individual Work)   1.31 0.21 .089
instructional_practice (Laboratory)   0.80 0.13 .169
instructional_practice (Lecture)   1.24 0.20 .184
instructional_practice (NaN)   0.82 0.38 .674
instructional_practice (Quiz and Test)   0.92 0.16 .654
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.006
ICCstud_ID   0.291
ICCteacher_ID   0.000
Observations   3613
Deviance   2640.616
m5 <- glmer(cluster_5 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice +
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m5_out <- broom::tidy(m5)
sjPlot::sjt.glmer(m5, show.re.var = F, show.ci = F, show.se = T)
    cluster_5
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.05 0.04 <.001
gender (Female)   1.74 0.43 .025
race (Other)   1.04 0.59 .946
race (Hispanic)   1.06 0.44 .888
race (White)   1.38 0.58 .450
grade (tenth)   0.45 0.16 .030
grade (eleventh)   0.53 0.18 .057
grade (twelth)   0.15 0.13 .033
gengrades   1.32 0.11 .001
interestsci   0.94 0.15 .701
expectwell   1.27 0.22 .173
instructional_practice (Individual Work)   1.34 0.26 .126
instructional_practice (Laboratory)   0.96 0.18 .833
instructional_practice (Lecture)   0.90 0.18 .593
instructional_practice (NaN)   2.25 1.28 .153
instructional_practice (Quiz and Test)   0.52 0.11 .002
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.055
ICCstud_ID   0.372
ICCteacher_ID   0.006
Observations   3613
Deviance   2491.594
m6 <- glmer(cluster_6 ~ 1 +
                gender + race + grade + gengrades + interestsci + expectwell + instructional_practice + 
                (1|teacher_ID) + (1|stud_ID) + (1|beep_ID),
            family = binomial(link = "logit"), 
            control=lme4::glmerControl(optimizer="bobyqa",
                                       optCtrl=list(maxfun=2e5)),
            data = df_j)

m6_out <- broom::tidy(m6)
sjPlot::sjt.glmer(m6, show.re.var = F, show.ci = F, show.se = T)
    cluster_6
    Odds Ratio std. Error p
Fixed Parts
(Intercept)   0.17 0.07 <.001
gender (Female)   0.99 0.15 .968
race (Other)   1.26 0.46 .526
race (Hispanic)   1.32 0.35 .294
race (White)   1.39 0.38 .225
grade (tenth)   1.46 0.32 .086
grade (eleventh)   1.46 0.30 .065
grade (twelth)   1.52 0.75 .399
gengrades   0.90 0.05 .051
interestsci   0.95 0.09 .604
expectwell   1.17 0.13 .155
instructional_practice (Individual Work)   0.59 0.10 .001
instructional_practice (Laboratory)   1.43 0.22 .021
instructional_practice (Lecture)   0.66 0.11 .013
instructional_practice (NaN)   0.86 0.40 .748
instructional_practice (Quiz and Test)   0.35 0.07 <.001
Random Parts
Nbeep_ID   244
Nstud_ID   203
Nteacher_ID   13
ICCbeep_ID   0.029
ICCstud_ID   0.166
ICCteacher_ID   0.000
Observations   3613
Deviance   2846.657

4B: ICCs

xm1 <- sjstats::icc(m1)
xm2 <- sjstats::icc(m2)
xm3 <- sjstats::icc(m3)
xm4 <- sjstats::icc(m4)
xm5 <- sjstats::icc(m5)
xm6 <- sjstats::icc(m6)

df1_with_ICCs <- data.frame(profile_1 = xm1, 
                            profile_2 = xm2, 
                            profile_3 = xm3,
                            profile_4 = xm4,
                            profile_5 = xm5, 
                            profile_6 = xm6)

row.names(df_with_ICCs) <- c("Momentary RE", "Student RE", "Teacher RE")

df_with_ICCs
##              profile_1  profile_2  profile_3   profile_4  profile_5
## Momentary RE 0.0000000 0.12122591 0.08742863 0.006808985 0.05680992
## Student RE   0.3951807 0.27105319 0.46877232 0.303127910 0.40472867
## Teacher RE   0.0000000 0.04996376 0.00000000 0.005911475 0.01024088
##               profile_6
## Momentary RE 0.04423836
## Student RE   0.17004291
## Teacher RE   0.01520584