Detritus Wiki

User Tools

Site Tools

Trace: • pgs
neuroimagen:pgs

Table of Contents

WMH pgs

  1. Se supone que un metodo dedicado a sacar las WMH es mejor
  2. De esta manera lo que saco es una mascara de las WMH y podria sacar los volumenes en cada tracto. 8-o
  3. Por lo menos no tendre que exlpicar a mas reviewers que es igual uno que otro para sacar relaciones lineales

HowTo

  1. https://pubmed.ncbi.nlm.nih.gov/33957235/
  2. https://hub.docker.com/r/wmhchallenge/pgs
  3. https://wmh.isi.uu.nl/wp-content/uploads/2019/08/pgs.pdf
[osotolongo@brick03 f5cehbi]$ singularity run --cleanenv -B /nas:/nas -B /old_nas/f5cehbi/output:/output -B /old_nas/f5cehbi/input:/input /nas/usr/local/opt/singularity/pgs.simg sh /WMHs_segmentation_PGS.sh sub-0020_T1w.nii.gz sub-0020_T2w.nii.gz sub-0020_WMH.nii.gz
Step 1  of 9: reading in images...Done! It took roughly 1 seconds
Step 2 of 9: registration...
  2a) Similarity transform...
  2b) Affine transform...
  2c) Resampling volume...
  Registration done! It took roughly 6 seconds
Step 3 of 9: rough bias field correction...Done! It took roughly 6 seconds
Step 4 of 9: calculating features...Done! It took roughly 28 seconds
Step 5 of 9: voxel classification...Done! It took roughly 17 seconds
Step 6 of 9: fitting the shape model...Done! It took roughly 6 seconds
Step 7 of 9: free deformation...Done! It took roughly 1 seconds
Step 8 of 9: building volume from mesh...Done! It took roughly 3 seconds
Step 9 of 9: warping back to original space...Done! It took roughly 0 seconds

 Output ready 

 The whole process took roughly 68 seconds

 Freeing memory and deleting temporary files...


Traceback (most recent call last):
  File "/WMHs_segmentation_PGS.py", line 89, in <module>
    FLAIR_bet_data = FLAIR_data*T1_bet_mask_data
ValueError: operands could not be broadcast together with shapes (176,252,256) (176,256,256) 
[osotolongo@brick03 f5cehbi]$ fslinfo output/sub-0020_T1w.nii.gz
data_type      INT16
dim1           176
dim2           256
dim3           256
dim4           1
datatype       4
pixdim1        1.000000
pixdim2        0.976562
pixdim3        0.976562
pixdim4        2.200000
cal_max        0.0000
cal_min        0.0000
file_type      NIFTI-1+
[osotolongo@brick03 f5cehbi]$ fslinfo output/sub-0020_T2w.nii.gz 
data_type      INT16
dim1           176
dim2           252
dim3           256
dim4           1
datatype       4
pixdim1        1.000000
pixdim2        1.015625
pixdim3        1.015625
pixdim4        5.000000
cal_max        0.0000
cal_min        0.0000
file_type      NIFTI-1+

Nota: Registrar T2 a T1 primero! 

[osotolongo@brick03 f5cehbi]$ antsRegistrationSyNQuick.sh -d 3 -t a -f input/pre/sub-0020_T1w.nii.gz -m input/pre/sub-0020_T2w.nii.gz -o input/pre/t2tot1_ | tee t2tot1Output.txt
...
...
...
[osotolongo@brick03 f5cehbi]$ antsApplyTransforms -d 3 -i input/pre/sub-0020_T2w.nii.gz -r input/pre/sub-0020_T1w.nii.gz -t input/pre/t2tot1_0GenericAffine.mat -o input/pre/sub-0020_T2w_resampled.nii.gz
[osotolongo@brick03 f5cehbi]$ ls input/pre/
sub-0020_T1w.nii.gz  sub-0020_T2w_resampled.nii.gz  t2tot1_InverseWarped.nii.gz
sub-0020_T2w.nii.gz  t2tot1_0GenericAffine.mat      t2tot1_Warped.nii.gz

