Постройте компоненты PCA, полученные из sklearn.decomposition.PCA, отдельно в трехмерном пространстве

Question

В своем проекте я работаю с трехмерными данными МРТ, где четвертое измерение представляет разные предметы (для этого я использую пакет nilearn). Я используюsklearn.decomposition.PCAчтобы извлечь из моих данных заданное количество основных компонентов. Теперь я хотел бы нанести компоненты отдельно на изображение мозга, то есть я хотел бы показать изображение мозга с моими извлеченными компонентами (в данном случае 2) в разных цветах.

Вот пример кода с использованием набора данных OASIS, который можно загрузить через API nilearn:

1. маскировка с использованием nilearn.input_data.NiftiMasker, который конвертирует мои четырехмерные данные в двухмерный массив (n_subjects x n_voxels).
2. стандартизация матрицы данных с использованием StandardScaler
3. запуск PCA с помощью sklearn.decomposition.PCA:

## set workspace
import numpy as np

from nilearn.datasets import fetch_oasis_vbm
from nilearn.input_data import NiftiMasker
from nilearn.image import index_img

from sklearn.preprocessing import StandardScaler
from sklearn.decomposition import PCA
from sklearn.pipeline import Pipeline

from nilearn import plotting

## Load Data  #################################################################

# take only first 30 subjects as example
oasis_dataset = fetch_oasis_vbm(n_subjects=30)
imgs = np.array(oasis_dataset['gray_matter_maps'])

## PIPELINE ###################################################################

# create a random number generator
rng = np.random.RandomState(42)

# Convert Images to 2D Data Array 
niftimasker = NiftiMasker(mask_strategy='template')

# z-standardize images
scaler = StandardScaler()

# Extract 2 Components
pca = PCA(n_components=2,
          svd_solver='full',
          random_state=rng)

# create pipeline
pipe = Pipeline([('niftimasker',niftimasker),
                 ('scaler',scaler),
                 ('pca',pca)])

# call fit_transform on pipeline
X = pipe.fit_transform(imgs)

Насколько я понимаю, что я получаю после запуска PCA, это загрузки PCA? К сожалению, я не понимаю, как получить из этого два изображения, каждое из которых содержит один компонент PCA.

Answer 1

Чтобы вернуть данные в формат изображения, вам нужно будет выполнить NiftiMasker.inverse_transform(). Для этого необходимо сохранить размеры в пространстве вокселей.

Итак, как сейчас работает конвейер, вы используете уменьшение размерности в пространстве вокселей. На всякий случай, если вы хотите уменьшить размерность в предметном пространстве, вы должны изменить следующее:

pipe = Pipeline([('niftimasker',niftimasker),
             ('scaler',scaler),
#                  ('pca',pca)
            ])

X = pipe.fit_transform(imgs)
X_reduced = pca.fit_transform(X.T).T

Затем вы применили бы обратное преобразование следующим образом:

component_image = niftimasker.inverse_transform(X_reduced)

Затем, чтобы получить изображение каждого отдельного предметного компонента, вы будете использовать index_image из nilearn.image. Например, это изображение для первого предметного компонента:

component1_image = index_img(component_image,0)

Однако я думаю, что вы заинтересованы в уменьшении размерности воксельного пространства. Следовательно, чтобы сохранить размеры вокселей для обратного преобразования, вам нужно будет получить индекс каждой функции вокселя, выбранной при уменьшении размерности PCA. Сохраните конвейер в том виде, в каком он был изначально, и сделайте следующее:

X = pipe.fit_transform(imgs)

components = pca.components_
#In your case 2, but replace range(2) with range(n_components)
most_important = [np.abs(components[i]).argmax() for i in range(2)]

Затем выложите массивы nan с объектами x и вокселями y: (в вашем случае 30 x 229007)

comp1, comp2 = np.tile(np.nan, [30,229007]), np.tile(np.nan, [30,229007])
for x,y in enumerate(X):
    comp1[x,most_important[0]] = y[0]
    comp2[x,most_important[1]] = y[1]

Затем примените обратное преобразование к каждому компоненту:

component1_image = niftimasker.inverse_transform(comp1)
component2_image = niftimasker.inverse_transform(comp2)

Теперь у вас будет 2 изображения, каждое с 30 объектами и 1 допустимым значением вокселя, представляющим выбранный компонент. Вам решать, как агрегировать компонентный воксель по 30 объектам, в этом случае я просто собираюсь использовать функцию среднего изображения из nilearn.image:

mean_component1_image = mean_img(component1_image)
mean_component2_image = mean_img(component2_image)

Наконец, в обоих случаях постройте соответствующее изображение. В версии с уменьшенным вокселем вы увидите небольшое изменение двух изображений в измерении X (вторая диаграмма), но вряд ли Y и Z. Я использую plot_glass_brain из nilearn.plotting:

plotting.plot_glass_brain(mean_component1_image)
plotting.plot_glass_brain(mean_component2_image)

Чтобы использовать наложения, настройте карты цветов, чтобы упростить визуализацию, и другие параметры печати обратитесь к этому и другим руководствам по построению графиков nilearn:

https://nilearn.github.io/plotting/index.html

Позвольте мне знать, если у вас есть еще вопросы.