Nilearn学习笔记2-从FMRI数据到时间序列

通过前面的学习了解，我们知道最原始的Fmri数据是4维的，包含三围的空间信息和一维的时间。在实际应用中，我们更多的是利用大脑图像时间序列做研究分析，因为无法直接使用fmri数据做相关研究。在我们使用数据之前，我们需要对原始数据做一些数据预处理和变换。
##1. mask
在所有的分析之中，我们第一步所做的事儿都是把四维fmri数据转换为二维矩阵，这个过程称为MASK。通俗的理解就是提取我们能利用的特征。通过mask得到的二维矩阵包含一维的时间和一维的特征，也就是将fmri数据中每一个时间片上的特征提取出来，再组在一起就是一个二维矩阵。如图所示：

在nilearn库中，提供了两个函数计算mask：

(1) nilearn.masking.compute_background_mask for brain images where the brain stands out of a constant background. This is typically the case when working on statistic maps output after a brain extraction
(2)nilearn.masking.compute_epi_mask for EPI images

(两个函数应用不同的数据，如果数据是EPI采样，就用第二个)
下面看例子：

from nilearn import masking
mask = masking.compute_background_mask(r'E:\home\bct_test\NC_01_0001\rs6_f8dGR_w3_rabrat_4D.nii')
print(mask.get_data().shape)
from nilearn.masking import apply_mask
masked_data = apply_mask(r'E:\home\bct_test\NC_01_0001\rs6_f8dGR_w3_rabrat_4D.nii', mask)
print(masked_data.shape)
# masked_data shape is (timepoints, voxels). We can plot the first 150
# timepoints from two voxels

# And now plot a few of these
import matplotlib.pyplot as plt
plt.figure(figsize=(7, 5))
plt.plot(masked_data[:230, 98:100])
plt.xlabel('Time [TRs]', fontsize=16)
plt.ylabel('Intensity', fontsize=16)
plt.xlim(0, 150)
plt.subplots_adjust(bottom=.12, top=.95, right=.95, left=.12)

plt.show()

得到的图像：

代码中可以看到，我选取了第98和99个特征，画出的图就是这样的。

##2. Timeseries
Mask之后的数据我们也可以做研究分析了，但是我们也发现特征实在是太多了，根本就无法下手，或者是计算两太大。一般而言我们需要找出随着时间在变化的量，只有通过变化的量我们才能得到有效的信息。最简单的一个例子我们想研究北京的雾霾与堵车之间是否有关，那么我们就要得到相应的数据后分析数据，可事实是北京每天都堵车，程度也还一样，没有什么变化，也就无法分析了（这个例子有点儿扯）。
在mask之后，接下来要做的任务就是提取时间序列。说到时间序列，不得不说的是图谱。前面mask之后的特征实在是太多，怎样将这些特征与我们已有的经过验证的图谱对应起来。这个时候就涉及到一个重采样。

#导入图谱，目前用的图谱有好几个，具体选择哪个没有定论。
from nilearn import datasets
dataset = datasets.fetch_atlas_harvard_oxford('cort-maxprob-thr25-2mm')
atlas_filename = dataset.maps
labels = dataset.labels

在nilearn中，重采样函数：resample_to_img

例如：from nilearn.image import resample_to_img
Atlas = resample_to_img(atlas_filename, mask, interpolation='nearest')

(这个是针对于已有参考图谱，最后一定要选最邻近的。如果没有参考图谱，有相应的函数)
重新采样之后就可以得到相应的时间序列了，在把时间序列转换为相关矩阵，并画出其图像。

from nilearn.input_data import NiftiLabelsMasker
masker = NiftiLabelsMasker(labels_img=Atlas, standardize=True,
                           memory='nilearn_cache', verbose=5)

time_series = masker.fit_transform\
    (r'E:\home\bct_test\NC_01_0001\rs6_f8dGR_w3_rabrat_4D.nii')

from nilearn.connectome import ConnectivityMeasure
correlation_measure = ConnectivityMeasure(kind='correlation')
correlation_matrix = correlation_measure.fit_transform([time_series])[0]

# Plot the correlation matrix
import numpy as np
from matplotlib import pyplot as plt
plt.figure(figsize=(10, 10))
# Mask the main diagonal for visualization:
np.fill_diagonal(correlation_matrix, 0)

plt.imshow(correlation_matrix, interpolation="nearest", cmap="RdBu_r",
           vmax=0.8, vmin=-0.8)

# Add labels and adjust margins
x_ticks = plt.xticks(range(len(labels) - 1), labels[1:], rotation=90)
y_ticks = plt.yticks(range(len(labels) - 1), labels[1:])
plt.gca().yaxis.tick_right()
plt.subplots_adjust(left=.01, bottom=.3, top=.99, right=.62)
plt.show()

图形：