User Tools

Site Tools


neuroimagen:altdti

This is an old revision of the document!


experimento DTI

pruebas

Siguiendo los pasos de la pipeline 3.0 procesamos las imagenes DTI en intentamos hallar una red que parta de una parte especifica del cortex.

El archivo dti_track.seed se fabrica con 42 regiones distintas,

grep ctx /usr/local/freesurfer/FreeSurferColorLUT.txt | grep 'parietal\|frontal' | awk {'print $1'} > /nas/data/facehbi/dti_track.seed 

Estas regiones son,

1003 ctx-lh-caudalmiddlefrontal
1008 ctx-lh-inferiorparietal
1012 ctx-lh-lateralorbitofrontal
1014 ctx-lh-medialorbitofrontal
1027 ctx-lh-rostralmiddlefrontal
1028 ctx-lh-superiorfrontal
1029 ctx-lh-superiorparietal
1032 ctx-lh-frontalpole
2003 ctx-rh-caudalmiddlefrontal
2008 ctx-rh-inferiorparietal
2012 ctx-rh-lateralorbitofrontal
2014 ctx-rh-medialorbitofrontal
2027 ctx-rh-rostralmiddlefrontal
2028 ctx-rh-superiorfrontal
2029 ctx-rh-superiorparietal
2032 ctx-rh-frontalpole
1106 ctx-lh-G_frontal_inf-Opercular_part
1107 ctx-lh-G_frontal_inf-Orbital_part
1108 ctx-lh-G_frontal_inf-Triangular_part
1109 ctx-lh-G_frontal_middle
1110 ctx-lh-G_frontal_superior
1122 ctx-lh-G_parietal_inferior-Angular_part
1123 ctx-lh-G_parietal_inferior-Supramarginal_part
1124 ctx-lh-G_parietal_superior
1154 ctx-lh-S_frontal_inferior
1155 ctx-lh-S_frontal_middle
1156 ctx-lh-S_frontal_superior
1159 ctx-lh-S_intraparietal-and_Parietal_transverse
1177 ctx-lh-S_subparietal
2106 ctx-rh-G_frontal_inf-Opercular_part
2107 ctx-rh-G_frontal_inf-Orbital_part
2108 ctx-rh-G_frontal_inf-Triangular_part
2109 ctx-rh-G_frontal_middle
2110 ctx-rh-G_frontal_superior
2122 ctx-rh-G_parietal_inferior-Angular_part
2123 ctx-rh-G_parietal_inferior-Supramarginal_part
2124 ctx-rh-G_parietal_superior
2154 ctx-rh-S_frontal_inferior
2155 ctx-rh-S_frontal_middle
2156 ctx-rh-S_frontal_superior
2159 ctx-rh-S_intraparietal-and_Parietal_transverse
2177 ctx-rh-S_subparietal

El resultado de probtrackx es una red bastante extensa,

Ahora voy a quitar los giros y demas y quedarme solo con las 16 primeras lineas.

$ cat dti_track.seed 
1003
1008
1012
1014
1027
1028
1029
1032
2003
2008
2012
2014
2027
2028
2029
2032

La red que obtengo es ahora bastante similar,

Aplicando esta red como mascara (25%) saco los datos de FA en la red y los puedo comparar con otros datos . Ejemplo, SUVR en la misma visita,

Aislando el cortex

Voy a intentar aislar todas las regiones del cortex para hacer una lista de lo que puedo parcelar

$ grep ctx /usr/local/freesurfer/FreeSurferColorLUT.txt | grep -v "G\|S\|_\|nknown\|#" | awk {'if($1<3000) print $1,$2'} | grep lh | sed 's/-lh//' | awk {'print $2'}

y me quedo con los nombres utiles para que un neurologo me edite la lista,

ctx-bankssts
ctx-caudalanteriorcingulate
ctx-caudalmiddlefrontal
ctx-corpuscallosum
ctx-cuneus
ctx-entorhinal
ctx-fusiform
ctx-inferiorparietal
ctx-inferiortemporal
ctx-isthmuscingulate
ctx-lateraloccipital
ctx-lateralorbitofrontal
ctx-lingual
ctx-medialorbitofrontal
ctx-middletemporal
ctx-parahippocampal
ctx-paracentral
ctx-parsopercularis
ctx-parsorbitalis
ctx-parstriangularis
ctx-pericalcarine
ctx-postcentral
ctx-posteriorcingulate
ctx-precentral
ctx-precuneus
ctx-rostralanteriorcingulate
ctx-rostralmiddlefrontal
ctx-superiorfrontal
ctx-superiorparietal
ctx-superiortemporal
ctx-supramarginal
ctx-frontalpole
ctx-temporalpole
ctx-transversetemporal
ctx-insula

