C:\JUAN\POST\2020-07-30-GAM-MODELS-FOR-NEUROSCIENCE-CODE-101
Modelos GAM para neurociencia - Código 101
2020-07-30
GAMS PARA ENFERMIDADES NEURODEXENERATIVAS
Vale, este é un resumo da primeira parte da tese ou, polo menos, da idea orixinal que tiñamos pensada e que agora está lixeiramente aparcada. En moi poucas palabras, a idea é crear a database inicial a partir dos PET e MRI xa normalizados con Matlab e SPM. Unha vez que esta sexa un data.frame estable con todos os nosos datos (como tamén será para os SCC), poderemos facer modelos de regresión non paramétricos nos que representar funcións que unan puntos coa mesma intensidade estimada de actividade cerebral.
O obxectivo inicial é obter imaxes (gam.vis e gam.check) similares ás que se expoñen debaixo. No código que presento aquí resúmese o proceso para obtelas, aínda que aínda queda unha etapa extra por levar a cabo: comparar estes resultados cos nativos obtidos con SPM, un software implementado con Matlab que é o estándar de mercado no mundiño do diagnóstico por imaxe médica (polo menos no caso da neuroimaxe). Mentres non poida demostrar de maneira fidedigna que a miña metodoloxía é mellor ca a de SPM, só serán imaxes moi bonitas.
#install.packages(c("cowplot","magick"))
library(magick);library(cowplot);library(ggplot2)
p1 <- ggdraw() + draw_image("GAMexample1.png")
p2 <- ggdraw() + draw_image("GAMexample2.png")
p3 <- ggdraw() + draw_image("GAMexample3.png")
p4 <- ggdraw() + draw_image("GAMexample4.png")
plot_grid(
p1, p2,p3, p4,
ncol = 2)

Claro que este é o mesmo problema ca o dos SCC, a comparabilidade a nivel matemático. Claro que nos SCC vou un paso por diante, xa que si podo obter imaxes que, a nivel visual, sexan practicamente idénticas no formato, e polo menos demostrar visualmente que os resultados son similares xa é un resultado positivo. Quizais non sexa o resultado final desexado, pero cando menos estaría demostrando que hai outra vía máis simple de obter resultados similares e desde un enfoque de datos funcionais.
Pero agora o que nos atinxe é obter eses modelos sobre diferentes cortes e ver a súa utilidade a nivel visual e que nos indican. Como comparalos entre eles e como comparalos con SPM é outra historia ben distinta.
PREÁMBULO:
#install.packages(c("mgcv","gamair","oro.nifti","memsic"))
library(mgcv);library(gamair);library(oro.nifti);library(memisc)
VISUALIZACIÓNS BÁSICAS:
Para a visualización de imaxes PET ou MRI eu non recomendaría en absoluto utilizar R. Existen varios paquetes para traballar con imaxes, pero está claro que o seu obxectivo é traballar cos datos subxacentes á imaxe, non tanto coa imaxe en si. Para traballar con imaxes a nivel visual recomendaría máis ben o uso de SPM (clásico), ou ben de MRIcro ou Mango. Calquera deles se parece máis a un editor de imaxes estándar, algo así como un Photoshop de cerebros. De todos os xeitos, ás veces pode ser necesario comprobar sen máis que todo está correcto, extraer algún valor puntual dun píxel, facer unha visualización ortogonal…
Aínda que o número de funcións dispoñibles é limitado, aquí van algunhas das máis útiles:
LECTURA DE IMAXES E DATOS XERAIS:
img_003_S_1059 <- readNIfTI("003_S_1059", verbose = FALSE, warn = -1, reorient = TRUE,
call = NULL, read_data = TRUE)
img_003_S_1059
## NIfTI-1 format
## Type : nifti
## Data Type : 16 (FLOAT32)
## Bits per Pixel : 32
## Slice Code : 0 (Unknown)
## Intent Code : 0 (None)
## Qform Code : 2 (Aligned_Anat)
## Sform Code : 2 (Aligned_Anat)
## Dimension : 79 x 95 x 79
## Pixel Dimension : 2 x 2 x 2
## Voxel Units : mm
## Time Units : sec
VISUALIZACIÓNS MÁIS ÚTILES
Unha lista das diferentes visualizacións que podemos obter, aínda que, como xa dixen, son moito máis completas e interactivas con outros programas.
Axial:
image(img_003_S_1059) #axial por defecto

Coronal:
image(img_003_S_1059,plane="coronal") #coronal

Sagital:
image(img_003_S_1059,plane="sagittal") #sagital

Plano ortográfico:
orthographic(img_003_S_1059,text="Plano Ortografico PET") #plano ortografico
## Warning in min(x, na.rm = na.rm): ningún argumento finito para min; retornando
## Inf
## Warning in max(x, na.rm = na.rm): ningun argumento finito para max; retornando -
## Inf

DATOS ASOCIADOS:
Estas imaxes, por suposto, teñen unha serie de datos asociados (se non, xa me dirás que clase de data science ía facer eu aquí). En xeral, estes datos só me serven de algo en bulk para sacar un data.frame útil. Non obstante, potencialmente podo pedir o valor dun dato concreto nunha coordenada concreta, e isto pode ser útil para comprobar en certos chanzos do desenvolvemento do código que todo vaia correctamente. Por se acaso, aquí queda:
Datos asociados á imaxe:
img_data=img_data(img_003_S_1059)
qplot(as.vector(img_data),
geom="histogram",
main = "Histogram for PET values",
xlab = "PET levels",
fill=I("blue"),
col=I("red"),
alpha=I(.2))

negatives=img_data[img_data<0]
length(negatives) # de hecho hay muchos valores inferiores a 0 por lo que lo más rápido será, en la database, igualarlos a 0 (MRIcro los identifica como 0's d ehecho)
## [1] 7590
Aquí, por exemplo, xa vemos no histograma que hai moitos valores próximos a cero ou mesmo negativos, polo que máis adiante haberá que telo en conta. Isto probablemente deriva dun proceso de normalización que reduciu valores nulos por debaixo de 0.
Coordenadas concretas na imaxe:
#Diferentes coordenadas de la imagen con valores negativos (comprobacion)
img_data[1,2,3]
## [1] 0
img_003_S_1059[45,45,16]
## [1] 0.5555722
img_003_S_1059[10,70,40];img_data[10,70,40] # valdría cualquiera de los dos comandos
## [1] 0.3873922
## [1] 0.3873922
img_data[66,34,24]
## [1] 0.8744069