Brain Atropy Index

train %>% count(col); test %>% count(col)

predictLDA <- function(model, df){
  df$LD1 <- predict(model, df)$x[,1]
  return(df)
}

mean_AD <- (predictLDA(model, train) %>% group_by(col) %>% summarise(LD1 = mean(LD1)))$LD1[1]
mean_CN <- (predictLDA(model, train) %>% group_by(col) %>% summarise(LD1 = mean(LD1)))$LD1[2]

ggplot(data = predictLDA(model, train)) + geom_density(aes(x=LD1, fill=col), alpha=0.5) + ggtitle("Density plot of train set") + geom_vline(xintercept = mean_AD, linetype="dotted") + geom_vline(xintercept = mean_CN, linetype="dotted") + scale_fill_manual(values = c("#FC4E07","#00AFBB"))

ggplot(data = predictLDA(model, test)) + geom_density(aes(x=LD1, fill=col), alpha=0.5) + ggtitle("Density plot of test set") + geom_vline(xintercept = mean_AD, linetype="dotted") + geom_vline(xintercept = mean_CN, linetype="dotted") + scale_fill_manual(values = c("#00AFBB","#E7B800","#FC4E07"))

ggplot(data = predictLDA(model, BAI_NRCD_SUM.FLW) %>% filter(!is.na(col))) + geom_density(aes(x=LD1, fill=col), alpha=0.5) + ggtitle("Density plot of follow-up set") + geom_vline(xintercept = mean_AD, linetype="dotted") + geom_vline(xintercept = mean_CN, linetype="dotted") + scale_fill_manual(values = c("#00AFBB","#E7B800","#FC4E07"))

ggplot(data = predictLDA(model, train)) + geom_point(aes(x=AGE_MRI, y= LD1, col=col)) + scale_color_manual(values = c("#FC4E07","#00AFBB")) + ggtitle("Distribution of LD1 score (Train set)")

ggplot(data = predictLDA(model, test) %>% filter(col!="MCI")) + geom_point(aes(x=AGE_MRI, y= LD1, col=col)) + scale_color_manual(values = c("#00AFBB","#FC4E07")) + ggtitle("Distribution of LD1 score (Test set)")

tmp <- predictLDA(model, train) %>% filter(col!="MCI") %>% mutate(col = as.character(col))
roc <- pROC::roc(tmp$col, tmp$LD1)

Setting levels: control = AD, case = CN
Setting direction: controls < cases

plot(roc, print.auc=TRUE, print.thres="best")
title("ROC curve of train set")

tmp <- predictLDA(model, test) %>% filter(col!="MCI") %>% mutate(col = as.character(col))
roc <- pROC::roc(tmp$col, tmp$LD1)

Setting levels: control = AD, case = CN
Setting direction: controls < cases

plot(roc, print.auc=TRUE, print.thres="best")
title("ROC curve of test set")

tmp <- predictLDA(model, BAI_NRCD_SUM.FLW) %>% filter(col!="MCI") %>% mutate(col = as.character(col))
roc <- pROC::roc(tmp$col, tmp$LD1)

Setting levels: control = AD, case = CN
Setting direction: controls < cases

plot(roc, print.auc=TRUE, print.thres="best")
title("ROC curve of follow-up set")

Others

tmp <- BAI_NRCD_SUM %>% filter(OID %in% (BAI_NRCD_SUM %>% count(OID, col) %>% count(OID) %>% filter(n==1))$OID)
ggplot(data = predictLDA(model, tmp %>% filter(OID %in% BAI_NRCD_SUM.FLW$OID, !is.na(col))), aes(x=AGE_MRI, y= LD1, col=col)) + geom_point() + geom_line(aes(group=OID),col="black") + scale_color_manual(values = c("#00AFBB","#E7B800","#FC4E07")) + facet_grid(.~col)

tmp <- BAI_NRCD_SUM %>% filter(OID %in% (BAI_NRCD_SUM %>% count(OID, col) %>% count(OID) %>% filter(n>1))$OID)

Factor `col` contains implicit NA, consider using `forcats::fct_explicit_na`

ggplot(data = predictLDA(model, tmp %>% filter(OID %in% BAI_NRCD_SUM.FLW$OID, !is.na(col))), aes(x=AGE_MRI, y= LD1, col=col)) + geom_point() + geom_line(aes(group=OID),col="black") + scale_color_manual(values = c("#00AFBB","#E7B800","#FC4E07"))