Mapa FP MB

Las regiones escogidas son,

ctx-caudalmiddlefrontal
ctx-inferiorparietal
ctx-middletemporal
ctx-parsopercularis
ctx-parstriangularis
ctx-postcentral
ctx-precentral
ctx-superiorfrontal
ctx-superiorparietal
ctx-superiortemporal
ctx-supramarginal

Voy a quedarme con lo que necesito ahora,

$ grep ctx /usr/local/freesurfer/FreeSurferColorLUT.txt | grep -v "G\|S\|_\|nknown\|#" | awk {'if($1<3000) print $1,$2'} > tofp.txt
$ cat tocut.txt
caudalmiddlefrontal
inferiorparietal
middletemporal
parsopercularis
parstriangularis
postcentral
precentral
superiorfrontal
superiorparietal
superiortemporal
supramarginal
$ grep -f tocut.txt tofp.txt | awk {'print $1'} > dti_track.seed 
$ cat dti_track.seed
1003
1008
1015
1018
1020
1022
1024
1028
1029
1030
1031
2003
2008
2015
2018
2020
2022
2024
2028
2029
2030
2031

con estas seeds ya puedo correr el experimento.

$ dti_track.pl facehbi
... 3dias ...
$ cd /nas/data/facehbi
$ for x in `ls -d working/*_probtrack_out`; do mv $x `echo $x | sed 's/out/FPCustom/'`;done
$ dti_metrics_alt.pl -path FPCustom facehbi
 
$ dti_track.pl v2MriPet
... 3dias ...
$ cd /nas/data/v2MriPet
$ for x in `ls -d working/*_probtrack_out`; do mv $x `echo $x | sed 's/out/FPCustom/'`;done
$ dti_metrics_alt.pl -path FPCustom v2MriPet

Comparando con DMN

Sorprendentemente, aunque estan relacionados, hay bastante diferencia entre la FA medida en la DMN y la medida en la nueva red,

> summary(v1m)
 
Call:
lm(formula = DMN_FA_v1 ~ FPCustom_FA_v1, data = dti_c)
 
Residuals:
      Min        1Q    Median        3Q       Max 
-0.031888 -0.006546  0.001095  0.006806  0.026302 
 
Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)     0.06472    0.01262   5.129 8.74e-07 ***
FPCustom_FA_v1  0.89746    0.04285  20.943  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.01099 on 152 degrees of freedom
Multiple R-squared:  0.7426,	Adjusted R-squared:  0.741 
F-statistic: 438.6 on 1 and 152 DF,  p-value: < 2.2e-16

> summary(v2m)
 
Call:
lm(formula = DMN_FA_v2 ~ FPCustom_FA_v2, data = dti_c)
 
Residuals:
      Min        1Q    Median        3Q       Max 
-0.036853 -0.005949  0.000129  0.006172  0.024941 
 
Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)     0.07965    0.01236   6.443 1.46e-09 ***
FPCustom_FA_v2  0.83691    0.04172  20.058  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.01042 on 152 degrees of freedom
Multiple R-squared:  0.7258,	Adjusted R-squared:  0.724 
F-statistic: 402.3 on 1 and 152 DF,  p-value: < 2.2e-16

Modelo mixto

Nada raro, la relacion con el SUVR va mas o menos igual.

> model.a <- lm(Custom_FA ~ SUVR , data=idata)
> summary(model.a)
 
Call:
lm(formula = Custom_FA ~ SUVR, data = idata)
 
Residuals:
      Min        1Q    Median        3Q       Max 
-0.055969 -0.015556 -0.001559  0.015067  0.055285 
 
Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  0.308468   0.007587  40.659   <2e-16 ***
SUVR        -0.011063   0.006002  -1.843   0.0662 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.02038 on 306 degrees of freedom
Multiple R-squared:  0.01098,	Adjusted R-squared:  0.00775 
F-statistic: 3.398 on 1 and 306 DF,  p-value: 0.06625
 
> model.c <- lme(Custom_FA ~ SUVR, random = ~ 1| Subject, data=idata)
> summary(model.c)
Linear mixed-effects model fit by REML
 Data: idata 
       AIC       BIC logLik
  -1566.58 -1551.686 787.29
 
Random effects:
 Formula: ~1 | Subject
        (Intercept)  Residual
StdDev:  0.01580869 0.0129124
 
Fixed effects: Custom_FA ~ SUVR 
                  Value   Std.Error  DF  t-value p-value
