rawdata = read_sav("data/AccessDB_T0T4Merge_7.12.19.1.sav")
names(rawdata)[1] = "SubjectID"
D <- rawdata
###check the number of participants at each timepoint
D <- D %>%
filter(T0SY1Y2C_T1C == 1 |T0SY1Y2C_T2C == 1 | T0SY1Y2C_T3C == 1 | T0SY1Y2C_T4C == 1) ### participants attend at least 1 timepionts
a <- nrow(D)
focal.var = c("T1AWMADMSS_N","T2AWMADMSS_N","T3AWMADMSS_N","T4AWMADMSS_N",
"T1AWMADMRS_N","T2AWMADMRS_N","T3AWMADMRS_N","T4AWMADMRS_N",
"T1GMPA", "T2GMPA", "T3GMPA", "T4GMPA",
"T1GPPA", "T2GPPA", "T3GPPA", "T4GPPA",
"T1PRWMPAE", "T2PRWMPAE", "T3PRWMPAE", "T4PRWMPAE",
"T1PMPA", "T2PMPA", "T3PMPA", "T4PMPA",
"T1NL10W_PAE","T2NL10W_PAE","T3NL10W_PAE","T4NL10W_PAE",
"T1NL100W_PAE", "T1NL1000W_PAE", "T2NL100W_PAE", "T2NL1000W_PAE",
"T3NL100W_PAE", "T3NL1000W_PAE", "T4NL100W_PAE", "T4NL1000W_PAE",
"T1PKC_PA", "T2PKC_PA", "T3PKC_PA", "T4PKC_PA",
"T1WJCW", "T2WJCW", "T3WJCW", "T4WJCW",
"T1WJCRS", "T2WJCRS", "T3WJCRS", "T4WJCRS",
"T1ASWmac", "T2ASWmac", "T3ASWmac", "T4ASWmac",
"T1LWIDWS", "T2LWIDWS", "T3LWIDWS", "T4LWIDWS")
D = D[rowSums(is.na(D[,focal.var])) != length(focal.var),]
D <- D %>%
dplyr::select("SubjectID", "T0SY1Y2G_Y1Grade", "T0SY1Y2G_Y2Grade",
"T0SY1Y2G_StartGrade", "T0SGENGender", "T0PDEMOMaxParentEducation",
"T0PDEMOIncomeCode", "T0PDEMOIncomeMid","T0PDEMORaceHispanic",
"T1AWMADMSS_N", "T2AWMADMSS_N", "T3AWMADMSS_N", "T4AWMADMSS_N",
"T1AWMADMRS_N","T2AWMADMRS_N","T3AWMADMRS_N","T4AWMADMRS_N",
"T1GMPA", "T2GMPA", "T3GMPA", "T4GMPA",
"T1GPPA", "T2GPPA", "T3GPPA", "T4GPPA",
"T1PRWMPAE", "T2PRWMPAE", "T3PRWMPAE", "T4PRWMPAE",
"T1PMPA", "T2PMPA", "T3PMPA", "T4PMPA",
"T1NL10W_PAE","T2NL10W_PAE","T3NL10W_PAE","T4NL10W_PAE",
"T1NL100W_PAE", "T1NL1000W_PAE", "T2NL100W_PAE", "T2NL1000W_PAE",
"T3NL100W_PAE", "T3NL1000W_PAE", "T4NL100W_PAE", "T4NL1000W_PAE",
"T1PKC_PA", "T2PKC_PA", "T3PKC_PA", "T4PKC_PA",
"T1WJCW", "T2WJCW", "T3WJCW", "T4WJCW",
"T1WJCRS", "T2WJCRS", "T3WJCRS", "T4WJCRS",
"T1ASWmac", "T2ASWmac", "T3ASWmac", "T4ASWmac",
"T1LWIDWS", "T2LWIDWS", "T3LWIDWS", "T4LWIDWS")
D[D>=(-999) & D<=(-50)] <- NA
##drop participants who repeated a grade
D1_dummy = D[c("SubjectID", "T0SY1Y2G_Y1Grade", "T0SY1Y2G_Y2Grade",
"T0SY1Y2G_StartGrade")]
D1_dummy$repeated = ifelse(D1_dummy$T0SY1Y2G_Y1Grade == D1_dummy$T0SY1Y2G_Y2Grade, 1,0)
repeatgrade = D1_dummy[which(D1_dummy$repeated==1),"SubjectID"]
D1 = D[which(D$SubjectID %!in% repeatgrade$SubjectID),] ##exclude 3 children who repeat grade
table(D1$T0SY1Y2G_StartGrade)
##
## -1 0 1 2 3
## 103 149 131 130 99
D1$GMPA_na = rowSums(is.na(D1[c("T1GMPA","T2GMPA","T3GMPA","T4GMPA")]))
D1$GPPA_na = rowSums(is.na(D1[c("T1GPPA","T2GPPA","T3GPPA","T4GPPA")]))
D1$AWMADMSS_na = rowSums(is.na(D1[c("T1AWMADMSS_N","T2AWMADMSS_N","T3AWMADMSS_N","T4AWMADMSS_N")]))
D1$PRWM_na = rowSums(is.na(D1[c("T1PRWMPAE","T2PRWMPAE","T3PRWMPAE","T4PRWMPAE")]))
D1$PKC_na = rowSums(is.na(D1[c("T1PKC_PA", "T2PKC_PA", "T3PKC_PA", "T4PKC_PA")]))
D1$WJCW_na = rowSums(is.na(D1[c("T1WJCW", "T2WJCW", "T3WJCW", "T4WJCW")]))
D1$ASW_na = rowSums(is.na(D1[c("T1ASWmac", "T2ASWmac", "T3ASWmac", "T4ASWmac")]))
D1$PMPA_na = rowSums(is.na(D1[c("T1PMPA", "T2PMPA", "T3PMPA", "T4PMPA")]))
D1$LWID_na = rowSums(is.na(D1[c("T1LWIDWS", "T2LWIDWS", "T3LWIDWS", "T4LWIDWS")]))
participant_summary = D1[c("SubjectID","T0SY1Y2G_Y1Grade","T0SY1Y2G_StartGrade","T0SY1Y2G_Y2Grade","GMPA_na","GPPA_na","AWMADMSS_na","PRWM_na","PKC_na","WJCW_na","ASW_na","PMPA_na", "LWID_na")]
table(participant_summary$T0SY1Y2G_StartGrade, participant_summary$T0SY1Y2G_Y1Grade)
##
## -1 0 1 2 3
## -1 88 0 0 0 0
## 0 0 143 0 0 0
## 1 0 0 126 0 0
## 2 0 0 0 116 0
## 3 0 0 0 0 95
## Mental rotation
participant_summary$mr = participant_summary$GMPA_na + participant_summary$GPPA_na
## Exact calculation
participant_summary$ec = participant_summary$PKC_na + participant_summary$WJCW_na
nrow(participant_summary)
## [1] 614
###exclude participants who missed more than two time points of the focal variables
participant_summary_2 <- participant_summary
participant_summary_2 = subset(participant_summary, AWMADMSS_na <3)
participant_summary_2 = subset(participant_summary_2, mr <7)
participant_summary_2 = subset(participant_summary_2, PMPA_na <3)
#participant_summary_2 = subset(participant_summary_2, mt <7)
participant_summary_2 = subset(participant_summary_2, ec <7)
participant_summary_2 = subset(participant_summary_2, ASW_na <3)
participant_summary_2 = subset(participant_summary_2, LWID_na <3)
participant_summary_2$SubjectID = as.factor(as.character(participant_summary_2$SubjectID))
###agg_spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt_spread is from "sf_recalculatepae.Rmd"
spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt<-read_csv("data/sf_proportionalreasoning_nltask_cleaned.csv")
names(spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt)[2] = "SubjectID"
spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt$SubjectID = as.factor(as.character(spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt$SubjectID))
participant_summary_2 = participant_summary_2 %>%
left_join(spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt[c("SubjectID","PRWMPAE_na")],"SubjectID")
participant_summary_2 = subset(participant_summary_2, PRWMPAE_na < 3)
nrow(participant_summary_2)
## [1] 421
#This is the amended proportional reasoning data (PAE was recalculated)
spatial_profiles_database = D1
spatial_profiles_database$SubjectID = as.factor(as.character(spatial_profiles_database$SubjectID))
spatial_profiles_database_atleast2tp = spatial_profiles_database[which(spatial_profiles_database$SubjectID %in% participant_summary_2$SubjectID),]
#This step adds proportional reasoning database
spatial_profiles_database_atleast2tp = spatial_profiles_database_atleast2tp %>%
left_join(spatial_proportional_data_nona_importantcols_names_alltimes_no20_noshortrt[c(2:6)], by = "SubjectID")
spatial_profiles_database_atleast2tp$T1PRWMPAE.asin = 2*asin(sqrt(spatial_profiles_database_atleast2tp$t1))
spatial_profiles_database_atleast2tp$T2PRWMPAE.asin = 2*asin(sqrt(spatial_profiles_database_atleast2tp$t2))
spatial_profiles_database_atleast2tp$T3PRWMPAE.asin = 2*asin(sqrt(spatial_profiles_database_atleast2tp$t3))
spatial_profiles_database_atleast2tp$T4PRWMPAE.asin = 2*asin(sqrt(spatial_profiles_database_atleast2tp$t4))
nrow(spatial_profiles_database_atleast2tp)
## [1] 421
var_zy1 = c("T1AWMADMSS_N","T2AWMADMSS_N",
"T1GMPA", "T2GMPA",
"T1GPPA", "T2GPPA",
"t1", "t2",
"T1PRWMPAE.asin", "T2PRWMPAE.asin",
"T1PMPA","T2PMPA",
"T1WJCW", "T2WJCW",
"T1PKC_PA", "T2PKC_PA",
"T1ASWmac", "T2ASWmac",
"T1LWIDWS", "T2LWIDWS")
var_zy2 = c("T3AWMADMSS_N","T4AWMADMSS_N",
"T3GPPA", "T4GPPA",
"T3GMPA", "T4GMPA",
"t3", "t4",
"T3PRWMPAE.asin", "T4PRWMPAE.asin",
"T3PMPA", "T4PMPA",
"T3WJCW", "T4WJCW",
"T3PKC_PA", "T4PKC_PA",
"T3ASWmac", "T4ASWmac",
"T3LWIDWS", "T4LWIDWS")
for (v in var_zy1) {
vs <- str_c(v, ".s")
spatial_profiles_database_atleast2tp[[vs]] <- ave(spatial_profiles_database_atleast2tp[[v]], spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade, FUN=scale)
}
for (v in var_zy2) {
vs <- str_c(v, ".s")
spatial_profiles_database_atleast2tp[[vs]] <- ave(spatial_profiles_database_atleast2tp[[v]], spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade, FUN=scale)
}
# reverse code NL & PM scores
spatial_profiles_database_atleast2tp$T1PRWMPAE.asinR = -spatial_profiles_database_atleast2tp$T1PRWMPAE.asin
spatial_profiles_database_atleast2tp$T2PRWMPAE.asinR = -spatial_profiles_database_atleast2tp$T2PRWMPAE.asin
spatial_profiles_database_atleast2tp$T3PRWMPAE.asinR = -spatial_profiles_database_atleast2tp$T3PRWMPAE.asin
spatial_profiles_database_atleast2tp$T4PRWMPAE.asinR = -spatial_profiles_database_atleast2tp$T4PRWMPAE.asin
spatial_profiles_database_atleast2tp$T1PRWMPAE.asin.sR = - spatial_profiles_database_atleast2tp$T1PRWMPAE.asin.s
spatial_profiles_database_atleast2tp$T2PRWMPAE.asin.sR = - spatial_profiles_database_atleast2tp$T2PRWMPAE.asin.s
spatial_profiles_database_atleast2tp$T3PRWMPAE.asin.sR = - spatial_profiles_database_atleast2tp$T3PRWMPAE.asin.s
spatial_profiles_database_atleast2tp$T4PRWMPAE.asin.sR = - spatial_profiles_database_atleast2tp$T4PRWMPAE.asin.s
spatial_profiles_database_atleast2tp$T1MR.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T1GPPA.s,spatial_profiles_database_atleast2tp$T1GMPA.s)
spatial_profiles_database_atleast2tp$T2MR.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T2GPPA.s,spatial_profiles_database_atleast2tp$T2GMPA.s)
spatial_profiles_database_atleast2tp$T3MR.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T3GPPA.s,spatial_profiles_database_atleast2tp$T3GMPA.s)
spatial_profiles_database_atleast2tp$T4MR.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T4GPPA.s,spatial_profiles_database_atleast2tp$T4GMPA.s)
# spatial_profiles_database_atleast2tp$T1MT.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T1CMTTPA.s,spatial_profiles_database_atleast2tp$T1THPA.s)
#
# spatial_profiles_database_atleast2tp$T2MT.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T2CMTTPA.s,spatial_profiles_database_atleast2tp$T2THPA.s)
#
# spatial_profiles_database_atleast2tp$T3MT.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T3CMTTPA.s,spatial_profiles_database_atleast2tp$T3THPA.s)
#
# spatial_profiles_database_atleast2tp$T4MT.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T4CMTTPA.s,spatial_profiles_database_atleast2tp$T4THPA.s)
#Exact Calculations
spatial_profiles_database_atleast2tp$T1EC.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T1PKC_PA.s,spatial_profiles_database_atleast2tp$T1WJCW.s)
spatial_profiles_database_atleast2tp$T2EC.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade <1, spatial_profiles_database_atleast2tp$T2PKC_PA.s,spatial_profiles_database_atleast2tp$T2WJCW.s)
spatial_profiles_database_atleast2tp$T3EC.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T3PKC_PA.s,spatial_profiles_database_atleast2tp$T3WJCW.s)
spatial_profiles_database_atleast2tp$T4EC.s = ifelse(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade <1, spatial_profiles_database_atleast2tp$T4PKC_PA.s,spatial_profiles_database_atleast2tp$T4WJCW.s)
nrow(spatial_profiles_database_atleast2tp)
## [1] 421
# write.csv(spatial_profiles_database_atleast2tp, "data/spatial_profiles_database_atleast2tp_062424_n421.csv")
#write.csv(spatial_profiles_database_atleast2tp, "data/spatial_profiles_database_atleast2tp_080524_n421.csv")
#table(spatial_profiles_database_atleast2tp$T0SY1Y2G_Y2Grade, spatial_profiles_database_atleast2tp$T0SY1Y2G_StartGrade)
D_info = read_sav("data/AccessDB_T0T4Merge_7.12.19.1.sav")
names(D_info)[1] = "SubjectID"
D_info <- D_info %>%
dplyr::select("SubjectID", "T0SY1Y2G_StartGrade", "T0SGENGender", "T0PDEMOMaxParentEducation", "T0PDEMOIncomeCode", "T0PDEMOIncomeMid", "T0PDEMORaceHispanic")
D_info = D_info[which(D_info$SubjectID %in% participant_summary_2$SubjectID),]
TAll_complete_info = D_info
nrow(TAll_complete_info)
## [1] 421
table(TAll_complete_info$T0SGENGender) #1=boys, 2=girls
##
## 1 2
## 230 190
table(TAll_complete_info$T0SY1Y2G_StartGrade)
##
## -1 0 1 2 3
## 63 105 90 96 67
table(TAll_complete_info$T0PDEMORaceHispanic)
##
## 1 2 3 4 5 6 7 10
## 1 13 185 1 82 1 39 53
#1 = American Indian/Alaskan Native
#2 = Asian/ Asian American,
#3 = Black/ African American,
#4 = Native Hawaiian/ Other Pacific Islander,
#5 = White,
#6 = Other,
#7 = Multirace,
#10 = Hispanic# ##
#frequency of race
race <- TAll_complete_info %>%
dplyr::group_by(T0PDEMORaceHispanic)%>%
dplyr::summarise(cnt = n())%>%
dplyr::mutate(freq = round (cnt/sum(cnt),3))%>%
arrange(desc(freq))
race
## # A tibble: 9 × 3
## T0PDEMORaceHispanic cnt freq
## <dbl+lbl> <int> <dbl>
## 1 3 [Black/ African American] 185 0.439
## 2 5 [White] 82 0.195
## 3 10 [Hispanic] 53 0.126
## 4 NA 46 0.109
## 5 7 [Multirace] 39 0.093
## 6 2 [Asian/ Asian American] 13 0.031
## 7 1 [American Indian/Alaskan Native] 1 0.002
## 8 4 [Native Hawaiian/ Other Pacific Islander] 1 0.002
## 9 6 [Other] 1 0.002
###parent educational level###
table(TAll_complete_info$T0PDEMOMaxParentEducation)
##
## 10 12 13 14 16 17 18
## 17 53 65 57 64 41 80
psych::describe(TAll_complete_info$T0PDEMOMaxParentEducation)
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 377 14.88 2.4 14 14.96 2.97 10 18 8 -0.16 -1.18 0.12
###family income###
table(TAll_complete_info$T0PDEMOIncomeMid)
##
## 7.5 25 42.5 62.5 87.5 100
## 49 90 61 58 33 66
des <- psych::describe(TAll_complete_info$T0PDEMOIncomeMid)
print(des,digits=3) #show the output, but round to 3 (trailing) digits
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 357 51.32 32.15 42.5 50.73 29.65 7.5 100 92.5 0.291 -1.29 1.702
#spatial_profiles_database_atleast2tp <- read_csv("data/spatial_profiles_database_atleast2tp_062424_n421.csv")
spatial_profiles_database_atleast2tp <- read_csv("data/spatial_profiles_database_atleast2tp_080524_n421.csv")
nrow(spatial_profiles_database_atleast2tp)
## [1] 421
T1 <- spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1AWMADMSS_N.s",
"T1MR.s",
"T1PRWMPAE.asin.sR",
"T1PMPA.s","T1EC.s","T1ASWmac.s", "T1LWIDWS.s")
T1 = as.data.frame(T1)
hist(T1[c(3:7)], pch = 19)
ggpairs(T1[c(3:7)], pch = 19)
T1 <- spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1AWMADMSS_N.s",
"T1MR.s",
"T1PRWMPAE.asin.sR",
"T1PMPA.s")
T1 = as.data.frame(T1)
T1_na = T1
T1_na$T1na = rowSums(is.na(T1_na[c( "T0SY1Y2G_StartGrade", "T1AWMADMSS_N.s","T1MR.s",
"T1PRWMPAE.asin.sR", "T1PMPA.s")]))
T1_complete = T1
T1_complete = subset(T1_complete, !is.na(T0SY1Y2G_StartGrade))
#Two data points in T1
T1_complete = T1_complete[rowSums(is.na(T1_complete))<2,]
# ###346 rows in T1_complete
# write.csv(T1_complete, "data/T1_complete_062424_n346.csv")
T1_complete = read.csv("data/T1_complete_062424_n346.csv")
set.seed(777)
vLPA_missing = T1_complete %>%
dplyr::select(-c(X,SubjectID, T0SY1Y2G_StartGrade))%>%
single_imputation(method = "missForest") %>%
estimate_profiles(1:10)
vLPA_missing
## tidyLPA analysis using mclust:
##
## Model Classes AIC BIC Entropy prob_min prob_max n_min n_max BLRT_p
## 1 1 3854.33 3885.10 1.00 1.00 1.00 1.00 1.00
## 1 2 3748.90 3798.91 0.75 0.92 0.92 0.46 0.54 0.01
## 1 3 3710.00 3779.23 0.76 0.73 0.93 0.09 0.46 0.01
## 1 4 3704.88 3793.35 0.79 0.68 0.92 0.06 0.46 0.05
## 1 5 3714.15 3821.85 0.68 0.47 0.94 0.06 0.45 0.44
## 1 6 3716.39 3843.33 0.66 0.41 0.89 0.02 0.35 0.30
## 1 7 3712.81 3858.97 0.66 0.53 0.86 0.01 0.38 0.14
## 1 8 3716.79 3882.18 0.68 0.59 0.87 0.04 0.35 0.55
## 1 9 3715.34 3899.97 0.69 0.60 0.90 0.01 0.34 0.10
## 1 10 3686.71 3890.57 0.76 0.55 0.89 0.01 0.34 0.01
compare_solutions(vLPA_missing, statistics = c("AIC","BIC"))
## Compare tidyLPA solutions:
##
## Model Classes AIC BIC Warnings
## 1 1 3854 3885
## 1 2 3749 3799
## 1 3 3710 3779
## 1 4 3705 3793
## 1 5 3714 3822
## 1 6 3716 3843
## 1 7 3713 3859 Warning
## 1 8 3717 3882
## 1 9 3715 3900 Warning
## 1 10 3687 3891 Warning
##
## Best model according to AIC is Model 1 with 10 classes.
## Best model according to BIC is Model 1 with 3 classes.
##
## An analytic hierarchy process, based on the fit indices AIC, AWE, BIC, CLC, and KIC (Akogul & Erisoglu, 2017), suggests the best solution is Model 1 with 3 classes.
set.seed(777)
vLPA4 = T1_complete %>%
dplyr::select(-c(X, SubjectID, T0SY1Y2G_StartGrade)) %>%
single_imputation(method = "missForest") %>%
estimate_profiles(3)
profile_data = as.data.frame(get_data(vLPA4))
profile_data = cbind(T1_complete, profile_data)
table( profile_data$T0SY1Y2G_StartGrade,profile_data$Class)
##
## 1 2 3
## -1 25 13 0
## 0 38 40 13
## 1 30 38 8
## 2 41 38 5
## 3 22 31 4
table(profile_data$Class)
##
## 1 2 3
## 156 160 30
profile_data$SubjectID=as.factor(as.character(profile_data$SubjectID))
TAll_complete_info$SubjectID = as.factor(as.character(TAll_complete_info$SubjectID))
T1_complete_class$Class_num = as.factor(as.numeric(T1_complete_class$Class))
T1_complete_class$Class = as.factor(as.numeric(T1_complete_class$Class))
T1_complete_class$Class = ifelse(T1_complete_class$Class==1, "lowmr",ifelse(T1_complete_class$Class==3, "low","high"))
T1_complete_class$Class = as.factor(as.character(T1_complete_class$Class))
T1_complete_class = T1_complete_class %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
spatial_profiles_database_atleast2tp$SubjectID = as.factor(as.character(spatial_profiles_database_atleast2tp$SubjectID))
T1_complete_class = T1_complete_class %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1LWIDWS")], by ="SubjectID")
T1_complete_class$Class = relevel(T1_complete_class$Class, ref = "lowmr")
fit_reflowmr = multinom(Class ~ T0PDEMOMaxParentEducation
+ T0SGENGender
+ T0SY1Y2G_StartGrade
+ T1LWIDWS,
data = T1_complete_class)
## # weights: 18 (10 variable)
## initial value 314.203115
## iter 10 value 253.353167
## final value 249.112738
## converged
# Print the model summary
summary(fit_reflowmr)
## Call:
## multinom(formula = Class ~ T0PDEMOMaxParentEducation + T0SGENGender +
## T0SY1Y2G_StartGrade + T1LWIDWS, data = T1_complete_class)
##
## Coefficients:
## (Intercept) T0PDEMOMaxParentEducation T0SGENGender T0SY1Y2G_StartGrade
## high -8.664 0.17575 0.2317 -0.3362
## low 2.126 0.01714 1.2594 0.6523
## T1LWIDWS
## high 0.01475
## low -0.01590
##
## Std. Errors:
## (Intercept) T0PDEMOMaxParentEducation T0SGENGender T0SY1Y2G_StartGrade
## high 1.9345 0.05594 0.2651 0.2084
## low 0.4652 0.08571 0.4450 0.2405
## T1LWIDWS
## high 0.004854
## low 0.003882
##
## Residual Deviance: 498.2
## AIC: 518.2
stargazer(fit_reflowmr, header=FALSE, type='text')
##
## ======================================================
## Dependent variable:
## ----------------------------
## high low
## (1) (2)
## ------------------------------------------------------
## T0PDEMOMaxParentEducation 0.176*** 0.017
## (0.056) (0.086)
##
## T0SGENGender 0.232 1.259***
## (0.265) (0.445)
##
## T0SY1Y2G_StartGrade -0.336 0.652***
## (0.208) (0.241)
##
## T1LWIDWS 0.015*** -0.016***
## (0.005) (0.004)
##
## Constant -8.664*** 2.126***
## (1.935) (0.465)
##
## ------------------------------------------------------
## Akaike Inf. Crit. 518.200 518.200
## ======================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
T1_complete_class$Class = relevel(T1_complete_class$Class, ref = "low")
fit_reflow = multinom(Class ~ T0PDEMOMaxParentEducation
+ T0SGENGender
+ T0SY1Y2G_StartGrade
+ T1LWIDWS,
data = T1_complete_class)
## # weights: 18 (10 variable)
## initial value 314.203115
## iter 10 value 254.678612
## final value 249.112737
## converged
# Print the model summary
summary(fit_reflow)
## Call:
## multinom(formula = Class ~ T0PDEMOMaxParentEducation + T0SGENGender +
## T0SY1Y2G_StartGrade + T1LWIDWS, data = T1_complete_class)
##
## Coefficients:
## (Intercept) T0PDEMOMaxParentEducation T0SGENGender T0SY1Y2G_StartGrade
## lowmr -2.119 -0.01719 -1.259 -0.6516
## high -10.786 0.15858 -1.028 -0.9880
## T1LWIDWS
## lowmr 0.01589
## high 0.03065
##
## Std. Errors:
## (Intercept) T0PDEMOMaxParentEducation T0SGENGender T0SY1Y2G_StartGrade
## lowmr 1.067 0.08590 0.4449 0.2534
## high 1.024 0.08689 0.4524 0.2492
## T1LWIDWS
## lowmr 0.004496
## high 0.004342
##
## Residual Deviance: 498.2
## AIC: 518.2
stargazer(fit_reflow, header=FALSE, type='text')
##
## ======================================================
## Dependent variable:
## ----------------------------
## lowmr high
## (1) (2)
## ------------------------------------------------------
## T0PDEMOMaxParentEducation -0.017 0.159*
## (0.086) (0.087)
##
## T0SGENGender -1.259*** -1.028**
## (0.445) (0.452)
##
## T0SY1Y2G_StartGrade -0.652** -0.988***
## (0.253) (0.249)
##
## T1LWIDWS 0.016*** 0.031***
## (0.004) (0.004)
##
## Constant -2.119** -10.790***
## (1.067) (1.024)
##
## ------------------------------------------------------
## Akaike Inf. Crit. 518.200 518.200
## ======================================================
## Note: *p<0.1; **p<0.05; ***p<0.01
T1_complete_class_gather = T1_complete_class_gather %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1LWIDWS")], by ="SubjectID")
T1_complete_class_gather = T1_complete_class_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
model_task_class <- lmer(zscore ~ Class * task + T1LWIDWS + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = T1_complete_class_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_task_class)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## zscore ~ Class * task + T1LWIDWS + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: T1_complete_class_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2570
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.958 -0.643 0.052 0.660 2.930
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.00593 0.077
## Residual 0.63431 0.796
## Number of obs: 1056, groups: SubjectID, 286
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -0.0683227 0.2838493 314.7142788 -0.24
## Class1 -0.5786159 0.1044230 1038.2358900 -5.54
## Class3 -1.5560530 0.1938704 1040.1122010 -8.03
## taskT1MR.s 0.5240379 0.0996535 794.5224605 5.26
## taskT1PMPA.s -0.0039582 0.1029171 821.4945545 -0.04
## taskT1PRWMPAE.asin.sR -0.1888127 0.1016285 810.6134985 -1.86
## T1LWIDWS 0.0000167 0.0004782 291.0511886 0.03
## T0SGENGender 0.0321752 0.0511097 286.6558268 0.63
## T0PDEMOMaxParentEducation 0.0267391 0.0106173 288.3450460 2.52
## Class1:taskT1MR.s -1.1528220 0.1432235 791.7955219 -8.05
## Class3:taskT1MR.s 0.0789968 0.2554240 841.8799580 0.31
## Class1:taskT1PMPA.s 0.1808449 0.1487490 822.6586510 1.22
## Class3:taskT1PMPA.s -0.8656883 0.2604256 862.6760249 -3.32
## Class1:taskT1PRWMPAE.asin.sR 0.3467853 0.1455110 805.1311786 2.38
## Class3:taskT1PRWMPAE.asin.sR 0.4391036 0.2562535 843.7229533 1.71
## Pr(>|t|)
## (Intercept) 0.80994
## Class1 0.0000000381018909 ***
## Class3 0.0000000000000027 ***
## taskT1MR.s 0.0000001867631674 ***
## taskT1PMPA.s 0.96933
## taskT1PRWMPAE.asin.sR 0.06355 .
## T1LWIDWS 0.97211
## T0SGENGender 0.52950
## T0PDEMOMaxParentEducation 0.01233 *
## Class1:taskT1MR.s 0.0000000000000030 ***
## Class3:taskT1MR.s 0.75719
## Class1:taskT1PMPA.s 0.22442
## Class3:taskT1PMPA.s 0.00092 ***
## Class1:taskT1PRWMPAE.asin.sR 0.01739 *
## Class3:taskT1PRWMPAE.asin.sR 0.08698 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Anova(model_task_class, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: zscore
## F Df Df.res Pr(>F)
## (Intercept) 0.06 1 303 0.810
## Class 38.42 2 1039 < 0.0000000000000002 ***
## task 19.25 3 797 0.0000000000047 ***
## T1LWIDWS 0.00 1 279 0.972
## T0SGENGender 0.40 1 275 0.530
## T0PDEMOMaxParentEducation 6.34 1 277 0.012 *
## Class:task 30.98 6 806 < 0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
means_model_task_class = emmeans(model_task_class, pairwise ~ Class|task, mult.name = "task")
summary(means_model_task_class, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "bonferroni")
## $emmeans
## task = T1AWMADMSS_N.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3846 0.0726 1040 0.210 0.5586 5.298 <.0001
## 1 -0.1940 0.0745 1039 -0.373 -0.0154 -2.605 0.0279
## 3 -1.1715 0.1794 1041 -1.602 -0.7412 -6.528 <.0001
##
## task = T1MR.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.9086 0.0699 1039 0.741 1.0762 13.001 <.0001
## 1 -0.8228 0.0727 1039 -0.997 -0.6484 -11.315 <.0001
## 3 -0.5684 0.1545 1040 -0.939 -0.1981 -3.680 0.0007
##
## task = T1PMPA.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3806 0.0744 1040 0.202 0.5590 5.117 <.0001
## 1 -0.0172 0.0790 1040 -0.207 0.1722 -0.217 1.0000
## 3 -2.0411 0.1607 1040 -2.426 -1.6558 -12.701 <.0001
##
## task = T1PRWMPAE.asin.sR:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.1958 0.0729 1040 0.021 0.3705 2.686 0.0220
## 1 -0.0361 0.0745 1039 -0.215 0.1424 -0.485 1.0000
## 3 -0.9212 0.1547 1040 -1.292 -0.5501 -5.953 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 3 estimates
## P value adjustment: bonferroni method for 3 tests
##
## $contrasts
## task = T1AWMADMSS_N.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.579 0.104 1038 0.328 0.829 5.540 <.0001
## Class2 - Class3 1.556 0.194 1040 1.091 2.021 8.022 <.0001
## Class1 - Class3 0.977 0.194 1040 0.511 1.443 5.031 <.0001
##
## task = T1MR.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 1.731 0.101 1038 1.489 1.974 17.095 <.0001
## Class2 - Class3 1.477 0.170 1039 1.069 1.885 8.687 <.0001
## Class1 - Class3 -0.254 0.171 1040 -0.664 0.155 -1.490 0.4094
##
## task = T1PMPA.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.398 0.109 1039 0.137 0.659 3.655 0.0008
## Class2 - Class3 2.422 0.178 1039 1.996 2.847 13.642 <.0001
## Class1 - Class3 2.024 0.179 1040 1.595 2.453 11.307 <.0001
##
## task = T1PRWMPAE.asin.sR:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.232 0.105 1038 -0.019 0.483 2.217 0.0806
## Class2 - Class3 1.117 0.172 1039 0.705 1.529 6.501 <.0001
## Class1 - Class3 0.885 0.172 1040 0.473 1.297 5.153 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
## Conf-level adjustment: bonferroni method for 3 estimates
## P value adjustment: bonferroni method for 3 tests
T1_complete_class = T1_complete_class %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1LWIDWS.s")], by ="SubjectID")
T1_complete_class_lw = T1_complete_class[c("SubjectID","T1LWIDWS.s","Class")]
T1_complete_class_lw = T1_complete_class_lw[rowSums(is.na(T1_complete_class_lw))==0,]
ezANOVA(T1_complete_class_lw, dv = T1LWIDWS.s, between = Class, wid = SubjectID)
## $ANOVA
## Effect DFn DFd F p p<.05 ges
## 1 Class 2 317 13.02 0.000003686 * 0.0759
##
## $`Levene's Test for Homogeneity of Variance`
## DFn DFd SSn SSd F p p<.05
## 1 2 317 0.3525 121 0.4619 0.6305
summarySE(T1_complete_class_lw, "T1LWIDWS.s","Class", na.rm = T)
## Class N T1LWIDWS.s sd se ci
## 1 low 29 -0.4380 1.0539 0.19570 0.4009
## 2 lowmr 146 -0.2200 0.9309 0.07705 0.1523
## 3 high 145 0.2839 0.9816 0.08151 0.1611
spatial_profiles_database_atleast2tp$SubjectID=as.factor(as.character(spatial_profiles_database_atleast2tp$SubjectID))
T1_complete_class = T1_complete_class %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1EC.s","T1ASWmac.s")], by ="SubjectID")
T1_complete_class_ec = T1_complete_class[c("SubjectID","T1EC.s","Class")]
T1_complete_class_ec = T1_complete_class_ec[rowSums(is.na(T1_complete_class_ec))==0,]
T1_complete_class_ec = T1_complete_class_ec %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1LWIDWS.s")], by ="SubjectID")
T1_complete_class_ec = subset(T1_complete_class_ec, !is.na(T1LWIDWS.s))
T1_complete_class_ec$Class = as.factor(T1_complete_class_ec$Class)
summary(lm(T1EC.s ~ Class + T1LWIDWS.s, T1_complete_class_ec))
##
## Call:
## lm(formula = T1EC.s ~ Class + T1LWIDWS.s, data = T1_complete_class_ec)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.5701 -0.5152 0.0275 0.4952 2.5192
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.6647 0.1633 -4.07 0.00006004 ***
## Classlowmr 0.4206 0.1781 2.36 0.019 *
## Classhigh 0.9715 0.1812 5.36 0.00000017 ***
## T1LWIDWS.s 0.2510 0.0508 4.94 0.00000131 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.856 on 297 degrees of freedom
## Multiple R-squared: 0.234, Adjusted R-squared: 0.226
## F-statistic: 30.2 on 3 and 297 DF, p-value: <0.0000000000000002
Anova(lm(T1EC.s ~ Class + T1LWIDWS.s, T1_complete_class_ec), test.statistic = "F", type= "II")
## Anova Table (Type II tests)
##
## Response: T1EC.s
## Sum Sq Df F value Pr(>F)
## Class 30.6 2 20.9 0.0000000033 ***
## T1LWIDWS.s 17.9 1 24.4 0.0000013070 ***
## Residuals 217.8 297
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
model_ec = lm(T1EC.s ~ Class + T1LWIDWS.s, T1_complete_class_ec)
means_model_ec = emmeans(model_ec, pairwise ~ Class)
summary(means_model_ec, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.667 0.1632 297 -0.988 -0.3457 -4.086 0.0001
## lowmr -0.246 0.0746 297 -0.393 -0.0995 -3.303 0.0011
## high 0.305 0.0745 297 0.158 0.4512 4.090 0.0001
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.421 0.178 297 -0.771 -0.0701 -2.362 0.0188
## low - high -0.972 0.181 297 -1.328 -0.6149 -5.362 <.0001
## lowmr - high -0.551 0.107 297 -0.761 -0.3403 -5.149 <.0001
##
## Confidence level used: 0.95
summarySE(T1_complete_class_ec, "T1EC.s","Class", na.rm = T)
## Class N T1EC.s sd se ci
## 1 low 28 -0.7714 1.0637 0.20102 0.4125
## 2 lowmr 135 -0.3017 0.8721 0.07506 0.1484
## 3 high 138 0.3800 0.8682 0.07390 0.1461
T1_complete_class_asw = T1_complete_class[c("SubjectID","T1ASWmac.s","Class")]
T1_complete_class_asw = T1_complete_class_asw[rowSums(is.na(T1_complete_class_asw))==0,]
T1_complete_class_asw = T1_complete_class_asw %>%
left_join(spatial_profiles_database_atleast2tp[c("SubjectID","T1LWIDWS.s")], by ="SubjectID")
T1_complete_class_asw = subset(T1_complete_class_asw, !is.na(T1LWIDWS.s))
T1_complete_class_asw$Class = as.factor(T1_complete_class_asw$Class)
summary(lm(T1ASWmac.s ~ Class + T1LWIDWS.s, T1_complete_class_asw))
##
## Call:
## lm(formula = T1ASWmac.s ~ Class + T1LWIDWS.s, data = T1_complete_class_asw)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.4113 -0.6827 0.0344 0.6862 2.0733
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.3387 0.1868 -1.81 0.0709 .
## Classlowmr 0.2449 0.2030 1.21 0.2286
## Classhigh 0.4994 0.2081 2.40 0.0171 *
## T1LWIDWS.s 0.1646 0.0597 2.76 0.0062 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.958 on 287 degrees of freedom
## Multiple R-squared: 0.0699, Adjusted R-squared: 0.0602
## F-statistic: 7.19 on 3 and 287 DF, p-value: 0.000114
model_aswmac = lm(T1ASWmac.s ~ Class + T1LWIDWS.s, T1_complete_class_asw)
means_model_aswmac = emmeans(model_aswmac, pairwise ~ Class)
summary(means_model_aswmac, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.340 0.1868 287 -0.7075 0.0277 -1.820 0.0698
## lowmr -0.095 0.0839 287 -0.2602 0.0701 -1.132 0.2584
## high 0.159 0.0857 287 -0.0092 0.3281 1.861 0.0638
##
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.245 0.203 287 -0.644 0.1546 -1.206 0.2286
## low - high -0.499 0.208 287 -0.909 -0.0897 -2.399 0.0171
## lowmr - high -0.255 0.122 287 -0.494 -0.0148 -2.090 0.0375
##
## Confidence level used: 0.95
Anova(model_aswmac, test.statistic = "F", type= "II")
## Anova Table (Type II tests)
##
## Response: T1ASWmac.s
## Sum Sq Df F value Pr(>F)
## Class 7 2 3.81 0.0233 *
## T1LWIDWS.s 7 1 7.61 0.0062 **
## Residuals 263 287
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
panamath_longitudinal_nonz_gather = panamath_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
panamath_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(panamath_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
panamath_longitudinal_nonz_gather$time_num = as.numeric(as.factor(panamath_longitudinal_nonz_gather$time))
model_long_panamath <- lmer(PMPA ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = panamath_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_panamath)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## PMPA ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: panamath_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -1490
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.446 -0.447 0.124 0.538 2.774
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.00586 0.0766
## Residual 0.00840 0.0916
## Number of obs: 995, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 0.49706 0.04443 468.75345 11.19
## T0SY1Y2G_StartGrade0 0.18980 0.02983 931.05323 6.36
## T0SY1Y2G_StartGrade1 0.24113 0.03129 937.80293 7.71
## T0SY1Y2G_StartGrade2 0.31827 0.03081 939.50409 10.33
## T0SY1Y2G_StartGrade3 0.30202 0.03404 946.89484 8.87
## time_num 0.06266 0.00861 768.64065 7.28
## T0SGENGender -0.01688 0.01060 305.47312 -1.59
## T0PDEMOMaxParentEducation 0.00567 0.00223 308.50105 2.54
## T0SY1Y2G_StartGrade0:time_num -0.03570 0.00986 754.18586 -3.62
## T0SY1Y2G_StartGrade1:time_num -0.03816 0.01031 754.63886 -3.70
## T0SY1Y2G_StartGrade2:time_num -0.05907 0.01023 765.12757 -5.77
## T0SY1Y2G_StartGrade3:time_num -0.05215 0.01121 757.38759 -4.65
## Pr(>|t|)
## (Intercept) < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade0 0.000000000311057 ***
## T0SY1Y2G_StartGrade1 0.000000000000033 ***
## T0SY1Y2G_StartGrade2 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade3 < 0.0000000000000002 ***
## time_num 0.000000000000854 ***
## T0SGENGender 0.11223
## T0PDEMOMaxParentEducation 0.01156 *
## T0SY1Y2G_StartGrade0:time_num 0.00031 ***
## T0SY1Y2G_StartGrade1:time_num 0.00023 ***
## T0SY1Y2G_StartGrade2:time_num 0.000000011273025 ***
## T0SY1Y2G_StartGrade3:time_num 0.000003915749444 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3
## T0SY1Y2G_StG0 -0.477
## T0SY1Y2G_StG1 -0.502 0.698
## T0SY1Y2G_StG2 -0.504 0.709 0.678
## T0SY1Y2G_StG3 -0.466 0.642 0.614 0.623
## time_num -0.452 0.672 0.640 0.650 0.589
## T0SGENGendr -0.317 0.032 0.027 0.012 0.039
## T0PDEMOMxPE -0.733 -0.035 0.031 0.031 0.027
## T0SY1Y2G_SG0: 0.394 -0.777 -0.559 -0.567 -0.514
## T0SY1Y2G_SG1: 0.382 -0.561 -0.782 -0.543 -0.492
## T0SY1Y2G_SG2: 0.372 -0.566 -0.538 -0.784 -0.495
## T0SY1Y2G_SG3: 0.344 -0.516 -0.491 -0.499 -0.789
## tim_nm T0SGEN T0PDEM T0SY1Y2G_SG0: T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG0
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr 0.021
## T0PDEMOMxPE -0.007 -0.063
## T0SY1Y2G_SG0: -0.874 -0.024 0.011
## T0SY1Y2G_SG1: -0.836 -0.019 0.002 0.730
## T0SY1Y2G_SG2: -0.842 -0.014 0.017 0.735 0.703
## T0SY1Y2G_SG3: -0.768 -0.021 0.013 0.671 0.642 0.647
Anova(model_long_panamath, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: PMPA
## F Df Df.res Pr(>F)
## (Intercept) 125.13 1 464 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade 31.07 4 936 < 0.0000000000000002 ***
## time_num 52.88 1 765 0.00000000000088 ***
## T0SGENGender 2.54 1 301 0.112
## T0PDEMOMaxParentEducation 6.45 1 304 0.012 *
## T0SY1Y2G_StartGrade:time_num 9.30 4 738 0.00000024763522 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_panamath, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio
## -1 0.06266 0.00862 765 0.04574 0.0796 7.272
## 0 0.02695 0.00480 705 0.01754 0.0364 5.620
## 1 0.02449 0.00566 718 0.01338 0.0356 4.326
## 2 0.00359 0.00553 753 -0.00726 0.0144 0.649
## 3 0.01051 0.00718 738 -0.00359 0.0246 1.463
## p.value
## <.0001
## <.0001
## <.0001
## 0.5162
## 0.1439
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade0 0.03570 0.00986 751 0.01634
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade1 0.03816 0.01031 751 0.01792
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade2 0.05907 0.01023 762 0.03897
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade3 0.05215 0.01122 754 0.03012
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade1 0.00246 0.00742 713 -0.01211
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade2 0.02337 0.00732 733 0.00900
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade3 0.01644 0.00864 728 -0.00051
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade2 0.02090 0.00791 735 0.00537
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade3 0.01398 0.00915 731 -0.00397
## T0SY1Y2G_StartGrade2 - T0SY1Y2G_StartGrade3 -0.00692 0.00906 744 -0.02471
## upper.CL t.ratio p.value
## 0.0551 3.620 0.0003
## 0.0584 3.702 0.0002
## 0.0792 5.771 <.0001
## 0.0742 4.648 <.0001
## 0.0170 0.332 0.7402
## 0.0377 3.194 0.0015
## 0.0334 1.904 0.0573
## 0.0364 2.642 0.0084
## 0.0319 1.529 0.1268
## 0.0109 -0.764 0.4452
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
means_model_task_class = emmeans(model_task_class, pairwise ~ Class|task, mult.name = "task")
summary(means_model_task_class, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## task = T1AWMADMSS_N.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3846 0.0726 1040 0.2421 0.5270 5.298 <.0001
## 1 -0.1940 0.0745 1039 -0.3402 -0.0479 -2.605 0.0093
## 3 -1.1715 0.1794 1041 -1.5236 -0.8193 -6.528 <.0001
##
## task = T1MR.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.9086 0.0699 1039 0.7715 1.0457 13.001 <.0001
## 1 -0.8228 0.0727 1039 -0.9655 -0.6801 -11.315 <.0001
## 3 -0.5684 0.1545 1040 -0.8715 -0.2654 -3.680 0.0002
##
## task = T1PMPA.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3806 0.0744 1040 0.2347 0.5266 5.117 <.0001
## 1 -0.0172 0.0790 1040 -0.1721 0.1378 -0.217 0.8281
## 3 -2.0411 0.1607 1040 -2.3565 -1.7258 -12.701 <.0001
##
## task = T1PRWMPAE.asin.sR:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.1958 0.0729 1040 0.0528 0.3388 2.686 0.0073
## 1 -0.0361 0.0745 1039 -0.1822 0.1100 -0.485 0.6281
## 3 -0.9212 0.1547 1040 -1.2248 -0.6175 -5.953 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## task = T1AWMADMSS_N.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.579 0.104 1038 0.3737 0.7835 5.540 <.0001
## Class2 - Class3 1.556 0.194 1040 1.1754 1.9367 8.022 <.0001
## Class1 - Class3 0.977 0.194 1040 0.5962 1.3587 5.031 <.0001
##
## task = T1MR.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 1.731 0.101 1038 1.5327 1.9302 17.095 <.0001
## Class2 - Class3 1.477 0.170 1039 1.1434 1.8107 8.687 <.0001
## Class1 - Class3 -0.254 0.171 1040 -0.5893 0.0806 -1.490 0.1365
##
## task = T1PMPA.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.398 0.109 1039 0.1842 0.6113 3.655 0.0003
## Class2 - Class3 2.422 0.178 1039 2.0734 2.7701 13.642 <.0001
## Class1 - Class3 2.024 0.179 1040 1.6727 2.3752 11.307 <.0001
##
## task = T1PRWMPAE.asin.sR:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.232 0.105 1038 0.0266 0.4371 2.217 0.0269
## Class2 - Class3 1.117 0.172 1039 0.7798 1.4541 6.501 <.0001
## Class1 - Class3 0.885 0.172 1040 0.5481 1.2222 5.153 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
mentalrotation_longitudinal_nonz_gather = mentalrotation_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
mentalrotation_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(mentalrotation_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
mentalrotation_longitudinal_nonz_gather$time_num = as.numeric(as.factor(mentalrotation_longitudinal_nonz_gather$time))
model_long_mr <- lmer(GMPA ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = mentalrotation_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_mr)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## GMPA ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: mentalrotation_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -93.1
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7097 -0.5277 0.0489 0.5224 2.6031
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.0457 0.214
## Residual 0.0282 0.168
## Number of obs: 859, groups: SubjectID, 271
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -0.27444 0.14234 695.81779 -1.93
## T0SY1Y2G_StartGrade1 0.42549 0.10867 751.50631 3.92
## T0SY1Y2G_StartGrade2 0.50868 0.10789 743.87733 4.71
## T0SY1Y2G_StartGrade3 0.64654 0.11127 772.10346 5.81
## time_num 0.14849 0.02809 602.56497 5.29
## T0SGENGender 0.04912 0.02896 262.79858 1.70
## T0PDEMOMaxParentEducation 0.00956 0.00617 262.56117 1.55
## T0SY1Y2G_StartGrade1:time_num -0.06012 0.02965 601.03989 -2.03
## T0SY1Y2G_StartGrade2:time_num -0.08576 0.02948 601.70868 -2.91
## T0SY1Y2G_StartGrade3:time_num -0.09120 0.03034 600.84724 -3.01
## Pr(>|t|)
## (Intercept) 0.0543 .
## T0SY1Y2G_StartGrade1 0.0000983614 ***
## T0SY1Y2G_StartGrade2 0.0000028885 ***
## T0SY1Y2G_StartGrade3 0.0000000091 ***
## time_num 0.0000001757 ***
## T0SGENGender 0.0910 .
## T0PDEMOMaxParentEducation 0.1227
## T0SY1Y2G_StartGrade1:time_num 0.0430 *
## T0SY1Y2G_StartGrade2:time_num 0.0038 **
## T0SY1Y2G_StartGrade3:time_num 0.0028 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3 tim_nm T0SGEN
## T0SY1Y2G_StG1 -0.702
## T0SY1Y2G_StG2 -0.707 0.893
## T0SY1Y2G_StG3 -0.691 0.866 0.872
## time_num -0.692 0.905 0.911 0.884
## T0SGENGendr -0.227 0.008 -0.006 0.014 0.006
## T0PDEMOMxPE -0.627 0.036 0.043 0.041 -0.001 -0.115
## T0SY1Y2G_SG1: 0.657 -0.926 -0.864 -0.838 -0.948 -0.007
## T0SY1Y2G_SG2: 0.658 -0.862 -0.930 -0.842 -0.953 -0.002
## T0SY1Y2G_SG3: 0.639 -0.838 -0.844 -0.913 -0.926 -0.005
## T0PDEM T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr
## T0PDEMOMxPE
## T0SY1Y2G_SG1: -0.001
## T0SY1Y2G_SG2: 0.002 0.903
## T0SY1Y2G_SG3: 0.003 0.877 0.882
Anova(model_long_mr, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: GMPA
## F Df Df.res Pr(>F)
## (Intercept) 3.72 1 697 0.0543 .
## T0SY1Y2G_StartGrade 13.05 3 715 0.000000026 ***
## time_num 27.93 1 604 0.000000176 ***
## T0SGENGender 2.88 1 264 0.0910 .
## T0PDEMOMaxParentEducation 2.40 1 264 0.1227
## T0SY1Y2G_StartGrade:time_num 4.35 3 596 0.0048 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_mr, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio
## 0 0.1485 0.02810 604 0.0933 0.2037 5.285
## 1 0.0884 0.00947 589 0.0698 0.1070 9.334
## 2 0.0627 0.00893 594 0.0452 0.0803 7.021
## 3 0.0573 0.01146 592 0.0348 0.0798 4.999
## p.value
## <.0001
## <.0001
## <.0001
## <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade1 0.06012 0.0297 602 0.001885
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade2 0.08576 0.0295 603 0.027861
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade3 0.09120 0.0303 602 0.031600
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade2 0.02565 0.0130 592 0.000081
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade3 0.03108 0.0149 591 0.001884
## T0SY1Y2G_StartGrade2 - T0SY1Y2G_StartGrade3 0.00543 0.0145 593 -0.023106
## upper.CL t.ratio p.value
## 0.1183 2.027 0.0431
## 0.1437 2.909 0.0038
## 0.1508 3.005 0.0028
## 0.0512 1.970 0.0493
## 0.0603 2.091 0.0370
## 0.0340 0.374 0.7087
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
means_model_task_class = emmeans(model_task_class, pairwise ~ Class|task, mult.name = "task")
summary(means_model_task_class, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## task = T1AWMADMSS_N.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3846 0.0726 1040 0.2421 0.5270 5.298 <.0001
## 1 -0.1940 0.0745 1039 -0.3402 -0.0479 -2.605 0.0093
## 3 -1.1715 0.1794 1041 -1.5236 -0.8193 -6.528 <.0001
##
## task = T1MR.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.9086 0.0699 1039 0.7715 1.0457 13.001 <.0001
## 1 -0.8228 0.0727 1039 -0.9655 -0.6801 -11.315 <.0001
## 3 -0.5684 0.1545 1040 -0.8715 -0.2654 -3.680 0.0002
##
## task = T1PMPA.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3806 0.0744 1040 0.2347 0.5266 5.117 <.0001
## 1 -0.0172 0.0790 1040 -0.1721 0.1378 -0.217 0.8281
## 3 -2.0411 0.1607 1040 -2.3565 -1.7258 -12.701 <.0001
##
## task = T1PRWMPAE.asin.sR:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.1958 0.0729 1040 0.0528 0.3388 2.686 0.0073
## 1 -0.0361 0.0745 1039 -0.1822 0.1100 -0.485 0.6281
## 3 -0.9212 0.1547 1040 -1.2248 -0.6175 -5.953 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## task = T1AWMADMSS_N.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.579 0.104 1038 0.3737 0.7835 5.540 <.0001
## Class2 - Class3 1.556 0.194 1040 1.1754 1.9367 8.022 <.0001
## Class1 - Class3 0.977 0.194 1040 0.5962 1.3587 5.031 <.0001
##
## task = T1MR.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 1.731 0.101 1038 1.5327 1.9302 17.095 <.0001
## Class2 - Class3 1.477 0.170 1039 1.1434 1.8107 8.687 <.0001
## Class1 - Class3 -0.254 0.171 1040 -0.5893 0.0806 -1.490 0.1365
##
## task = T1PMPA.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.398 0.109 1039 0.1842 0.6113 3.655 0.0003
## Class2 - Class3 2.422 0.178 1039 2.0734 2.7701 13.642 <.0001
## Class1 - Class3 2.024 0.179 1040 1.6727 2.3752 11.307 <.0001
##
## task = T1PRWMPAE.asin.sR:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.232 0.105 1038 0.0266 0.4371 2.217 0.0269
## Class2 - Class3 1.117 0.172 1039 0.7798 1.4541 6.501 <.0001
## Class1 - Class3 0.885 0.172 1040 0.5481 1.2222 5.153 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
mentalrotation2_longitudinal_nonz_gather = mentalrotation2_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
mentalrotation2_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(mentalrotation2_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
mentalrotation2_longitudinal_nonz_gather$time_num = as.numeric(as.factor(mentalrotation2_longitudinal_nonz_gather$time))
model_long_mr2 <- lmer(GPPA ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = mentalrotation2_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_mr2)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## GPPA ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: mentalrotation2_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -219.9
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7856 -0.6321 0.0488 0.6292 3.0732
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.00785 0.0886
## Residual 0.01806 0.1344
## Number of obs: 287, groups: SubjectID, 118
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 0.15356 0.08621 131.43852 1.78
## T0SY1Y2G_StartGrade0 0.17949 0.04747 267.96725 3.78
## time_num 0.06270 0.01086 182.46448 5.77
## T0SGENGender 0.04158 0.02316 112.95406 1.80
## T0PDEMOMaxParentEducation 0.01643 0.00468 109.34072 3.51
## T0SY1Y2G_StartGrade0:time_num -0.03283 0.02382 176.84746 -1.38
## Pr(>|t|)
## (Intercept) 0.07717 .
## T0SY1Y2G_StartGrade0 0.00019 ***
## time_num 0.000000033 ***
## T0SGENGender 0.07521 .
## T0PDEMOMaxParentEducation 0.00065 ***
## T0SY1Y2G_StartGrade0:time_num 0.16978
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 tim_nm T0SGEN T0PDEM
## T0SY1Y2G_StG0 -0.235
## time_num -0.319 0.541
## T0SGENGendr -0.428 0.010 0.002
## T0PDEMOMxPE -0.837 -0.027 0.025 0.035
## T0SY1Y2G_SG0: 0.144 -0.842 -0.456 0.007 -0.013
Anova(model_long_mr2, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: GPPA
## F Df Df.res Pr(>F)
## (Intercept) 3.17 1 129 0.07740 .
## T0SY1Y2G_StartGrade 14.29 1 268 0.00019 ***
## time_num 33.26 1 180 0.000000034 ***
## T0SGENGender 3.22 1 111 0.07548 .
## T0PDEMOMaxParentEducation 12.29 1 107 0.00067 ***
## T0SY1Y2G_StartGrade:time_num 1.90 1 175 0.16994
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_mr2, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio
## -1 0.0627 0.0109 180 0.0413 0.0842 5.767
## 0 0.0299 0.0212 173 -0.0120 0.0717 1.409
## p.value
## <.0001
## 0.1606
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade0 0.0328 0.0238 175 -0.0142
## upper.CL t.ratio p.value
## 0.0799 1.378 0.1699
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
means_model_task_class = emmeans(model_task_class, pairwise ~ Class|task, mult.name = "task")
summary(means_model_task_class, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## task = T1AWMADMSS_N.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3846 0.0726 1040 0.2421 0.5270 5.298 <.0001
## 1 -0.1940 0.0745 1039 -0.3402 -0.0479 -2.605 0.0093
## 3 -1.1715 0.1794 1041 -1.5236 -0.8193 -6.528 <.0001
##
## task = T1MR.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.9086 0.0699 1039 0.7715 1.0457 13.001 <.0001
## 1 -0.8228 0.0727 1039 -0.9655 -0.6801 -11.315 <.0001
## 3 -0.5684 0.1545 1040 -0.8715 -0.2654 -3.680 0.0002
##
## task = T1PMPA.s:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.3806 0.0744 1040 0.2347 0.5266 5.117 <.0001
## 1 -0.0172 0.0790 1040 -0.1721 0.1378 -0.217 0.8281
## 3 -2.0411 0.1607 1040 -2.3565 -1.7258 -12.701 <.0001
##
## task = T1PRWMPAE.asin.sR:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## 2 0.1958 0.0729 1040 0.0528 0.3388 2.686 0.0073
## 1 -0.0361 0.0745 1039 -0.1822 0.1100 -0.485 0.6281
## 3 -0.9212 0.1547 1040 -1.2248 -0.6175 -5.953 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## task = T1AWMADMSS_N.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.579 0.104 1038 0.3737 0.7835 5.540 <.0001
## Class2 - Class3 1.556 0.194 1040 1.1754 1.9367 8.022 <.0001
## Class1 - Class3 0.977 0.194 1040 0.5962 1.3587 5.031 <.0001
##
## task = T1MR.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 1.731 0.101 1038 1.5327 1.9302 17.095 <.0001
## Class2 - Class3 1.477 0.170 1039 1.1434 1.8107 8.687 <.0001
## Class1 - Class3 -0.254 0.171 1040 -0.5893 0.0806 -1.490 0.1365
##
## task = T1PMPA.s:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.398 0.109 1039 0.1842 0.6113 3.655 0.0003
## Class2 - Class3 2.422 0.178 1039 2.0734 2.7701 13.642 <.0001
## Class1 - Class3 2.024 0.179 1040 1.6727 2.3752 11.307 <.0001
##
## task = T1PRWMPAE.asin.sR:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## Class2 - Class1 0.232 0.105 1038 0.0266 0.4371 2.217 0.0269
## Class2 - Class3 1.117 0.172 1039 0.7798 1.4541 6.501 <.0001
## Class1 - Class3 0.885 0.172 1040 0.5481 1.2222 5.153 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
workingmemory_longitudinal_nonz_gather = workingmemory_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
workingmemory_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(workingmemory_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
workingmemory_longitudinal_nonz_gather$time_num = as.numeric(as.factor(workingmemory_longitudinal_nonz_gather$time))
model_long_wm <- lmer(AWMADMRS_N ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = workingmemory_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_wm)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: AWMADMRS_N ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender +
## T0PDEMOMaxParentEducation + (1 | SubjectID)
## Data: workingmemory_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 6059
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.696 -0.573 0.039 0.562 2.680
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 7.71 2.78
## Residual 10.52 3.24
## Number of obs: 1092, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 2.3774 1.5378 438.3529 1.55
## T0SY1Y2G_StartGrade0 2.0542 1.0183 969.9144 2.02
## T0SY1Y2G_StartGrade1 5.5941 1.0530 969.7183 5.31
## T0SY1Y2G_StartGrade2 8.0849 1.0274 968.5050 7.87
## T0SY1Y2G_StartGrade3 11.2573 1.1265 972.3782 9.99
## time_num 1.5725 0.2689 817.9148 5.85
## T0SGENGender -0.1685 0.3753 301.8362 -0.45
## T0PDEMOMaxParentEducation 0.3499 0.0789 301.4215 4.44
## T0SY1Y2G_StartGrade0:time_num 0.1792 0.3203 816.0332 0.56
## T0SY1Y2G_StartGrade1:time_num -0.2161 0.3295 810.7887 -0.66
## T0SY1Y2G_StartGrade2:time_num -0.1977 0.3219 814.2002 -0.61
## T0SY1Y2G_StartGrade3:time_num -0.7027 0.3518 812.0247 -2.00
## Pr(>|t|)
## (Intercept) 0.123
## T0SY1Y2G_StartGrade0 0.044 *
## T0SY1Y2G_StartGrade1 0.0000001340539488 ***
## T0SY1Y2G_StartGrade2 0.0000000000000095 ***
## T0SY1Y2G_StartGrade3 < 0.0000000000000002 ***
## time_num 0.0000000071660447 ***
## T0SGENGender 0.654
## T0PDEMOMaxParentEducation 0.0000127782125618 ***
## T0SY1Y2G_StartGrade0:time_num 0.576
## T0SY1Y2G_StartGrade1:time_num 0.512
## T0SY1Y2G_StartGrade2:time_num 0.539
## T0SY1Y2G_StartGrade3:time_num 0.046 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3
## T0SY1Y2G_StG0 -0.443
## T0SY1Y2G_StG1 -0.482 0.671
## T0SY1Y2G_StG2 -0.493 0.688 0.668
## T0SY1Y2G_StG3 -0.455 0.627 0.610 0.625
## time_num -0.429 0.628 0.608 0.623 0.569
## T0SGENGendr -0.307 0.013 0.022 -0.005 0.016
## T0PDEMOMxPE -0.750 -0.028 0.037 0.049 0.042
## T0SY1Y2G_SG0: 0.352 -0.757 -0.510 -0.523 -0.477
## T0SY1Y2G_SG1: 0.354 -0.512 -0.758 -0.509 -0.464
## T0SY1Y2G_SG2: 0.353 -0.524 -0.508 -0.756 -0.475
## T0SY1Y2G_SG3: 0.320 -0.480 -0.464 -0.476 -0.758
## tim_nm T0SGEN T0PDEM T0SY1Y2G_SG0: T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG0
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr 0.001
## T0PDEMOMxPE 0.016 -0.075
## T0SY1Y2G_SG0: -0.839 -0.001 -0.003
## T0SY1Y2G_SG1: -0.816 -0.011 -0.013 0.685
## T0SY1Y2G_SG2: -0.835 0.004 -0.009 0.701 0.681
## T0SY1Y2G_SG3: -0.764 0.006 -0.005 0.641 0.623 0.638
Anova(model_long_wm, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: AWMADMRS_N
## F Df Df.res Pr(>F)
## (Intercept) 2.39 1 439 0.123
## T0SY1Y2G_StartGrade 39.76 4 976 < 0.0000000000000002 ***
## time_num 34.19 1 819 0.0000000072 ***
## T0SGENGender 0.20 1 303 0.654
## T0PDEMOMaxParentEducation 19.69 1 302 0.0000127757 ***
## T0SY1Y2G_StartGrade:time_num 2.50 4 808 0.041 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
prop_longitudinal_nonz_gather = prop_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
prop_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(prop_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
prop_longitudinal_nonz_gather$time_num = as.numeric(as.factor(prop_longitudinal_nonz_gather$time))
model_long_propreasoning <- lmer(PRWMPAE.asinR ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = prop_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_propreasoning)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: PRWMPAE.asinR ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender +
## T0PDEMOMaxParentEducation + (1 | SubjectID)
## Data: prop_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -284.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.219 -0.423 0.106 0.561 2.707
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.0218 0.148
## Residual 0.0299 0.173
## Number of obs: 1114, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -1.58523 0.08303 458.17760 -19.09
## T0SY1Y2G_StartGrade0 0.12969 0.05528 996.00294 2.35
## T0SY1Y2G_StartGrade1 0.19948 0.05699 993.31368 3.50
## T0SY1Y2G_StartGrade2 0.39792 0.05568 992.34524 7.15
## T0SY1Y2G_StartGrade3 0.49856 0.06011 983.40230 8.29
## time_num 0.09687 0.01498 863.25565 6.47
## T0SGENGender 0.00352 0.01991 298.72289 0.18
## T0PDEMOMaxParentEducation 0.01980 0.00420 301.67952 4.72
## T0SY1Y2G_StartGrade0:time_num -0.00987 0.01758 853.59774 -0.56
## T0SY1Y2G_StartGrade1:time_num 0.00128 0.01804 848.25646 0.07
## T0SY1Y2G_StartGrade2:time_num -0.03266 0.01763 852.14221 -1.85
## T0SY1Y2G_StartGrade3:time_num -0.06806 0.01892 842.73014 -3.60
## Pr(>|t|)
## (Intercept) < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade0 0.01916 *
## T0SY1Y2G_StartGrade1 0.00049 ***
## T0SY1Y2G_StartGrade2 0.00000000000171915 ***
## T0SY1Y2G_StartGrade3 0.00000000000000036 ***
## time_num 0.00000000016596876 ***
## T0SGENGender 0.85957
## T0PDEMOMaxParentEducation 0.00000363212220121 ***
## T0SY1Y2G_StartGrade0:time_num 0.57453
## T0SY1Y2G_StartGrade1:time_num 0.94345
## T0SY1Y2G_StartGrade2:time_num 0.06427 .
## T0SY1Y2G_StartGrade3:time_num 0.00034 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3
## T0SY1Y2G_StG0 -0.463
## T0SY1Y2G_StG1 -0.498 0.688
## T0SY1Y2G_StG2 -0.514 0.704 0.687
## T0SY1Y2G_StG3 -0.483 0.652 0.636 0.652
## time_num -0.449 0.646 0.628 0.644 0.596
## T0SGENGendr -0.299 0.007 0.004 -0.014 0.017
## T0PDEMOMxPE -0.749 -0.016 0.049 0.065 0.054
## T0SY1Y2G_SG0: 0.377 -0.766 -0.535 -0.548 -0.507
## T0SY1Y2G_SG1: 0.371 -0.536 -0.765 -0.534 -0.494
## T0SY1Y2G_SG2: 0.380 -0.549 -0.534 -0.764 -0.506
## T0SY1Y2G_SG3: 0.352 -0.511 -0.497 -0.509 -0.760
## tim_nm T0SGEN T0PDEM T0SY1Y2G_SG0: T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG0
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr -0.009
## T0PDEMOMxPE 0.028 -0.071
## T0SY1Y2G_SG0: -0.852 0.004 -0.015
## T0SY1Y2G_SG1: -0.830 0.009 -0.022 0.707
## T0SY1Y2G_SG2: -0.850 0.014 -0.025 0.724 0.705
## T0SY1Y2G_SG3: -0.791 0.006 -0.017 0.674 0.657 0.672
Anova(model_long_propreasoning, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: PRWMPAE.asinR
## F Df Df.res Pr(>F)
## (Intercept) 364.49 1 463 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade 28.29 4 983 < 0.0000000000000002 ***
## time_num 41.81 1 866 0.00000000017 ***
## T0SGENGender 0.03 1 302 0.86
## T0PDEMOMaxParentEducation 22.27 1 305 0.00000361770 ***
## T0SY1Y2G_StartGrade:time_num 6.60 4 831 0.00003126854 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_propreasoning, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio
## -1 0.0969 0.01498 866 0.06747 0.1263 6.466
## 0 0.0870 0.00921 830 0.06892 0.1051 9.445
## 1 0.0982 0.01006 818 0.07840 0.1179 9.755
## 2 0.0642 0.00930 827 0.04596 0.0825 6.905
## 3 0.0288 0.01157 811 0.00609 0.0515 2.489
## p.value
## <.0001
## <.0001
## <.0001
## <.0001
## 0.0130
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade0 0.00987 0.0176 857 -0.02464
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade1 -0.00128 0.0180 851 -0.03670
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade2 0.03266 0.0176 855 -0.00195
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade3 0.06806 0.0189 846 0.03091
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade1 -0.01115 0.0136 823 -0.03793
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade2 0.02279 0.0131 828 -0.00290
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade3 0.05819 0.0148 818 0.02916
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade2 0.03394 0.0137 822 0.00705
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade3 0.06934 0.0153 814 0.03924
## T0SY1Y2G_StartGrade2 - T0SY1Y2G_StartGrade3 0.03540 0.0148 817 0.00626
## upper.CL t.ratio p.value
## 0.0444 0.561 0.5746
## 0.0341 -0.071 0.9435
## 0.0673 1.852 0.0643
## 0.1052 3.596 0.0003
## 0.0156 -0.818 0.4139
## 0.0485 1.741 0.0820
## 0.0872 3.935 0.0001
## 0.0608 2.477 0.0134
## 0.0994 4.522 <.0001
## 0.0645 2.385 0.0173
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
multiplot(graph_workingmemory_longitudinal_nonz, graph_prop_longitudinal_nonz,graph_mentalrotation2_longitudinal_nonz, graph_panamath_longitudinal_nz, graph_mentalrotation_longitudinal_nonz, cols = 3)
EC_longitudinal_nonz = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
#"T1WJCW","T2WJCW","T3WJCW","T4WJCW")
"T1WJCRS","T2WJCRS","T3WJCRS","T4WJCRS")
EC_longitudinal_nonz = EC_longitudinal_nonz[which(EC_longitudinal_nonz$SubjectID %in% T1_complete_class$SubjectID),]
EC_longitudinal_nonz_gather = gather(EC_longitudinal_nonz, "time", "WJCRS", -c("SubjectID","T0SY1Y2G_StartGrade"))
EC_longitudinal_nonz_gather_summary = summarySE(EC_longitudinal_nonz_gather, "WJCRS", c("time","T0SY1Y2G_StartGrade"), na.rm = T)
EC_longitudinal_nonz_gather_summary$time = as.factor(EC_longitudinal_nonz_gather_summary$time)
EC_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade = as.factor(EC_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade)
EC_longitudinal_nonz_gather_summary$task = "Exact Calculation for Older Children"
graph_EC_longitudinal_nonz = ggplot(EC_longitudinal_nonz_gather_summary, aes(x = time, y = WJCRS, group = T0SY1Y2G_StartGrade)) +
#geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (WJCRS-se), ymax = (WJCRS+se), fill = T0SY1Y2G_StartGrade),
alpha = .15, color = NA)+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
scale_color_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
scale_fill_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
facet_grid(.~task)+
# scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Exact Calculation for Older Children")+
xlab("Time")+
#scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="bottom",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_EC_longitudinal_nonz
EC_longitudinal_nonz_gather = EC_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
EC_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(EC_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
EC_longitudinal_nonz_gather$time_num = as.numeric(as.factor(EC_longitudinal_nonz_gather$time))
model_long_ecolder <- lmer(WJCRS ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = EC_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_ecolder)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## WJCRS ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: EC_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 3977
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -6.336 -0.473 0.017 0.543 2.401
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 7.60 2.76
## Residual 5.86 2.42
## Number of obs: 775, groups: SubjectID, 271
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -0.76118 2.06141 708.72964 -0.37
## T0SY1Y2G_StartGrade1 6.10779 1.69688 658.31422 3.60
## T0SY1Y2G_StartGrade2 10.51143 1.69248 652.69883 6.21
## T0SY1Y2G_StartGrade3 15.84014 1.73107 672.30390 9.15
## time_num 2.86043 0.45114 569.19170 6.34
## T0SGENGender 0.00705 0.38674 253.22165 0.02
## T0PDEMOMaxParentEducation 0.18666 0.08263 254.74192 2.26
## T0SY1Y2G_StartGrade1:time_num -0.13167 0.47312 563.51681 -0.28
## T0SY1Y2G_StartGrade2:time_num 0.12989 0.47188 564.25345 0.28
## T0SY1Y2G_StartGrade3:time_num -0.12253 0.48346 561.83598 -0.25
## Pr(>|t|)
## (Intercept) 0.71205
## T0SY1Y2G_StartGrade1 0.00034 ***
## T0SY1Y2G_StartGrade2 0.00000000094 ***
## T0SY1Y2G_StartGrade3 < 0.0000000000000002 ***
## time_num 0.00000000047 ***
## T0SGENGender 0.98546
## T0PDEMOMaxParentEducation 0.02473 *
## T0SY1Y2G_StartGrade1:time_num 0.78088
## T0SY1Y2G_StartGrade2:time_num 0.78322
## T0SY1Y2G_StartGrade3:time_num 0.80002
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3 tim_nm T0SGEN
## T0SY1Y2G_StG1 -0.751
## T0SY1Y2G_StG2 -0.751 0.913
## T0SY1Y2G_StG3 -0.742 0.892 0.895
## time_num -0.743 0.927 0.929 0.908
## T0SGENGendr -0.203 -0.004 -0.017 0.005 -0.009
## T0PDEMOMxPE -0.559 0.004 0.007 0.010 -0.030 -0.109
## T0SY1Y2G_SG1: 0.710 -0.945 -0.886 -0.866 -0.953 0.006
## T0SY1Y2G_SG2: 0.707 -0.886 -0.948 -0.868 -0.956 0.012
## T0SY1Y2G_SG3: 0.695 -0.865 -0.867 -0.935 -0.933 0.004
## T0PDEM T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr
## T0PDEMOMxPE
## T0SY1Y2G_SG1: 0.027
## T0SY1Y2G_SG2: 0.032 0.912
## T0SY1Y2G_SG3: 0.027 0.890 0.892
Anova(model_long_ecolder, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: WJCRS
## F Df Df.res Pr(>F)
## (Intercept) 0.14 1 711 0.712
## T0SY1Y2G_StartGrade 62.79 3 701 < 0.0000000000000002 ***
## time_num 40.16 1 576 0.00000000047 ***
## T0SGENGender 0.00 1 262 0.985
## T0PDEMOMaxParentEducation 5.10 1 263 0.025 *
## T0SY1Y2G_StartGrade:time_num 0.71 3 537 0.546
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
# em=emtrends(model_long_propreasoning, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
# summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
EC2_longitudinal_nonz = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1PKC_PA","T2PKC_PA","T3PKC_PA","T4PKC_PA")
EC2_longitudinal_nonz = EC2_longitudinal_nonz[which(EC2_longitudinal_nonz$SubjectID %in% T1_complete_class$SubjectID),]
EC2_longitudinal_nonz_gather = gather(EC2_longitudinal_nonz, "time", "PKC_PA", -c("SubjectID","T0SY1Y2G_StartGrade"))
EC2_longitudinal_nonz_gather_summary = summarySE(EC2_longitudinal_nonz_gather, "PKC_PA", c("time","T0SY1Y2G_StartGrade"), na.rm = T)
EC2_longitudinal_nonz_gather_summary$time = as.factor(EC2_longitudinal_nonz_gather_summary$time)
EC2_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade = as.factor(EC2_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade)
EC2_longitudinal_nonz_gather_summary$task = "Exact Calculation for Younger Children"
graph_EC2_longitudinal_nonz = ggplot(EC2_longitudinal_nonz_gather_summary, aes(x = time, y = PKC_PA, group = T0SY1Y2G_StartGrade)) +
#geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (PKC_PA-se), ymax = (PKC_PA+se), fill = T0SY1Y2G_StartGrade),
alpha = .15, color = NA)+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
scale_color_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
scale_fill_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
facet_grid(.~task)+
scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Exact Calculation for Younger Children")+
xlab("Time")+
#scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="bottom",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_EC2_longitudinal_nonz
multiplot(graph_EC2_longitudinal_nonz, graph_EC_longitudinal_nonz, cols =2)
EC2_longitudinal_nonz_gather = EC2_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
EC2_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(EC2_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
EC2_longitudinal_nonz_gather$time_num = as.numeric(as.factor(EC2_longitudinal_nonz_gather$time))
model_long_ecyounger <- lmer(PKC_PA ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = EC2_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_ecyounger)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## PKC_PA ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: EC2_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -11.4
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2859 -0.5355 -0.0862 0.4835 2.8767
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.0477 0.218
## Residual 0.0247 0.157
## Number of obs: 272, groups: SubjectID, 118
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -0.33840 0.16367 114.41809 -2.07
## T0SY1Y2G_StartGrade0 0.15107 0.06919 257.49519 2.18
## time_num 0.07694 0.01344 153.66516 5.72
## T0SGENGender 0.00162 0.04506 108.37372 0.04
## T0PDEMOMaxParentEducation 0.03023 0.00915 106.44731 3.30
## T0SY1Y2G_StartGrade0:time_num 0.03883 0.02901 153.51151 1.34
## Pr(>|t|)
## (Intercept) 0.0409 *
## T0SY1Y2G_StartGrade0 0.0299 *
## time_num 0.000000053 ***
## T0SGENGender 0.9714
## T0PDEMOMaxParentEducation 0.0013 **
## T0SY1Y2G_StartGrade0:time_num 0.1827
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 tim_nm T0SGEN T0PDEM
## T0SY1Y2G_StG0 -0.194
## time_num -0.199 0.456
## T0SGENGendr -0.441 0.021 0.011
## T0PDEMOMxPE -0.857 -0.051 0.002 0.035
## T0SY1Y2G_SG0: 0.087 -0.699 -0.463 -0.004 0.005
Anova(model_long_ecyounger, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: PKC_PA
## F Df Df.res Pr(>F)
## (Intercept) 4.27 1 119 0.0409 *
## T0SY1Y2G_StartGrade 4.77 1 258 0.0299 *
## time_num 32.73 1 158 0.000000051 ***
## T0SGENGender 0.00 1 113 0.9714
## T0PDEMOMaxParentEducation 10.92 1 111 0.0013 **
## T0SY1Y2G_StartGrade:time_num 1.79 1 158 0.1829
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
AC_longitudinal_nonz = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1ASWmac","T2ASWmac","T3ASWmac","T4ASWmac")
AC_longitudinal_nonz = AC_longitudinal_nonz[which(AC_longitudinal_nonz$SubjectID %in% T1_complete_class$SubjectID),]
AC_longitudinal_nonz_gather = gather(AC_longitudinal_nonz, "time", "ASWmac", -c("SubjectID","T0SY1Y2G_StartGrade"))
AC_longitudinal_nonz_gather_summary = summarySE(AC_longitudinal_nonz_gather, "ASWmac", c("time","T0SY1Y2G_StartGrade"), na.rm = T)
AC_longitudinal_nonz_gather_summary$time = as.factor(AC_longitudinal_nonz_gather_summary$time)
AC_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade = as.factor(AC_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade)
AC_longitudinal_nonz_gather_summary$task = "Approximate Calculation"
graph_AC_longitudinal_nonz = ggplot(AC_longitudinal_nonz_gather_summary, aes(x = time, y = ASWmac, group = T0SY1Y2G_StartGrade)) +
#geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (ASWmac-se), ymax = (ASWmac+se), fill = T0SY1Y2G_StartGrade),
alpha = .15, color = NA)+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
scale_color_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
scale_fill_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
facet_grid(.~task)+
# scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Approximate Calculation")+
xlab("Time")+
#scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="bottom",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_AC_longitudinal_nonz
AC_longitudinal_nonz_gather = AC_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
AC_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(AC_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
AC_longitudinal_nonz_gather$time_num = as.numeric(as.factor(AC_longitudinal_nonz_gather$time))
model_long_ac <- lmer(ASWmac ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = AC_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_ac)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## ASWmac ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: AC_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: -1282
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6508 -0.6112 0.0126 0.6269 2.7489
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.00503 0.071
## Residual 0.01280 0.113
## Number of obs: 1069, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 0.36127 0.04721 527.33765 7.65
## T0SY1Y2G_StartGrade0 0.12744 0.03465 1044.11352 3.68
## T0SY1Y2G_StartGrade1 0.07647 0.03536 1041.60629 2.16
## T0SY1Y2G_StartGrade2 0.20878 0.03479 1042.96511 6.00
## T0SY1Y2G_StartGrade3 0.16817 0.03789 1043.10016 4.44
## time_num 0.04089 0.01013 864.75824 4.04
## T0SGENGender 0.01738 0.01075 295.21608 1.62
## T0PDEMOMaxParentEducation 0.00858 0.00227 299.31535 3.78
## T0SY1Y2G_StartGrade0:time_num -0.02936 0.01190 848.89667 -2.47
## T0SY1Y2G_StartGrade1:time_num 0.00222 0.01207 838.91950 0.18
## T0SY1Y2G_StartGrade2:time_num -0.04938 0.01188 843.72200 -4.16
## T0SY1Y2G_StartGrade3:time_num -0.02850 0.01284 834.43334 -2.22
## Pr(>|t|)
## (Intercept) 0.000000000000094 ***
## T0SY1Y2G_StartGrade0 0.00025 ***
## T0SY1Y2G_StartGrade1 0.03082 *
## T0SY1Y2G_StartGrade2 0.000000002693389 ***
## T0SY1Y2G_StartGrade3 0.000010039430217 ***
## time_num 0.000058943166806 ***
## T0SGENGender 0.10714
## T0PDEMOMaxParentEducation 0.00019 ***
## T0SY1Y2G_StartGrade0:time_num 0.01383 *
## T0SY1Y2G_StartGrade1:time_num 0.85403
## T0SY1Y2G_StartGrade2:time_num 0.000035789501215 ***
## T0SY1Y2G_StartGrade3:time_num 0.02674 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3
## T0SY1Y2G_StG0 -0.516
## T0SY1Y2G_StG1 -0.549 0.705
## T0SY1Y2G_StG2 -0.555 0.717 0.705
## T0SY1Y2G_StG3 -0.512 0.658 0.647 0.658
## time_num -0.526 0.706 0.693 0.704 0.646
## T0SGENGendr -0.283 0.006 0.006 -0.009 0.007
## T0PDEMOMxPE -0.707 -0.021 0.039 0.042 0.035
## T0SY1Y2G_SG0: 0.442 -0.830 -0.589 -0.599 -0.550
## T0SY1Y2G_SG1: 0.446 -0.592 -0.827 -0.591 -0.542
## T0SY1Y2G_SG2: 0.443 -0.602 -0.590 -0.827 -0.551
## T0SY1Y2G_SG3: 0.405 -0.557 -0.546 -0.555 -0.828
## tim_nm T0SGEN T0PDEM T0SY1Y2G_SG0: T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG0
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr -0.007
## T0PDEMOMxPE 0.013 -0.070
## T0SY1Y2G_SG0: -0.851 0.001 -0.001
## T0SY1Y2G_SG1: -0.839 0.002 -0.017 0.714
## T0SY1Y2G_SG2: -0.852 0.008 -0.006 0.725 0.715
## T0SY1Y2G_SG3: -0.789 0.013 -0.001 0.671 0.661 0.672
Anova(model_long_ac, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: ASWmac
## F Df Df.res Pr(>F)
## (Intercept) 58.54 1 534 0.000000000000094 ***
## T0SY1Y2G_StartGrade 11.91 4 1041 0.000000001836696 ***
## time_num 16.28 1 868 0.000059434699282 ***
## T0SGENGender 2.61 1 300 0.10716
## T0PDEMOMaxParentEducation 14.30 1 304 0.00019 ***
## T0SY1Y2G_StartGrade:time_num 9.81 4 809 0.000000094180766 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_ac, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio
## -1 0.04089 0.01013 868 0.021000 0.06078 4.035
## 0 0.01153 0.00626 810 -0.000756 0.02381 1.842
## 1 0.04311 0.00657 782 0.030212 0.05601 6.560
## 2 -0.00849 0.00622 791 -0.020702 0.00371 -1.366
## 3 0.01239 0.00790 788 -0.003117 0.02789 1.568
## p.value
## 0.0001
## 0.0658
## <.0001
## 0.1724
## 0.1172
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade0 0.02936 0.01191 853 0.00599
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade1 -0.00222 0.01208 843 -0.02593
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade2 0.04938 0.01189 848 0.02605
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade3 0.02850 0.01285 838 0.00328
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade1 -0.03159 0.00907 795 -0.04940
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade2 0.02002 0.00882 801 0.00270
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade3 -0.00086 0.01008 797 -0.02064
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade2 0.05161 0.00905 786 0.03384
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade3 0.03072 0.01028 785 0.01055
## T0SY1Y2G_StartGrade2 - T0SY1Y2G_StartGrade3 -0.02088 0.01005 790 -0.04061
## upper.CL t.ratio p.value
## 0.05274 2.466 0.0139
## 0.02149 -0.184 0.8541
## 0.07272 4.154 <.0001
## 0.05372 2.218 0.0268
## -0.01377 -3.480 0.0005
## 0.03734 2.270 0.0235
## 0.01892 -0.085 0.9319
## 0.06937 5.703 <.0001
## 0.05090 2.990 0.0029
## -0.00115 -2.077 0.0381
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
LetterWordID_longitudinal_nonz = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1LWIDWS","T2LWIDWS","T3LWIDWS","T4LWIDWS")
LetterWordID_longitudinal_nonz = LetterWordID_longitudinal_nonz[which(LetterWordID_longitudinal_nonz$SubjectID %in% T1_complete_class$SubjectID),]
LetterWordID_longitudinal_nonz_gather = gather(LetterWordID_longitudinal_nonz, "time", "LWIDWS", -c("SubjectID","T0SY1Y2G_StartGrade"))
LetterWordID_longitudinal_nonz_gather_summary = summarySE(LetterWordID_longitudinal_nonz_gather, "LWIDWS", c("time","T0SY1Y2G_StartGrade"), na.rm = T)
LetterWordID_longitudinal_nonz_gather_summary$time = as.factor(LetterWordID_longitudinal_nonz_gather_summary$time)
LetterWordID_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade = as.factor(LetterWordID_longitudinal_nonz_gather_summary$T0SY1Y2G_StartGrade)
LetterWordID_longitudinal_nonz_gather_summary$task = "Letter Word ID"
graph_LetterWordID_longitudinal_nonz = ggplot(LetterWordID_longitudinal_nonz_gather_summary, aes(x = time, y = LWIDWS, group = T0SY1Y2G_StartGrade)) +
#geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (LWIDWS-se), ymax = (LWIDWS+se), fill = T0SY1Y2G_StartGrade),
alpha = .15, color = NA)+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
scale_color_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
scale_fill_manual(values = c("#74a9cf","#3690c0","#0570b0","#045a8d","#023858"))+
geom_line(aes(color = T0SY1Y2G_StartGrade)) +
facet_grid(.~task)+
# scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Letter Word ID")+
xlab("Time")+
#scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="bottom",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_LetterWordID_longitudinal_nonz
LetterWordID_longitudinal_nonz_gather = LetterWordID_longitudinal_nonz_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
LetterWordID_longitudinal_nonz_gather$T0SY1Y2G_StartGrade = as.factor(as.numeric(LetterWordID_longitudinal_nonz_gather$T0SY1Y2G_StartGrade))
LetterWordID_longitudinal_nonz_gather$time_num = as.numeric(as.factor(LetterWordID_longitudinal_nonz_gather$time))
model_long_lwid <- lmer(LWIDWS ~ T0SY1Y2G_StartGrade*time_num + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = LetterWordID_longitudinal_nonz_gather, control = lmerControl(optimizer = "bobyqa"))
summary(model_long_lwid)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## LWIDWS ~ T0SY1Y2G_StartGrade * time_num + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: LetterWordID_longitudinal_nonz_gather
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 9040
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -5.192 -0.505 -0.008 0.511 3.094
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 566 23.8
## Residual 135 11.6
## Number of obs: 1064, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) 290.966 10.719 337.356 27.15
## T0SY1Y2G_StartGrade0 29.103 5.717 509.004 5.09
## T0SY1Y2G_StartGrade1 75.337 5.927 511.054 12.71
## T0SY1Y2G_StartGrade2 123.973 5.796 513.980 21.39
## T0SY1Y2G_StartGrade3 148.932 6.328 512.242 23.54
## time_num 13.599 1.001 762.255 13.59
## T0SGENGender 1.904 2.816 305.180 0.68
## T0PDEMOMaxParentEducation 2.121 0.592 304.929 3.58
## T0SY1Y2G_StartGrade0:time_num 6.961 1.181 760.951 5.89
## T0SY1Y2G_StartGrade1:time_num 3.310 1.217 760.184 2.72
## T0SY1Y2G_StartGrade2:time_num -4.939 1.200 761.384 -4.12
## T0SY1Y2G_StartGrade3:time_num -8.307 1.303 760.906 -6.38
## Pr(>|t|)
## (Intercept) < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade0 0.00000050240 ***
## T0SY1Y2G_StartGrade1 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade2 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade3 < 0.0000000000000002 ***
## time_num < 0.0000000000000002 ***
## T0SGENGender 0.49944
## T0PDEMOMaxParentEducation 0.00039 ***
## T0SY1Y2G_StartGrade0:time_num 0.00000000573 ***
## T0SY1Y2G_StartGrade1:time_num 0.00666 **
## T0SY1Y2G_StartGrade2:time_num 0.00004281987 ***
## T0SY1Y2G_StartGrade3:time_num 0.00000000031 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) T0SY1Y2G_StG0 T0SY1Y2G_StG1 T0SY1Y2G_StG2 T0SY1Y2G_StG3
## T0SY1Y2G_StG0 -0.351
## T0SY1Y2G_StG1 -0.410 0.678
## T0SY1Y2G_StG2 -0.419 0.693 0.674
## T0SY1Y2G_StG3 -0.396 0.635 0.618 0.631
## time_num -0.224 0.412 0.398 0.407 0.373
## T0SGENGendr -0.329 0.017 0.018 -0.007 0.027
## T0PDEMOMxPE -0.805 -0.038 0.048 0.061 0.054
## T0SY1Y2G_SG0: 0.187 -0.490 -0.337 -0.345 -0.316
## T0SY1Y2G_SG1: 0.189 -0.339 -0.493 -0.335 -0.307
## T0SY1Y2G_SG2: 0.185 -0.344 -0.332 -0.494 -0.311
## T0SY1Y2G_SG3: 0.169 -0.317 -0.306 -0.313 -0.492
## tim_nm T0SGEN T0PDEM T0SY1Y2G_SG0: T0SY1Y2G_SG1: T0SY1Y2G_SG2:
## T0SY1Y2G_StG0
## T0SY1Y2G_StG1
## T0SY1Y2G_StG2
## T0SY1Y2G_StG3
## time_num
## T0SGENGendr 0.001
## T0PDEMOMxPE 0.004 -0.073
## T0SY1Y2G_SG0: -0.847 -0.002 0.000
## T0SY1Y2G_SG1: -0.823 -0.006 -0.006 0.697
## T0SY1Y2G_SG2: -0.834 0.003 -0.003 0.707 0.686
## T0SY1Y2G_SG3: -0.768 -0.001 0.000 0.651 0.632 0.641
Anova(model_long_lwid, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: LWIDWS
## F Df Df.res Pr(>F)
## (Intercept) 736.86 1 336 < 0.0000000000000002 ***
## T0SY1Y2G_StartGrade 249.75 4 509 < 0.0000000000000002 ***
## time_num 184.55 1 761 < 0.0000000000000002 ***
## T0SGENGender 0.46 1 304 0.49945
## T0PDEMOMaxParentEducation 12.85 1 304 0.00039 ***
## T0SY1Y2G_StartGrade:time_num 75.15 4 758 < 0.0000000000000002 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(model_long_lwid, pairwise ~ T0SY1Y2G_StartGrade, var="time_num", mult.name = "T0SY1Y2G_StartGrade")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## T0SY1Y2G_StartGrade time_num.trend SE df lower.CL upper.CL t.ratio p.value
## -1 13.60 1.001 761 11.63 15.56 13.590 <.0001
## 0 20.56 0.628 756 19.33 21.79 32.760 <.0001
## 1 16.91 0.692 755 15.55 18.27 24.450 <.0001
## 2 8.66 0.662 758 7.36 9.96 13.080 <.0001
## 3 5.29 0.834 758 3.65 6.93 6.340 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade0 -6.96 1.181 760 -9.28
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade1 -3.31 1.217 759 -5.70
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade2 4.94 1.200 760 2.58
## (T0SY1Y2G_StartGrade-1) - T0SY1Y2G_StartGrade3 8.31 1.303 760 5.75
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade1 3.65 0.934 756 1.82
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade2 11.90 0.912 757 10.11
## T0SY1Y2G_StartGrade0 - T0SY1Y2G_StartGrade3 15.27 1.044 757 13.22
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade2 8.25 0.957 756 6.37
## T0SY1Y2G_StartGrade1 - T0SY1Y2G_StartGrade3 11.62 1.083 757 9.49
## T0SY1Y2G_StartGrade2 - T0SY1Y2G_StartGrade3 3.37 1.065 758 1.28
## upper.CL t.ratio p.value
## -4.641 -5.891 <.0001
## -0.922 -2.721 0.0067
## 7.295 4.115 <.0001
## 10.865 6.376 <.0001
## 5.484 3.909 0.0001
## 13.690 13.045 <.0001
## 17.317 14.628 <.0001
## 10.128 8.616 <.0001
## 13.744 10.723 <.0001
## 5.459 3.164 0.0016
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
names(spatial_profiles_database_atleast2tp)
## [1] "...1" "...2"
## [3] "SubjectID" "T0SY1Y2G_Y1Grade"
## [5] "T0SY1Y2G_Y2Grade" "T0SY1Y2G_StartGrade"
## [7] "T0SGENGender" "T0PDEMOMaxParentEducation"
## [9] "T0PDEMOIncomeCode" "T0PDEMOIncomeMid"
## [11] "T0PDEMORaceHispanic" "T1AWMADMSS_N"
## [13] "T2AWMADMSS_N" "T3AWMADMSS_N"
## [15] "T4AWMADMSS_N" "T1AWMADMRS_N"
## [17] "T2AWMADMRS_N" "T3AWMADMRS_N"
## [19] "T4AWMADMRS_N" "T1GMPA"
## [21] "T2GMPA" "T3GMPA"
## [23] "T4GMPA" "T1GPPA"
## [25] "T2GPPA" "T3GPPA"
## [27] "T4GPPA" "T1PRWMPAE"
## [29] "T2PRWMPAE" "T3PRWMPAE"
## [31] "T4PRWMPAE" "T1PMPA"
## [33] "T2PMPA" "T3PMPA"
## [35] "T4PMPA" "T1NL10W_PAE"
## [37] "T2NL10W_PAE" "T3NL10W_PAE"
## [39] "T4NL10W_PAE" "T1NL100W_PAE"
## [41] "T1NL1000W_PAE" "T2NL100W_PAE"
## [43] "T2NL1000W_PAE" "T3NL100W_PAE"
## [45] "T3NL1000W_PAE" "T4NL100W_PAE"
## [47] "T4NL1000W_PAE" "T1PKC_PA"
## [49] "T2PKC_PA" "T3PKC_PA"
## [51] "T4PKC_PA" "T1WJCW"
## [53] "T2WJCW" "T3WJCW"
## [55] "T4WJCW" "T1WJCRS"
## [57] "T2WJCRS" "T3WJCRS"
## [59] "T4WJCRS" "T1ASWmac"
## [61] "T2ASWmac" "T3ASWmac"
## [63] "T4ASWmac" "T1LWIDWS"
## [65] "T2LWIDWS" "T3LWIDWS"
## [67] "T4LWIDWS" "GMPA_na"
## [69] "GPPA_na" "AWMADMSS_na"
## [71] "PRWM_na" "PKC_na"
## [73] "WJCW_na" "ASW_na"
## [75] "PMPA_na" "LWID_na"
## [77] "t1" "t2"
## [79] "t3" "t4"
## [81] "T1PRWMPAE.asin" "T2PRWMPAE.asin"
## [83] "T3PRWMPAE.asin" "T4PRWMPAE.asin"
## [85] "T1AWMADMSS_N.s" "T2AWMADMSS_N.s"
## [87] "T1GMPA.s" "T2GMPA.s"
## [89] "T1GPPA.s" "T2GPPA.s"
## [91] "t1.s" "t2.s"
## [93] "T1PRWMPAE.asin.s" "T2PRWMPAE.asin.s"
## [95] "T1PMPA.s" "T2PMPA.s"
## [97] "T1WJCW.s" "T2WJCW.s"
## [99] "T1PKC_PA.s" "T2PKC_PA.s"
## [101] "T1ASWmac.s" "T2ASWmac.s"
## [103] "T1LWIDWS.s" "T2LWIDWS.s"
## [105] "T3AWMADMSS_N.s" "T4AWMADMSS_N.s"
## [107] "T3GPPA.s" "T4GPPA.s"
## [109] "T3GMPA.s" "T4GMPA.s"
## [111] "t3.s" "t4.s"
## [113] "T3PRWMPAE.asin.s" "T4PRWMPAE.asin.s"
## [115] "T3PMPA.s" "T4PMPA.s"
## [117] "T3WJCW.s" "T4WJCW.s"
## [119] "T3PKC_PA.s" "T4PKC_PA.s"
## [121] "T3ASWmac.s" "T4ASWmac.s"
## [123] "T3LWIDWS.s" "T4LWIDWS.s"
## [125] "T1PRWMPAE.asinR" "T2PRWMPAE.asinR"
## [127] "T3PRWMPAE.asinR" "T4PRWMPAE.asinR"
## [129] "T1PRWMPAE.asin.sR" "T2PRWMPAE.asin.sR"
## [131] "T3PRWMPAE.asin.sR" "T4PRWMPAE.asin.sR"
## [133] "T1MR.s" "T2MR.s"
## [135] "T3MR.s" "T4MR.s"
## [137] "T1EC.s" "T2EC.s"
## [139] "T3EC.s" "T4EC.s"
panamath_longitudinal_profiles = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1PMPA.s","T2PMPA.s","T3PMPA.s","T4PMPA.s")
panamath_longitudinal_profiles_gather = gather(panamath_longitudinal_profiles, "time", "PMPA.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
panamath_longitudinal_profiles_gather$SubjectID = as.factor(as.character(panamath_longitudinal_profiles_gather$SubjectID))
panamath_longitudinal_profiles_gather = panamath_longitudinal_profiles_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
panamath_longitudinal_gather_summary = summarySE(panamath_longitudinal_profiles_gather, "PMPA.s", c("time","Class"), na.rm = T)
panamath_longitudinal_gather_summary$time = as.factor(panamath_longitudinal_gather_summary$time)
panamath_longitudinal_gather_summary$Class = as.factor(panamath_longitudinal_gather_summary$Class)
panamath_longitudinal_gather_summary$task = "Panamath"
panamath_longitudinal_gather1 = subset(panamath_longitudinal_profiles_gather, !is.na(PMPA.s))
panamath_longitudinal_gather1$timenum = as.numeric(as.factor(panamath_longitudinal_gather1$time))
panamath_longitudinal_gather1$Class = as.factor(panamath_longitudinal_gather1$Class)
panamath_longitudinal_gather1 = panamath_longitudinal_gather1 %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
panamath_model <- lmer(PMPA.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
(1 + timenum|SubjectID), data = subset(panamath_longitudinal_gather1, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(panamath_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## PMPA.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 + timenum | SubjectID)
## Data: subset(panamath_longitudinal_gather1, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2460
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.943 -0.475 0.121 0.584 2.488
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## SubjectID (Intercept) 0.1241 0.352
## timenum 0.0134 0.116 0.53
## Residual 0.4679 0.684
## Number of obs: 995, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -2.3130 0.3177 374.4580 -7.28
## Classlowmr 2.0792 0.2038 261.2691 10.20
## Classhigh 2.5788 0.2013 257.2196 12.81
## timenum 0.5031 0.0654 241.2779 7.70
## T0SGENGender -0.0654 0.0770 295.2227 -0.85
## T0PDEMOMaxParentEducation 0.0185 0.0160 295.2603 1.15
## Classlowmr:timenum -0.4888 0.0730 248.2053 -6.70
## Classhigh:timenum -0.5774 0.0722 246.0346 -8.00
## Pr(>|t|)
## (Intercept) 0.00000000000197 ***
## Classlowmr < 0.0000000000000002 ***
## Classhigh < 0.0000000000000002 ***
## timenum 0.00000000000035 ***
## T0SGENGender 0.40
## T0PDEMOMaxParentEducation 0.25
## Classlowmr:timenum 0.00000000014238 ***
## Classhigh:timenum 0.00000000000005 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Clsslw Clsshg timenm T0SGEN T0PDEM Clssl:
## Classlowmr -0.560
## Classhigh -0.506 0.821
## timenum -0.431 0.663 0.670
## T0SGENGendr -0.362 0.108 0.081 -0.006
## T0PDEMOMxPE -0.711 0.001 -0.065 0.011 -0.060
## Clsslwmr:tm 0.381 -0.737 -0.600 -0.895 0.010 -0.005
## Clsshgh:tmn 0.388 -0.600 -0.738 -0.905 0.002 -0.006 0.810
Anova(panamath_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: PMPA.s
## F Df Df.res Pr(>F)
## (Intercept) 52.72 1 391 0.00000000000209 ***
## Class 82.09 2 274 < 0.0000000000000002 ***
## timenum 59.18 1 241 0.00000000000037 ***
## T0SGENGender 0.72 1 302 0.40
## T0PDEMOMaxParentEducation 1.32 1 302 0.25
## Class:timenum 31.99 2 258 0.00000000000039 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(ggpredict(panamath_model, terms=c("timenum","Class")))
emeans=emmeans(panamath_model, pairwise ~ Class|timenum, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## timenum = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -1.6333 0.1421 281 -1.9131 -1.3535 -11.490 <.0001
## lowmr -0.0429 0.0678 297 -0.1762 0.0905 -0.633 0.5274
## high 0.3681 0.0630 283 0.2441 0.4921 5.843 <.0001
##
## timenum = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -1.1301 0.1239 295 -1.3740 -0.8862 -9.120 <.0001
## lowmr -0.0285 0.0589 305 -0.1444 0.0874 -0.484 0.6290
## high 0.2938 0.0548 299 0.1861 0.4016 5.366 <.0001
##
## timenum = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.6270 0.1381 282 -0.8987 -0.3552 -4.541 <.0001
## lowmr -0.0141 0.0669 308 -0.1457 0.1174 -0.211 0.8329
## high 0.2195 0.0626 299 0.0963 0.3427 3.506 0.0005
##
## timenum = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.1238 0.1770 268 -0.4723 0.2246 -0.700 0.4847
## lowmr 0.0003 0.0871 297 -0.1712 0.1717 0.003 0.9977
## high 0.1452 0.0821 286 -0.0163 0.3067 1.769 0.0779
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## timenum = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -1.590 0.1580 286 -1.901 -1.2794 -10.065 <.0001
## low - high -2.001 0.1559 282 -2.308 -1.6944 -12.835 <.0001
## lowmr - high -0.411 0.0928 291 -0.594 -0.2283 -4.427 <.0001
##
## timenum = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -1.102 0.1379 298 -1.373 -0.8304 -7.992 <.0001
## low - high -1.424 0.1360 295 -1.692 -1.1563 -10.471 <.0001
## lowmr - high -0.322 0.0808 301 -0.481 -0.1633 -3.988 0.0001
##
## timenum = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.613 0.1540 289 -0.916 -0.3098 -3.980 0.0001
## low - high -0.847 0.1520 285 -1.146 -0.5472 -5.567 <.0001
## lowmr - high -0.234 0.0920 305 -0.415 -0.0527 -2.541 0.0116
##
## timenum = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.124 0.1977 276 -0.513 0.2651 -0.628 0.5307
## low - high -0.269 0.1954 272 -0.654 0.1157 -1.377 0.1697
## lowmr - high -0.145 0.1200 294 -0.381 0.0912 -1.208 0.2280
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
em=emtrends(panamath_model, pairwise ~ Class, var="timenum", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class timenum.trend SE df lower.CL upper.CL t.ratio p.value
## low 0.5031 0.0654 241 0.3743 0.6320 7.693 <.0001
## lowmr 0.0144 0.0326 278 -0.0497 0.0785 0.441 0.6592
## high -0.0743 0.0308 268 -0.1349 -0.0137 -2.415 0.0164
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.4888 0.0731 247 0.3449 0.633 6.690 <.0001
## low - high 0.5774 0.0723 245 0.4351 0.720 7.990 <.0001
## lowmr - high 0.0887 0.0448 273 0.0005 0.177 1.980 0.0488
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
graph_panamath_longitudinal = ggplot(subset(panamath_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = PMPA.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (PMPA.s-se), ymax = (PMPA.s+se), fill = Class),
alpha = .15, color = NA)+
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme_bw()+
ylab("Panamath")+
xlab("Time")+
facet_grid(.~task)+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_panamath_longitudinal
###mental rotation longitudinal
mentalrotation_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1MR.s","T2MR.s","T3MR.s","T4MR.s")
mentalrotation_longitudinal_gather = gather(mentalrotation_longitudinal, "time", "MR.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
mentalrotation_longitudinal_gather$SubjectID = as.factor(mentalrotation_longitudinal_gather$SubjectID)
mentalrotation_longitudinal_gather = mentalrotation_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
mentalrotation_longitudinal_gather_summary = summarySE(mentalrotation_longitudinal_gather, "MR.s", c("time","Class"), na.rm = T)
mentalrotation_longitudinal_gather_summary$time = as.factor(mentalrotation_longitudinal_gather_summary$time)
mentalrotation_longitudinal_gather_summary$Class = as.factor(mentalrotation_longitudinal_gather_summary$Class)
mentalrotation_longitudinal_gather_summary$task = "Mental Rotation"
mentalrotation_longitudinal_gather1 = subset(mentalrotation_longitudinal_gather, !is.na(MR.s))
mentalrotation_longitudinal_gather1$timenum = as.numeric(as.factor(mentalrotation_longitudinal_gather1$time))
mentalrotation_longitudinal_gather1$Class = as.factor(mentalrotation_longitudinal_gather1$Class)
mentalrotation_longitudinal_gather1 = mentalrotation_longitudinal_gather1 %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
mentalrotation_model <- lmer(MR.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = subset(mentalrotation_longitudinal_gather1, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(mentalrotation_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: MR.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: subset(mentalrotation_longitudinal_gather1, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2632
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.569 -0.588 0.085 0.569 3.967
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.246 0.496
## Residual 0.416 0.645
## Number of obs: 1146, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value
## (Intercept) -1.00072 0.29226 488.19968 -3.42
## Classlowmr -0.17541 0.19597 1033.59285 -0.90
## Classhigh 1.63852 0.19432 1035.86588 8.43
## timenum 0.10223 0.05544 844.68455 1.84
## T0SGENGender 0.14673 0.06964 305.94887 2.11
## T0PDEMOMaxParentEducation 0.00602 0.01446 305.38620 0.42
## Classlowmr:timenum 0.05128 0.06122 845.65744 0.84
## Classhigh:timenum -0.23109 0.06079 847.62581 -3.80
## Pr(>|t|)
## (Intercept) 0.00067 ***
## Classlowmr 0.37094
## Classhigh < 0.0000000000000002 ***
## timenum 0.06555 .
## T0SGENGender 0.03594 *
## T0PDEMOMaxParentEducation 0.67771
## Classlowmr:timenum 0.40247
## Classhigh:timenum 0.00015 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Clsslw Clsshg timenm T0SGEN T0PDEM Clssl:
## Classlowmr -0.587
## Classhigh -0.530 0.828
## timenum -0.476 0.701 0.707
## T0SGENGendr -0.351 0.105 0.079 0.000
## T0PDEMOMxPE -0.688 -0.001 -0.072 0.008 -0.071
## Clsslwmr:tm 0.430 -0.772 -0.640 -0.906 0.004 -0.008
## Clsshgh:tmn 0.429 -0.640 -0.773 -0.912 -0.003 0.002 0.826
Anova(mentalrotation_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: MR.s
## F Df Df.res Pr(>F)
## (Intercept) 11.72 1 487 0.00067 ***
## Class 134.48 2 1032 < 0.0000000000000002 ***
## timenum 3.40 1 844 0.06556 .
## T0SGENGender 4.44 1 305 0.03595 *
## T0PDEMOMaxParentEducation 0.17 1 305 0.67772
## Class:timenum 32.08 2 850 0.000000000000037 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(ggpredict(mentalrotation_model, terms=c("timenum","Class")))
emeans=emmeans(mentalrotation_model, pairwise ~ Class|timenum, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## timenum = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.589 0.1389 642 -0.862 -0.317 -4.243 <.0001
## lowmr -0.714 0.0648 628 -0.841 -0.586 -11.008 <.0001
## high 0.818 0.0614 642 0.698 0.939 13.318 <.0001
##
## timenum = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.487 0.1152 336 -0.714 -0.260 -4.227 <.0001
## lowmr -0.560 0.0541 333 -0.666 -0.454 -10.348 <.0001
## high 0.689 0.0510 338 0.589 0.790 13.501 <.0001
##
## timenum = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.385 0.1158 342 -0.613 -0.157 -3.323 0.0010
## lowmr -0.406 0.0548 347 -0.514 -0.299 -7.415 <.0001
## high 0.560 0.0518 353 0.459 0.662 10.815 <.0001
##
## timenum = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.283 0.1404 660 -0.558 -0.007 -2.014 0.0444
## lowmr -0.253 0.0665 665 -0.384 -0.122 -3.801 0.0002
## high 0.431 0.0633 680 0.307 0.556 6.815 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## timenum = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.1241 0.1537 636 -0.178 0.426 0.808 0.4196
## low - high -1.4074 0.1523 640 -1.706 -1.108 -9.243 <.0001
## lowmr - high -1.5316 0.0896 631 -1.708 -1.355 -17.085 <.0001
##
## timenum = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.0729 0.1278 335 -0.179 0.324 0.570 0.5692
## low - high -1.1763 0.1265 336 -1.425 -0.927 -9.298 <.0001
## lowmr - high -1.2492 0.0748 334 -1.396 -1.102 -16.706 <.0001
##
## timenum = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.0216 0.1287 343 -0.232 0.275 0.168 0.8670
## low - high -0.9453 0.1273 344 -1.196 -0.695 -7.423 <.0001
## lowmr - high -0.9668 0.0758 349 -1.116 -0.818 -12.760 <.0001
##
## timenum = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.0297 0.1558 658 -0.336 0.276 -0.191 0.8488
## low - high -0.7142 0.1543 660 -1.017 -0.411 -4.627 <.0001
## lowmr - high -0.6845 0.0921 668 -0.865 -0.504 -7.431 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
em=emtrends(mentalrotation_model, pairwise ~ Class, var="timenum", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class timenum.trend SE df lower.CL upper.CL t.ratio p.value
## low 0.102 0.0554 844 -0.0066 0.2111 1.844 0.0656
## lowmr 0.154 0.0260 850 0.1026 0.2045 5.912 <.0001
## high -0.129 0.0249 861 -0.1778 -0.0799 -5.164 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.0513 0.0612 845 -0.171 0.0689 -0.838 0.4025
## low - high 0.2311 0.0608 847 0.112 0.3504 3.801 0.0002
## lowmr - high 0.2824 0.0360 855 0.212 0.3531 7.841 <.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
graph_mentalrotation_longitudinal = ggplot(subset(mentalrotation_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = MR.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (MR.s-se), ymax = (MR.s+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
#scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Mental Rotation")+
xlab("Time")+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_mentalrotation_longitudinal
workingmemory_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1AWMADMSS_N.s","T2AWMADMSS_N.s","T3AWMADMSS_N.s","T4AWMADMSS_N.s")
workingmemory_longitudinal_gather = gather(workingmemory_longitudinal, "time", "AWMADMSS_N.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
workingmemory_longitudinal_gather$SubjectID = as.factor(workingmemory_longitudinal_gather$SubjectID)
workingmemory_longitudinal_gather = workingmemory_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
workingmemory_longitudinal_gather_summary = summarySE(workingmemory_longitudinal_gather, "AWMADMSS_N.s", c("time","Class"), na.rm = T)
workingmemory_longitudinal_gather_summary$time = as.factor(workingmemory_longitudinal_gather_summary$time)
workingmemory_longitudinal_gather_summary$Class = as.factor(workingmemory_longitudinal_gather_summary$Class)
workingmemory_longitudinal_gather_summary$task = "Working Memory"
workingmemory_longitudinal_gather_summary = summarySE(workingmemory_longitudinal_gather, "AWMADMSS_N.s", c("time","Class"), na.rm = T)
workingmemory_longitudinal_gather_summary$time = as.factor(workingmemory_longitudinal_gather_summary$time)
workingmemory_longitudinal_gather_summary$Class = as.factor(workingmemory_longitudinal_gather_summary$Class)
workingmemory_longitudinal_gather_summary$task = "Working Memory"
workingmemory_longitudinal_gather1 = subset(workingmemory_longitudinal_gather, !is.na(AWMADMSS_N.s))
workingmemory_longitudinal_gather1$timenum = as.numeric(as.factor(workingmemory_longitudinal_gather1$time))
workingmemory_longitudinal_gather1$Class = as.factor(workingmemory_longitudinal_gather1$Class)
workingmemory_longitudinal_gather1 = workingmemory_longitudinal_gather1 %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
workingmemory_model <- lmer(AWMADMSS_N.s ~ Class * (timenum)+ T0SGENGender + T0PDEMOMaxParentEducation +
(1+timenum|SubjectID), data = subset(workingmemory_longitudinal_gather1, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(workingmemory_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## AWMADMSS_N.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 + timenum | SubjectID)
## Data: subset(workingmemory_longitudinal_gather1, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2834
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.6290 -0.5823 -0.0095 0.6021 2.8126
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## SubjectID (Intercept) 0.3615 0.601
## timenum 0.0131 0.115 -0.29
## Residual 0.5357 0.732
## Number of obs: 1091, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.91451 0.35471 439.63020 -5.40 0.00000011 ***
## Classlowmr 0.78131 0.24597 323.99833 3.18 0.00163 **
## Classhigh 1.33601 0.24428 325.27170 5.47 0.00000009 ***
## timenum 0.18904 0.07229 326.09389 2.62 0.00933 **
## T0SGENGender -0.00705 0.08199 301.30862 -0.09 0.93156
## T0PDEMOMaxParentEducation 0.06425 0.01703 300.90502 3.77 0.00019 ***
## Classlowmr:timenum -0.17012 0.07903 318.15008 -2.15 0.03210 *
## Classhigh:timenum -0.22904 0.07861 319.77759 -2.91 0.00382 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Clsslw Clsshg timenm T0SGEN T0PDEM Clssl:
## Classlowmr -0.615
## Classhigh -0.564 0.847
## timenum -0.516 0.738 0.743
## T0SGENGendr -0.345 0.101 0.077 0.003
## T0PDEMOMxPE -0.667 -0.003 -0.070 0.005 -0.069
## Clsslwmr:tm 0.469 -0.799 -0.679 -0.915 0.002 -0.004
## Clsshgh:tmn 0.468 -0.679 -0.800 -0.919 -0.004 0.006 0.841
Anova(workingmemory_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: AWMADMSS_N.s
## F Df Df.res Pr(>F)
## (Intercept) 28.99 1 450 0.000000118 ***
## Class 18.67 2 315 0.000000022 ***
## timenum 6.83 1 323 0.0094 **
## T0SGENGender 0.01 1 304 0.9318
## T0PDEMOMaxParentEducation 14.13 1 304 0.0002 ***
## Class:timenum 4.39 2 299 0.0132 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(ggpredict(workingmemory_model, terms=c("timenum","Class")))
emeans=emmeans(workingmemory_model, pairwise ~ Class|timenum, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## timenum = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.785 0.1726 345 -1.125 -0.4458 -4.551 <.0001
## lowmr -0.174 0.0760 296 -0.324 -0.0246 -2.291 0.0226
## high 0.322 0.0722 299 0.180 0.4637 4.458 <.0001
##
## timenum = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.596 0.1386 330 -0.869 -0.3237 -4.303 <.0001
## lowmr -0.155 0.0633 300 -0.280 -0.0305 -2.450 0.0148
## high 0.282 0.0597 303 0.164 0.3992 4.718 <.0001
##
## timenum = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.407 0.1381 295 -0.679 -0.1353 -2.948 0.0035
## lowmr -0.136 0.0655 301 -0.265 -0.0073 -2.080 0.0384
## high 0.242 0.0619 301 0.120 0.3636 3.903 0.0001
##
## timenum = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.218 0.1716 292 -0.556 0.1196 -1.271 0.2046
## lowmr -0.117 0.0813 294 -0.277 0.0427 -1.443 0.1502
## high 0.202 0.0776 291 0.049 0.3544 2.600 0.0098
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## timenum = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.611 0.1891 337 -0.983 -0.2393 -3.232 0.0013
## low - high -1.107 0.1875 338 -1.476 -0.7381 -5.904 <.0001
## lowmr - high -0.496 0.1052 300 -0.703 -0.2887 -4.712 <.0001
##
## timenum = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.441 0.1530 324 -0.742 -0.1401 -2.883 0.0042
## low - high -0.878 0.1514 325 -1.176 -0.5801 -5.798 <.0001
## lowmr - high -0.437 0.0875 302 -0.609 -0.2646 -4.992 <.0001
##
## timenum = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.271 0.1535 297 -0.573 0.0312 -1.765 0.0787
## low - high -0.649 0.1519 296 -0.948 -0.3499 -4.272 <.0001
## lowmr - high -0.378 0.0906 301 -0.556 -0.1997 -4.173 <.0001
##
## timenum = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.101 0.1904 294 -0.476 0.2739 -0.530 0.5968
## low - high -0.420 0.1887 292 -0.791 -0.0485 -2.225 0.0268
## lowmr - high -0.319 0.1127 293 -0.541 -0.0972 -2.831 0.0050
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
em=emtrends(workingmemory_model, pairwise ~ Class, var="timenum", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class timenum.trend SE df lower.CL upper.CL t.ratio p.value
## low 0.1890 0.0724 322 0.0467 0.3314 2.613 0.0094
## lowmr 0.0189 0.0320 277 -0.0440 0.0818 0.592 0.5543
## high -0.0400 0.0309 284 -0.1009 0.0209 -1.293 0.1969
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.1701 0.0791 315 0.0145 0.326 2.151 0.0323
## low - high 0.2290 0.0787 316 0.0742 0.384 2.911 0.0039
## lowmr - high 0.0589 0.0445 280 -0.0286 0.146 1.325 0.1862
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
graph_workingmemory_longitudinal = ggplot(subset(workingmemory_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = AWMADMSS_N.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (AWMADMSS_N.s-se), ymax = (AWMADMSS_N.s+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
#scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("working memory")+
xlab("Time")+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_workingmemory_longitudinal
prop_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1PRWMPAE.asin.sR","T2PRWMPAE.asin.sR","T3PRWMPAE.asin.sR","T4PRWMPAE.asin.sR")
prop_longitudinal_gather = gather(prop_longitudinal, "time", "PRWMPAE.asin.sR", -c("SubjectID","T0SY1Y2G_StartGrade"))
prop_longitudinal_gather$SubjectID = as.factor(prop_longitudinal_gather$SubjectID)
prop_longitudinal_gather = prop_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
prop_longitudinal_gather_summary = summarySE(prop_longitudinal_gather, "PRWMPAE.asin.sR", c("time","Class"), na.rm = T)
prop_longitudinal_gather_summary$time = as.factor(prop_longitudinal_gather_summary$time)
prop_longitudinal_gather_summary$Class = as.factor(prop_longitudinal_gather_summary$Class)
prop_longitudinal_gather_summary$task = "Proportional Reasoning"
prop_longitudinal_gather_summary = summarySE(prop_longitudinal_gather, "PRWMPAE.asin.sR", c("time","Class"), na.rm = T)
prop_longitudinal_gather_summary$time = as.factor(prop_longitudinal_gather_summary$time)
prop_longitudinal_gather_summary$Class = as.factor(prop_longitudinal_gather_summary$Class)
prop_longitudinal_gather_summary$task = "Proportional Reasoning"
prop_longitudinal_gather1 = subset(prop_longitudinal_gather, !is.na(PRWMPAE.asin.sR))
prop_longitudinal_gather1$timenum = as.numeric(as.factor(prop_longitudinal_gather1$time))
prop_longitudinal_gather1$Class = as.factor(prop_longitudinal_gather1$Class)
prop_longitudinal_gather1 = prop_longitudinal_gather1 %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
prop_model <- lmer(PRWMPAE.asin.sR ~ Class * (timenum)+T0SGENGender+T0PDEMOMaxParentEducation+
(1+timenum|SubjectID), data = subset(prop_longitudinal_gather1, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(prop_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## PRWMPAE.asin.sR ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 + timenum | SubjectID)
## Data: subset(prop_longitudinal_gather1, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2817
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.173 -0.454 0.108 0.576 2.319
##
## Random effects:
## Groups Name Variance Std.Dev. Corr
## SubjectID (Intercept) 0.3552 0.596
## timenum 0.0163 0.128 -0.34
## Residual 0.4944 0.703
## Number of obs: 1114, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.9899 0.3312 413.8560 -6.01 0.0000000041 ***
## Classlowmr 0.8570 0.2217 289.7756 3.87 0.00014 ***
## Classhigh 1.0881 0.2198 290.8124 4.95 0.0000012509 ***
## timenum 0.0766 0.0645 269.3125 1.19 0.23567
## T0SGENGender 0.0819 0.0793 301.6468 1.03 0.30256
## T0PDEMOMaxParentEducation 0.0644 0.0165 303.4853 3.90 0.00012 ***
## Classlowmr:timenum -0.0680 0.0717 272.1583 -0.95 0.34344
## Classhigh:timenum -0.0588 0.0712 273.7389 -0.83 0.40922
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Clsslw Clsshg timenm T0SGEN T0PDEM Clssl:
## Classlowmr -0.574
## Classhigh -0.519 0.823
## timenum -0.463 0.696 0.702
## T0SGENGendr -0.362 0.113 0.087 0.006
## T0PDEMOMxPE -0.688 -0.015 -0.085 -0.007 -0.065
## Clsslwmr:tm 0.412 -0.771 -0.632 -0.900 -0.001 0.009
## Clsshgh:tmn 0.413 -0.631 -0.772 -0.906 -0.009 0.017 0.815
Anova(prop_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: PRWMPAE.asin.sR
## F Df Df.res Pr(>F)
## (Intercept) 35.93 1 415 0.0000000044 ***
## Class 12.32 2 298 0.0000072406 ***
## timenum 1.41 1 268 0.23577
## T0SGENGender 1.06 1 305 0.30424
## T0PDEMOMaxParentEducation 15.12 1 307 0.00012 ***
## Class:timenum 0.45 2 279 0.63536
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
plot(ggpredict(prop_model, terms=c("timenum","Class")))
emeans=emmeans(prop_model, pairwise ~ Class|timenum, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## timenum = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.8362 0.1555 299 -1.1424 -0.5301 -5.376 <.0001
## lowmr -0.0473 0.0737 303 -0.1923 0.0978 -0.641 0.5219
## high 0.1930 0.0695 305 0.0563 0.3297 2.779 0.0058
##
## timenum = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.7596 0.1295 300 -1.0145 -0.5047 -5.865 <.0001
## lowmr -0.0386 0.0613 304 -0.1593 0.0820 -0.630 0.5291
## high 0.2108 0.0576 306 0.0974 0.3242 3.659 0.0003
##
## timenum = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.6830 0.1329 287 -0.9446 -0.4214 -5.138 <.0001
## lowmr -0.0300 0.0636 301 -0.1551 0.0951 -0.472 0.6373
## high 0.2286 0.0603 304 0.1100 0.3472 3.792 0.0002
##
## timenum = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.6063 0.1639 278 -0.9291 -0.2836 -3.699 0.0003
## lowmr -0.0214 0.0793 294 -0.1774 0.1346 -0.270 0.7875
## high 0.2464 0.0759 296 0.0970 0.3959 3.245 0.0013
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## timenum = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.789 0.1727 301 -1.129 -0.4492 -4.569 <.0001
## low - high -1.029 0.1709 302 -1.366 -0.6929 -6.021 <.0001
## lowmr - high -0.240 0.1017 306 -0.440 -0.0402 -2.363 0.0188
##
## timenum = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.721 0.1440 301 -1.004 -0.4377 -5.008 <.0001
## low - high -0.970 0.1423 301 -1.250 -0.6903 -6.817 <.0001
## lowmr - high -0.249 0.0846 305 -0.416 -0.0830 -2.948 0.0034
##
## timenum = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.653 0.1480 291 -0.944 -0.3617 -4.412 <.0001
## low - high -0.912 0.1464 290 -1.200 -0.6234 -6.225 <.0001
## lowmr - high -0.259 0.0881 303 -0.432 -0.0854 -2.937 0.0036
##
## timenum = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.585 0.1827 283 -0.945 -0.2254 -3.202 0.0015
## low - high -0.853 0.1810 281 -1.209 -0.4965 -4.711 <.0001
## lowmr - high -0.268 0.1101 297 -0.484 -0.0511 -2.432 0.0156
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
emeans=emmeans(prop_model, pairwise ~ Class, mult.name = "Class")
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.7234 0.1271 294 -0.974 -0.4732 -5.690 <.0001
## lowmr -0.0346 0.0604 303 -0.153 0.0843 -0.572 0.5677
## high 0.2192 0.0569 307 0.107 0.3313 3.850 0.0001
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.689 0.1414 296 -0.967 -0.4105 -4.870 <.0001
## low - high -0.943 0.1399 296 -1.218 -0.6674 -6.740 <.0001
## lowmr - high -0.254 0.0835 305 -0.418 -0.0895 -3.040 0.0026
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
em=emtrends(prop_model, pairwise ~ Class, var="timenum", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class timenum.trend SE df lower.CL upper.CL t.ratio p.value
## low 0.07663 0.0645 268 -0.0503 0.2036 1.188 0.2358
## lowmr 0.00862 0.0313 283 -0.0530 0.0702 0.276 0.7830
## high 0.01780 0.0301 293 -0.0415 0.0771 0.590 0.5554
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.06801 0.0717 271 -0.0731 0.2091 0.949 0.3436
## low - high 0.05883 0.0712 272 -0.0813 0.1990 0.826 0.4094
## lowmr - high -0.00918 0.0435 288 -0.0947 0.0764 -0.211 0.8329
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
graph_prop_longitudinal = ggplot(subset(prop_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = PRWMPAE.asin.sR, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (PRWMPAE.asin.sR-se), ymax = (PRWMPAE.asin.sR+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
#scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Proportional Reasoning")+
xlab("Time")+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_prop_longitudinal
exactcalculation_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1EC.s","T2EC.s","T3EC.s","T4EC.s")
exactcalculation_longitudinal_gather = gather(exactcalculation_longitudinal, "time", "EC.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
exactcalculation_longitudinal_gather$SubjectID = as.factor(as.character(exactcalculation_longitudinal_gather$SubjectID))
exactcalculation_longitudinal_gather = exactcalculation_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
exactcalculation_longitudinal_gather_summary = summarySE(exactcalculation_longitudinal_gather, "EC.s", c("time","Class"), na.rm = T)
exactcalculation_longitudinal_gather_summary$time = as.factor(exactcalculation_longitudinal_gather_summary$time)
exactcalculation_longitudinal_gather_summary$Class = as.factor(exactcalculation_longitudinal_gather_summary$Class)
exactcalculation_longitudinal_gather_summary$task = "Exact Calculation"
exactcalculation_longitudinal_gather1 = subset(exactcalculation_longitudinal_gather, !is.na(EC.s))
exactcalculation_longitudinal_gather1$timenum = as.numeric(as.factor(exactcalculation_longitudinal_gather1$time))
exactcalculation_longitudinal_gather1$Class = as.factor(exactcalculation_longitudinal_gather1$Class)
exactcalculation_longitudinal_gather1 = exactcalculation_longitudinal_gather1 %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
exactcalculation_model <- lmer(EC.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = subset(exactcalculation_longitudinal_gather1, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(exactcalculation_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: EC.s ~ Class * (timenum) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: subset(exactcalculation_longitudinal_gather1, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2561
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -6.327 -0.548 0.022 0.567 2.478
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.373 0.610
## Residual 0.445 0.667
## Number of obs: 1047, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -1.60223 0.33881 440.58831 -4.73 0.0000030417
## Classlowmr 0.60901 0.21758 902.81981 2.80 0.0052
## Classhigh 1.31038 0.21554 904.55415 6.08 0.0000000018
## timenum 0.05225 0.05812 737.60493 0.90 0.3690
## T0SGENGender -0.00674 0.08281 302.07329 -0.08 0.9352
## T0PDEMOMaxParentEducation 0.04807 0.01720 301.71370 2.79 0.0055
## Classlowmr:timenum -0.00392 0.06473 741.80721 -0.06 0.9517
## Classhigh:timenum -0.10362 0.06423 743.33157 -1.61 0.1071
##
## (Intercept) ***
## Classlowmr **
## Classhigh ***
## timenum
## T0SGENGender
## T0PDEMOMaxParentEducation **
## Classlowmr:timenum
## Classhigh:timenum
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Clsslw Clsshg timenm T0SGEN T0PDEM Clssl:
## Classlowmr -0.559
## Classhigh -0.501 0.822
## timenum -0.423 0.661 0.668
## T0SGENGendr -0.369 0.118 0.095 0.006
## T0PDEMOMxPE -0.704 -0.008 -0.082 -0.005 -0.064
## Clsslwmr:tm 0.377 -0.734 -0.599 -0.898 0.001 0.004
## Clsshgh:tmn 0.382 -0.600 -0.735 -0.905 -0.016 0.011 0.813
Anova(exactcalculation_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: EC.s
## F Df Df.res Pr(>F)
## (Intercept) 22.36 1 445 0.000003036058 ***
## Class 25.93 2 910 0.000000000011 ***
## timenum 0.81 1 741 0.3690
## T0SGENGender 0.01 1 306 0.9352
## T0PDEMOMaxParentEducation 7.81 1 305 0.0055 **
## Class:timenum 3.60 2 756 0.0277 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(exactcalculation_model, pairwise ~ Class, var="timenum", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class timenum.trend SE df lower.CL upper.CL t.ratio p.value
## low 0.0522 0.0581 741 -0.0619 0.16636 0.899 0.3690
## lowmr 0.0483 0.0285 762 -0.0076 0.10425 1.696 0.0902
## high -0.0514 0.0273 774 -0.1050 0.00227 -1.880 0.0605
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr 0.00392 0.0647 745 -0.1232 0.131 0.061 0.9517
## low - high 0.10362 0.0642 747 -0.0225 0.230 1.613 0.1071
## lowmr - high 0.09970 0.0395 767 0.0222 0.177 2.525 0.0118
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
emeans=emmeans(exactcalculation_model, pairwise ~ Class|timenum, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## timenum = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.8473 0.1582 559 -1.1581 -0.5365 -5.355 <.0001
## lowmr -0.2422 0.0751 569 -0.3897 -0.0947 -3.225 0.0013
## high 0.3594 0.0706 573 0.2207 0.4982 5.088 <.0001
##
## timenum = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.7951 0.1358 324 -1.0623 -0.5278 -5.854 <.0001
## lowmr -0.1939 0.0642 328 -0.3201 -0.0676 -3.021 0.0027
## high 0.3081 0.0603 329 0.1895 0.4266 5.112 <.0001
##
## timenum = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.7428 0.1365 330 -1.0112 -0.4744 -5.443 <.0001
## lowmr -0.1456 0.0650 342 -0.2733 -0.0178 -2.241 0.0257
## high 0.2567 0.0614 349 0.1360 0.3774 4.182 <.0001
##
## timenum = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## low -0.6906 0.1598 576 -1.0045 -0.3766 -4.320 <.0001
## lowmr -0.0972 0.0771 608 -0.2487 0.0542 -1.261 0.2078
## high 0.2053 0.0735 627 0.0611 0.3496 2.795 0.0053
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## timenum = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.605 0.1757 559 -0.950 -0.2600 -3.444 0.0006
## low - high -1.207 0.1739 560 -1.548 -0.8652 -6.940 <.0001
## lowmr - high -0.602 0.1035 567 -0.805 -0.3984 -5.814 <.0001
##
## timenum = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.601 0.1509 325 -0.898 -0.3043 -3.984 0.0001
## low - high -1.103 0.1492 325 -1.397 -0.8096 -7.395 <.0001
## lowmr - high -0.502 0.0885 328 -0.676 -0.3278 -5.671 <.0001
##
## timenum = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.597 0.1519 333 -0.896 -0.2985 -3.933 0.0001
## low - high -1.000 0.1501 332 -1.295 -0.7042 -6.659 <.0001
## lowmr - high -0.402 0.0899 345 -0.579 -0.2255 -4.476 <.0001
##
## timenum = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## low - lowmr -0.593 0.1781 580 -0.943 -0.2435 -3.331 0.0009
## low - high -0.896 0.1762 581 -1.242 -0.5498 -5.084 <.0001
## lowmr - high -0.303 0.1070 614 -0.513 -0.0925 -2.829 0.0048
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
plot(ggpredict(exactcalculation_model, terms=c("timenum","Class")))
graph_exactcalculation_longitudinal = ggplot(subset(exactcalculation_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = EC.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (EC.s-se), ymax = (EC.s+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
#scale_y_continuous(breaks = seq(0,1,.25), limits = c(0,1.1))+
theme_bw()+
ylab("Exact Calculation")+
xlab("Time")+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_exactcalculation_longitudinal
approx_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1ASWmac.s","T2ASWmac.s","T3ASWmac.s","T4ASWmac.s")
approx_longitudinal_gather = gather(approx_longitudinal, "time", "ASWmac.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
approx_longitudinal_gather$SubjectID = as.factor(as.character(approx_longitudinal_gather$SubjectID))
approx_longitudinal_gather = approx_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
approx_longitudinal_gather = approx_longitudinal_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
approx_longitudinal_gather_summary = summarySE(approx_longitudinal_gather, "ASWmac.s", c("time","Class"), na.rm = T)
approx_longitudinal_gather_summary$time = as.factor(approx_longitudinal_gather_summary$time)
approx_longitudinal_gather_summary$Class = as.factor(approx_longitudinal_gather_summary$Class)
approx_longitudinal_gather_summary$task = "Approximate"
approx_longitudinal_gather$time_num = as.numeric(as.factor(approx_longitudinal_gather$time))
approx_longitudinal_gather$Class = as.factor(as.character(approx_longitudinal_gather$Class))
approx_model <- lmer(ASWmac.s ~ Class * (time_num) + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = subset(approx_longitudinal_gather, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(approx_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## ASWmac.s ~ Class * (time_num) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: subset(approx_longitudinal_gather, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2856
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.7524 -0.6190 0.0013 0.6842 3.0194
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.225 0.475
## Residual 0.662 0.813
## Number of obs: 1069, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.7156 0.2770 365.7109 -2.58 0.01018 *
## Classlow -0.7970 0.2365 1050.3675 -3.37 0.00078 ***
## Classlowmr -0.3003 0.1406 1050.3358 -2.14 0.03293 *
## time_num -0.0103 0.0333 811.3583 -0.31 0.75815
## T0SGENGender 0.1559 0.0754 300.4603 2.07 0.03940 *
## T0PDEMOMaxParentEducation 0.0479 0.0157 305.6013 3.05 0.00251 **
## Classlow:time_num 0.0132 0.0804 801.9007 0.16 0.86956
## Classlowmr:time_num 0.0124 0.0477 804.5459 0.26 0.79531
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Classlw Clsslwm tim_nm T0SGEN T0PDEM Clsslw:t_
## Classlow -0.167
## Classlowmr -0.346 0.283
## time_num -0.305 0.341 0.572
## T0SGENGendr -0.343 -0.077 0.039 -0.013
## T0PDEMOMxPE -0.848 0.064 0.104 0.023 -0.066
## Clsslw:tm_n 0.116 -0.826 -0.236 -0.414 0.009 0.001
## Clsslwmr:t_ 0.209 -0.238 -0.827 -0.698 0.012 -0.014 0.289
Anova(approx_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: ASWmac.s
## F Df Df.res Pr(>F)
## (Intercept) 6.67 1 368 0.0102 *
## Class 6.44 2 1050 0.0017 **
## time_num 0.09 1 813 0.7582
## T0SGENGender 4.28 1 302 0.0394 *
## T0PDEMOMaxParentEducation 9.28 1 307 0.0025 **
## Class:time_num 0.04 2 804 0.9627
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(approx_model, pairwise ~ Class, var="time_num", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class time_num.trend SE df lower.CL upper.CL t.ratio p.value
## high -0.01025 0.0333 813 -0.0756 0.0551 -0.308 0.7582
## low 0.00296 0.0732 801 -0.1408 0.1467 0.040 0.9678
## lowmr 0.00213 0.0342 800 -0.0650 0.0692 0.062 0.9504
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low -0.013210 0.0804 804 -0.171 0.1447 -0.164 0.8696
## high - lowmr -0.012379 0.0477 806 -0.106 0.0813 -0.259 0.7954
## low - lowmr 0.000831 0.0808 801 -0.158 0.1595 0.010 0.9918
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
emeans=emmeans(approx_model, pairwise ~ Class|time_num, mult.name = "Class", at=list(timenum=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## time_num = 2.45:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.2039 0.0540 307 0.0976 0.3103 3.774 0.0002
## low -0.5607 0.1215 299 -0.7997 -0.3217 -4.616 <.0001
## lowmr -0.0661 0.0571 297 -0.1785 0.0464 -1.157 0.2484
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## time_num = 2.45:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.765 0.1335 300 0.502 1.027 5.728 <.0001
## high - lowmr 0.270 0.0791 301 0.114 0.426 3.415 0.0007
## low - lowmr -0.495 0.1348 299 -0.760 -0.229 -3.669 0.0003
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
plot(ggpredict(approx_model, terms=c("time_num","Class")))
graph_approx_longitudinal = ggplot(subset(approx_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = ASWmac.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (ASWmac.s-se), ymax = (ASWmac.s+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme_bw()+
ylab("Approximate")+
xlab("Time")+
theme(legend.position="none",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_approx_longitudinal
lw_longitudinal = spatial_profiles_database_atleast2tp %>%
dplyr::select("SubjectID","T0SY1Y2G_StartGrade",
"T1LWIDWS.s","T2LWIDWS.s","T3LWIDWS.s","T4LWIDWS.s")
lw_longitudinal_gather = gather(lw_longitudinal, "time", "LWIDWS.s", -c("SubjectID","T0SY1Y2G_StartGrade"))
lw_longitudinal_gather$SubjectID = as.factor(as.character(lw_longitudinal_gather$SubjectID))
lw_longitudinal_gather = lw_longitudinal_gather %>%
left_join(T1_complete_class[c("SubjectID","Class")], by ="SubjectID")
lw_longitudinal_gather = lw_longitudinal_gather %>%
left_join(TAll_complete_info[c("SubjectID","T0SGENGender","T0PDEMOMaxParentEducation")], by ="SubjectID")
lw_longitudinal_gather_summary = summarySE(lw_longitudinal_gather, "LWIDWS.s", c("time","Class"), na.rm = T)
lw_longitudinal_gather_summary$time = as.factor(lw_longitudinal_gather_summary$time)
lw_longitudinal_gather_summary$Class = as.factor(lw_longitudinal_gather_summary$Class)
lw_longitudinal_gather_summary$task = "LetterWordID"
lw_longitudinal_gather$time_num = as.numeric(as.factor(lw_longitudinal_gather$time))
lw_longitudinal_gather$Class = as.factor(as.character(lw_longitudinal_gather$Class))
lw_model <- lmer(LWIDWS.s ~ Class * (time_num) + T0SGENGender + T0PDEMOMaxParentEducation +
(1|SubjectID), data = subset(lw_longitudinal_gather, !is.na(Class)), control = lmerControl(optimizer = "bobyqa"))
summary(lw_model)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula:
## LWIDWS.s ~ Class * (time_num) + T0SGENGender + T0PDEMOMaxParentEducation +
## (1 | SubjectID)
## Data: subset(lw_longitudinal_gather, !is.na(Class))
## Control: lmerControl(optimizer = "bobyqa")
##
## REML criterion at convergence: 2023
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -4.050 -0.518 -0.013 0.514 3.758
##
## Random effects:
## Groups Name Variance Std.Dev.
## SubjectID (Intercept) 0.710 0.843
## Residual 0.175 0.418
## Number of obs: 1064, groups: SubjectID, 311
##
## Fixed effects:
## Estimate Std. Error df t value Pr(>|t|)
## (Intercept) -0.81081 0.35639 314.92610 -2.28 0.0236 *
## Classlow -0.60069 0.20678 511.52982 -2.90 0.0038 **
## Classlowmr -0.37738 0.12222 510.33904 -3.09 0.0021 **
## time_num 0.00954 0.01712 761.30190 0.56 0.5775
## T0SGENGender 0.11502 0.10128 305.81639 1.14 0.2570
## T0PDEMOMaxParentEducation 0.05676 0.02104 305.40975 2.70 0.0074 **
## Classlow:time_num -0.08666 0.04202 759.13558 -2.06 0.0395 *
## Classlowmr:time_num 0.00761 0.02469 758.81931 0.31 0.7580
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) Classlw Clsslwm tim_nm T0SGEN T0PDEM Clsslw:t_
## Classlow -0.136
## Classlowmr -0.329 0.279
## time_num -0.127 0.204 0.344
## T0SGENGendr -0.361 -0.107 0.071 -0.006
## T0PDEMOMxPE -0.877 0.097 0.156 0.013 -0.066
## Clsslw:tm_n 0.049 -0.494 -0.140 -0.407 -0.001 -0.001
## Clsslwmr:t_ 0.090 -0.142 -0.493 -0.694 0.007 -0.012 0.282
Anova(lw_model, test.statistic = "F", type= "III")
## Analysis of Deviance Table (Type III Wald F tests with Kenward-Roger df)
##
## Response: LWIDWS.s
## F Df Df.res Pr(>F)
## (Intercept) 5.18 1 315 0.02357 *
## Class 7.03 2 511 0.00097 ***
## time_num 0.31 1 762 0.57751
## T0SGENGender 1.29 1 306 0.25697
## T0PDEMOMaxParentEducation 7.28 1 306 0.00737 **
## Class:time_num 2.56 2 759 0.07819 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
em=emtrends(lw_model, pairwise ~ Class, var="time_num", mult.name = "Class")
summary(em, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emtrends
## Class time_num.trend SE df lower.CL upper.CL t.ratio p.value
## high 0.00954 0.0171 762 -0.0241 0.04315 0.557 0.5775
## low -0.07712 0.0384 759 -0.1525 -0.00177 -2.009 0.0449
## lowmr 0.01715 0.0178 757 -0.0178 0.05206 0.964 0.3352
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.08666 0.0420 760 0.00416 0.1692 2.062 0.0395
## high - lowmr -0.00761 0.0247 759 -0.05607 0.0409 -0.308 0.7580
## low - lowmr -0.09427 0.0423 759 -0.17732 -0.0112 -2.228 0.0261
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
emeans=emmeans(lw_model, pairwise ~ Class|time_num, mult.name = "Class", at=list(time_num=c(1,2,3,4)))
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## time_num = 1:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.213 0.0766 380 0.0621 0.3631 2.777 0.0058
## low -0.475 0.1727 376 -0.8143 -0.1352 -2.749 0.0063
## lowmr -0.157 0.0811 373 -0.3167 0.0024 -1.937 0.0535
##
## time_num = 2:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.222 0.0728 313 0.0789 0.3654 3.051 0.0025
## low -0.552 0.1646 312 -0.8757 -0.2280 -3.353 0.0009
## lowmr -0.140 0.0775 311 -0.2924 0.0124 -1.808 0.0716
##
## time_num = 3:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.232 0.0730 316 0.0881 0.3753 3.174 0.0016
## low -0.629 0.1653 318 -0.9542 -0.3038 -3.805 0.0002
## lowmr -0.123 0.0778 316 -0.2759 0.0301 -1.580 0.1151
##
## time_num = 4:
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.241 0.0771 390 0.0897 0.3928 3.130 0.0019
## low -0.706 0.1746 392 -1.0494 -0.3628 -4.043 0.0001
## lowmr -0.106 0.0820 388 -0.2670 0.0556 -1.289 0.1983
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## time_num = 1:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.687 0.190 376 0.315 1.0601 3.625 0.0003
## high - lowmr 0.370 0.112 376 0.149 0.5902 3.298 0.0011
## low - lowmr -0.318 0.192 375 -0.694 0.0591 -1.658 0.0982
##
## time_num = 2:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.774 0.181 313 0.418 1.1296 4.283 <.0001
## high - lowmr 0.362 0.107 312 0.152 0.5725 3.388 0.0008
## low - lowmr -0.412 0.183 312 -0.771 -0.0523 -2.254 0.0249
##
## time_num = 3:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.861 0.181 317 0.504 1.2176 4.744 <.0001
## high - lowmr 0.354 0.107 316 0.144 0.5655 3.307 0.0010
## low - lowmr -0.506 0.183 318 -0.867 -0.1451 -2.758 0.0062
##
## time_num = 4:
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.947 0.192 391 0.571 1.3240 4.945 <.0001
## high - lowmr 0.347 0.113 388 0.125 0.5693 3.068 0.0023
## low - lowmr -0.600 0.194 391 -0.981 -0.2195 -3.099 0.0021
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
emeans=emmeans(lw_model, pairwise ~ Class, mult.name = "Class")
summary(emeans, infer=c(TRUE,TRUE), null=0, type = "response", adjust = "none")
## $emmeans
## Class emmean SE df lower.CL upper.CL t.ratio p.value
## high 0.227 0.0724 306 0.084 0.3689 3.129 0.0019
## low -0.587 0.1638 307 -0.909 -0.2645 -3.583 0.0004
## lowmr -0.132 0.0771 306 -0.284 0.0195 -1.715 0.0873
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
##
## $contrasts
## contrast estimate SE df lower.CL upper.CL t.ratio p.value
## high - low 0.813 0.180 306 0.460 1.1672 4.523 <.0001
## high - lowmr 0.359 0.106 306 0.149 0.5679 3.373 0.0008
## low - lowmr -0.455 0.182 307 -0.812 -0.0968 -2.500 0.0129
##
## Results are averaged over the levels of: T0SGENGender
## Degrees-of-freedom method: kenward-roger
## Confidence level used: 0.95
plot(ggpredict(lw_model, terms=c("time_num","Class")))
graph_lw_longitudinal = ggplot(subset(lw_longitudinal_gather_summary,!is.na(Class)), aes(x = time, y = LWIDWS.s, group = Class)) +
geom_hline(yintercept = 0, linetype = "dashed")+
#scale_fill_manual(values = c("#000000","#228B22"))+
#scale_color_manual(values = c("#000000","#228B22"))+
geom_ribbon(aes(ymin = (LWIDWS.s-se), ymax = (LWIDWS.s+se), fill = Class),
alpha = .15, color = NA)+
geom_line(aes(color = Class)) +
scale_color_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
scale_fill_manual(values = c("#e41a1c","#377eb8","#4daf4a","#984ea3"))+
geom_line(aes(color = Class)) +
facet_grid(.~task)+
scale_y_continuous(breaks = seq(-2.5,1,.5), limits = c(-2.7,1.0))+
theme_bw()+
ylab("LetterWordID")+
xlab("Time")+
theme(legend.position="bottom",
axis.title.x = element_text(size=size_text),
axis.text.x = element_text(size=size_text),
#axis.title.x = element_blank(),
panel.grid.major = element_blank(), panel.grid.minor = element_blank(),
panel.background = element_rect(fill = "white", colour = "grey50"),
strip.background =element_rect(fill="#f0f0f0"),
strip.text = element_text(size = size_text),
axis.text.y = element_text(size=size_text),
axis.title.y = element_text(size=size_text),
legend.text=element_text(size=size_text))
graph_lw_longitudinal
multiplot(graph_panamath_longitudinal,graph_mentalrotation_longitudinal,graph_prop_longitudinal,graph_exactcalculation_longitudinal,graph_workingmemory_longitudinal,graph_approx_longitudinal,graph_lw_longitudinal, cols = 4)