Paulina J. Villaseñor
Laboratorio de Conectividad Cerebral C13
LANIREM
---
En este tutorial veremos paso a paso como exportar tus imagenes DWI desde el Bruker, como convertirlas en formatos Nifti (.nii) y finalmente como pre-procesarlas.
Debido a que usaremos la terminal y bash para ejecutar estos pasos, si no estas tan familiarizado en como se usa este lenguaje. [Aquí](https://seankross.com/the-unix-workbench/) te dejo un manual muy sencillo y didactico de como empezar y abrazar con cariño el sistema Unix/Linux/Bash.
### EXPORTACIÓN DE DWI DESDE EL BRUKER
Primero hay que abrir la terminal (ctrl + T) e irnos al directorio donde las imagenes del Bruker se almacenan:
```
cd /misc/bruker7/data02/user/conchalab
```
> Ojo, aqui pongo como ejemplo el `user` como `conchalab` pero hay que ingresar al usuario del lab en el que se inicio sesion en la consola del Bruker a la hora de escanear.
Lo siguiente es localizar los archivos que desear convertir. Puedes buscarlos al usar el comando `ls` o mas facil, buscarlo utilizando un `*` si sabes el nombre de la rata que buscas.
```
ls *irm150d_rata64A*
```
> Al hacer este filtro, yo estoy buscando especificamente por la rata 64A y el archivo que me encontro es el siguiente:
>
> 20220104_085643_INB_C13_hluna_irm150d_rata64A_INB_C13_hluna_1_1
>
> Y es el que voy a utilizar de ahora en adelante como ejemplo.
Si nosotros enlistamos (`ls`) esta carpeta para ver que hay adentro, veremos que hay numeros al inicio, estos corresponden a cada adquisicion en el orden en el que fueron tomadas y son las que vamos a ir convirtiendo una a una.
```
ls 20220104_085643_INB_C13_hluna_irm150d_rata64A_INB_C13_hluna_1_1/
1 2 3 4 5 6 7 8 AdjResult AdjStatePerStudy Mapshim ResultState ScanProgram.scanProgram subject
```
Bien, hasta aqui ya sabemos como acceder a tus imágenes del Bruker, siguiente paso es exportarlas en formato Nifti.
Paso numero uno es cargar el modulo de Bruker (gracias a Ricardo Rios que nos hizo la vida mas facil al crear los modulos, si aun no te familiarizas con ellos, da click [aquí](https://github.com/c13inb/c13inb.github.io/wiki/Modules) y aprende mas a como usarlos).
```
module load brkraw/0.3.11
```
Una vez cargado el modulo estas listo para utilizarlo. Primero queremos saber cual es la información detallada de cada una de nuestras adquisiciones. Para esto utilizamos el comando `brkraw info` y nos despliega la siguiente información:
```
brkraw info 20220104_085643_INB_C13_hluna_irm150d_rata64A_INB_C13_hluna_1_1/
```
```
Paravision 7.0.0
----------------
UserAccount: conchalab
Date: 2022-01-04
Researcher: rata64A
Subject ID: INB_C13_hluna_irm150d_rata64A
Session ID: INB_C13_hluna_irm150d_rata64A
Study ID: 1
Date of Birth: 07 Aug 2021
Sex: male
Weight: 0.433 kg
Subject Type: Quadruped
Position: Prone Entry: HeadFirst
[ScanID] Sequence::Protocol::[Parameters]
[001] Bruker:FLASH::1_Localizer::1_Localizer (E1)
[ TR: 100 ms, TE: 2.50 ms, pixelBW: 159.22 Hz, FlipAngle: 30 degree]
[01] dim: 2D, matrix_size: 256 x 256 x 3, fov_size: 50 x 50 (unit:mm)
spatial_resol: 0.195 x 0.195 x 2.000 (unit:mm), temporal_resol: 12800.000 (unit:msec)
[002] Bruker:FLASH::1_Localizer::1_Localizer (E2)
[ TR: 100 ms, TE: 2.50 ms, pixelBW: 159.22 Hz, FlipAngle: 30 degree]
[01] dim: 2D, matrix_size: 256 x 256 x 3, fov_size: 50 x 50 (unit:mm)
spatial_resol: 0.195 x 0.195 x 2.000 (unit:mm), temporal_resol: 12800.000 (unit:msec)
[003] Bruker:FLASH::1_Localizer::1_Localizer (E3)
[ TR: 100 ms, TE: 2.50 ms, pixelBW: 159.22 Hz, FlipAngle: 30 degree]
[01] dim: 2D, matrix_size: 256 x 256 x 3, fov_size: 50 x 50 (unit:mm)
spatial_resol: 0.195 x 0.195 x 2.000 (unit:mm), temporal_resol: 12800.000 (unit:msec)
[004] Bruker:FLASH::T1_FLASH::T1_FLASH (E4)
[ TR: 201.57 ms, TE: 3.50 ms, pixelBW: 98.64 Hz, FlipAngle: 30 degree]
[01] dim: 2D, matrix_size: 384 x 384 x 13, fov_size: 25.6 x 25.6 (unit:mm)
spatial_resol: 0.067 x 0.067 x 1.100 (unit:mm), temporal_resol: 309614.466 (unit:msec)
[005] Bruker:FieldMap::B0Map-ADJ_B0MAP::T1_FLASH
[ TR: 20 ms, TE: 0 ms, pixelBW: 1860.12 Hz, FlipAngle: 30 degree]
[01] dim: 3D, matrix_size: 64 x 64 x 64, fov_size: 45 x 45 x 45 (unit:mm)
spatial_resol: 0.703 x 0.703 x 0.703 (unit:mm), temporal_resol: 81920.000 (unit:msec)
[006] Bruker:DtiEpi::DTI_EPI_30dir::DWIzoom (E6)
[ TR: 2000 ms, TE: 22.86 ms, pixelBW: 2289.38 Hz, FlipAngle: 90 degree]
[01] dim: 2D, matrix_size: 126 x 86 x 25 x 285, fov_size: 22 x 15 (unit:mm)
spatial_resol: 0.175 x 0.174 x 1.250 (unit:mm), temporal_resol: 4000.000 (unit:msec)
[02] dim: 2D, matrix_size: 126 x 86 x 22 x 25, fov_size: 22 x 15 (unit:mm)
spatial_resol: 0.175 x 0.174 x 0.006 (unit:mm), temporal_resol: 0.000 (unit:msec)
[007] Bruker:DtiEpi::DTI_EPI_30dir::DWI-IVIM-zoom(E11) (E7)
[ TR: 2000 ms, TE: 22.86 ms, pixelBW: 2289.38 Hz, FlipAngle: 90 degree]
[01] dim: 2D, matrix_size: 126 x 86 x 25 x 63, fov_size: 22 x 15 (unit:mm)
spatial_resol: 0.175 x 0.174 x 1.250 (unit:mm), temporal_resol: 4000.000 (unit:msec)
[02] dim: 2D, matrix_size: 126 x 86 x 22 x 25, fov_size: 22 x 15 (unit:mm)
spatial_resol: 0.175 x 0.174 x 0.006 (unit:mm), temporal_resol: 0.000 (unit:msec)
[008] Bruker:RARE::T2_TurboRARE::T2_TurboRARE (E8)
[ TR: 4212.78 ms, TE: 33 ms, pixelBW: 140.85 Hz, FlipAngle: 141.72 degree]
[01] dim: 2D, matrix_size: 256 x 256 x 26, fov_size: 30 x 30 (unit:mm)
spatial_resol: 0.117 x 0.117 x 1.200 (unit:mm), temporal_resol: 269617.981 (unit:msec)
```
Podría parecer mucha información al inicio, pero al final no es mas que cada adquisición enumerada del `001` al `007`, donde yo sé que la imagén DWI es la numero `006`, por la siguiente informacion:
#
`[006] Bruker:DtiEpi::DTI_EPI_30dir::DWIzoom (E6)` Justo aqui menciona que es una adquisicion EPI (Echo Planar Image)
`[ TR: 2000 ms, TE: 22.86 ms, pixelBW: 2289.38 Hz, FlipAngle: 90 degree]` Aqui estan los parametros de adquisicion
` [01] dim: 2D, matrix_size: 126 x 86 x 25 x 285, fov_size: 22 x 15 (unit:mm)` y aqui justo vemos la dimensionalidad de la imagen, donde 285 son los volumenes de difusion
` spatial_resol: 0.175 x 0.174 x 1.250 (unit:mm), temporal_resol: 4000.000 (unit:msec)` Aqui especifica el tamaño del voxel en los planos x, y, z
#
Una vez que sabemos que la imagén DWI que queremos convertir a Nifti es la numero 6, podemos ejecutar el siguiente comando:
```
brkraw tonii 20220104_085643_INB_C13_hluna_irm150d_rata64A_INB_C13_hluna_1_1/ -o ./64A_dwi -r 1 -s 6
```
En otras palabras:
`tonii` es el comando que convierte de Bruker a Nifti.
```-o``` es el output de como quieres que se llame tu imagen y en donde quieres guardarla, en este caso 64A_dwi es el nombre que yo le doy y `./` hago referencia de que se guarde en el directorio actual.
```-r``` es la reconstruccion que queremos, en este caso es la primera y por eso ponemos 1
```-s``` es la imagen que queremos convertir, en este caso es la numero 6 porque es la DWI
Vamos a ver que despues de la conversion tendremos **tres** outputs:
```
ls
64A_dwi.bvec
64A_dwi.bval
64A_dwi.nii.gz
```
Los archivos `.bvec` contiene información acerca de los autovectores, mientras que el archivo `.bval` son los autovalores. Fundamentales para las DWI. Y por ultimo tenemos la imagen `.nii.gz` que es la que podremos visualizar.
Para ver que tus imagenes se convirtieron exitosamente en formato Nifti, vamos a visualizarlas utilizando `mrview` del software `mrtrix`. Para esto, no olvides cargar tu modulo: `module load mrtrix/3.0.4`
```
mrview 64A_dwi.nii.gz
```
Y el resultado es esto:

Ahora, las imágenes DWI así crudas son feísimas! y con mucho ruido. Pero para eso esta el pre-procesamiento y aquí te diré que es y como se hace.
---
### PRE-PROCESAMIENTO
Las imágenes sensibles a difusón deben de pasar por el pre-procesamiento antes de someterlas a cualquier análisis.
Existen varios pasos en el pre-procesamiento, que van desde quitar el ruido de la señal hasta remover ciertos tipos de artefactos. Aquí vamos a introducir primero paso por paso en como se hace y al final encontraras un script que pertenece al repositorio del INB que hace todos los pasos.
> La gran mayoria de los comandos son parte del repositorio de MRtrix3, si quieres aprender mas acerca de como se usa da click [aqui](https://mrtrix.readthedocs.io/en/latest/#).
1. **Denoising**
Este paso es fundamental y normalmente el primer paso antes de cualquier otro. Consiste en remover el ruido proveniente de la señal. Aquí puedes utilizar el comando `dwidenoise`:
`dwidenoise 64A_dwi.nii.gz denoised_64A_dwi.nii.gz -noise noise_64A_dwi.nii.gz`
> Donde `dwidenoise` es el comando, despues viene el `input` (DWI cruda), seguido del `output` (mi nueva imagen con denoise) y por ultimo `-noise` y su correspondiente `output` para el ruido estimado.
> Si quieres saber mas en como funciona haz clic aquí [aquí](https:/mrtrix.readthedocs.io/en/dev/dwi_preprocessing/denoising.html/)
Ten paciencia que el denoising es tardadito... pero una vez completado puedes ver tu nueva imagen:
```
mrview denoised_64A_dwi.nii.gz
```

_
2. **Eddy Correction**
Las corrientes de Eddy (tambien conocidas como Foucault) son el resultado del cambio rapído en los gradientes del campo magnético y estos inducen ciertos tipos de artefactos en nuestras imágenes DWI. Los desarrolladores de FSL crearon una herramienta llamada `eddy_cuda10.2` que ejecuta esta corrección y mucho más.
Antes de correr eddy, necesitamos hacer una serie de primeros pasos para preparar los datos de acuerdo a como lo pide el software. En su [pagina web](https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/eddy/UsersGuide/) tienen toda la información detallada de como hacerlo. Aquí lo resumiré con el ejemplo de la rata 64A.
1) Primero necesitamos sacar una máscara binaria del cerebro de la rata
```
dwi2mask -fslgrad 64A_dwi.bvec 64A_dwi.bval 64A_dwi.nii.gz mascara_64A_dwi.nii.gz
```
2) Ahora necesitamos un archivo que describa los parametros de la adquisición de cada imágen.
```
topup= 0.05
echo "0 -1 0" $topup > acqp_64A_dwi.txt
```
```
cat acqp_64A_dwi.txt
0 -1 0 0.05
```
Vemos que en el output tenemos `0 -1 0` que no es nada mas que la codificación en fase y `0.05` es la multiplicación entre el factor EPI y los ms de espacio entre ecos. Toda esta información al final son los parámetros de adquisición. Más información [aquí](https://fsl.fmrib.ox.ac.uk/fsl/fslwiki/eddy/Faq#How_do_I_know_what_to_put_into_my_--acqp_file)
3) Hay que crear un archivo índice que ayude a indicar que volúmenes (aquí 285) de DWI fueron tomadas con ciertos parametros de acuerdo al archivo acqp_64A_dwi.txt. En este caso, todos los volúemenes fueron adquiridos de igual manera.
```
indx=""
for ((i=1; i<=285; i+=1)); do indx="$indx 1"; done
echo $indx > indice_64A_dwi.txt
```
```
echo $indx
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
```
Ahora sí cargamos nuestro modulo:
```
module load fsl/6.0.7.4
```
Despues usaremos el siguiente código ya con todos los pasos anteriores:
```
eddy_cuda10.2 /
--imain=denoised_64A_dwi.nii.gz /
--mask=mascara_64A_dwi.nii.gz /
--index=indice_64A_dwi.txt /
--acqp=acqp_64A_dwi.txt /
--bvecs=64A_dwi.bvec /
--bvals=64A_d
--out 64A_dwi
```
No olvides checar tus outputs!
```
ls eddy*
64A_dwi.nii.gz (con denoise y eddy)
64A_dwi_eddy_parameters
64A_dwi_eddy_command_txt
64A_dwi_eddy_post_eddy_shell_alignment_parameters
64A_dwi_eddy_movement_rms
64A_dwi_eddy_post_eddy_shell_PE_translation_parameters
64A_dwi_eddy_outlier_map
64A_dwi_eddy_restricted_movement_rms
64A_dwi_eddy_outlier_n_sqr_stdev_map
64A_dwi_eddy_rotated_bvecs
64A_dwi_eddy_outlier_n_stdev_map
64A_dwi_eddy_values_of_all_input_parameters
64A_dwi_eddy_outlier_report
```
Y voalá! tenemos nuestra nueva imágen!
```
mrview eddy_denoised_65A_dwi.nii.gz
```

Quizás necesites realizar otras correcciones a tus imágenes DWI como Gibb's ringing removal. Pero vuelvo a lo mismo, eso dependerá de si tus imágenes lo necesitan. Si no, con la corrección de denoising, eddy e inhomogeneidades tienes suficiente para empezar a analaizarlas!
:::success
:thumbsup: Success
:::
---
Siempre es bueno aprender a procesar tus imágenes paso por paso para entender el proceso y porque no, crear tu propio código de pre-procesamiento. Sin embargo, el profesor Dr. Luis Concha (Lab C-13) nos hizo la vida mucho mas fácil y creo un script que hace tooooooodo en una sola exhibición!
El primer paso es cargar el modulo `inb_tools`, aun que este modulo debería de estar ya cargado automaticamente.
El script lo puedes mandar a llamar con solo escribir en la terminal `inb_dwi_bruker_preproc.sh` y al dar `enter` podemos ver un manual de que es lo que hace y que opciones tiene. Vemos que utiliza basicamente los mismos pasos que vimos antes, incluyendo el bias field correction:
```
inb_dwi_bruker_preproc.sh
inb_dwi_bruker_preproc.sh <-i dwi.nii.gz> [-i dwi2.nii.gz] <-o outbase>
Take one or more 2D-EPI DWI acquisitions and preprocess them according to:
0. Concatenate the input DWIs if there is more than one input.
1. dwidenoise (mrtrix, Exp2 estimator - Cordero-Grande 2019).
2. eddy (fsl), including eddy_quad for quality check
3. bias-field correction (N4BiasFieldCorrection). Parameters set for rat imaging.
```
Vemos que primero pide un `-i` input (imágen DWI cruda) y despues un `-o` output (tu nueva imágen)
También el script viene con una serie de opciones de acuerdo a tus necesidades. Ya sea el permutar los axes, re-escalar el voxel, corregir el movimiento (muy recomendado) y/o voltear alguno de los vectores. Este ultimo es necesario ya que al convertir desde Bruker, uno de los vectores sale volteado! Hay que corroborar cual es de acuerdo a tus imágenes.
```
Options:
-p Permute axes to 0,2,1,3 (don't do it)
-s <factor> Scale the image voxel dimensions by some factor (e.g. 2, or 10).
Useful for eddy, as it is expecting human data, not from rodents.
-m Perform motion correction (mcflirt) before running eddy.
This is useful for removing image drift during acquisition.
Flip diffusion gradient vector components:
You can use none, one or any combination of the following.
This is useful if your conversion from bruker data messes up the gradients.
-x Flip x component of diffusion gradient direction
-y Flip y component of diffusion gradient direction
-z Flip z component of diffusion gradient direction
-t Keep temporary directory.
```
Listo, una vez que sabemos que hace el script, lo podemos correr! (spoiler, tarda unos minutos)
```
module load ANTs/2.4.4
module load fsl/6.0.7
module load mrtrix/3.0.4
inb_dwi_bruker_preproc.sh -i 64A_dwi.nii.gz -o inb_64A_dwi -m -s 10 -z
```
Veamos nuestros outputs:
```
ls inb*
inb_64A_dwi_d.bval
inb_64A_dwi_d.bvec
inb_64A_dwi_deb.bval
inb_64A_dwi_deb.bvec
inb_64A_dwi_deb.nii.gz
inb_64A_dwi_de.bval
inb_64A_dwi_de.bvec
inb_64A_dwi_de.nii.gz
inb_64A_dwi_d_mask.nii.gz
inb_64A_dwi_d.nii.gz
inb_64A_dwi_de.files
```
Ahora, vas a notar que hay tres archivos `.nii.gz`, `.bvec` y `bval`, pero cada uno tiene le antecede ya sea`d`, `de` y `deb`. ¿Que significa esto? Esto no es nada mas qué los outputs deribados de cada parte del pre-procesamiento y que el script los nombra asi como guía para saber que datos pertenecen a cada paso del pre-procesamiento:
denoising:
```
inb_64A_dwi_d.bval
inb_64A_dwi_d.bvec
inb_64A_dwi_d.nii.gz
```
denoising + eddy:
```
inb_64A_dwi_de.bval
inb_64A_dwi_de.bvec
inb_64A_dwi_de.nii.gz
```
denoising + eddy + bias field correction:
```
inb_64A_dwi_deb.bval
inb_64A_dwi_deb.bvec
inb_64A_dwi_deb.nii.gz
```
...y todos los archivos deribados del eddy:
```
inb_64A_dwi_de.files
```
Y nuestra nueva imágen!
```
mrview inb_64A_dwi_deb.nii.gz
```

Y al final, esta es la imágen que utilizarás para comenzar tus análisis.