(Intercept)  0.30835579 0.008911474 153 34.60211  0.0000
SUVR        -0.01097326 0.007035937 153 -1.55960  0.1209
 Correlation: 
     (Intr)
SUVR -0.986
 
Standardized Within-Group Residuals:
        Min          Q1         Med          Q3         Max 
-2.10991932 -0.50342918 -0.01151992  0.50847071  1.97546782 
 
Number of Observations: 308
Number of Groups: 154
 
> anova(model.c, model.a)
        Model df       AIC       BIC   logLik   Test  L.Ratio p-value
model.c     1  4 -1566.580 -1551.686 787.2900                        
model.a     2  3 -1500.089 -1488.918 753.0443 1 vs 2 68.49128  <.0001
 
> model.b <- lmer(Custom_FA ~ SUVR + (1| Subject), data=idata)
> summary(model.b)
Linear mixed model fit by REML ['lmerMod']
Formula: Custom_FA ~ SUVR + (1 | Subject)
   Data: idata
 
REML criterion at convergence: -1574.6
 
Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.10992 -0.50343 -0.01152  0.50847  1.97547 
 
Random effects:
 Groups   Name        Variance  Std.Dev.
 Subject  (Intercept) 0.0002499 0.01581 
 Residual             0.0001667 0.01291 
Number of obs: 308, groups:  Subject, 154
 
Fixed effects:
             Estimate Std. Error t value
(Intercept)  0.308356   0.008911   34.60
SUVR        -0.010973   0.007036   -1.56
 
Correlation of Fixed Effects:
     (Intr)
SUVR -0.986
 
> anova(model.b, model.a)
refitting model(s) with ML (instead of REML)
Data: idata
Models:
model.a: Custom_FA ~ SUVR
model.b: Custom_FA ~ SUVR + (1 | Subject)
        Df     AIC   BIC logLik deviance  Chisq Chi Df Pr(>Chisq)    
model.a  3 -1520.2 -1509 763.08  -1526.2                             
model.b  4 -1585.9 -1571 796.94  -1593.9 67.704      1  < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Modelos

Estratificando por APOE

Tenemos los genotipos de APOE en un CSV,

> head updateAPOE_FACEHBI_051218.csv
code_facehbi;Interno;APOE
F079;20120457;e3e3
F103;20121210;e3e3
F080;20121207;e2e2
F097;20130137;e3e3
F018;20130447;e3e4
F096;20130432;e3e4
F002;20131084;e3e3
F113;20131063;e3e4
F027;20131385;e3e4

Limpiamos un poco esto,

> awk -F";" {'print $1,$3'} updateAPOE_FACEHBI_051218.csv | sed 's/F//; s/ /,/; s/code_facehbi/Subject/' | sort -n > facehbi_apoe.csv
> head facehbi_apoe.csv 
Subject,APOE
001,e3e4
002,e3e3
003,e3e3
004,e3e3
005,e2e3
006,e3e4
007,e3e3
008,e3e3
009,e3e3

Lo separamos para estratificar,

> sed 's/e2e2\|e2e3/0/; s/e2e4\|e3e3/1/; s/e3e4\|e4e4/2/' facehbi_apoe.csv > facehbi_apoe_strats.csv
> head facehbi_apoe_strats.csv
Subject,APOE
001,2
002,1
003,1
004,1
005,0
006,2
007,1
008,1
009,1

Esto lo importo en R (mas o menos),

> idatawnp <- read.csv("facehbi_dmn_fa_suvr_np_tmp_v1.csv", sep=";", header=TRUE)
> iapoe <-read.csv("facehbi_apoe_strats.csv", sep=";", header= TRUE)
> idatacustom <- read.csv("facehbi_suvr_dmnfa_customfa_v12.csv", sep=";", header=TRUE)
> idatatmp <-merge(idatawnp, idatacustom, by="Subject")
> idata <-merge(idatatmp, iapoe, by="Subject")
> idata_apoe0 <- idata[idata$APOE == "0",]

Un monton de data aqui pero nos centramos en la visita basal,

> m0 <- lm(idata_apoe0$FPCustom_FA_v1 ~ idata_apoe0$Global_v1.x)
> summary(m0)
 
Call:
lm(formula = idata_apoe0$FPCustom_FA_v1 ~ idata_apoe0$Global_v1.x)
 
Residuals:
      Min        1Q    Median        3Q       Max 
-0.051447 -0.010969  0.001284  0.009302  0.040698 
 
Coefficients:
                        Estimate Std. Error t value Pr(>|t|)  
