PyMVPA_manual 4.2&4.3基础知识和处理流程

 

所有环境和包全部准备完毕后,来看mvpa的数据结构

所有volume的数据全部都保存在 dataset中,所有的操作都离不开dataset,所以需要了解pymvpa中的dataset的数据结构

4.2数据集基础和概念 

pymvpa中的大多数数据集被表示为二维数组,其中第一轴是样本轴,第二轴代表样本的特征。在最简单的情况下,数据集只包含数值矩阵的数据。

在上面的示例中,数据矩阵中的每一行向量都成为数据集中的一个观察值,一个样本samples,每个列向量表示一个单独的变量,一个特征feature。样本和特征的概念对于数据集是必不可少的,因此我们将进行更深入的研究。

数据集假定数据的第一个轴将用于定义单个样本。

但是如果数据集只有一维,那么将认为这一维数据集的每个元素是一个样本,只有一个特征,但是这样一维的情况很少见。

另一方面,如果数据集是从多维数据创建的,则只有其第二轴表示特征。

在本例中,我们有一个包含三个样本和四个特征的数据集,其中每个特征都是一个2x3矩阵。

 

4.2.1属性(Attributes

在MVPA上下文中,我们通常需要更多地了解每个示例,而不仅仅是它的特性的价值。例如,为了训练一种有监督的学习算法来区分两类样本,我们需要每个样本的目标值来标记每个样本与其各自的类别。然后可以使用这些信息,例如,将数据集分割成特定的样本组。对于这种类型的辅助信息,DataSet还可以包含三种属性的集合:样本属性、特征属性和数据集属性。

1.样本属性

数据集中的每个样本都可以具有任意数量的附加属性。它们作为与集合中样本数相同长度的向量存储,并可通过sa属性访问。集合类似于标准python的dict字典,因此添加示例属性就像向字典中添加元素一样。

但是,样本属性不是直接存储在普通数据中,而是出于各种原因作为所谓的可收集的Collectable,然后将numpy数组嵌入实际属性:

但除了基本合理性检查之外,不存在对样本属性值的进一步限制。只要属性向量的长度与数据集中的采样数匹配,并且属性值可以存储在numpy数组中,则允许使用任何值。因此,甚至可以将n维数组(而不仅仅是向量)作为属性,只要它们的第一个轴与数据集中的采样数相匹配即可。此外,存储文字(非数字)属性是完全可能和支持的。还应注意,每个属性可以具有其自己的单个数据类型,因此可以在同一数据集中具有文字和数字属性

2.特征属性

特性属性几乎与样本属性相同,唯一的区别是,与每个样本拥有一个属性值不同,特性属性每个feature有一个特征。此外,它们存储在DataSet中的一个单独的集合中,称为fa:

意思是说有几列就有几个feature的值,ds.nfeature的数值是几,那么feature矩阵中的feature的个数就是几。

3.数据集属性

最后,数据集也可以有属性,不是针对每个样本,也不是针对每个特征,而是作为一个整体的数据集:所谓的数据集属性。分配这些属性并在以后访问这些属性的方式与其他两种类型的属性完全相同,但DataSet属性存储在它们自己的集合中,可以通过DataSet的a属性访问它们。然而,与Sample和Feature的Attribute不同的是,没有对类型或大小施加任何约束-任何东西都可以存储。让我们存储一个包含当前目录中所有文件名称的列表,因为我们可以

4.2.2切片、重采样、特征选择

关于特征选择,参看知乎:https://www.zhihu.com/question/28641663

此时,我们已经可以从简单数组构建数据集,并通过任意数量的附加属性丰富数据集。但是仅仅拥有一个数据集是不够的。我们通常需要能够选择数据集的子集进行进一步的处理。需要用到切片处理

切片数据集(即选择特定的子集)与分割numpy数组非常相似。实际上,它的工作原理几乎是一样的。DataSet支持python的切片语法,但也支持布尔掩码和索引的选择。通过始终选择数据集中的每个其他样本,以下三种切片操作将产生等效的输出数据集:

所有三种切片样式都同样适用于数据集中特征子集的选择。记住,特征是在数据集的第二个轴上表示的。

通过将按索引选择的方法应用于第二轴,我们可以很容易地获得示例数据集的最后两个特性。请注意:是为第一轴切片提供的。这是python表示沿此轴获取所有内容的方法,因此包括所有示例。

您可以猜到,同时选择样本和特征子集也是可能的。

    

                            

note:如果你以前有过NumPy的经验,你现在可能会感到困惑。你可能期望的是:  

上面的代码将同样的切片直接应用于.Sample的numpy数组,结果是根本不同的。对于numpy数组,这种类型的切片允许在数组的每个轴上根据它们的索引选择特定的元素。对于pymvpa的数据集,这种模式并不十分有用,相反,我们通常希望选择行和列,即由它们的索引提供的样本和特性。

通过对数据集进行切割,他们相应的属性值也会跟着被切割,同样,如果根据属性值对他们进行切割,数据也会跟着变,也就是说,无论对数据集还是数据属性进行切割,都奏效,都是粘连在一起的。

4.2.3加载fMRI数据

理论足够了,让我们看个例子。pymvpa有几个辅助函数来从专门的格式加载数据,还有一个用于fmri数据的函数fmri_dataset().

最简单的情况下,我们现在利用 fmri_dataset() 完成它的工作,只需将它指向fMRI数据文件。数据存储为一个Nifti文件,该文件将一个实验的所有volume连接到一个文件中。

 data_niifile = os.path.join(data_path,waexample.nii)    #将。nii文件导入
 data_maskfile = os.path.join(data_path, 'mask', example.nii) # 将掩码mask导入

dataper= fmri_dataset(data_niifile,  mask = data_maskfile) #此时的dataper才是真的有掩码,有很多volume的文件

如果不加mask

我们得到了一个两维的数据集,这个数据集的第一个维数是 samples的个数,代表了volume的个数,另外一个维度是feature的个数,这个feature是一个volume中voxel的个数。体素被表示为一维向量,它们似乎已经失去了与三维体素空间的关联。然而,情况并非如此,我们稍后会看到。pYMVPA以这种简单的格式表示数据,使它与大量的通用算法兼容,这些算法预期数据是一个简单的矩阵。

为何要加mask

我们从Nifti文件中加载了所有数据,但通常我们只对一个子集感兴趣,fmri_dataset()能够执行数据屏蔽,我们只需要指定一个掩码图像,这样的掩模图像在相当多的fMRI分析中产生,  它可以是在颅骨剥离过程中计算的全脑掩模,或者是来自功能定位器的激活图。这里的mask是我们提取的感兴趣的roi 

加了mask后,samples的个数没变,但是feature的个数变了,因为加了mask后的voxels的个数就少了,剩下的voxel对应于掩模图像中的非零元素的体素. 现在,让我们进一步研究这个数据集.

 

除了采样之外,数据集还提供了一些属性,这些属性通过在“nifti图像头文件”中存在的信息来增强数据。每个样本都有关于其在时间序列中的体积指数和实际获取时间的信息(相对于文件的开始)。此外,每个特征的原始体素索引(有时称为ijk)也是可用的。最后,数据集还包含有关输入卷的维数、体素大小以及任何其他特定于Nifti的信息,因为它还包括完整的Nifti图像头的转储。

一个体素的大小是3毫米× 3毫米× 3毫米

数据集还带有一个关键的附加属性:映射器。映射器mapper是PyMVPA中的一个重要概念,因此它有自己的教程章节

拥有所有这些属性作为数据集的一部分通常是一件有用的事情,但在某些情况下(例如,涉及到效率和/或非常大的数据集),人们可能希望拥有一个更精简的数据集,其中包含真正必要的信息。实现这一目标的一种方法是去掉所有不需要的属性。DataSet类的Copy()方法可以提供帮助。

我们可以看到除了 time_coords 之外的所有属性都已被过滤掉。将deep 参数设置为false会导致复制功能重新使用源数据集中的数据生成新的剥离数据,而不会复制内存中的所有数据-这意味着两个数据集现在共享样本数据,对 数据集 进行的任何更改也会影响stripped

4.2.4中间存储

用处不大就不写了……

请注意,这种类型的数据集存储不适合数据的长期归档,因为它依赖于一个稳定的软件环境。对于长期存储,请使用其他格式。

4.3获取形状数据

在教程部分,数据集的基础和概念中,我们已经发现了数据集的魔法成分:映射器Mapper。映射器可能是PyMVPA中最强大的概念,没有它们就没有什么可以做的。

通常,映射器是一种变换数据的算法。该变换可以简单地选择数据的子集,或者与多级预处理流水线一样复杂。有些转换是可逆的,一些转换是不可逆的。有些是简单一步计算,有些是迭代算法,必须在使用之前对数据进行训练。在PyMVPA中,所有这些变换都是映射器。

让我们创建一个虚拟数据集(5个示例,12个特性)。这一次,我们将使用一个新的方法来创建数据集,DataSet_wizard。在这里,它完全等价于常规的构造函数调用(即DataSet),但我们很快就会看到一些方便的方面。

有些数据集(例如带有掩码的fMRI_DataSet())包含作为数据集属性.a.mapper.的映射器,但是,并非每个数据集实际上都具有映射器。例如,我们刚创建的简单的虚拟数据集没有

现在让我们来看看一个非常相似的数据集,它只在一个很小但非常重要的细节上有所不同:

我们看到,生成的数据集看起来与上面的数据集相同,但这一次它是从一个3D样本数组中创建的(即五个样本,其中每个样本是一个4x3矩阵)。不知何故,这个3D数组被转换成数据集中的2d样本数组。这种神奇行为是通过观察数据集的映射器是一个FlattenMapper

这个映射器的目的正是我们刚才观察到的:将数据数组整形为2d。它是通过保持第一个轴(在pymvpa数据集中,这是分离样本的轴)并将所有其他轴连接到第二个轴来实现的。

因为映射器代表特定的转换,所以它们也可以被看作是已经完成的工作的协议。如果我们查看数据集,我们就知道它在从原始数据集到DataSet中的样本数组的过程中已经被扁平化了。如果处理变得更加复杂,这个特性就会变得非常有用。让我们看下一个可能的步骤-选择一个有趣特性的子集:

现在情况发生了变化:数据集中出现了两个新的映射器-一个ChainMapper 链表映射器和一个StaticFeatureSelection静态特征重新选择。后者描述(并实际执行)我们刚刚进行的切片操作,而前者将两个映射器封装到一个处理管道中。我们可以看到,映射链表示数据集的处理历史,就像面包屑轨迹一样。

下面这一块没看懂

如前所述,映射器不仅可以转换单个数据集,还可以输入其他数据(只要它与映射器兼容)。

尽管图中的 subds 具有比我们输入的数据更少的特性,但正向映射 forward mapping 应用于DataSet本身的转换也适用于我们的测试4x3数组。该过程生成的特征向量与subds中的特征向量的形状相同。通过查看forward-mapped前向映射的数据,我们可以验证是否选择了正确的特性。

4.3.1加载实际数据

我们几乎所有的部分都可以开始第一次分析。我们知道如何从时间序列图像中加载fMRI数据,我们知道如何在数据集中添加和访问属性,我们知道如何对数据集进行切片,我们知道我们可以使用映射器来操作数据集。

现在,我们的目标是将所有这些小片段合并到生成数据集的代码中,就像Haxby等人在开创性工作中使用的数据集一样。(2001)-在块设计实验中,参与者被动地观看八个对象类别的灰度图像。从原始的粗BOLD时间序列中,我们对第一个被试进行了12次完整的记录,他们计算出:在每一半的数据中每个刺激类别的激活模式pattern of activation(由奇数和偶数运算分割;即16个样本),包括对数据进行交叉验证的分类分析所必需的相关样本属性sample attributes

我们已经看到了如何从Nifti图像加载fMRI数据,但这一次我们需要的不仅仅是epi图像。对于分类分析,我们还需要将每个样本与相应的实验条件相关联,即类标签,有时也称为目标值target。此外,对于交叉验证过程,我们还需要将整个数据集划分为独立的块chunks。独立性是实现对分类器泛化性能的无偏估计的关键,即它在预测新数据的正确类别标签方面的准确性,这在训练过程中是看不见的。那么,我们从哪里得到这些信息呢?

target目标值和块chunks都由实验的设计定义。在最简单的情况下,fmri volume 的samples 样本的target目标值是在获取volume时已经存在/激活的实验条件。但是,我们将在以后查看更复杂的场景。独立数据块对应于假定的fmri volume是独立的。mri采集过程的特性导致随后采集的体积非常相似,因此它们不能被认为是独立的。理想情况下,将实验分成几个采集sessions,其中的sessions定义相应的数据块。

将此信息导入pymvpa有许多方法。最简单的方法是创建一个两列文本文件,其中第一列中有目标值target,第二列中有块chunk标识符,在nifti图像中每卷有一行。

这里是指:target就是我们实验的条件,而chunks的数量与一个run的volume数目是一样的

SampleAttributes允许我们加载这种类型的文件,并访问其内容。我们有121个标签和块值,每个volume一个。另外,我们发现有九个不同的条件,所有的样本都与同一块相关联。不同扫描/运行的属性文件将增加块值。

attr_fname = os.path.join(data_path, 'sub001',

... 'BOLD', 'task001_run001', 'attributes.txt')

attr = SampleAttributes(attr_fname,header=['chunks','targets'])       

'attributes.txt' 是需要提前生成的,根据扫描的时间,分别经历了多少volume,时间上的target(扫volume对应的条件),每个条件的具体信息,是一个run中的第几个trail等信息。  行的数量是一个run中的volume的数量,列数是信息(target,trial等)的数量。

现在,我们可以加载fMRI数据,就像我们以前所做的那样-只加载对应于颞叶腹侧皮层面的体素,并将样本属性分配给数据集。fmri_data()允许我们直接传递它们:

我们从最后一部分中得到了我们已经知道的数据集,但这次也有关于块和目标的信息。

4.3.2增加结构,减少工作重复

虽然可以为每次fMRI扫描创建单独的属性文件,但这样做的效果并不理想。通常情况下,刺激与fMRI容积采样率不同步,因此会丢失定时信息。此外,关于刺激的信息,或者一般的实验设计,很可能已经以不同的形式或形状提供了。

为了方便使用广泛的数据集,pymvpa提供了专用的数据集支持,遵循openfmri.org数据共享平台使用的规范。这些是有关文件名约定和设计规范的简单准则,可以方便地为您自己的数据采用。

通过处理程序访问此类数据集,只需知道数据集存储在磁盘上的位置即可。此处理程序提供对基本信息的方便访问,例如主题数、任务描述和其他属性。

更重要的是,它支持访问实验设计信息:

如您所见,刺激设计信息可以在每个事件的标准python字典列表中获得。这包括刺激的开始和持续时间,以及文字条件标签和任务描述。

使用实用程序函数,将这样的事件列表转换为样本属性数组是很简单的,就像我们以前从文件中加载的那样

events2sample_attr()使用存储在DataSet的样本属性time_coords 中的样本获取时间信息来匹配刺激事件和数据样本。

请注意,将刺激事件转换为属性数组是标记fmri数据的一种相当粗糙的方法,它只适用于块设计类实验。我们稍后将在本教程中看到其他方法。

除了实验设计信息外,数据集处理程序还提供了对实际BoldfMRI数据的方便访问:

get_bold_run_dataset()方法的工作方式与我们以前看到的fMRI_DataSet()相同,并且支持相同的参数。但是,不提供自定义文件名,BOLD数据由被试、任务和获取运行ID标识定义。

多session的数据

许多fMRI实验涉及多次运行。加载这样的数据最好在循环中完成。下面的代码片段从DataSet中的示例主题加载对象查看任务的所有可用运行。

现在是获得数据集的summary()的很好的时间:基本统计、每个chunk目标之间样本数量的平衡等。比如:

在使用openfmri.org数据时,您可以查看如何以更紧凑的方式执行上述数据准备的示例。

下一步是从我们感兴趣的数据集中提取激活模式。但等等!我们知道fMRI数据通常受到许多噪声的污染,或者实际上是我们不感兴趣的信息。例如,在数据中存在时间漂移(当扫描器升温时信号趋于增加)。我们也知道信号在整个大脑中并不是完全均匀的。

所有这些工件都带有很大的方差,即(希望)与实验设计无关,我们应该尝试删除它,以尽可能干净的信号呈现分类器。有无数的方法可以对数据进行预处理以试图实现这个目标。一些关键词有:高/低/带通滤波、去尖峰、运动校正、强度归一化等。在本教程中,我们保持简单。我们刚加载的数据已经进行了运动校正。对于每一个长于几分钟的实验,就像在这种情况下一样,时间趋势去除或去趋势是至关重要的。

4.3.3基本预处理

消除长期趋势Detrending

pymvpa提供了从数据中删除多项式趋势的功能(其他方法也可用),这意味着多项式被拟合到时间序列中,并且只有它们没有解释的东西留在数据集中。在线性去趋势的情况下,这意味着通过线性回归与每个体素的时间序列拟合一条直线,并以残差作为新的特征值。Detrending可以看作是数据转换的一种类型,因此在pymvpa中它被实现为映射器。

我们刚刚创建的是一个映射器,它将执行块向线性(一阶多项式)的去趋势。由于我们的数据来源于12个不同的run,而假设在所有运行过程中都存在一个连续的线性趋势是不合适的,因此,按块的方式减少趋势是可取的。映射程序将使用 chunks 属性来标识数据集中的块。

我们已经看到,我们可以简单地用这个映射器对我们的数据集进行前向映射。但是,如果我们想让映射器出现在处理历史记录的数据集中,我们可以使用它的 get_mapped() 方法。此方法将导致DataSet使用给定的映射程序映射自身的浅副本,并返回该副本。让我们试试:

detrended_fds很容易被识别为已被扁平、切片和线性去趋势的数据集。

注意,detrending 并不总是一个明确的步骤。例如,如果您计划在一般线性模型中使用血流动力学响应函数对数据建模(如使用openfmri.org数据处理nipy),多项式去趋势可以作为建模的一部分同时进行。

 标准化Normalization

虽然这有望解决数据中的时间漂移问题,但我们仍然存在一些非均匀的体素强度,这对于某些机器学习算法来说是一个问题。这对于一些机器学习算法来说可能是一个问题。在本教程中,我们再次遵循一个简单的方法来解决此问题,并执行对数据的逐个特征的、按块的Z-分数转化。这具有许多优点。首先,它将把所有的特征扩展到大约相同的范围内,并去除它们的平均值。后者非常重要,因为一些分类器在处理具有大偏移量的数据时会受到损害。然而,我们不打算执行一个非常简单的z分数去除全局平均值,而是使用数据集的rest条件(基线)样本来估计平均值和标准差。使用这些参数来缩放数据集特性会产生一个分数,该分数对应于每个时间点的体素强度与其余平均值的差值。

这种类型的数据标准化是,您已经猜到了,也被实现为映射器:

此映射器配置实现了按chunk块划分(默认)z评分,同时在相应的数据块中使用“REST”估计与样本目标的均值和标准差。

请记住,所有映射器都返回只具有已修改内容副本的新数据集。然而,去趋势和z评分都有或将修改样本本身.这意味着内存消耗将增加三倍!我们将获得原始数据、去趋势数据和z分数数据,但通常我们只对最后的处理阶段感兴趣。为了减少内存占用,两个映射器都有执行相同处理的兄弟节点,但不复制数据。对于PolyDetrendMapper,这是poly_detrend(),而对于ZScoreMapper,则是zscore()。以下调用将与上面创建的映射器相同,但使用较少的内存:

所得到的数据集现在都是去趋势化和标准化的。信息在映射器中显示得很好。从此,我们不再对REST类别的样本使用,因此我们从数据集中删除这些样本:

计算激活模式Patterns Of Activation

最后一个预处理步骤是计算实际激活模式。在原始研究中,haxby和同事分别对奇数和偶数运行数据进行了glm分析,并使用相应的对比统计(刺激类别与rest)作为分类器输入。在本教程中,我们将使用简单得多的快捷方式,只计算奇数和偶数运行的每个条件的平均样本。

为此,我们首先添加一个新的SAMPLE属性,为DataSet中的每个示例分配相应的标签,该标签指示它属于哪两种运行类型:

其余的都是琐碎的,对于这种情况-应用函数(如 mean)到一组样本组(stimulus category刺激类别和run type运行类型的所有组合)。PyMVPA中的FxMapper 它有许多方便的功能,我们需要的是 mean_group_sample()。它获取样本属性列表,确定其唯一值的所有可能组合,选择与这些组合相对应的数据集样本,并对其进行平均。最后,由于这也是映射器,因此返回了具有平均采样的新数据集:

现在可以了!我们现在有一个完全预处理的数据集:加了mask、去趋势、归一化,每个刺激条件下一个采样分别是奇数和偶数的平均值。现在我们可以做一些正式的分类,这将在Classififiers – All Alike, Yet Different这章显示出来,但是还有一个重要的方面是mappers,我们必须先看看。

4.3.4一遍又一遍-一个Mapper's的故事

让我们回顾一下从教程部件开始的简单数据集。

映射器的一个非常重要的特性是,如果可能的话,它们允许逆转转换。在简单数据集的情况下,我们可以要求映射器撤销扁平,并将我们的样本返回到原来的3D形状。

在交互式脚本编写会话中,这将是一个相对庞大的输入命令,尽管它可能会被相当频繁地使用。为了减少手指的痛苦,有一条同样的捷径:

重要的是要认识到,反向映射不仅适用于单个映射程序,而且还适用于ChainMapper链式映射器。回到我们的演示数据集,我们可以看到它是如何工作的:

通过映射器链的单个采样(一维特征向量)的反向映射创建了一个4x3阵列,该阵列对应于原始数据空间中的样本的尺寸。此外,我们还看到,每个特征值都精确地放置在与以前的数据集切片操作中选择的特征对应的位置上。

但是现在让我们再来看看我们的fMRI数据集。在这里,映射链稍微复杂一点:

Initial flflattening followed by mask, detrending, Z-scoring and fifinally averaging. 在这种情况下,我们会进行反向映射吗?让我们来测试:

所发生的正是我们所期望的:初始的一维向量通过映射链向后传递。恢复基于组的平均对于单个向量没有多大意义,因此它被忽略了。z得分和时差趋势也是如此.但是,对于所有剩余的映射程序,转换都是相反的。首先 unmasked,将mask去掉,然后重塑到原始的维度-大脑体积。 (40,64,64)

我们只需通过链中的前两个映射器进行反向映射,并比较结果,就可以验证是否确实如此:

如果您想知道:ChainMapper的行为类似于一个常规python列表。我们刚刚选择了列表中的前两个映射器作为另一个ChainMapper链表映射器,并将其用于反向映射。

回到Nifti

在反向映射的上下文中的最后一个有趣的方面:每当需要从pymvpa导出数据(如结果)时,DataSet映射器 mappers也起着至关重要的作用。例如,我们可以很容易地将revtest向量导出到NIfTI脑体积图像中。这是可能的,因为映射器可以将其放回3D空间,而且数据集还存储有关原始源nifti图像的信息.

pymvpa提供了map2nifti(),这是一个将这两件事结合起来并将任何向量转换成相应的Nifti映像的函数:

这个映像现在可以存储为一个文件(例如,nimg.to_filename('mytest.nii.gz') )。在这种格式下,它现在与绝大多数神经成像软件兼容。

在pymvpa中有更多的映射程序,比我们在教程部分所能介绍的要多得多。一些更多的将用于其他部分,但更多的是可以找到映射器mappers模块。尽管它们都实现了不同的转换,但它们都可以相同的方式使用,并且都可以组合成一个链。

 

 

  • 1
    点赞
  • 9
    收藏
    觉得还不错? 一键收藏
  • 0
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值