tmp <- predictLDA(model, BAI_NRCD_SUM %>% filter(col=="CN") %>% getFLW() %>% filter(FLW=="FIRST")) %>% select(LD1, AGE_MRI)

Joining, by = c("OID", "ID_MRI", "data", "col", "PET", "colPET", "SEX", "AGE_MRI", "c4", "c10", "c11", "c12", "c13", "c17", "c18", "c26", "c28", "c1002", "c1003", "c1005", "c1006", "c1007", "c1008", "c1009", "c1010", "c1011", "c1012", "c1013", "c1014", "c1015", "c1016", "c1017", "c1018", "c1019", "c1020", "c1021", "c1022", "c1023", "c1024", "c1025", "c1026", "c1027", "c1028", "c1029", "c1030", "c1031", "c1034", "c1035", "THICK", "ICV", "WB_V", "GM_V", "WM_V", "CB_V", "DGM_V", "CSF_V", "BS_V", "WB_P", "GM_P", "WM_P", "CB_P", "DGM_P", "CSF_P", "BS_P")

model_lms <- lms(LD1, AGE_MRI, data = tmp)

*** Initial  fit***  
GAMLSS-RS iteration 1: Global Deviance = 983.8082 
GAMLSS-RS iteration 2: Global Deviance = 983.8082 
*** Fitting BCCGo *** 
BCCGo  failed 
*** Fitting BCPEo *** 
BCPEo  failed 
*** Fitting BCTo *** 
BCTo  failed 
*** Refitting NO *** 
GAMLSS-RS iteration 1: Global Deviance = 982.474 
GAMLSS-RS iteration 2: Global Deviance = 982.4813 
% of cases below  0.4798274 centile is  0.5347594 
% of cases below  2.371295 centile is  2.406417 
% of cases below  8.095366 centile is  9.358289 
% of cases below  23.97418 centile is  25.40107 
% of cases below  52.82418 centile is  50 
% of cases below  77.16275 centile is  74.59893 
% of cases below  88.99815 centile is  90.64171 
% of cases below  96.65526 centile is  97.59358 
% of cases below  99.32433 centile is  99.46524 

tmp <- predictLDA(model_noAge, BAI_NRCD_SUM %>% filter(col=="CN") %>% getFLW() %>% filter(FLW=="FIRST")) %>% select(LD1, AGE_MRI)

Joining, by = c("OID", "ID_MRI", "data", "col", "PET", "colPET", "SEX", "AGE_MRI", "c4", "c10", "c11", "c12", "c13", "c17", "c18", "c26", "c28", "c1002", "c1003", "c1005", "c1006", "c1007", "c1008", "c1009", "c1010", "c1011", "c1012", "c1013", "c1014", "c1015", "c1016", "c1017", "c1018", "c1019", "c1020", "c1021", "c1022", "c1023", "c1024", "c1025", "c1026", "c1027", "c1028", "c1029", "c1030", "c1031", "c1034", "c1035", "THICK", "ICV", "WB_V", "GM_V", "WM_V", "CB_V", "DGM_V", "CSF_V", "BS_V", "WB_P", "GM_P", "WM_P", "CB_P", "DGM_P", "CSF_P", "BS_P")

model_lms <- lms(LD1, AGE_MRI, data = tmp)

*** Initial  fit***  
GAMLSS-RS iteration 1: Global Deviance = 970.4225 
GAMLSS-RS iteration 2: Global Deviance = 970.4225 
*** Fitting BCCGo *** 
BCCGo  failed 
*** Fitting BCPEo *** 
BCPEo  failed 
*** Fitting BCTo *** 
BCTo  failed 
*** Refitting NO *** 

additive.fit convergence not obtained in  30  iterations

GAMLSS-RS iteration 1: Global Deviance = 969.7177 

additive.fit convergence not obtained in  30  iterations

GAMLSS-RS iteration 2: Global Deviance = 969.7413 

additive.fit convergence not obtained in  30  iterations

GAMLSS-RS iteration 3: Global Deviance = 969.7476 
% of cases below  0.4257131 centile is  0.5347594 
% of cases below  2.048408 centile is  2.406417 
% of cases below  7.958279 centile is  9.358289 
% of cases below  26.40779 centile is  25.40107 
% of cases below  53.53071 centile is  50 
% of cases below  75.71152 centile is  74.59893 
% of cases below  88.06334 centile is  90.64171 
% of cases below  96.82558 centile is  97.59358 
% of cases below  99.46934 centile is  99.46524