(Intercept)              0.18674    0.07971   2.343   0.0278 *
idata_apoe0$Global_v1.x  0.09175    0.06785   1.352   0.1889  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.02192 on 24 degrees of freedom
Multiple R-squared:  0.07078,	Adjusted R-squared:  0.03207 
F-statistic: 1.828 on 1 and 24 DF,  p-value: 0.1889

OK. ahi no hay nada pero esto solo era para ver si funcionaba la cosa. Vamos a organizar un poco la tabla,

> colnames(idata)
 [1] "Subject"                        "DMN_FA_v1.x"                    "Global_v1.x"                   
 [4] "Escolaridad"                    "male"                           "LL_Namingtotal_NP"             
 [7] "LL_total_NP"                    "total_NP"                       "Piramides_Y_Palmeras_Palab_FAC"
[10] "Piramides_Y_Plameras_Imag_FAC"  "Kissing_Dancing_Imagenes_FAC"   "Kissing_Dancing_Palabras_FAC"  
[13] "FNAME_TOTAL_FAC"                "Action_Namimg_Libre_FAC"        "Boston_Libre_FAC"              
[16] "Boston_Total_FAC"               "Global_v1.y"                    "Global_v2"                     
[19] "DMN_FA_v1.y"                    "DMN_FA_v2"                      "FPCustom_FA_v1"                
[22] "FPCustom_FA_v2"                 "APOE" 

Con esos nombres me voy a equivocar seguro, asi que

> colnames(idata)[colnames(idata) == "Piramides_Y_Plameras_Imag_FAC"] <- "PyP_i"
> colnames(idata)[colnames(idata) == "Piramides_Y_Palmeras_Palab_FAC"] <- "PyP_p"
> colnames(idata)[colnames(idata) == "Kissing_Dancing_Imagenes_FAC"] <- "KD_i"
> colnames(idata)[colnames(idata) == "Kissing_Dancing_Palabras_FAC"] <- "KD_p"
> colnames(idata)[colnames(idata) == "LL_Namingtotal_NP"] <- "LL_Naming"
> colnames(idata)[colnames(idata) == "LL_total_NP"] <- "LL"
> colnames(idata)[colnames(idata) == "total_NP"] <- "NP"
> colnames(idata)[colnames(idata) == "FNAME_TOTAL_FAC"] <- "FName"
> colnames(idata)[colnames(idata) == "Action_Namimg_Libre_FAC"] <- "ANaming"
> colnames(idata)[colnames(idata) == "Boston_Libre_FAC"] <- "Boston_Libre"
> colnames(idata)[colnames(idata) == "Boston_Total_FAC"] <- "Boston"
> colnames(idata)
 [1] "Subject"        "DMN_FA_v1.x"    "Global_v1.x"    "Escolaridad"    "male"           "LL_Naming"     
 [7] "LL"             "NP"             "PyP_p"          "PyP_i"          "KD_i"           "KD_p"          
[13] "FName"          "ANaming"        "Boston_Libre"   "Boston"         "Global_v1.y"    "Global_v2"     
[19] "DMN_FA_v1.y"    "DMN_FA_v2"      "FPCustom_FA_v1" "FPCustom_FA_v2" "APOE" 

y tambien,

> colnames(idata)[colnames(idata) == "DMN_FA_v1.x"] <- "DMN"
> colnames(idata)[colnames(idata) == "Global_v1.x"] <- "SUVR"
> colnames(idata)[colnames(idata) == "FPCustom_FA_v1"] <- "FPCustom"
> drops <- c("Global_v1.y", "Global_v2", "DMN_FA_v1.y", "DMN_FA_v2", "FPCustom_FA_v2")
> idata[ , !(names(idata) %in% drops)]
> okdata <- idata[ , !(names(idata) %in% drops)]

con lo cual queda mucho mas potable esto,

> colnames(okdata)
 [1] "Subject"      "DMN"          "SUVR"         "Escolaridad"  "male"         "LL_Naming"    "LL"          
 [8] "NP"           "PyP_p"        "PyP_i"        "KD_i"         "KD_p"         "FName"        "ANaming"     
[15] "Boston_Libre" "Boston"       "FPCustom"     "APOE"  

