Differences

This shows you the differences between two versions of the page.

--- neuroimagen:bioface_atn [2021/08/30 08:49] – [Otras ROI en ADNI DB] osotolongo
+++ neuroimagen:bioface_atn [2022/06/14 08:22] (current) – [evaluando WMH] osotolongo
@@ Line 86: / Line 86: @@
 <code>
 > library("e1071")
-> library("caTools")
+> library("caret")
 > library("caTools")
 > library("ADNIMERGE")
@@ Line 173: / Line 173: @@
 **Nice!**
+===== What we got so far =====
+A partir de la base de datos de ADNI, tomo los sujetos diagnosticados como //Dementia// y //CN// en el baseline, así como su edad y ciertos valores calculados a partir de la segmentacion. Tomo el diagnostico de //Dementia// como positivo en neurodegeneracion y negativo en caso contrario.
+<code>
+library("e1071")
+library("caret")
+library("caTools")
+library("ADNIMERGE")
+pop <- adnimerge[adnimerge$VISCODE=="bl" & (adnimerge$DX=="Dementia" | adnimerge$DX=="CN"),];
+pop$ND = as.factor(ifelse(pop$DX == "Dementia", 1, 0))
+xt <- pop[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "ND")]
+</code>
+Con estos datos, entreno el clasificador,
+<code>
+classifier_cl <- naiveBayes(ND ~ ., data = xt)
+</code>
+Ahora, tomo un proyecto y reuno los datos de segmentacion y demograficos. No es neceario todo pero la edad y la segmentacion basica si. Esto puede cambiar de proyeto en proyecto pero es mas o menos asi, (lo repito por tenerlo todo junto)
+<code bash>
+xnat_pullfs.pl -s aseg -x bioface19 -o bf_base_aseg.csv
+xnat_pullfs.pl -s aparc -x bioface19 -o bf_base_aparc.csv
+join -t, bf_base_aseg.csv bf_base_aparc.csv > bf_base.csv
+awk -F"," '{print $3","$5}' bioface_np.csv | sed 's/Subjecte/Subject_ID/; s/Edad/AGE/' > bf_age.csv
+join -t, bf_base.csv bf_age.csv > bf_data.csv
+</code>
+Cargo estos datos y construyo las variables identicas a la de los datos de ADNIMERGE.
+<code>
+setwd("/home/data/bioface")
+base <- read.csv("bf_data.csv")
+base$Hippocampus = base$Left.Hippocampus + base$Right.Hippocampus
+base$Entorhinal = base$lh.entorhinal.GrayVol + base$rh.entorhinal.GrayVol
+base$Ventricles <- base$Left.Inf.Lat.Vent + base$Right.Inf.Lat.Vent + base$Left.Lateral.Ventricle + base$Right.Lateral.Ventricle
+base$MidTemp = base$lh.middletemporal.GrayVol + base$rh.middletemporal.GrayVol
+base$ICV = base$eTIV
+</code>
+Ahora puedo calcular el valor de //(N)//,
+<code>
+base$ND <- predict(classifier_cl, newdata = base)
+</code>
+o la probabilidad de que exista neurodegeneración
+<code>
+base$post <- predict(classifier_cl, newdata = base, type="raw")
+</code>
+^ Clasificación ^ Probabilidades de N+ ^
+| {{ :neuroimagen:bioface_nplus_classifier.png?600 | my pretty cool neurodegeneration classifier output }} | {{ :neuroimagen:nplus_probabilities_histogram.png?600 | how the probabilities are distributed? }} |
+Y esto voy a hacerlo para FACEHBI también para compararlo,
+<code>
+setwd("/home/data/facehbi")
+read.csv("face_data.csv") -> face
+face$Hippocampus = face$Left.Hippocampus + face$Right.Hippocampus
+face$Entorhinal = face$lh.entorhinal.GrayVol + face$rh.entorhinal.GrayVol
+face$Ventricles <- face$Left.Inf.Lat.Vent + face$Right.Inf.Lat.Vent + face$Left.Lateral.Ventricle + face$Right.Lateral.Ventricle
+face$MidTemp = face$lh.middletemporal.GrayVol + face$rh.middletemporal.GrayVol
+face$ICV = face$eTIV
+face$ND <- predict(classifier_cl, newdata = face)
+face$post <- predict(classifier_cl, newdata = face, type="raw")
+</code>
+^ Clasificación ^ Probabilidades de N+ ^
+| {{ :neuroimagen:facehbi_nbayes_classifier.png?600 | my pretty cool neurodegeneration classifier output (FACEHBI)}} | {{ :neuroimagen:facehbi_nplus_prob_histogram.png?600 |how the probabilities are distributed?}} |
+==== Métodos no lineales ====
+Vamos a intentar con métodos no lineales. Ojo, estos  métodos no tienen interpretación probabilistica. Son construcciones matemáticas que funcionan.
+=== random Forest ===
+** BIOFACE **
+<code>
+rf_cl <- randomForest(ND ~ ., data = xt, na.action = na.omit)
+base$ND <- predict(rf_cl, newdata = base)
+base$post <- predict(rf_cl, newdata = base, type="prob")
+plot(base$AGE, base$Hippocampus, main = "Hippocampus Volume versus Age", xlab="Age", ylab="HV", pch=19, col=ifelse(base$ND==1,"red","green"))
+hist(base$post[,2], col="red", density = 50, probability=TRUE, breaks = 10, main="Probabilities profile", xlab = "Probability of N+")
+</code>
+^ Clasificación ^ Probabilidades de N+ ^
+| {{ :neuroimagen:rf_classifier_bioface.png?600 |}} | {{ :neuroimagen:rf_prob_profile_bioface.png?600 |}} |
+**FACEHBI**
+<code>
+face$ND <- predict(rf_cl, newdata = face)
+face$post <- predict(rf_cl, newdata = face, type = "prob")
+plot(face$AGE, face$Hippocampus, main = "Hippocampus Volume versus Age", xlab="Age", ylab="HV", pch=19, col=ifelse(face$ND==1,"red","green"))
+hist(face$post[,2], col="red", density = 50, probability=TRUE, breaks = 10, main="Probabilities profile", xlab = "Probability of N+")
+</code>
+^ Clasificación ^ Probabilidades de N+ ^
+| {{ :neuroimagen:rf_classifier_facehbi.png?600 |}} | {{ :neuroimagen:rf_prob_profile_facehbi.png?600 |}} |
+=== support-vector machine ===
+**BIOFACE**
+<code>
+svm_cl <- svm(ND ~ ., data = xt, na.action = na.omit)
+base$ND <- predict(svm_cl, newdata = base)
+plot(base$AGE, base$Hippocampus, main = "Hippocampus Volume versus Age", xlab="Age", ylab="HV", pch=19, col=ifelse(base$ND==1,"red","green"))
+</code>
+{{ :neuroimagen:svm_classifier_bioface.png?600 |}}
+** FACEHBI **
+<code>
+face$ND <- predict(svm_cl, newdata = face)
+plot(face$AGE, face$Hippocampus, main = "Hippocampus Volume versus Age", xlab="Age", ylab="HV", pch=19, col=ifelse(face$ND==1,"red","green"))
+</code>
+{{ :neuroimagen:svm_classifier_facehbi.png?600 |}}
+===== Introduciendo WMH en AT(N)=====
+¿Son las afecciones vasculares parte de la neurodegeneracion? ¿Debe ser considerado el aumento en el volumen de WMH una parte fundamental de la evaluacion de N?
+Pues si queremos hacerlo de esta manera, hay basicamente dos maneras. Utilizando las WMH de Freesurfer (T1w) o las calculadas por un soft independiente (T2w). Si tenemos un FLAIR 3D podemos hacerlo de ambas maneras pero en caso contrario.
+==== evaluando WMH ====
+Las WMH vienen ya calculadas en el procedimiento //standard// asi que esto deberia ser //painless//. Pero en la DB de ADNIMERGE, no vienen integradas en la tabla de resumen que estamos usando, asi que hay que incluirlas. Pero revisando el procedimiento de integracion (adnimerger) veo que hay que arreglar un poco las //VISCODE//.
+<code RMarkdown>
+ucsffsx -> ucsffsx_t
+ucsffsx_t$VISCODE <- ifelse(ucsffsx_t$VISCODE == "sc", "bl", ucsffsx_t$VISCODE)
+a1 <- merge(pop, ucsffsx_t, by=c("RID", "VISCODE"))
+a1t <- a1[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "ST128SV", "ND")]
+...
+...
+axt <- rbind(a1t, a2t, a3t)
+xt <- rename(axt, "WMH" = "ST128SV")
+</code>
+y esto hay que hacerlo para cada una de las tablas y despues concatenar los resutados. Ahora habria tambien que introducir las WMH en la DB del test,
+<code>
+base$WMH = base$WM.hypointensities
+</code>
+++++ que se dice un poco mas rapido de lo que se hace, tambien es cierto, |
+<code Rmarkdown>
+library("e1071")
+library("caret")
+library("caTools")
+library("tidyverse")
+library("ADNIMERGE")
+input_file="input_data.csv"
+output_file="classifier_output_w.csv"
+output_fig="classifier_output_hippocampus_w.ps"
+pop <- adnimerge[(adnimerge$DX=="Dementia" | adnimerge$DX=="CN"),];
+pop$ND = as.factor(ifelse(pop$DX == "Dementia", 1, 0))
+ucsffsx -> ucsffsx_t
+ucsffsx_t$VISCODE <- ifelse(ucsffsx_t$VISCODE == "sc", "bl", ucsffsx_t$VISCODE)
+a1 <- merge(pop, ucsffsx_t, by=c("RID", "VISCODE"))
+ucsffsx51 -> ucsffsx_t
+ucsffsx_t$VISCODE <- ifelse(ucsffsx_t$VISCODE == "sc", "bl", ucsffsx_t$VISCODE)
+a2 <- merge(pop, ucsffsx_t, by=c("RID", "VISCODE"))
+ucsffsx6 -> ucsffsx_t
+ucsffsx_t$VISCODE <- ifelse(ucsffsx_t$VISCODE == "sc", "bl", ucsffsx_t$VISCODE)
+a3 <- merge(pop, ucsffsx_t, by=c("RID", "VISCODE"))
+a1t <- a1[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "ST128SV", "ND")]
+a2t <- a2[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "ST128SV", "ND")]
+a3t <- a3[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "ST128SV", "ND")]
+axt <- rbind(a1t, a2t, a3t)
+xt <- rename(axt, "WMH" = "ST128SV")
+classifier_cl <- naiveBayes(ND ~ ., data = xt)
+base <- read.csv(input_file)
+base$Hippocampus = base$Left.Hippocampus + base$Right.Hippocampus
+base$Entorhinal = base$lh.entorhinal.GrayVol + base$rh.entorhinal.GrayVol
+base$Ventricles <- base$Left.Inf.Lat.Vent + base$Right.Inf.Lat.Vent + base$Left.Lateral.Ventricle + base$Right.Lateral.Ventricle
+base$MidTemp = base$lh.middletemporal.GrayVol + base$rh.middletemporal.GrayVol
+base$ICV = base$eTIV
+base$WMH = base$WM.hypointensities
+base$ND <- predict(classifier_cl, newdata = base)
+base2e = base[, c("Subject_ID", "ND")]
+write.csv(base2e, file=output_file, row.names=FALSE, quote=FALSE)
+a <- lm(base$Hippocampus ~ base$ICV)
+base$aHV = base$Hippocampus - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aHV, main = "Hippocampus volume versus Age", xlab="Age", ylab="adjusted HV", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+a <- lm(base$MidTemp ~ base$ICV)
+base$aMidTemp = base$MidTemp - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+output_fig="classifier_output_middletemporal_w.ps"
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aMidTemp, main = "Middle temporal cortex volume versus Age", xlab="Age", ylab="adjusted MidTemp Volume", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+a <- lm(base$Entorhinal ~ base$ICV)
+base$aEntorhinal = base$Entorhinal - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+output_fig="classifier_output_entorhinal_w.ps"
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aEntorhinal, main = "Entorhinal cortex volume versus Age", xlab="Age", ylab="adjusted Entorhinal Volume", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+</code>
+++++
+OK, ya puestos también podemos hacerlo con las WMH calculadas desde T2W. Yo uso PGS, para las WMH de T2w, asi que tengo que pegarlas a las tablas,
+(Ejemplo con facehbi)
+<code>
+[osotolongo@brick03 facehbi]$ sed 's/;/,/' facehbi_mri.csv > codes.csv
+[osotolongo@brick03 facehbi]$ sed -i '1iSubject,Subject_ID' codes.csv
+[osotolongo@brick03 facehbi]$ join -t, codes.csv facehbi_wmh_metrics.csv | awk -F',' '{print $2","$3}' > facehbi_wmh.csv
+[osotolongo@brick03 facehbi]$ ./xnat_get_age.pl -x facehbi
+[osotolongo@brick03 facehbi]$ xnat_pullfs.pl -s aseg -p facehbi -o base_aseg.csv
+[osotolongo@brick03 facehbi]$ xnat_pullfs.pl -s aparc -p facehbi -o base_aparc.csv
+[osotolongo@brick03 facehbi]$ join -t, base_aseg.csv base_aparc.csv > base_full.csv
+[osotolongo@brick03 facehbi]$ (head -n 1 facehbi_age_data.csv && tail -n +2 facehbi_age_data.csv | sort -t, -n)> facehbi_age_data_sorted.csv
+[osotolongo@brick03 facehbi]$ join -t, facehbi_age_data_sorted.csv facehbi_wmh.csv > facehbi_age_wmh.csv
+[osotolongo@brick03 facehbi]$ join -t, base_full.csv facehbi_age_wmh.csv > input_data.csv
+</code>
+Claro que aqui hay que usar la tabla de UCD WMH,
+<code>
+pop <- adnimerge[(adnimerge$DX=="Dementia" | adnimerge$DX=="CN"),];
+ucd_wmh -> ucd_t
+ucd_t$VISCODE <- ifelse(ucd_t$VISCODE == "scmri", "bl", ucd_t$VISCODE)
+ucd_t$VISCODE <- ifelse(ucd_t$VISCODE == "sc", "bl", ucd_t$VISCODE)
+pop <- merge(pop, ucd_t, by=c("RID", "VISCODE"))
+</code>
+y despues modificar un poco las de nuestro proyecto para que cuadren las unidades,
+<code>
+base$TOTAL_WMH = base$WMH*0.001
+</code>
+++++ OK, ya puestos también podemos hacerlo con las WMH calculadas desde T2W |
+<code Rmarkdown>
+library("e1071")
+library("caret")
+library("caTools")
+library("ADNIMERGE")
+input_file="input_dataw.csv"
+output_file="classifier_output_w2.csv"
+output_fig="classifier_output_hippocampus_w2.ps"
+pop <- adnimerge[(adnimerge$DX=="Dementia" | adnimerge$DX=="CN"),];
+ucd_wmh -> ucd_t
+ucd_t$VISCODE <- ifelse(ucd_t$VISCODE == "scmri", "bl", ucd_t$VISCODE)
+ucd_t$VISCODE <- ifelse(ucd_t$VISCODE == "sc", "bl", ucd_t$VISCODE)
+pop <- merge(pop, ucd_t, by=c("RID", "VISCODE"))
+pop$ND = as.factor(ifelse(pop$DX == "Dementia", 1, 0))
+xt <- pop[, c("Hippocampus", "Entorhinal", "Ventricles", "MidTemp", "AGE", "ICV", "TOTAL_WMH", "ND")]
+classifier_cl <- naiveBayes(ND ~ ., data = xt)
+base <- read.csv(input_file)
+base$Hippocampus = base$Left.Hippocampus + base$Right.Hippocampus
+base$Entorhinal = base$lh.entorhinal.GrayVol + base$rh.entorhinal.GrayVol
+base$Ventricles <- base$Left.Inf.Lat.Vent + base$Right.Inf.Lat.Vent + base$Left.Lateral.Ventricle + base$Right.Lateral.Ventricle
+base$MidTemp = base$lh.middletemporal.GrayVol + base$rh.middletemporal.GrayVol
+base$ICV = base$eTIV
+base$TOTAL_WMH = base$WMH*0.001
+base$ND <- predict(classifier_cl, newdata = base)
+base$posterior <- predict(classifier_cl, newdata = base, type="raw")
+base2e = base[, c("Subject_ID", "ND", "posterior")]
+write.csv(base2e, file=output_file, row.names=FALSE, quote=FALSE)
+a <- lm(base$Hippocampus ~ base$ICV)
+base$aHV = base$Hippocampus - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aHV, main = "Hippocampus volume versus Age", xlab="Age", ylab="adjusted HV", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+a <- lm(base$MidTemp ~ base$ICV)
+base$aMidTemp = base$MidTemp - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+output_fig="classifier_output_middletemporal_w2.ps"
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aMidTemp, main = "Middle temporal cortex volume versus Age", xlab="Age", ylab="adjusted MidTemp", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+a <- lm(base$Entorhinal ~ base$ICV)
+base$aEntorhinal = base$Entorhinal - a$coefficients[[2]]*(base$ICV - mean(base$ICV, na.rm=TRUE))
+output_fig="classifier_output_entorhinal_w2.ps"
+postscript(output_fig, width=1024, height=600, bg="white")
+plot(base$AGE, base$aEntorhinal, main = "Entorhinal cortex volume versus Age", xlab="Age", ylab="adjusted Entorhinal", pch=19, col=ifelse(base$ND==1,"red","green"))
+dev.off()
+</code>
+++++ (**NO WAY**, explicar como se hace que hay que meter el ucd_wmh de ADNI y sacar wl WMH con PGS)
+Con BIOFACE:
+^ no WMH ^ T1w WMH ^ T2w WMH ^
+| {{ :neuroimagen:bioface_classifier_output_hippocampus.png?400 |}} | {{ :neuroimagen:bioface_classifier_output_hippocampus_w.png?400 |}} | {{ :neuroimagen:bioface_classifier_output_hippocampus_w2.png?400 |}} |
+^ ^ T1w WMH ^ T2w WMH ^
+| -2+ | 0 | 1 (1%) |
+| +2- | 9 (10%) | 4 (4%) |
+Con FACEHBI:
+^ no WMH ^ T1w WMH ^ T2w WMH ^
+| {{ :neuroimagen:facehbi_classifier_output_hippocampus.png?400 |}} | {{ :neuroimagen:facehbi_classifier_output_hippocampus_w.png?400 |}} | {{ :neuroimagen:facehbi_classifier_output_hippocampus_w2.png?400 |}} |
+^ ^ T1w WMH ^ T2w WMH ^
+| -2+ | 0 | 0 |
+| +2- | 21 (4%) | 10 (2%) |