neuroimagen:bioface_sbm
Table of Contents
Analisis SBM con Freesurfer (a.k.a. FSGA) para BIOFACE
Haciendo los archivos
Para construir los archivos .fsgd a partir de la base de datos han de hacerse una serie de pasos. Voy a poner el ejemplo con la variable RLPP. Para cualquier otra variable simplemente se cambia.
[osotolongo@brick03 fsga2]$ sed 's/;/,/g' ../bioface_mri.csv > bioface_join.csv [osotolongo@brick03 fsga2]$ awk -F"," '{print $3","$5","$6","$7","$21}' ../bioface_np.csv > bioface_rlpp.csv [osotolongo@brick03 fsga2]$ sed -i 's/Subjecte/PSubject/g' bioface_rlpp.csv [osotolongo@brick03 fsga2]$ join -t"," -1 2 -2 1 bioface_join.csv bioface_rlpp.csv > bioface_data_rlpp.csv [osotolongo@brick03 fsga2]$ awk -F"," '{print "bioface_"$2","$3","$4","$5","$6}' bioface_data_rlpp.csv | sed 's/bioface_Subject/Variables/;s/Edad/Age/;s/Sexo/Gender/;s/Años_Estudios/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > rlpp_body.csv [osotolongo@brick03 fsga2]$ head rlpp_body.csv Variables Age Gender Education RLPP Input bioface_0001 Main 64 1 11 2 Input bioface_0002 Main 53 2 9 10 Input bioface_0003 Main 60 2 8 3 Input bioface_0004 Main 64 2 7 3 Input bioface_0005 Main 62 2 16 7 Input bioface_0006 Main 64 2 5 7 Input bioface_0007 Main 62 2 11 5 Input bioface_0008 Main 58 2 11 8 Input bioface_0009 Main 51 2 11 7 [osotolongo@brick03 fsga2]$ cat headers_rlpp.txt GroupDescriptorFile 1 Title BIOFACE_rlpp Class Main [osotolongo@brick03 fsga2]$ cat headers_rlpp.txt rlpp_body.csv > rlpp.fsgd
y por otra parte el modelo,
[osotolongo@brick03 fsga2]$ cat rlpp.mtx 0 0 0 0 1
Preparando la ejecucion
Ahora voy a hacer un script que corra todo el proceso, para no saltarme ningun paso
[osotolongo@brick03 fsga2]$ cat fsga_rlpp.sh #!/bin/bash mris_preproc --fsgd rlpp.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.rlpp.thickness.10.mgh mris_preproc --fsgd rlpp.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.rlpp.thickness.10.mgh mri_glmfit --y lh.rlpp.thickness.10.mgh --fsgd rlpp.fsgd --C rlpp.mtx --surf fsaverage lh --glmdir lh.rlpp.glmdir mri_glmfit --y rh.rlpp.thickness.10.mgh --fsgd rlpp.fsgd --C rlpp.mtx --surf fsaverage rh --glmdir rh.rlpp.glmdir mri_glmfit-sim --glmdir lh.rlpp.glmdir --cache 2 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.rlpp.glmdir --cache 2 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.rlpp.glmdir --cache 2 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.rlpp.glmdir --cache 2 pos --cwp 0.05 --2spaces
y se lanza todo con,
[osotolongo@brick03 fsga2]$ ./fsga_rlpp.sh
Buscando resultados
Para encontrar los clusters hacemos,
[osotolongo@brick03 fsga2]$ for x in `find ./ -name "*.summary"`; do echo ${x}; grep -v "^#" ${x}; done ./lh.rlpn.glmdir/rlpn/cache.th20.neg.sig.cluster.summary ./lh.rlpn.glmdir/rlpn/cache.th20.pos.sig.cluster.summary ./lh.rlpp.glmdir/rlpp/cache.th20.neg.sig.cluster.summary ./lh.rlpp.glmdir/rlpp/cache.th20.pos.sig.cluster.summary 1 3.519 112413 3477.26 -18.5 -75.6 28.2 0.00020 0.00000 0.00040 6595 15036.07 precuneus 2 3.947 157071 2138.53 -58.6 -51.1 1.7 0.00020 0.00000 0.00040 4377 10588.38 bankssts 3 4.314 99360 1135.62 -44.3 7.6 5.7 0.00878 0.00719 0.01057 2749 7526.91 parsopercularis ./rh.rlpn.glmdir/rlpn/cache.th20.neg.sig.cluster.summary ./rh.rlpn.glmdir/rlpn/cache.th20.pos.sig.cluster.summary ./rh.rlpp.glmdir/rlpp/cache.th20.neg.sig.cluster.summary ./rh.rlpp.glmdir/rlpp/cache.th20.pos.sig.cluster.summary 1 4.982 104537 7071.55 46.5 -57.5 23.1 0.00020 0.00000 0.00040 14015 35733.48 inferiorparietal 2 3.430 29029 2647.45 27.8 19.3 42.7 0.00020 0.00000 0.00040 5370 11735.45 caudalmiddlefrontal 3 4.229 45003 1848.86 12.5 -52.8 30.8 0.00020 0.00000 0.00040 4269 10182.11 precuneus 4 2.481 143807 1116.95 30.5 -49.7 -13.8 0.01276 0.01077 0.01475 2172 4473.95 fusiform
Ahora, en el hemisferio/modelo que hemos encontrado resultados podemos identificar estos clusters con,
$ freeview -f $SUBJECTS_DIR/fsaverage/surf/rh.inflated:overlay=rh.rlpp.glmdir/rlpp/cache.th20.pos.sig.cluster.mgh:overlay_threshold=0.1,5:annot=rh.rlpp.glmdir/rlpp/cache.th20.pos.sig.ocn.annot -viewport 3d &
Para superponer el atlas DK,
$ freeview -f $SUBJECTS_DIR/fsaverage/surf/rh.inflated:overlay=rh.rlpp.glmdir/rlpp/cache.th20.pos.sig.cluster.mgh:overlay_threshold=0.1,5:annot=aparc.annot:annot_outline=1 -viewport 3d &
y para ver las areas de Brodman afectadas,
$ freeview -f $SUBJECTS_DIR/fsaverage/surf/rh.inflated:overlay=rh.rlpp.glmdir/rlpp/cache.th20.pos.sig.cluster.mgh:overlay_threshold=0.1,5:annot=PALS_B12_Brodmann.annot:annot_outline=1 -viewport 3d &
Cambiando a Volumen
Para area y/o otra variable es los mismo, con un solo ejemplo ya la idea la tenemos. Hay que hacer las 9 combinaciones pero no deberia tomar demasiado.
Extraer ICV
[osotolongo@brick03 fsga]$ awk -F"," '{print $1","$67}' ../fsrecon/aseg_stats.csv > bioface_icv.csv [osotolongo@brick03 fsga]$ head bioface_icv.csv Subject,EstimatedTotalIntraCranialVol 0001,1486505.98285 0002,1249118.11926 0003,1333969.51364 0004,1344500.90389 0005,1408317.34379 0006,1321218.04232 0007,1501903.59123 0008,1480617.98715 0009,1368664.47062 [osotolongo@brick03 fsga]$ cat convert_icv.r setwd("/nas/data/bioface/fsga") read.csv("bioface_icv.csv") -> x x$zICV <- (x$EstimatedTotalIntraCranialVol - mean(x$EstimatedTotalIntraCranialVol))/sd(x$EstimatedTotalIntraCranialVol) write.csv(x, file="bioface_zic_tempv.csv", row.names=FALSE) [osotolongo@brick03 fsga]$ Rscript convert_icv.r [osotolongo@brick03 fsga]$ head bioface_zic_tempv.csv "Subject","EstimatedTotalIntraCranialVol","zICV" 1,1486505.98285,0.41844790218558 2,1249118.11926,-1.20476544891787 3,1333969.51364,-0.624567648851558 4,1344500.90389,-0.552555992480593 5,1408317.34379,-0.116191236089841 6,1321218.04232,-0.711759800041669 7,1501903.59123,0.523733838540968 8,1480617.98715,0.378186900266852 9,1368664.47062,-0.387330089465778 [osotolongo@brick03 fsga]$ awk -F"," '{$1 = sprintf("%04d",$1); print $1","$3}' bioface_zic_tempv.csv | sed 's/"//g;s/0000/Subject/' > bioface_zicv.csv [osotolongo@brick03 fsga]$ head bioface_zicv.csv Subject,zICV 0001,0.41844790218558 0002,-1.20476544891787 0003,-0.624567648851558 0004,-0.552555992480593 0005,-0.116191236089841 0006,-0.711759800041669 0007,0.523733838540968 0008,0.378186900266852 0009,-0.387330089465778
Integrando con los datos
[osotolongo@brick03 fsga]$ join -t"," -1 2 -2 1 bioface_data.csv bioface_zicv.csv > bioface_vdata.csv [osotolongo@brick03 fsga]$ head bioface_vdata.csv Subject,PSubject,Edad,Sexo,Estudios,FACEmem,zICV 0001,B001,64,1,11,13,0.41844790218558 0002,B002,53,2,9,48,-1.20476544891787 0003,B003,60,2,8,10,-0.624567648851558 0004,B004,64,2,7,21,-0.552555992480593 0005,B005,62,2,16,30,-0.116191236089841 0006,B006,64,2,5,25,-0.711759800041669 0007,B007,62,2,11,39,0.523733838540968 0008,B008,58,2,11,45,0.378186900266852 0009,B009,51,2,11,28,-0.387330089465778 [osotolongo@brick03 fsga]$ awk -F"," '{print "bioface_"$1","$3","$4","$5","$6","$7}' bioface_vdata.csv | sed 's/bioface_Subject/Variables/;s/Edad/Age/;s/Sexo/Gender/;s/Estudios/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > memv_body.csv [osotolongo@brick03 fsga]$ head memv_body.csv Variables Age Gender Education FACEmem zICV Input bioface_0001 Main 64 1 11 13 0.41844790218558 Input bioface_0002 Main 53 2 9 48 -1.20476544891787 Input bioface_0003 Main 60 2 8 10 -0.624567648851558 Input bioface_0004 Main 64 2 7 21 -0.552555992480593 Input bioface_0005 Main 62 2 16 30 -0.116191236089841 Input bioface_0006 Main 64 2 5 25 -0.711759800041669 Input bioface_0007 Main 62 2 11 39 0.523733838540968 Input bioface_0008 Main 58 2 11 45 0.378186900266852 Input bioface_0009 Main 51 2 11 28 -0.387330089465778 [osotolongo@brick03 fsga]$ cat headers.txt GroupDescriptorFile 1 Title BIOFACE_mem Class Main [osotolongo@brick03 fsga]$ cat headers.txt memv_body.csv > memv.fsgd [osotolongo@brick03 fsga]$ head memv.fsgd GroupDescriptorFile 1 Title BIOFACE_mem Class Main Variables Age Gender Education FACEmem zICV Input bioface_0001 Main 64 1 11 13 0.41844790218558 Input bioface_0002 Main 53 2 9 48 -1.20476544891787 Input bioface_0003 Main 60 2 8 10 -0.624567648851558 Input bioface_0004 Main 64 2 7 21 -0.552555992480593 Input bioface_0005 Main 62 2 16 30 -0.116191236089841 Input bioface_0006 Main 64 2 5 25 -0.711759800041669
Preparando el run
[osotolongo@brick03 fsga]$ cat memv.mtx 0 0 0 0 1 0
- fsgav.sh
#!/bin/bash mris_preproc --fsgd memv.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.bioface.volume.10.mgh mris_preproc --fsgd memv.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.bioface.volume.10.mgh mri_glmfit --y lh.bioface.volume.10.mgh --fsgd memv.fsgd --C memv.mtx --surf fsaverage lh --glmdir lh.bioface.vglmdir mri_glmfit --y rh.bioface.volume.10.mgh --fsgd memv.fsgd --C memv.mtx --surf fsaverage rh --glmdir rh.bioface.vglmdir mri_glmfit-sim --glmdir lh.bioface.vglmdir --cache 2 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.bioface.vglmdir --cache 2 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.bioface.vglmdir --cache 2 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.bioface.vglmdir --cache 2 pos --cwp 0.05 --2spaces
Sacar resultados
Ahora,
$ ./fsgav.sh
y buscamos los resultados,
[osotolongo@brick03 fsga]$ for x in `find *vglmdir/memv/ -name "*.th20.*.summary"`; do echo ${x}; grep -v "^#" ${x}; done lh.bioface.vglmdir/memv/cache.th20.neg.sig.cluster.summary lh.bioface.vglmdir/memv/cache.th20.pos.sig.cluster.summary 1 3.959 79787 828.17 -32.6 -76.7 -11.8 0.00180 0.00100 0.00260 1178 3049.31 fusiform 2 2.513 56552 625.22 -47.9 -26.0 -10.1 0.01574 0.01355 0.01792 1292 2610.65 middletemporal rh.bioface.vglmdir/memv/cache.th20.neg.sig.cluster.summary rh.bioface.vglmdir/memv/cache.th20.pos.sig.cluster.summary 1 5.991 47819 954.75 45.3 -50.6 24.1 0.00060 0.00020 0.00100 2214 6368.82 inferiorparietal 2 3.719 5713 674.59 54.4 3.4 -15.3 0.01157 0.00958 0.01355 1372 3261.33 superiortemporal
SBM vs NBACE
Antes de ejecutar el analisis para las variabes de NBACE, vamos a mirar que hay en la base de datos,
[osotolongo@brick03 fsga_nbace]$ head -n 1 bioface_full_db.csv | sed 's/,/\n/g' | cat -n 1 SUJETO 2 Interno 3 INV_Edinburgh 4 ISIS 5 Suenyo_Pitsb 6 P_Orofonatorias 7 R_Pseudopalabras 8 R_Palabras 9 R_Frases 10 C_Palabras 11 C_Ordenes 12 Proces_Sintactico 13 L_Numeros 14 L_Logotomos 15 L_Textos 16 Escritura 17 TAP 18 Secuencias_Posturas_D 19 Secuencias_Posturas_I 20 Grafestesia_D 21 Grafestesia_I 22 Recon_Digital_D 23 Recon_Digital_I 24 Orientacion_D_I 25 Probl_aritmeticos 26 Num_LCR 27 PLASMA_Exosomes 28 PrimarioÚltimo 29 Primera_Prueba 30 Prueba 31 EIN1 32 EIP1 33 EIN2 34 EIP2 35 RCPN 36 RCPP 37 RF 38 RLPN 39 RLPP 40 NEM 41 PEM 42 Resultado 43 RLPTOTAL 44 TABLETAPTOTAL 45 Fecha_NPS_adicional 46 FCSRT_Total_Free_Recall 47 FCSRT_Total_Recall 48 FCSRT_Delayed_Free_Recall 49 FCSRT_Delayed_Total_Recall 50 ROCFT_Delayed 51 DMS_48_SET1 52 DMS_48_SET2 53 TMT_A_Time 54 TMT_A_errors 55 TMT_B_Time 56 TMT_B_errors 57 Fluency_M 58 Fluency_R 59 Fluency_P 60 Stroop_P 61 Stroop_C 62 Stroop_PC 63 Stroop_Formula 64 Stroop_Interferencia 65 ROCFT_Copy 66 ROCFT_Time 67 VOSP_Incomp_Letters 68 VOSP_Numb_Locat 69 JLO 70 Boston_free 71 Boston_free_CS 72 Boston_free_CS_CF 73 Piram_Palm_Imagenes 74 Piram_Palm_palabras 75 ATN_grupo 76 PrimaryLast 77 Fecha_Nacimiento 78 Sexo 79 Diagnostico_Sindromico 80 Diagnost_secund_agrupado 81 Diagnostico_secundario 82 Diagnostico_NPS 83 Anyos_escolaridad 84 Nivel_escolaridad 85 Situacion_laboral 86 Convivencia 87 Estado_civil 88 Lengua_materna 89 Bilingue 90 Lateralidad 91 Fecha_NRL 92 Fecha_Inicio_sintomas 93 Edad_inicio 94 Edad_visita 95 Sintoma_Inicio1 96 Sintoma_Inicio2 97 Sintoma_Inicio3 98 Sintoma_Inicio4 99 Sintoma_Inicio5 100 Sintoma_Inicio6 101 Sintoma_Inicio7 102 Progresion 103 Fluctuacion 104 HTA 105 DL 106 DM 107 EAP 108 Ictus_isq 109 Ictus_hem 110 AIT 111 Cardiopatia_isq 112 Tabaquismo 113 SAHS 114 EPOC 115 ASMA 116 TAnsiedad_generalizada 117 TDepresivo_persistente 118 EnfPsiOTRAS 119 Sind_AnsioDepre 120 Fibromialgia 121 Fatiga_cronica 122 SQM 123 Enf_autoinmunes 124 RBD 125 APNEASSUEnyO 126 Roncopatia 127 Bruxismo 128 MovVigorSuenyos 129 Sominiloquios 130 ParaSuenyViole 131 Insomnio 132 SomnolecenciaDiurna 133 Parasomnias 134 SPI 135 Migranya_aura 136 Migranya_sin_aura 137 Migranya_cronica 138 Migranya_Episodica 139 IRC 140 hipotiroidismo 141 Hipo_B12 142 Consumo_riesgo_OH 143 Consumo_otras_sustancias 144 m1.Antiagregantes 145 m2.Antihipertensivos 146 m3.Antidiabeticos1 147 m4.Hipolipemientes2 148 m5.Antidepresivos 149 m6.Ansioliticos 150 m7.Antiepilecticos 151 m8.Antiparkinsonianos 152 m9.Neurolecticos 153 m10.Antiinflamatorios 154 m11.Analgesicos 155 m12.Opiaceos 156 m13.Souvenaid 157 m14.GinkgoBiloba 158 m.15HormonasTiroideas 159 m16.Hinopticos 160 m17.B12 161 M18.Antihistaminicos 162 Antecedentes_familiares_1G 163 Antecedentes_familiares_presenil 164 IPAQ 165 UPDRSS_III 166 MMSE 167 CDR 168 Talla 169 Peso 170 IMC 171 Perimetro_abdominal 172 TA_S 173 TA_D 174 Obesidad_abdominal 175 Obesidad 176 Atrofia_global 177 Atrofia_parietal_D 178 Atrofia_parietal_I 179 Atrofia_temporal_D 180 Atrofia_temporal_I 181 Fazekas_profundo_D 182 Fazekas_profundo_I 183 Fazekas_periventricular 184 ARWMC_Frontal_D 185 ARWMC_Frontal_I 186 ARWMC_Parieto_occipital_D 187 ARWMC_Parieto_occipital_I 188 ARWMC_Temporal_D 189 ARWMC_Temporal_I 190 ARWMC_Basal_ganglia 191 ARWMC_infratentorial 192 ARWMC_total 193 EspaciosPV_centrosemioval 194 EspaciosPV_gangliosbasales 195 Infartos_lacunares 196 Microhemorragias_menor5_profundas 197 Microhemorragias_menor5_corticales 198 Microhemorragias_5_10_profundas 199 Microhemorragias_5_10_corticales 200 Siderosis_superficial 201 APOE 202 E4_positivo 203 PL 204 Fecha_PL 205 ABETA42_LUMI 206 AB42resAB40 207 LumipTau 208 LumihTau 209 Lumi_Abeta_ratio_dicot 210 Lumi_T_tau_dicot 211 Lumi_P_tau_dicot 212 Lumi_ATNr 213 FechaNeuropsicologia 214 O_Total_NP 215 O_T_NP 216 O_E_NP 217 O_P_NP 218 M_numparaules_NP 219 M_WMS_total_NP 220 M_ret_NP 221 M_recon_NP 222 M_f_rec_NP 223 M_digtspan_direct_NP 224 M_digtspan_invers_NP 225 M_digttotal_direct_NP 226 M_digttotal_invers_NP 227 LL_Namingtotal_NP 228 LL_comprensio_NP 229 LL_R_total_NP 230 LL_escriptura_NP 231 G_pop_total_NP 232 G_15obj_NP 233 G_15obj_nf_NP 234 G_Luria_NP 235 P_Ecopraxiatotal_NP 236 P_ideo_total_NP 237 P_constr_total_NP 238 FE_SKTtemps_NP 239 FE_SKTerrors_NP 240 FE_Pflu_NP 241 FE_anflu_NP 242 Fluidesa_Verbal_Accio_NP 243 FE_R_abstracte_NP 244 FE_T_rellotge_NP 245 HAD_ansiedad_NP 246 HAD_depresion_NP
Entonces voy primero a enlazarlos ID de neuroimagen del proyecto (Subject) con los ID del proyecto (PSubject). Y luego voy a extraer de la DB las variables que necesito. Voy a empezar con M_ret_NP y por supuesto las covariables. Entonces voy a unir los dos archivos y esta sera mi base de trabajo.
[osotolongo@brick03 fsga_nbace]$ sed 's/;/,/g' ../bioface_mri.csv | sort -t, -k 2 > bioface_join.csv [osotolongo@brick03 fsga_nbace]$ sed -i '1iSubject,PSubject' bioface_join.csv [osotolongo@brick03 fsga_nbace]$ head bioface_join.csv Subject,PSubject 0001,B001 0002,B002 0003,B003 0004,B004 0005,B005 0006,B006 0007,B007 0008,B008 0009,B009 [osotolongo@brick03 fsga_nbace]$ awk -F"," '{print $1","$94","$78","$83","$220}' bioface_full_db.csv > fsga_220/bioface_220.csv [osotolongo@brick03 fsga_nbace]$ sed -i 's/SUJETO/PSubject/g' fsga_220/bioface_220.csv [osotolongo@brick03 fsga_nbace]$ join -t"," -1 2 -2 1 bioface_join.csv fsga_220/bioface_220.csv > fsga_220/bioface_data_220.csv
A partir de la DB de trabajo tengo que construir el archivo FSGD,
[osotolongo@brick03 fsga_220]$ (head -n 1 bioface_data_220.csv; tail -n +2 bioface_data_220.csv | sort -t, -k 2) > bioface_data_sorted_220.csv [osotolongo@brick03 fsga_220]$ awk -F"," '{print "bioface_"$2","$3","$4","$5","$6}' bioface_data_220.csv | sed 's/bioface_Subject/Variables/;s/Edad_visita/Age/;s/Sexo/Gender/;s/Anyos_escolaridad/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > 220_body.csv [osotolongo@brick03 fsga_220]$ head 220_body.csv Variables Age Gender Education M_ret_NP Input bioface_0001 Main 64 0 13 6 Input bioface_0002 Main 53 1 10 5 Input bioface_0003 Main 60 1 8 5 Input bioface_0004 Main 64 1 7 5 Input bioface_0005 Main 62 1 16 4 Input bioface_0006 Main 64 1 6 7 Input bioface_0007 Main 62 1 9 5 Input bioface_0008 Main 58 1 15 5 Input bioface_0009 Main 51 1 9 6 [osotolongo@brick03 fsga_220]$ cat headers_220.txt GroupDescriptorFile 1 Title BIOFACE_220 Class Main [osotolongo@brick03 fsga_220]$ cat headers_220.txt 220_body.csv > 220.fsgd
Y ahora me ejecuto un script que haga todos los analisis de la variable correspondiente,
- fsga_220.sh
#!/bin/bash mris_preproc --fsgd 220.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.220.thickness.10.mgh mris_preproc --fsgd 220.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.220.thickness.10.mgh mri_glmfit --y lh.220.thickness.10.mgh --fsgd 220.fsgd --C 220.mtx --surf fsaverage lh --glmdir lh.220.glmdir mri_glmfit --y rh.220.thickness.10.mgh --fsgd 220.fsgd --C 220.mtx --surf fsaverage rh --glmdir rh.220.glmdir mri_glmfit-sim --glmdir lh.220.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.220.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220.glmdir --cache 3 pos --cwp 0.05 --2spaces
Ahora bien, esto me vale para el thickness pero para area y volume necesito el valor del intracraneal volume, pasado a z-scores para evitar ill-posed problems.
awk -F"," '{print $1","$65}' ../../fsrecon/aseg_stats.csv > bioface_icv.csv Rscript convert_icv.r awk -F"," '{$1 = sprintf("%04d",$1); print $1","$3}' bioface_zic_tempv.csv | sed 's/"//g;s/0000/Subject/' > bioface_zicv.csv
Hago el join con los datos que ya tengo y construyo un nuevo FSGD con los zICV incluidos.
join -t"," -1 2 -2 1 bioface_data_sorted_220.csv bioface_zicv.csv > bioface_vdata.csv
awk -F"," '{print "bioface_"$1","$3","$4","$5","$6","$7}' bioface_vdata.csv | sed 's/bioface_Subject/Variables/;s/Edad_visita/Age/;s/Sexo/Gender/;s/Anyos_escolaridad/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > 220v_body.csv
cat headers_220.txt 220v_body.csv > 220v.fsgd
Por supuesto que esto lo he de añadir al modelo (mtx),
[osotolongo@brick03 fsga_220]$ cat 220a.mtx 0 0 0 0 1 0 y los nuevos scripts de //area// y //volume// tienen que usar estos nuevos archivos, [osotolongo@brick03 fsga_220]$ cat fsga_220a.sh #!/bin/bash mris_preproc --fsgd 220v.fsgd --cache-in area.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.220.area.10.mgh mris_preproc --fsgd 220v.fsgd --cache-in area.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.220.area.10.mgh mri_glmfit --y lh.220.area.10.mgh --fsgd 220v.fsgd --C 220a.mtx --surf fsaverage lh --glmdir lh.220a.glmdir mri_glmfit --y rh.220.area.10.mgh --fsgd 220v.fsgd --C 220a.mtx --surf fsaverage rh --glmdir rh.220a.glmdir mri_glmfit-sim --glmdir lh.220a.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.220a.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220a.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220a.glmdir --cache 3 pos --cwp 0.05 --2spaces [osotolongo@brick03 fsga_220]$ cat fsga_220v.sh #!/bin/bash mris_preproc --fsgd 220v.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.220.volume.10.mgh mris_preproc --fsgd 220v.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.220.volume.10.mgh mri_glmfit --y lh.220.volume.10.mgh --fsgd 220v.fsgd --C 220a.mtx --surf fsaverage lh --glmdir lh.220v.glmdir mri_glmfit --y rh.220.volume.10.mgh --fsgd 220v.fsgd --C 220a.mtx --surf fsaverage rh --glmdir rh.220v.glmdir mri_glmfit-sim --glmdir lh.220v.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.220v.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220v.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.220v.glmdir --cache 3 pos --cwp 0.05 --2spaces
En cualquier caso, para encontrar rapidamente si hay resultados basta con hacer un grep rapido,
[osotolongo@brick03 fsga_220]$ grep -v "^#" *glmdir/220*/*sig.cluster.summary lh.220.glmdir/220/cache.th30.pos.sig.cluster.summary: 1 4.7802 70529 384.91 -29.6 -35.2 -15.6 0.00200 0.00120 0.00280 899 3173.92 parahippocampal rh.220v.glmdir/220a/cache.th30.pos.sig.cluster.summary: 1 3.5435 158360 418.06 45.4 -66.0 29.0 0.00400 0.00280 0.00519 690 2060.05 inferiorparietal
Automatizando
Bueno este proceso es mas o menos el mismo para todas las variables, asi que como que se puede automatizar con un poco de maña.
El primer procedimiento es extraer la info de la DB, respecto a la variable a estudiar. Para esto voy a hacer un plantilla y en tonces solo hay que sustituir el codigo de la variables en la plantilla, ejecutar el script resultante y todo deberia quedar en su sitio. Para cada variable hare un subdirectorio distinto, asi no se mezclan la cosas.
El segundo paso es ejecutar las ordenes propiamente del analisis. Es decir que hago otra plantilla donde escribir los codigo de variable y ya esta. Igual, que antes, escribo, ejecuto y quedan los analisis en su sitio.
[osotolongo@brick03 fsga_nbace]$ cat templates/organize.sh #!/bin/bash mkdir fsga_<code> awk -F"," '{print $1","$94","$78","$83","$<code>}' bioface_full_db.csv > fsga_<code>/bioface_<code>.csv sed -i 's/SUJETO/PSubject/g' fsga_<code>/bioface_<code>.csv join -t"," -1 2 -2 1 bioface_join.csv fsga_<code>/bioface_<code>.csv > fsga_<code>/bioface_data_<code>.csv (head -n 1 fsga_<code>/bioface_data_<code>.csv; tail -n +2 fsga_<code>/bioface_data_<code>.csv | sort -t, -k 2) > fsga_<code>/bioface_data_sorted_<code>.csv awk -F"," '{print "bioface_"$2","$3","$4","$5","$6}' fsga_<code>/bioface_data_sorted_<code>.csv | sed 's/bioface_Subject/Variables/;s/Edad_visita/Age/;s/Sexo/Gender/;s/Anyos_escolaridad/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > fsga_<code>/<code>_body.csv echo 'GroupDescriptorFile 1' > fsga_<code>/headers_<code>.txt echo 'Title BIOFACE_<code>' >> fsga_<code>/headers_<code>.txt echo 'Class Main' >> fsga_<code>/headers_<code>.txt echo '0 0 0 0 1' > fsga_<code>/<code>.mtx cat fsga_<code>/headers_<code>.txt fsga_<code>/<code>_body.csv > fsga_<code>/<code>.fsgd join -t"," -1 2 -2 1 fsga_<code>/bioface_data_sorted_<code>.csv bioface_zicv.csv > fsga_<code>/bioface_vdata.csv awk -F"," '{print "bioface_"$1","$3","$4","$5","$6","$7}' fsga_<code>/bioface_vdata.csv | sed 's/bioface_Subject/Variables/;s/Edad_visita/Age/;s/Sexo/Gender/;s/Anyos_escolaridad/Education/' | sed 's/bioface_\([^,]*\),/Input bioface_\1 Main /; s/,/ /g' > fsga_<code>/<code>v_body.csv cat fsga_<code>/headers_<code>.txt fsga_<code>/<code>v_body.csv > fsga_<code>/<code>v.fsgd echo '0 0 0 0 1 0' > fsga_<code>/<code>a.mtx [osotolongo@brick03 fsga_nbace]$ cat templates/runner.sh #!/bin/bash ##thickness mris_preproc --fsgd <code>v.fsgd --cache-in area.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.<code>.area.10.mgh mris_preproc --fsgd <code>v.fsgd --cache-in area.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.<code>.area.10.mgh mri_glmfit --y lh.<code>.area.10.mgh --fsgd <code>v.fsgd --C <code>a.mtx --surf fsaverage lh --glmdir lh.<code>a.glmdir mri_glmfit --y rh.<code>.area.10.mgh --fsgd <code>v.fsgd --C <code>a.mtx --surf fsaverage rh --glmdir rh.<code>a.glmdir mri_glmfit-sim --glmdir lh.<code>a.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.<code>a.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>a.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>a.glmdir --cache 3 pos --cwp 0.05 --2spaces ##area mris_preproc --fsgd <code>.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.<code>.thickness.10.mgh mris_preproc --fsgd <code>.fsgd --cache-in thickness.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.<code>.thickness.10.mgh mri_glmfit --y lh.<code>.thickness.10.mgh --fsgd <code>.fsgd --C <code>.mtx --surf fsaverage lh --glmdir lh.<code>.glmdir mri_glmfit --y rh.<code>.thickness.10.mgh --fsgd <code>.fsgd --C <code>.mtx --surf fsaverage rh --glmdir rh.<code>.glmdir mri_glmfit-sim --glmdir lh.<code>.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.<code>.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>.glmdir --cache 3 pos --cwp 0.05 --2spaces ##volume mris_preproc --fsgd <code>v.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi lh --out lh.<code>.volume.10.mgh mris_preproc --fsgd <code>v.fsgd --cache-in volume.fwhm10.fsaverage --target fsaverage --hemi rh --out rh.<code>.volume.10.mgh mri_glmfit --y lh.<code>.volume.10.mgh --fsgd <code>v.fsgd --C <code>a.mtx --surf fsaverage lh --glmdir lh.<code>v.glmdir mri_glmfit --y rh.<code>.volume.10.mgh --fsgd <code>v.fsgd --C <code>a.mtx --surf fsaverage rh --glmdir rh.<code>v.glmdir mri_glmfit-sim --glmdir lh.<code>v.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir lh.<code>v.glmdir --cache 3 pos --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>v.glmdir --cache 3 neg --cwp 0.05 --2spaces mri_glmfit-sim --glmdir rh.<code>v.glmdir --cache 3 pos --cwp 0.05 --2spaces
y ahora debo hacer un script para llenar las plantillas y ejecutarlas. Esto es cosa bastante simple.
[osotolongo@brick03 fsga_nbace]$ cat run_all.sh #!/bin/bash code=$1 shift sed "s/<code>/${code}/g" templates/organize.sh > org${code}.sh chmod +x org${code}.sh ./org${code}.sh rm org${code}.sh sed "s/<code>/${code}/g" templates/runner.sh > fsga_${code}/run.sh cd fsga_${code}/ chmod +x run.sh ./run.sh cd .. grep -v "^#" fsga_${code}/*glmdir/${code}*/*sig.cluster.summary
El grep final devuelve los resultados, en caso de haber alguno. y ahora solo tendria que ir ejecutando, por ejemplo,
$ ./run_all 240
y cambiando el codigo para cada caso que me interese.
Resultados
Despues de ejecutar todo, puedo mirar los resultados globales con facilidad,
[osotolongo@brick03 fsga_nbace]$ ls -d fsga_2* fsga_219 fsga_220 fsga_221 fsga_227 fsga_232 fsga_238 fsga_240 fsga_241 fsga_242 [osotolongo@brick03 fsga_nbace]$ grep -v "^#" fsga_2*/*glmdir/2*/*sig.cluster.summary fsga_220/lh.220.glmdir/220/cache.th30.pos.sig.cluster.summary: 1 4.7802 70529 384.91 -29.6 -35.2 -15.6 0.00200 0.00120 0.00280 899 3173.92 parahippocampal fsga_220/rh.220v.glmdir/220a/cache.th30.pos.sig.cluster.summary: 1 3.5435 158360 418.06 45.4 -66.0 29.0 0.00400 0.00280 0.00519 690 2060.05 inferiorparietal fsga_238/lh.238.glmdir/238/cache.th30.neg.sig.cluster.summary: 1 -4.9777 159437 405.48 -8.4 -58.6 13.1 0.00020 0.00000 0.00040 949 -3320.19 precuneus fsga_238/lh.238.glmdir/238/cache.th30.neg.sig.cluster.summary: 2 -5.4753 28501 389.38 -51.8 -53.9 -16.1 0.00100 0.00040 0.00160 555 -2051.28 inferiortemporal fsga_238/lh.238.glmdir/238/cache.th30.neg.sig.cluster.summary: 3 -3.9238 23271 347.97 -58.7 -49.8 14.8 0.00140 0.00080 0.00200 717 -2300.38 superiortemporal fsga_238/lh.238.glmdir/238/cache.th30.neg.sig.cluster.summary: 4 -4.1765 125637 231.02 -55.4 -32.9 24.6 0.02524 0.02247 0.02800 409 -1345.49 supramarginal fsga_240/lh.240a.glmdir/240a/cache.th30.pos.sig.cluster.summary: 1 3.6678 29871 377.44 -34.8 44.8 19.1 0.01950 0.01693 0.02208 566 1720.94 rostralmiddlefrontal fsga_241/rh.241a.glmdir/241a/cache.th30.pos.sig.cluster.summary: 1 4.3157 1024 596.17 6.6 -19.1 62.8 0.00140 0.00080 0.00200 1482 4717.49 paracentral
Ahora para mirar especificamente el cluster que nos interesa hagamos algo como,
freeview -f $SUBJECTS_DIR/fsaverage/surf/lh.inflated:annot=aparc.annot:annot_outline=1:overlay=fsga_220/lh.220.glmdir/220/cache.th30.pos.sig.masked.mgh:overlay_threshold=2.3,5 -viewport 3d
que nos muestra algo asi,
Pero manipulando un poco la configuracion del overlay y moviendo un poc la figura, nos queda algo mas enseñable,
Aqui estamos usando aparc.annot como annotation y podemos ver las regiones del cortex afectadas por el cluster. Pero tambien podemos cargar otro archivo, por ejemplo lh.PALS_B12_Brodmann.annot y esto nos identifica las regiones de Brodman afectadas.
neuroimagen/bioface_sbm.txt · Last modified: 2021/12/14 10:33 by osotolongo