Ahora, el objetivo es estudiar las variables neurosicologicas como funcion de la FA en las redes DTI (DMN y FPCustom), estratificando por APOE y usando como covariables SUVR, Genero, Escolaridad y Edad. m( joer q me falta la edad.

> awk -F";" '{print $2,$8}' faceHBI_matriuREF_14-1-19-1.csv | sed 's/ /;/; s/edat/Edad/; s/subject/Subject/' > edad.csv
> scp -P 20022 edad.csv detritus.fundacioace.org:facehbi/dti_model/
edad.csv                                                      100% 1312   200.6KB/s   00:00  . 

Vale,la meto

> edad_dlc <-read.csv("edad.csv", sep=";", header=TRUE)
> okdata <- merge(okdata, edad_dlc, by = "Subject")
> colnames(okdata)
 [1] "Subject"      "DMN"          "SUVR"         "Escolaridad"  "male"         "LL_Naming"    "LL"          
 [8] "NP"           "PyP_p"        "PyP_i"        "KD_i"         "KD_p"         "FName"        "ANaming"     
[15] "Boston_Libre" "Boston"       "FPCustom"     "APOE"         "Edad" 

Ahora si. Pero sigue siendo una puñeta. Mejor me hago un script que corra a traves de los modelos ;-). Venga, primero a lo bruto, just in case,

library(QuantPsyc)
x<-read.csv("facehbi_dti_np.csv")
Color=c("red","blue")
scan("npvars.names", what = character())->np
scan("nivars.names", what = character())->ni
sink(file = "facehbi_dti_np_models.txt", append = TRUE, type = "output", split = TRUE)
 
for(i in 1:length(np)){
        for(j in 1:length(ni)){
                y.data <- x[c(ni[j], np[i], "male", "Edad", "Escolaridad", "SUVR")]
                y.data <- y.data[complete.cases(y.data),]
                a <- lm( paste ('y.data$', np[i], ' ~ y.data$', ni[j], ' + y.data$male + y.data$Edad + y.data$Escolaridad + y.data$SUVR'))
                writeLines(paste("NP: ", np[i], " NI: ", ni[j]))
                writeLines(paste("R2: ", summary(a)$adj.r.squared, " p-value: ", summary(a)$coef[2,4]))
                beta <- lm.beta(a)
                for(k in 1:length(beta)){
                        writeLines(paste(names(beta[k]), ": ", beta[k]))
                }
                writeLines(paste("-------"))
        }
}
sink()

a ver,

> write.csv(okdata, "facehbi_dti_np.csv")
> source("get_lms.r")
Read 11 items
Read 2 items

Un vistazo a la salida, y effectivamente, no hay nada que hacer aqui. el mas prometedor es este,

NP:  NP  NI:  FPCustom
R2:  0.245411609394026  p-value:  0.0332379182003353

Puaf. Vamos a estratificar aver,

 okdata0 <- okdata[okdata$APOE == "0",]
> write.csv(okdata0, "facehbi_dti_np.csv")
> source("get_lms.r")
Read 11 items
Read 2 items

Poca cosa aqui,

NP:  NP  NI:  DMN
R2:  0.283699960037635  p-value:  0.624010265461206
NP:  FName  NI:  DMN
R2:  0.260050114900356  p-value:  0.141503076825646

Seguimos,

 okdata1 <- okdata[okdata$APOE == "1",]
> write.csv(okdata1, "facehbi_dti_np.csv")
> source("get_lms.r")
Read 11 items
Read 2 items

Nop.

 okdata2 <- okdata[okdata$APOE == "2",]
> write.csv(okdata2, "facehbi_dti_np.csv")
> source("get_lms.r")
Read 11 items
Read 2 items

grrrrrrr…

NP:  LL_Naming  NI:  DMN
R2:  0.338308455461783  p-value:  0.00172399799381871
NP:  NP  NI:  FPCustom
R2:  0.423468586516873  p-value:  0.344571425163692
NP:  KD_p  NI:  FPCustom
R2:  0.300292939478652  p-value:  0.382224137819204
NP:  KD_i  NI:  FPCustom
R2:  0.3278128869582  p-value:  0.511816038951647
NP:  FName  NI:  FPCustom
R2:  0.377022589868983  p-value:  0.457838263496454
NP:  Boston  NI:  FPCustom
R2:  0.305029202245965  p-value:  0.251505671812935
NP:  Boston_Libre  NI:  FPCustom
R2:  0.349670089711  p-value:  0.127585742382962

Odio estas cosas, todo porqueria que hay que mirar.

> m0 <- lm(okdata2$LL_Naming ~ okdata2$DMN + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR +okdata2$male)
> summary(m0)
 
Call:
lm(formula = okdata2$LL_Naming ~ okdata2$DMN + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)
 
Residuals:
    Min      1Q  Median      3Q     Max 
