Differences

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

--- neuroimagen:md_nph [2021/10/15 10:29]
osotolongo [Uso]
+++ neuroimagen:md_nph [2021/10/19 09:04] (current)
osotolongo [Uso]
@@ Line 1: / Line 1: @@
 ====== MD-NPH analisis ======
-//Esto es una idea de Ester Esteban, asi que todo es culpa suya. //
+//Esto es una idea de Ester Esteban, asi que toda la culpa es suya. //
@@ Line 45: / Line 45: @@
 asi que todo deberia ser sencillo.
-Primero leo la DB del proyecto y por cada unos de las imagenes MD que he llegado a procesar ejecuto el procedimiento anterior,
+Ok. Pero no taaaan sencillo. :-P
+Antes que nada, no voy a sacar el valor medio de toda la imagen MD sino de la porcion correspondiente a la WM. Es decir, que tengo que empezar por hacer una mascara de la WM. Es to no deberia ser dificil si reutilizo el //aseg// de FS,
+<code bash get_wm_mask.sh>
+#!/bin/sh
+subject=$1
+shift
+tmp_dir=$1
+shift
+if [ ! -f ${tmp_dir}/register.dat ]; then
+        if [  ! -d ${tmp_dir} ]; then mkdir -p ${tmp_dir}; fi;
+        tkregister2 --mov $SUBJECTS_DIR/${subject}/mri/rawavg.mgz --noedit --s ${subject} --regheader --reg ${tmp_dir}/register.dat;
+fi
+if [ ! -f ${tmp_dir}/all_aseg.nii.gz ]; then
+        mri_label2vol --seg $SUBJECTS_DIR/${subject}/mri/aseg.mgz --temp $SUBJECTS_DIR/${subject}/mri/rawavg.mgz --o ${tmp_dir}/all_aseg.nii.gz --reg ${tmp_dir}/register.dat;
+fi
+${FSLDIR}/bin/fslmaths ${tmp_dir}/all_aseg.nii.gz -uthr 2  -thr 2 -div 2 ${tmp_dir}/lhwm
+${FSLDIR}/bin/fslmaths ${tmp_dir}/all_aseg.nii.gz -uthr 41  -thr 41 -div 41 ${tmp_dir}/rhwm
+${FSLDIR}/bin/fslmaths ${tmp_dir}/lhwm -add ${tmp_dir}/rhwm -bin ${tmp_dir}/wm_mask.nii.gz
+~
+</code>
+{{ :neuroimagen:wm_mask.png?400 |}}
+Pero ojo, que esa mascara esta en espacio nativo T1, asi que tendria que hacer algo como,
+<code bash>
+antsApplyTransforms -d 3 -i ${tmp_dir}/wm_mask.nii.gz -r ${tmp_dir}/hifi_b0.nii.gz -t ${tmp_dir}/$subject_epi_reg_ANTS.mat -n GenericLabel -o ${tmp_dir}/wm_tmp.nii.gz
+</code>
+La mascara resultante hay que aplicarla a la imagen MD,
+{{ :neuroimagen:wm_mask_indti.png?400 |}}
+y a esta ultima imagen,
+{{ :neuroimagen:dti_md_masked.png?400 |}}
+es a la que hay que sacarle el histograma.
 <code perl>
-foreach my $subject (@dtis){
-        my $dti_md = $w_dir.'/'.$subject.'_dti_MD.nii.gz';
+                my $order = $ENV{'FSLDIR'}.'/bin/fslstats '.$dti_md_masked.' -l 0 -u '.$up_threshold.' -h '.$nbins;
-        if($subject and -e $dti_md){
-                my $order = $ENV{'FSLDIR'}.'/bin/fslstats '.$dti_md.' -l 0 -u '.$up_threshold.' -h '.$nbins;
                 my @shist = qx/$order/;
-                chomp @shist;
-                push @allsums, sum(@shist);
-                $nph{$subject} = max(@shist);
-        }
-}
 </code>
-ahora, el valor para normalizar es el medio de todas las areas,
+ahora, el valor para normalizar es el numero de voxels de la mascara de WM
 <code perl>
-my $norm = sum(@allsums)/scalar(@allsums);
+                $order = $ENV{'FSLDIR'}.'/bin/fslstats '.$dti_md_masked.' -V';
+                my $oout =  qx/$order/;
+                my @gout = split ' ', $oout;
+                my $norm = $gout[0];
 </code>
-y saco todo por //STDOUT//,
+y todo esto se guarda,
 <code perl>
-print "Subject,NPH\n";
+                $nph{$subject}{'sum'} = sum(@shist);
+                $nph{$subject}{'norm'} = $norm;
+                $nph{$subject}{'max'} = max(@shist);
+</code>
+para escribirlo luego,
+<code perl>
+open ODF, ">$ofile";
+print ODF "Subject,NPH\n";
 foreach my $subject (sort keys %nph){
-        my $nphv = $nph{$subject}/$norm;
+        my $nphv = $nph{$subject}{'max'}/$nph{$subject}{'norm'};
-        print "$subject,$nphv\n";
+        print ODF "$subject,$nphv\n";
 }
+close ODF;
 </code>
@@ Line 83: / Line 131: @@
 <code bash>
-[osotolongo@brick03 nph]$ ./nph.pl bioface > nph.csv
+[osotolongo@brick03 nph]$ ./nph.pl bioface
+[osotolongo@brick03 nph]$ head ../bioface_normalised_peak_height.csv
+Subject,NPH
+,0.142459141831309
+,0.101291664281762
+,0.128800261523374
+,0.131614509246088
+,0.116625745352508
+,0.142582104338399
+,0.140059946569362
+,0.129072904095975
+,0.12497310864109
 [osotolongo@brick03 nph]$ sed 's/;/,/;1iSubject,Subject_ID' ../bioface_mri.csv > bioface_codes.csv
-[osotolongo@brick03 nph]$ join -t, -j 1 bioface_codes.csv nph.csv | awk -F"," '{print $2","$3}' | csvsort -c 1 > normalised_peak_height.csv
+[osotolongo@brick03 nph]$ join -t, -j 1 bioface_codes.csv ../bioface_normalised_peak_height.csv | awk -F"," '{print $2","$3}' | csvsort -c 1 > nph_sorted.csv
-[osotolongo@brick03 nph]$ head normalised_peak_height.csv
+[osotolongo@brick03 nph]$ head nph_sorted.csv
 Subject_ID,NPH
-B002,0.0731249323687369
+B002,0.142459141831309
-B003,0.0732767676009448
+B003,0.101291664281762
-B004,0.0745456763272544
+B004,0.128800261523374
-B005,0.0746812434988687
+B005,0.131614509246088
-B006,0.0697899799470258
+B006,0.116625745352508
-B007,0.0788675577583173
+B007,0.142582104338399
-B008,0.0932810594443464
+B008,0.140059946569362
-B009,0.0787428359604322
+B009,0.129072904095975
-B010,0.0681414831401963
+B010,0.12497310864109
 </code>
 **Nota:** Para simplificar he usado //csvsort// de [[https://csvkit.readthedocs.io/en/latest/ | csvkit ]] pero tambien se puede dar una vuelta y usar //gnu sort//.

Detritus Wiki

User Tools

Site Tools

Differences

Page Tools