Ahora si, 

[osotolongo@brick03 f5cehbi]$ singularity run --cleanenv -B /nas:/nas -B /old_nas/f5cehbi/output:/output -B /old_nas/f5cehbi/input:/input /nas/usr/local/opt/singularity/pgs.simg sh /WMHs_segmentation_PGS.sh sub-0020_T1w.nii.gz sub-0020_T2w_resampled.nii.gz sub-0020_WMH.nii.gz

This is a highly massive computation, :-o Dont take it lightly!

but results are amazing, 

[osotolongo@brick05 f5cehbi]$ fslstats output/sub-0020_WMH.nii.gz -V 
2353 2243.995605

Scripting

Primero hay que sacar las iamgenes de WMH para el proyecto,

https://github.com/asqwerty666/acenip/blob/main/bin/wmh.pl

y luego las metricas de cada imagen,

https://github.com/asqwerty666/acenip/blob/main/bin/wmh_metrics.pl

comparando

Para cada proyecto, 

[osotolongo@brick03 facehbi]$ xnatapic list_subjects --project_id facehbi --label > xnat_subjects.list
[osotolongo@brick03 facehbi]$ for x in `awk -F"," '{print $2}' xnat_subjects.list`; do e=$(xnatapic list_experiments --project_id facehbi --subject_id ${x} --modality MRI --label); echo "${x},${e}"; done > xnat_sub_exp.list
[osotolongo@brick03 facehbi]$ for l in `cat xnat_sub_exp.list`; do s=$(echo ${l} | awk -F"," '{print $1}'); e=$(echo ${l} | awk -F"," '{print $2}'); mkdir -p fsresults/${s}; xnatapic get_fsresults --experiment_id ${e} --stats aseg fsresults/${s}/; done
[osotolongo@brick03 facehbi]$ for x in `ls fsresults/*/aseg.stats`; do s=$(echo ${x} | awk -F"/" '{print $2}'); v=$(grep " WM-hypointensities" ${x} | awk '{print $4}'); echo "${s},${v}"; done | sed '1iPSubject,WM-hypointensities' > fswmh_results.csv
[osotolongo@brick03 facehbi]$ sed 's/;/,/g' facehbi_mri.csv | sed '1iSubject,PSubject' > joiner.csv
[osotolongo@brick03 facehbi]$ wmh_metrics.pl facehbi
[osotolongo@brick03 facehbi]$ join -t"," joiner.csv facehbi_wmh_metrics.csv | awk -F"," '{print $2","$3}' > pgs_wmh.csv 
[osotolongo@brick03 facehbi]$ (head -n 1 pgs_wmh.csv; tail -n +2 pgs_wmh.csv | sort -t",") > pgs_wmh_sorted.csv
[osotolongo@brick03 facehbi]$ join -t"," fswmh_results.csv pgs_wmh_sorted.csv > compare.csv
[osotolongo@brick03 facehbi]$ head compare.csv

y despues se pegan, 

[osotolongo@brick03 f5cehbi]$ cat compare.csv /old_nas/facehbi/compare.csv | grep -v PSubject | sed '1iPSubject,WM-hypointensities,WMH' > compare_v05.csv

y el script en R hace el plot,

myplot.r
w <- read.csv("compare_v05.csv")
wl <- lm(w$WMH ~ w$WM.hypointensities)
pwl <- summary(wl)
postscript("wmh_comparing.ps", width=1024, height=600, bg="white")
plot(w$WM.hypointensities, w$WMH, main="WMH comparison pgs VS FS", xlab="Freesurfer WMH", ylab="pgs WMH", pch=15, col="blue")
abline(wl$coefficients[1],wl$coefficients[2],col="red",lwd=3)
text(2500,25000, bquote(R^2 == .(round(pwl$adj.r.squared, digits=3))), cex=2)
dev.off()

Poniendo las tres visitas de FACEHBI ($n \sim 550$) juntas tenemos algo asi,