-1.0747 -0.2115  0.1003  0.2248  0.4808 
 
Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          20.55919    1.58086  13.005 2.56e-14 ***
okdata2$DMN         -13.72502    4.01223  -3.421  0.00172 ** 
okdata2$Escolaridad   0.04295    0.01474   2.913  0.00648 ** 
okdata2$Edad         -0.03229    0.00971  -3.326  0.00222 ** 
okdata2$SUVR          0.16477    0.26430   0.623  0.53743    
okdata2$male          0.02292    0.13024   0.176  0.86142    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.3551 on 32 degrees of freedom
Multiple R-squared:  0.4277,	Adjusted R-squared:  0.3383 
F-statistic: 4.783 on 5 and 32 DF,  p-value: 0.002238
 
> m0 <- lm(okdata2$LL_Naming ~ okdata2$DMN + okdata2$Escolaridad + okdata2$Edad)
> summary(m0)
 
Call:
lm(formula = okdata2$LL_Naming ~ okdata2$DMN + okdata2$Escolaridad + 
    okdata2$Edad)
 
Residuals:
     Min       1Q   Median       3Q      Max 
-1.13733 -0.21729  0.09543  0.23242  0.42227 
 
Coefficients:
                      Estimate Std. Error t value Pr(>|t|)    
(Intercept)          20.529979   1.521191  13.496 3.25e-15 ***
okdata2$DMN         -13.286680   3.840520  -3.460  0.00148 ** 
okdata2$Escolaridad   0.041280   0.014099   2.928  0.00605 ** 
okdata2$Edad         -0.030158   0.008949  -3.370  0.00188 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
 
Residual standard error: 0.3469 on 34 degrees of freedom
Multiple R-squared:  0.4198,	Adjusted R-squared:  0.3686 
F-statistic: 8.201 on 3 and 34 DF,  p-value: 0.0003054

Baaaah, pero he escrito papers con menos. Esta es la variable que se llama en la tabla LLNaming_total_NP, a saber que sera.

Seguimos, la segunda asociacion no dice nada,

> m0 <- lm(okdata2$NP ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$NP ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
    Min      1Q  Median      3Q     Max 
-26.514  -7.122   1.120   5.021  32.711 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)         176.2271    44.0144   4.004 0.000346 ***
okdata2$FPCustom    121.2372   126.3732   0.959 0.344571    
okdata2$Escolaridad   1.2925     0.5144   2.513 0.017215 *  
okdata2$Edad         -0.7723     0.3146  -2.455 0.019701 *  
okdata2$SUVR         -5.1581     9.0583  -0.569 0.573036    
okdata2$male         -7.5469     4.4865  -1.682 0.102280    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 12.25 on 32 degrees of freedom
Multiple R-squared:  0.5014,	Adjusted R-squared:  0.4235 
F-statistic: 6.435 on 5 and 32 DF,  p-value: 0.0003044

> m0 <- lm(okdata2$NP ~ okdata2$Escolaridad + okdata2$Edad)
> summary(m0)

Call:
lm(formula = okdata2$NP ~ okdata2$Escolaridad + okdata2$Edad)

Residuals:
    Min      1Q  Median      3Q     Max 
-27.939  -8.188  -0.075   6.659  37.356 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)         209.0838    23.1923   9.015 1.19e-10 ***
okdata2$Escolaridad   1.4459     0.4692   3.082  0.00400 ** 
okdata2$Edad         -0.9078     0.3010  -3.016  0.00475 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 12.41 on 35 degrees of freedom
Multiple R-squared:  0.4404,	Adjusted R-squared:  0.4084 
F-statistic: 13.77 on 2 and 35 DF,  p-value: 3.87e-05

la variable NP esta asociada casi exclusivamente, en este grupo a la edad y la escolaridad del sujeto. Lomismo vale para el resto,

> m0 <- lm(okdata2$KD_p ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$KD_p ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.0610 -0.7722  0.1195  0.6817  2.1740 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          59.01916    8.95276   6.592 6.19e-06 ***
okdata2$FPCustom    -23.02432   25.62433  -0.899   0.3822    
okdata2$Escolaridad   0.27985    0.09696   2.886   0.0107 *  
okdata2$Edad         -0.08073    0.05345  -1.511   0.1504    
okdata2$SUVR         -0.74099    1.30206  -0.569   0.5772    
okdata2$male         -0.11501    0.77994  -0.147   0.8846    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.562 on 16 degrees of freedom
  (16 observations deleted due to missingness)
Multiple R-squared:  0.4669,	Adjusted R-squared:  0.3003 
F-statistic: 2.803 on 5 and 16 DF,  p-value: 0.05284

> m0 <- lm(okdata2$KD_i ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$KD_i ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.9060 -0.8693  0.0712  1.3062  3.7389 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          67.92106   14.52237   4.677 0.000217 ***
okdata2$FPCustom    -27.68934   41.32332  -0.670 0.511816    
okdata2$Escolaridad   0.31847    0.15604   2.041 0.057094 .  
okdata2$Edad         -0.25900    0.08683  -2.983 0.008356 ** 
okdata2$SUVR          1.23577    2.07456   0.596 0.559237    
okdata2$male          1.14538    1.26034   0.909 0.376169    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.539 on 17 degrees of freedom
  (15 observations deleted due to missingness)
Multiple R-squared:  0.4806,	Adjusted R-squared:  0.3278 
F-statistic: 3.146 on 5 and 17 DF,  p-value: 0.03432

> m0 <- lm(okdata2$FName ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$FName ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
    Min      1Q  Median      3Q     Max 
-44.333  -7.319  -1.934   6.399  33.950 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)   
(Intercept)          143.2128    55.3000   2.590  0.01434 * 
okdata2$FPCustom    -119.3217   158.7761  -0.752  0.45784   
okdata2$Escolaridad    1.3655     0.6463   2.113  0.04251 * 
okdata2$Edad          -1.1730     0.3952  -2.968  0.00564 **
okdata2$SUVR         -11.9182    11.3809  -1.047  0.30285   
okdata2$male          -7.9011     5.6369  -1.402  0.17064   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 15.39 on 32 degrees of freedom
Multiple R-squared:  0.4612,	Adjusted R-squared:  0.377 
F-statistic: 5.478 on 5 and 32 DF,  p-value: 0.0009417

> m0 <- lm(okdata2$Boston ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$Boston ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
     Min       1Q   Median       3Q      Max 
-15.5752  -1.7398   0.0496   2.4670   7.4385 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)          72.9727    15.7856   4.623 6.31e-05 ***
okdata2$FPCustom    -52.9013    45.2722  -1.169  0.25151    
okdata2$Escolaridad   0.5936     0.1842   3.223  0.00298 ** 
okdata2$Edad         -0.2660     0.1129  -2.357  0.02493 *  
okdata2$SUVR          2.6940     3.2501   0.829  0.41350    
okdata2$male         -0.4492     1.6135  -0.278  0.78257    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.386 on 31 degrees of freedom
  (1 observation deleted due to missingness)
Multiple R-squared:  0.4016,	Adjusted R-squared:  0.305 
F-statistic:  4.16 on 5 and 31 DF,  p-value: 0.005237

> m0 <- lm(okdata2$Boston_Libre ~ okdata2$FPCustom + okdata2$Escolaridad + okdata2$Edad+ okdata2$SUVR + okdata2$male)
> summary(m0)

Call:
lm(formula = okdata2$Boston_Libre ~ okdata2$FPCustom + okdata2$Escolaridad + 
    okdata2$Edad + okdata2$SUVR + okdata2$male)

Residuals:
    Min      1Q  Median      3Q     Max 
-9.6225 -1.9833 -0.2538  1.9214  6.9440 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          76.69651   12.78225   6.000 1.08e-06 ***
okdata2$FPCustom    -57.40968   36.70011  -1.564  0.12759    
okdata2$Escolaridad   0.49401    0.14938   3.307  0.00234 ** 
okdata2$Edad         -0.26234    0.09136  -2.872  0.00719 ** 
okdata2$SUVR          1.92228    2.63062   0.731  0.47026    
okdata2$male         -0.14670    1.30294  -0.113  0.91106    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 3.558 on 32 degrees of freedom
Multiple R-squared:  0.4376,	Adjusted R-squared:  0.3497 
F-statistic: 4.979 on 5 and 32 DF,  p-value: 0.001748

Resumiendo, En los sujetos con el alelo $\epsilon$-4 presente, hay una ligera asociacion entre la variable LLNaming_total_NP y la fraccion de anisotropia en la Default Mode Network, covariada por edady escolaridad de los sujetos de estudio. Todos los demas efectos que puedan observarse son debido al efecto de asociacion de los resultados de los test con la edad y escolaridad de los sujetos. :-P

Cosas raras

Si hacemos un modelo con todo observamos una cosa curiosa,

> m1 <- lm(okdata$LL_Naming ~ okdata$FPCustom + okdata$SUVR + okdata$Escolaridad + okdata$male +okdata$Edad + okdata$APOE)
> summary(m1)

Call:
lm(formula = okdata$LL_Naming ~ okdata$FPCustom + okdata$SUVR + 
    okdata$Escolaridad + okdata$male + okdata$Edad + okdata$APOE)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.75246  0.02246  0.08062  0.14430  0.30217 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)    
(Intercept)        16.118960   0.507667  31.751   <2e-16 ***
okdata$FPCustom    -2.409148   1.269032  -1.898   0.0596 .  
okdata$SUVR        -0.173373   0.179192  -0.968   0.3349    
okdata$Escolaridad  0.005906   0.005947   0.993   0.3223    
okdata$male        -0.006072   0.056222  -0.108   0.9141    
okdata$Edad        -0.005047   0.003865  -1.306   0.1936    
okdata$APOE        -0.054227   0.045219  -1.199   0.2324    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.3216 on 147 degrees of freedom
Multiple R-squared:  0.06926,	Adjusted R-squared:  0.03127 
F-statistic: 1.823 on 6 and 147 DF,  p-value: 0.09836

El modelo es una porqueria pero los residuales indican que hay un estratificacion para la variable LL_Naming,

Y tal y como indica el ajuste,

> m1 <- lm(okdata$LL_Naming ~ okdata$FPCustom)
> summary(m1)

Call:
lm(formula = okdata$LL_Naming ~ okdata$FPCustom)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.88482  0.05739  0.08721  0.12847  0.22546 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)    
(Intercept)      15.5248     0.3728  41.647   <2e-16 ***
okdata$FPCustom  -2.1407     1.2661  -1.691   0.0929 .  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.3248 on 152 degrees of freedom
Multiple R-squared:  0.01846,	Adjusted R-squared:  0.012 
F-statistic: 2.859 on 1 and 152 DF,  p-value: 0.09294