Guardando WMH en XNAT (deprecated)

Para guardar los valores calculados de WMH en XNAT necesito un archivo con el formato, 

PSubject,WMH
B014,1240.071167
B021,15298.660156
B030,906.359924
B035,966.427979
B037,3140.890869
B038,2161.114502
...
...
...

Ahora voy a aprovechar dos funciones de interaccion con XNAT; xput_res y xget_res. (ver XNATACE lib)

  • xput_res toma la referencia a un hash y lo sube como un .json. Asi que para cada sujeto descrito en el archivo de datos habria que armar un pequeño hash e invocar la funcion,
foreach my $sbj (sort keys %wmhs){
        my $experiment = xget_mri($xconf{'HOST'}, $xconf{'JSESSION'}, $xprj, $sbj);
        my %wmh_data;
        $wmh_data{'WMH'} = $wmhs{$sbj};
        xcreate_res($xconf{'HOST'}, $xconf{'JSESSION'}, $experiment, 'data');
        xput_res($xconf{'HOST'}, $xconf{'JSESSION'}, $experiment, 'data', 'wmh.json', \%wmh_data);
}
  • xget_res toma a localizacion de un archivo .json y devuelve un hash. Asi que por cada sujeto invoco la funcion y escribo el contenido del hash,
foreach my $sbj (sort keys %subjects){
        my $experiment = xget_mri($xconf{'HOST'}, $xconf{'JSESSION'}, $xprj, $sbj);
        my $label = xget_sbj_data($xconf{'HOST'}, $xconf{'JSESSION'}, $sbj, 'label');
        my %wmh_data = xget_res($xconf{'HOST'}, $xconf{'JSESSION'}, $experiment, 'data', 'wmh.json');
        print "$label,$wmh_data{'WMH'}\n";
}

Los scripts xnat_up_wmh.pl y xnat_get_wmh.pl son entonces sencillos de usar para subir y bajar el contenido del archivo CSV.

Calculando WMH en XNAT

Hemos conseguido configurar el pipeline pgsWMH para que calcule el volumen de WMH directamente en XNAT. (Lo archivos xml aqui) ¿Porque si ya disponiamos de un metodo valido?

  1. Al integrar el calculo de WMH en un pipeline de XNAT nos ahorramos lo pasos mas costosos (en esfuerzo humano y tiempo) que son localizar, extraer y organizar los DICOM, así como convertir los T1w y T2w a formato nifti.
  2. Si el pipeline esta corectamente configurado los datos se guardan silenciosamente en XNAT y estan disponibles automaticamente para consultarse en cualquier momento.
  3. La localizacion de los sujeto para los que aun no se ha calculado el volumen de WMH es ahora trivial y su calculo por separado aun mas simple.
  4. Al dismminuir el nivel de interaccion necesario para este calculo, los posibles errores de procesamiento intermedio se disminuyen a un minimo.

Así que, ahora para calcular el WMH de un proyecto completo basta hacer algo como, 

[osotolongo@brick03 f5cehbi]$ for x in `xnatapic list_experiments --project_id f5cehbi --modality MRI`; do xnatapic run_pipeline --project_id f5cehbi --pipeline pgsWMH --experiment_id ${x}; done

Dependiendo lo bien escritos que esten los TAGs de las imagenes, puede haber algun error, pero llegara un email anunciandolo.

Ahora, con un pequeño cambio en el script de bajar los datos de WMH, podemos obtener los volumenes en mm3. 

        my %wmh_data = xget_res_data($xconf{'HOST'}, $xconf{'JSESSION'}, $experiment, 'WMH', 'wmh.json');
        if (exists($wmh_data{'WMH mm3'}) and $wmh_data{'WMH mm3'}){
                print "$label,$wmh_data{'WMH mm3'}\n";
        }else{
                print "$label,NA\n";
        }

Esto es gracias a que la funcion xget_res_data es capaz de bajar los datos de cualquier JSON y almacenarlos en un hash. Solo necesitamos saber la estructura del JSON y preguntar por la variable que se quiera.