Nota: Estas variables estan muy sesgadas por los valores dados en la clinica, de ahi este comportamiento raro.

Composites

Vamos a intentar el procedimiento utilizando los composites de NP. Las variables de los composites son,

funcioExecutiva_fluencia
funcioExecutiva_velocprocess_IM
funcioExecutiva_atencio
memoria_fnameProf
memoria_fnameNom
memoria_wms
memoria_rbans
gnosia
praxia

Saco las variables que necesito,

> awk -F";" '{print $2,$8,$9,$11,$338,$340,$341,$342,$343,$344,$345,$346,$347,$350,$412}' faceHBI_matriuREF_14-1-19-1.csv | sed 's/ /;/g; s/edat/Edad/; s/subject/Subject/; s/Anyos_Escolaridad_FAC/Escolaridad/; s/Sex_1H_0M/female/; s/Global_v1/SUVR/' > facehbi_data.csv
> awk -F";" {'print $1,$3'} updateAPOE_FACEHBI_051218.csv | sed 's/F//; s/ /,/; s/code_facehbi/Subject/' | sort -n > facehbi_apoe.csv
> sed 's/e2e2\|e2e3/0/; s/e2e4\|e3e3/1/; s/e3e4\|e4e4/2/; s/,/;/' facehbi_apoe.csv > facehbi_apoe_strats.csv
> scp -P 20022 facehbi_apoe_strats.csv detritus.fundacioace.org:facehbi/dti_model/
facehbi_apoe_strats.csv             
> scp -P 20022 facehbi_data.csv detritus.fundacioace.org:facehbi/dti_model/
............
facehbi_[osotolongo@detritus dti_model]$ awk -F";" '{print $1,$4,$6}' facehbi_suvr_dmnfa_customfa_v12.csv | sed 's/ /;/g; s/DMN_FA_v1/DMN/; s/FPCustom_FA_v1/FPCustom/' > facehbi_dti.csv
data.csv

La importo en R,

> fdata <-read.csv("facehbi_data.csv", sep = ";", header=TRUE)
> fapoe <-read.csv("facehbi_apoe_strats.csv", sep = ";", header=TRUE)
> fdti <- read.csv("facehbi_dti.csv", sep=";", header=TRUE)
> okdata <- merge(fdata, fapoe, by = "Subject")
> okdata <- merge(okdata, fdti, by = "Subject")

preparo la lista de variables,

[osotolongo@detritus dti_model]$ cat npvars.names 
funcioExecutiva_fluencia
funcioExecutiva_velocprocess_IM
funcioExecutiva_atencio
memoria_fnameProf
memoria_fnameNom
memoria_wms
memoria_rbans
gnosia
praxia
[osotolongo@detritus dti_model]$ cat nivars.names 
DMN
FPCustom

y hacemos la primera prueba,

> write.csv(okdata, file="facehbi_dti_np.csv")
> source("get_lms.r")
Read 10 items
Read 2 items

No muy bien hasta ahora,

NP:  funcioExecutiva_velocprocess_IM  NI:  FPCustom
R2:  0.315511475678145  p-value:  0.240278609939977
neuroimagen/altdti.1552561541.txt.gz · Last modified: 2020/08/04 10:46 (external edit)