Entonces, el resultado se baja como, 

[osotolongo@brick03 f5cehbi]$ xnat_get_wmh.pl -x f5cehbi
Subject_ID,WMH
F204,5162.239258
F246,663.757324
F230,3814.697266
F227,379.562378
F013,5701.064941
F132,387.191772
F033,1088.142334
F073,4213.333008
F027,7450.103516
F074,794.410706
F071,3265.380859
...
...
...

o si lo deseamos, especificando un archivo de output.

Ahora tenemso el problema de ir ejecutando por lotes, pero utilizando la API de XNAT esto deberia ser simple. Algo como, 

[osotolongo@brick03 mri_face]$ for x in `xnatapic list_experiments --project_id mriface --modality MRI`; do if [[ -z `curl -f -X GET -b "JSESSIONID=542119A78B8CD068BACB5862961601F3" "http://detritus.fundacioace.com:8088/data/experiments/${x}/files?format=json" 2>/dev/null| jq '.ResultSet.Result[].Name' | grep wmh.json` ]]; then echo ${x}; fi; done

Debria darnos la lista de los experimentos que NO tienen el archivo wmh.json. Asi que bastaria con enviar el pipeline a ejecutarse sobre estos experimentos.

Por ejemplo, 

[osotolongo@brick03 mri_face]$ for x in `xnatapic list_experiments --project_id mriface --modality MRI`; do if [[ -z `curl -f -X GET -b "JSESSIONID=FE98FA6B6AB5520A3A329757C3BC1127" "http://detritus.fundacioace.com:8088/data/experiments/${x}/files?format=json" 2>/dev/null| jq '.ResultSet.Result[].Name' | grep wmh.json` ]]; then echo ${x}; fi; done
XNAT05_E00040
XNAT05_E00066
XNAT_E00937
XNAT_E00938
XNAT_E00939
XNAT_E00940
[osotolongo@brick03 mri_face]$ for x in $(seq 37 40); do xnatapic run_pipeline --project_id mriface --pipeline pgsWMH --experiment_id XNAT_E009${x}; done
[osotolongo@brick03 mri_face]$ queue
   JOBID PARTITION                                 NAME     USER ST       TIME  NODES NODELIST(REASON)
  323751      fast                   pgsWMH.XNAT_E00937     xnat  R       2:55      1 brick01
  323752      fast                   pgsWMH.XNAT_E00938     xnat  R       2:55      1 brick02
  323753      fast                   pgsWMH.XNAT_E00939     xnat  R       2:51      1 brick05
  323754      fast                   pgsWMH.XNAT_E00940     xnat  R       2:51      1 brick04

Y si, las dos primeras las he excluido porque ya se que han fallado antes por no tener bien el T2w. Entonces lo que voy a hacer es excluir estos experimentos. Voy a hacer un archivo wmh.json para los fallidos, 

[osotolongo@brick03 mri_face]$ cat wmh.json
{"ResultSet":{"Result":[{"ID":"NA","WMH voxels":"NA","WMH mm3":"NA"}]}}

y entonces, 

[osotolongo@brick03 mri_face]$ xnatapic get_jsession
FE98FA6B6AB5520A3A329757C3BC1127
[osotolongo@brick03 mri_face]$ curl -f -X PUT -b "JSESSIONID=FE98FA6B6AB5520A3A329757C3BC1127" "http://detritus.fundacioace.com:8088/data/experiments/XNAT05_E00040/resources/WMH/files/wmh.json?overwrite=true" -F file="@wmh.json"
[osotolongo@brick03 mri_face]$ curl -f -X PUT -b "JSESSIONID=FE98FA6B6AB5520A3A329757C3BC1127" "http://detritus.fundacioace.com:8088/data/experiments/XNAT05_E00066/resources/WMH/files/wmh.json?overwrite=true" -F file="@wmh.json"

y cuando haga la proxima busqueda estos resultados no apareceran.

neuroimagen/pgs.txt · Last modified: 2022/12/07 10:55 by osotolongo