4.4分类器-All Alike, Yet Different
在本章中,我们将继续研究形状数据,以便复制Haxby等人(2001年)的研究。对于本教程,有一个小帮手函数,可以在以下情况下手动生成我们所生成的数据集:
最初的研究使用了一个所谓的1个最近的相邻分类器,使用相关性作为距离测量。在pymvpa中,这种类型的分类器由knn类提供,这使得可以指定所需的参数。
k最近邻分类器基于训练数据集中的每个样本的样本相似性来执行分类。k的值指定要导出预测的邻居数,dfx设置确定邻居的距离度量,表决从分配给这些邻居的目标集中选择单个标签的策略。
现在我们有了一个分类器实例,通过将DataSet传递给它的train()方法,可以轻松地对其进行训练。
经过训练的分类器随后可用于对未标记的样本进行分类。分类性能可以通过将这些预测与目标标签进行比较来评估。
>>> predictions = clf.predict(ds.samples)
>>> np.mean(predictions == ds.sa.targets)
1.0
我们看到,分类器在我们的数据集上表现得非常好-它甚至连一个预测错误都没有。然而,大多数情况下,我们不会特别感兴趣的一个分类器的预测精度,在同一数据集,它被训练。
相反,我们感兴趣的是分类器对新的、看不见的数据的通用性。原则上,这将允许我们使用它为未标记的数据分配标签。因为我们只有一个数据集,因此需要将其拆分为(至少)两个部分来实现这一点。在最初的研究中,haxby和他的同事将数据集分为奇数和偶数运行的**模式。我们的DataSet在runtype SAMPLE属性中具有以下信息:
>>> print ds.sa.runtype
['even' 'even' 'even' 'even' 'even' 'even' 'even' 'even' 'odd' 'odd' 'odd'
'odd' 'odd' 'odd' 'odd' 'odd']
使用此属性,我们现在可以很容易地将数据集分割为一半。pYMVPA数据集可以与Numpy的ND阵列类似的方式进行切片。以下调用选择样本子集(即数据集的行),其中RunType属性的值是字符串“偶数”或“奇数”。
>>> ds_split1 = ds[ds.sa.runtype == 'odd']
>>> len(ds_split1)
8
>>> ds_split2 = ds[ds.sa.runtype == 'even']
>>> len(ds_split2)
8
现在我们可以重复上述步骤:调用一个数据集一半并与另一个数据集进行预测()并手动计算预测精度的步骤。然而,更方便的方法是让分类器为我们这样做。PyMVPA中的许多对象支持后处理步骤,我们可以使用该步骤从实际结果中计算内容。下面的示例计算分类器预测和存储在数据集中的目标值之间的失配误差。为了进行这项工作,我们不再调用分类器的预测()方法,而是直接“呼叫”分类器直接与测试数据集联系起来。在PyMVPA中,这是一个非常常见的使用模式,我们在本教程的过程中可以看到很多。请注意,我们现在计算错误,因此较低的值表示更准确的分类。
>>> clf.set_postproc(BinaryFxNode(mean_mismatch_error, 'targets'))
>>> clf.train(ds_split2)
>>> err = clf(ds_split1)
>>> print np.asscalar(err)
0.125
在这种情况下,我们选择数据集的哪一半用于培训,哪一半用于测试完全是任意的,因此我们也可以在交换角色后估计传输错误:
>>> clf.train(ds_split1)
>>> err = clf(ds_split2)
>>> print np.asscalar(err)
0.0
我们认为,平均分类器误差非常低,并且我们获得了与原始研究中报告的结果相当的精度水平。
4.4.1 交叉验证
我们刚刚做的是手动将数据集分为培训和测试数据集的组合,给定特定的示例属性-在这种情况下,**模式还是样本来自偶数运行还是奇数运行。我们对每个分割进行分类分析,以估计分类器模型的性能。一般来说,这种方法称为交叉验证,涉及将数据集划分为多对子集,根据某些标准选择样本组,并通过在分割中的第一个数据集上对分类器性能进行培训,并根据同一拆分测试第二数据集来估计分类器性能
PyMVPA提供了一种方法,以允许完全的交叉验证过程完全自动运行,而无需手动拆分数据集。使用CrossValidation类,通过指定应在每个数据集拆分上计算什么度量以及应如何生成数据集拆分,来设置交叉验证。通常计算的度量是我们在上一节中已经查看过的传输错误。第二个元素是一个用于数据集的生成器,是pYMVPA中的另一个非常常见的工具。以下示例使用半分区器,生成器在使用数据集调用时,将其关联的所有样本标记为数据集的第一个或第二个半部分。这基于指定的示例属性的值-在这种情况下运行类型-非常类似于我们之前执行的手动数据集拆分。半分区器将确保随后将样本分配给两部分,即第一生成数据集中前半部分的样本将位于第二生成数据集的后半部分。使用这两种技术,我们可以轻松地复制手工交叉验证-重用现有的分类器,但不需要定制的后处理步骤。
>>> # disable post-processing again
>>> clf.set_postproc(None)
>>> # dataset generator
>>> hpart = HalfPartitioner(attr='runtype')
>>> # complete cross-validation facility
>>> cv = CrossValidation(clf, hpart)
一旦创建了cv对象,就可以使用DataSet调用它,就像我们以前对分类器所做的那样。它将在内部执行所有数据集分区,将每个生成的数据集拆分为训练和测试集(基于分区),并对分类器进行反复训练和测试。最后,它将返回所有交叉验证折叠的结果。
>>> cv_results = cv(ds)
>>> np.mean(cv_results)
0.0625
实际上,交叉验证结果是作为另一个数据集返回的,该数据集每倍有一个样本,而一个特征具有计算出的每倍传输误差。
>>> len(cv_results)
2
>>> cv_results.samples
array([[ 0. ],
[ 0.125]])
4.4.2任何分类器,真的
迄今为止开发的用于分析的所有代码的简短总结如下:
>>> clf = kNN(k=1, dfx=one_minus_correlation, voting='majority')
>>> cvte = CrossValidation(clf, HalfPartitioner(attr='runtype'))
>>> cv_results = cvte(ds)
>>> np.mean(cv_results)
0.0625
看看这个小代码片段,我们可以很好地看到交叉验证的分类分析的逻辑部分。
1.加载数据;
2.选择分类器;
3.设置错误函数
4.在交叉验证过程中评估错误
5.检查结果
我们以前选择的分类器是以复制Haxby等人的意图为指导的。(2001),但是如果我们想尝试一种不同的算法呢?在这种情况下,pymvpa的另一个好特性开始发挥作用。所有分类器都实现了一个公共接口,使得它们可以很容易地互换,而不需要调整分析代码的任何其他部分。例如,如果我们想在我们的示例数据集中尝试流行的svm(支持向量机),它如下所示:
>>> clf = LinearCSVMC()
>>> cvte = CrossValidation(clf, HalfPartitioner(attr='runtype'))
>>> cv_results = cvte(ds)
>>> np.mean(cv_results)
0.1875
我们使用流行的libSVM库(请注意,pYMVPA提供了额外的SVM实现),而不是k-最近的邻居,而是创建了一个线性SVM分类器。其余代码保持相同。SVM及其默认设置似乎比简单的KNN分类器要稍微差一些。我们很快就会回到分类器了。让我们首先看看这个分析的剩余部分。
我们已经知道CrossValidation可用于计算错误。迄今为止,我们只使用实际目标与分类器预测之间的平均失配数作为误差函数(缺省值为默认值)。但是,PyMVPA在ErrorRFx模块中提供了许多可选功能,但是指定自定义函数也是微不足道的。例如,如果我们不想报告错误,而是精度,我们可以执行以下操作:
>>> cvte = CrossValidation(clf, HalfPartitioner(attr='runtype'),
... errorfx=lambda p, t: np.mean(p == t))
>>> cv_results = cvte(ds)
>>> np.mean(cv_results)
0.8125
此示例重用我们以前创建的svm分类器,并从前面的结果中得到我们期望的结果。
交叉验证过程的细节也是高度可定制的。我们已经看到,分区器用于为每个交叉验证折叠生成培训和测试数据集。到目前为止,我们只使用半分区器将数据集划分为奇数和偶数运行(基于自定义的SAMPLE属性runtype)。但是,通常情况下,执行所谓的“留出一出交叉验证”leave-one-out cross-validation更常见,其中选择数据集的一个独立部分作为测试数据集,而其他部分构成培训数据集。此过程将被重复,直到所有部分都被用作测试数据集一次。在我们数据集的情况下,我们可以认为12个运行中的每一个作为独立的测量结果(fMRI数据不允许我们考虑在时间上相邻的数据被认为是独立的)。
要运行这样的分析,我们首先需要重做我们的数据集预处理,因为在当前的数据集预处理中,对于奇数和偶数运行,每个刺激类别只有一个样本。要获取每个运行的每个刺激类别一个样本的数据集,我们需要修改平均步骤。使用我们从上一个教程部分学到的内容,以下代码段应该是合理的:
>>> # directory that contains the data files
>>> datapath = os.path.join(tutorial_data_path, 'haxby2001')
>>> # load the raw data
>>> ds = load_tutorial_data(roi='vt')
>>> # pre-process
>>> poly_detrend(ds, polyord=1, chunks_attr='chunks')
>>> zscore(ds, param_est=('targets', ['rest']))
>>> ds = ds[ds.sa.targets != 'rest']
>>> # average
>>> run_averager = mean_group_sample(['targets', 'chunks'])
>>> ds = ds.get_mapped(run_averager)
>>> ds.shape
(96, 577)
在整个数据集中,每个类别没有两个样本,现在我们每个类别都有一个样本,每个实验运行,因此在整个数据集中有96个样本。为了建立一个12倍离开一个运行交叉验证,我们可以使用nfoldparater。默认情况下,将一次从一个块中选择示例:
>>> cvte = CrossValidation(clf, NFoldPartitioner(),
... errorfx=lambda p, t: np.mean(p == t))
>>> cv_results = cvte(ds)
>>> np.mean(cv_results)
0.78125
我们获得了几乎相同的预测精度(重用svm分类器和自定义错误函数)。请注意,这一次我们对更多的样本进行了分析,每个样本都是从几个fMRI卷(每个大约9个)中计算出来的。
到目前为止,我们只研究了平均精度或误差。让我们进一步调查交叉验证分析的结果。
返回的值实际上是具有所有交叉验证折叠的结果的数据集。由于我们的错误函数仅计算每个折叠的单个标量值,因此数据集仅包含单个特征(在这种情况下,精度),并且每个折叠的采样数
4.4.3我们需要仔细观察
现在,我们已经使用两种不同的分类器和不同的预处理策略进行了一些交叉验证分析。在所有这些情况下,我们只研究了泛化性能或误差。然而,错误率隐藏了许多有趣的信息,这对于结果的解释是非常重要的。在我们的案例中,我们分析了具有八个不同类别的数据集。平均误分类率没有告诉我们每个类别对预测误差的贡献。可能是每个类别的一半样本被错误分类,但是相同的平均误差可能是由于所有类别的样本被完全错误地分类,而来自其余类别的样本的预测精度是完美的。这两个结果必须以完全不同的方式来解释,尽管相同的平均错误率。
在心理学研究中,这类结果通常被表示为一个contingency table或预期结果与经验结果的cross tabulation 。信号检测理论提供了一系列的技术来表征这样的结果。从这个角度看,分类分析和心理实验几乎没有什么不同,在心理实验中,人类观察者执行检测任务,因此同样的分析程序也适用于这里。
Pymvpa提供了对混淆矩阵的方便访问,即目标的应急表与实际的预调。但是,为防止浪费CPU时间和内存,默认情况下不会计算它们,而是必须显式启用。类似于数据集的SA和FA的条件属性的专用集合中可以使用这样的可选分析结果,它命名为CA。让我们看看它的工作方式:
>>> cvte = CrossValidation(clf, NFoldPartitioner(),
... errorfx=lambda p, t: np.mean(p == t),
... enable_ca=['stats'])
>>> cv_results = cvte(ds)
通过Enable_ca参数,我们为所有交叉验证折叠触发了计算混乱表,但否则代码中没有任何更改。之后,CA集合中提供了整个交叉验证过程的聚合混乱。让我们看一看(注意,在印刷手册中,由于页面宽度的限制,输出被截断-请参考基于html的版本,完整的矩阵)。
此输出是对所执行的分析的全面总结。我们可以看到,混淆矩阵具有很强的对角线,混淆主要发生在小对象之间。除了简单的应急表外,还有一些有用的汇总统计数据,包括平均准确性。
特别是对于多类数据集,矩阵很快变得难以理解。对于这些情况,混淆矩阵也可以通过它的plot()方法来绘制。如果将混淆用作进一步处理的输入,也可以纯矩阵格式访问:
>>> print cvte.ca.stats.matrix
[[ 6 0 3 0 0 5 0 1]
[ 0 10 0 0 0 0 0 0]
[ 0 0 7 0 0 0 0 0]
[ 0 2 0 12 0 0 0 0]
[ 0 0 0 0 12 0 0 0]
[ 2 0 1 0 0 6 0 0]
[ 2 0 1 0 0 0 12 1]
[ 2 0 0 0 0 1 0 10]]
分类器融合仅仅是由PyMVPA中的许多对象支持的条件属性的一般机制的示例。
4.5看这里和那里-探照灯Searchlights
在分类器中,我们看到了我们如何实现分类分析,但是我们仍然不知道我们感兴趣的信号在大脑(或我们所选择的玫瑰)中的位置。尽管我们使用不同的分类器对数据进行了反复的分析,但是我们还是研究了错误率和混淆矩阵。那我们能做什么呢?
理想地,我们希望有某种方式来估计每个特征的分数,该分数指示特定特征(体素的大部分时间)在特定分类任务的上下文中的重要性。在PyMVPA中存在这样的每个特征分数的向量的各种可能性。我们可以简单地计算每个特征的ANOVAF-得分,得到分数,这将告诉我们,在我们的数据集中的任何类别之间,哪些特征显著地变化。
在我们看一下实现细节之前,让我们先重新创建我们的预处理演示数据集。该代码非常类似于4.4类器-所有相同,但不同,不应提出任何问题。我们得到一个数据集,每次运行都有一个样本。
>>> from mvpa2.tutorial_suite import *
>>> ds = get_haxby2001_data(roi='vt')
>>> ds.shape
(16, 577)
4.5.1措施
现在我们有了数据集,计算期望的anovaf分数相对简单:
>>> aov = OneWayAnova()
>>> f = aov(ds)
>>> print f
<Dataset: [email protected], <fa: fprob>>
如果上面的代码片段并不奇怪,那么您可能就有了基本的想法。我们创建了一个对象实例AOV作为一个onewayanova。随后使用数据集调用此实例,并生成打包到数据集中的f-得分。我们以前在哪见过这个?右(边),正确的!这个方法与对交叉验证的调用没有什么不同。这两个对象都是实例化的对象(可能带有一些自定义参数),并在使用输入数据集调用时在DataSet中生成结果。这被称为处理对象,是pymvpa中的一个常见概念。
但是,两个处理对象之间存在差异。交叉验证返回具有单个特征的数据集-准确性或错误率,而一次性验证则返回每个特征值为一个的向量。后者被称为特征-智慧。但是,除了返回的数据集中的特性数量之外,没有太大的差别。例如,pymvpa中的所有措施都支持可选的后处理步骤.在实例化度量值期间,可以指定在内部调用任意映射器,以便在返回结果之前对结果进行前向映射。如果由于某种原因,f-得分需要缩放到间隔[0,1],则可以使用fxmapper来实现这一点。
>>> aov = OneWayAnova(
... postproc=FxMapper('features',
... lambda x: x / x.max(),
... attrfx=None))
>>> f = aov(ds)
>>> print f.samples.max()
1.0
现在我们知道如何计算分数的特征智慧,我们可以开始担心它们。我们最初的目标是解密在脑**的多变量模式中编码的信息。但现在我们正在使用一个单变量量度来定位重要的体素?一定还有别的东西--还有!
4.5.2搜索.
Kriegeskorte et al. (2006) 提出了一种算法,该算法占用脑体素的小的球形邻域,并计算多元度量,以量化在其模式中编码的信息量。在此之后,探照灯方法被扩展到在大脑的每一个可能的区域中运行一个完整的分类器交叉验证。自那时以来,许多研究都采用了这种方法,以当地限制的方式本地化有关信息。
我们几乎都知道在PYMVPA中实施探照灯分析。我们可以加载和预处理数据集,我们可以建立一个交叉验证过程。
>>> clf = kNN(k=1, dfx=one_minus_correlation, voting='majority')
>>> cv = CrossValidation(clf, HalfPartitioner())
唯一要做的是,我们必须把数据集分割成所有可能与大脑相交的球体邻里。为此,我们可以使用sphere_searchlight():
>>> sl = sphere_searchlight(cv, radius=3, postproc=mean_sample())
此单行配置一个探照灯分析,在数据集中的每个可能的范围内运行完全交叉验证。每个球体的半径为三个体素。该算法使用存储在输入数据集的特征属性中的坐标(缺省的voxelindices)来确定本地邻域。从postproc参数中,您可能已经猜到这个对象也是一个度量-而且您是正确的。此度量返回由基本度量(这里是交叉验证)计算的任何值,并将其分配给表示输出数据集中球体中心的特性。对于这个初始示例,我们不感兴趣的是完整的交叉验证输出(每个折叠的错误),而只对平均错误感兴趣,因此我们使用适当的映射器进行后处理。与任何其他处理对象一样,我们必须使用DataSet调用它来运行实际分析:
>>> res = sl(ds)
>>> print res
<Dataset: [email protected], <sa: cvfolds>, <fa: center_ids>, <a: mapper>>
就这样了。然而,这只是一个玩具的例子,只有我们的腹颞ROI。现在让我们在更大的卷上运行它,这样我们就可以真正地本地化一些东西(即使加载和预处理也需要几秒钟)。我们将重用相同的探照灯设置,并在此数据上运行它。由于数据的大小,根据系统中CPU的数量,计算结果可能需要几分钟时间。
>>> ds = get_haxby2001_data_alternative(roi=0)
>>> print ds.nfeatures
34888
>>> res = sl(ds)
现在让我们看看我们得到了什么。由于具有35K元素的矢量有点难以理解,所以我们不得不求助于一些统计信息。
>>> sphere_errors = res.samples[0]
>>> res_mean = np.mean(res)
>>> res_std = np.std(res)
>>> # we deal with errors here, hence 1.0 minus
>>> chance_level = 1.0 - (1.0 / len(ds.uniquetargets))
正如你将看到的,平均经验误差略低于机会水平。然而,无论如何,我们不会期望在所有领域都有一个完美的分类性能的信号。让我们看看,对于多少个球体,误差越大,两个标准偏差就越小,而不是机会。
>>> frac_lower = np.round(np.mean(sphere_errors < chance_level - 2 * res_std), 3)
因此,在几乎10%的所有球体中,误差远远低于我们预期的对分类器的随机猜测,即超过3000个球体!
4.5.3真的!
现在我们对探照灯分析的可能发生了什么,让我们再做一个,但这次是更熟悉的ROI-全大脑。
您现在已经执行了一些探照灯分析,调查了结果,可能试图解释这些结果。在神经科学方面,您从这些分析中得出了哪些结论?你在大脑中的对象表示方面学到了什么?在这种情况下,我们已经运行了8路分类分析,并查看了脑部数千球形ROI中所有条件的平均错误率。在一些球体中,分类器可以很好地执行,即它可以同样地预测所有样本。然而,这只适用于我们测试过的少数超过30K个球体,并不显示该分类器是否能够对所有条件进行分类,或仅对一些条件进行分类。在绝大多数情况下,我们观察到理论概率水平与零之间的误差,我们不知道导致误差减小的原因。我们甚至不知道哪些样本被错误地分类。
从分类器-所有相同,但不同的,我们知道,有办法摆脱这一困境。我们可以查看分类器的混淆矩阵,以获得更多其他隐藏的信息。然而,我们不能合理地对成千上万的探照灯球这样做(请注意,这不是完全正确的)。看 Connolly etal., 2012 对于一些创造性的使用-案例探照灯)。很明显,探照灯分析可能不是数据探索的结束,而是一次粗略的起飞,因为它提出的问题多于回答的问题。
此外,探照灯不能检测到信号延伸到一个小的地方附近。该属性有效地限制了可以使用此策略的分析范围。一项着眼于全球脑回路的研究很难将分析限制在几立方毫米的脑组织上。正如我们以前看到的,探照灯还有另一个不好的方面。虽然他们为我们提供了多变量定位措施,但也继承了单变量fMRI数据分析的– multiple comparisons.。PYMVPA提供了一种算法,可帮助解决组分析中的问题:组簇阈值GroupClusterThreshold. 。
尽管有这些限制,如果使用得当,探照灯分析可能是一个有价值的探索性工具。pymvpa的探照灯实现功能超出了我们在本教程中看到的内容。它不仅可以运行空间探照灯,而且可以同时考虑多个空间。这一点在多维探照灯中得到了进一步的说明。
4.6做更多的分类器-元分类器
在分类器中-它们都是相同的,但不同的-我们看到,可以将整个交叉验证分析封装到一个单独的对象中,该对象可以与任何数据集一起调用以生成所需的结果。我们还看到,尽管进行了这种封装,我们仍然可以获得大量有关所执行的分析的信息。但是,如果我们想进一步处理交叉验证分析中的数据,会发生什么呢?这似乎很困难,因为我们将整个数据集输入到分析中,并且只有在内部才能将其分割成各个部分。
当然,还有一个解决这个问题的解决方案-一个元分类器。这是一个分类器,它本身没有实现分类算法,而是使用另一个分类器来完成实际的工作。此外,元分类器添加另一处理步骤,该处理步骤在实际的基础分类器看到数据之前执行。
这种元分类器的一个例子是MappedClassifier。它的目的是简单的:在将数据传送到内部基础分类器之前,将映射器应用于训练和测试数据。通过这种技术,可以在交叉验证分析中实现任意的预处理。
在我们开始之前,让我们从分类器中复制数据集,所有这些都是一样的,但不同的:
现在,假设我们想对体素强度本身进行分类,而是在训练数据的奇异向量所跨越的空间中的相同样本上,它将如下所示:
首先,我们注意到代码和结果中几乎没有什么变化-错误略有减少,但仍然具有可比性。关键线是第二条,在这里,我们从svm分类器实例创建mppd分类器,并创建作为映射器实现奇异值分解的svdmapper。
我们知道映射器可以组合成复杂的处理管道,而且由于映射分类器以任何映射器为参数,所以可以在交叉验证过程中实现任意的预处理步骤。假设我们听到传言说,只有SvD向量所覆盖的空间的前两个维度覆盖了“有趣”的方差,其余的是噪声。我们可以很容易地用适当的映射器来检查:
很明显,丢弃的元件不能只是噪音,因为误差大幅度增加了。但可能是分类器无法处理数据。由于这段代码中没有任何特定于实际分类算法的内容,所以我们可以很容易地回到过去为我们服务过的knn分类器。
哦,更糟的是。我们得仔细看一下数据,找出这里发生的事情。
在本教程中,我们研究了分类器。我们已经看到,不管实际算法是什么,所有的分类器都在实现相同的接口。因此,它们可以被另一个分类器替换,而不必更改分析代码的任何其他部分。此外,我们已经看到,可以启用和访问处理管道的特定部分提供的可选信息。
4.7分类模型参数.灵敏度分析
在这里和那里-探照灯我们第一次尝试在大脑中定位与特定分类分析相关的信息。虽然我们比较成功,但我们遇到了一些问题,也不得不等一下。在这里,我们希望查看另一种本地化方法。在开始之前,我们对数据进行预处理,因为我们以前做过并执行体积平均,以获得每个刺激类别和原始实验会话的单个样本。
在这个数据集上的探照灯分析看起来就像我们在这里和那里看到的一样-探照灯,但是由于样本数量的增加,这需要花费更长的时间。作为探照灯分析结果的误差图只提供了近似定位。第一,它被重叠的球体涂抹,第二,球状的玫瑰花可能不能反映大脑中功能亚区的真实形状和范围。因此,它混合和匹配可能不属于一起的事物。通过使用更聪明的探照灯算法(参见基于fmri数据的表面探照灯),这可以在一定程度上得到缓解。但是,如果我们能够获得一个每个特征度量,其中每个值都可以真正地归因于各自的特征,而不仅仅是它周围的一个区域,那就更好了。
4.7.1这是一种魔法
获得这种测量的一种方法是检查分类器本身。每个分类器创建模型以从训练数据映射到相应的目标值上。在该模型中,分类器通常将某种权重与作为对分类器决策的影响的指示的每个特征相关联。如何从分类器获取此信息将成为本教程的主题。
但是,如果我们要检查一个经过训练的分类器,我们首先要训练一个分类器。但是,我们有一个完整的大脑数据集,里面有将近40K的特征。我们能做到吗?那么,让我们试一试(并希望您正在运行的机器上仍然有保修.)。
>>> clf = LinearCSVMC()
>>> cvte = CrossValidation(clf, NFoldPartitioner(),
... enable_ca=['stats'])
Ready, set, go!
>>> results = cvte(ds)
太快了,不是吗?但这有什么好处吗?
嗯,准确度并不准确,但混淆矩阵似乎没有任何显著的对角线。虽然我们可以很容易地在完整的脑数据集中训练支持向量机,但是它不能构建可靠的预测模型。至少我们在幸运的情况下已经知道数据中有一些信号,因此我们可以将这种失败归因于分类器。在大多数情况下,数据中实际上没有信号.
人们常常声称,分类性能随着特征选择而提高。如果我们能够将数据集简化为重要的数据集,分类器就不再需要处理所有的噪声了。一种简单的方法是计算一个全脑变体,并且只与显示类别间一定程度差异的体素一起使用。从这里到那里的探照灯,我们知道如何计算想要的f-分数,我们可以用它们来手动选择一些阈值的特征。然而,pymvpa提供了一种更方便的方法-特性选择器:
上面的代码片段配置了这样的选择器。它使用一个阿诺瓦度量选择500个特征与最高的f分。有更多的方式来执行选择,但我们将采用这一个目前。基于敏感性的特性选择实例是另一个处理对象,可以用DataSet调用该对象来执行特征选择:
这是我们想要的数据集,因此我们可以重新运行交叉验证并查看它是否有帮助。但首先,后退一步,再次查看此代码段。有一个对象使用数据集调用并返回数据集。不能防止注意到pymvpa或mapper中的度量之间的显著相似。是的,要素选择过程也是处理对象和工作,就像度量或贴图器一样。现在返回分析:
是的!Yes!我们做到了。对于8路分类,几乎80%的正确分类和混淆矩阵具有强烈的对角线。显然,ANOVA选择的特征是正确的。
让我们在数据的一个子集上重复这一分析。我们只选择瓶和鞋样品。在我们刚才所做的分析中,他们比较容易被分类器混淆。让我们看看整个大脑svm是如何处理这个二进制问题的。
>>> bin_demo = ds[np.array([i in ['bottle', 'shoe'] for i in ds.sa.targets])]
>>> results = cvte(bin_demo)
>>> print np.round(cvte.ca.stats.stats['ACC%'], 1)
62.5
不多,但这并不令人惊讶。让我们看看基于anova的特征选择有什么影响。
>>> fsel.train(bin_demo)
>>> bin_demo_p = fsel(bin_demo)
>>> results = cvte(bin_demo_p)
>>> print cvte.ca.stats.stats["ACC%"]
100.0
哇,那是一次跳跃。完美的分类性能,即使相同的类别不能被同一分类器区分,当对所有八个类别进行训练时。我想,很明显,我们选择功能的方式有些可疑-如果不是非法的话。ANOVA度量使用完整的数据集计算f-分数,从而确定哪些特征在整个数据集中显示类别差异,包括我们假定的独立测试数据。一旦我们发现了这些差异,我们就试图用分类器重新发现它们。能够做到这一点并不令人惊讶,而且准确地构成了双浸过程。因此,所获得的预测精度和建立的模型都是完全没有意义的。
4.7.2谢谢你的鱼
为了正确地实现基于anova的特征选择,我们必须只在训练数据集上这样做。使用pymvpa的方法是通过一个特性选择分类器FeatureSelectionClassifier:
>>> fclf = FeatureSelectionClassifier(clf, fsel)
这是一个元分类器,它只需要两样东西:一个基本分类器来做实际的分类工作和一个特征选择对象。我们可以简单地重用我们已经拥有的对象实例。现在我们有了一个元分类器,它可以和其他分类器一样使用。最重要的是,我们可以将其插入到交叉验证过程中(几乎与我们在开始时所使用的过程相同)。
>>> cvte = CrossValidation(fclf, NFoldPartitioner(),
... enable_ca=['stats'])
>>> results = cvte(bin_demo)
>>> print np.round(cvte.ca.stats.stats['ACC%'], 1)
70.8
这是一个更糟糕和更接近真相-或所谓的无偏估计分类器模型的泛化性。现在,我们还可以在原始的8类数据集上运行这个改进的过程。
>>> results = cvte(ds)
>>> print np.round(cvte.ca.stats.stats['ACC%'], 1)
78.1
>>> print cvte.ca.stats.matrix
[[ 5 0 2 0 0 4 0 2]
[ 0 10 0 0 0 0 0 0]
[ 0 0 8 0 0 1 0 0]
[ 2 2 0 12 0 0 0 0]
[ 0 0 0 0 12 0 0 0]
[ 1 0 1 0 0 7 0 0]
[ 0 0 1 0 0 0 12 1]
[ 4 0 0 0 0 0 0 9]]
这仍然是一个值得尊敬的准确性8路分类和混淆表也证实了这一点。
4.7.3解剖分类器
但是现在回到我们最初的目标:获取分类器对数据集中特征的重要性的意见。采用上述方法,对分类器进行了500个特征的训练。我们只能有它对这些的看法。虽然这仅仅是一个典型的探照灯球体的几倍大,但我们已经解除了探照灯的空间限制-这些特征可以来自所有的ROI。
然而,我们仍然希望考虑更多的特性,因此我们正在改变特征选择以保留更多的功能。
精度下降8%,约为特征数量的4倍。这次我们要求获得5%的分数。
但我们怎样才能最终得到重量呢?在pymvpa中,许多分类器都伴随着所谓的灵敏度分析器。这是一个知道如何从特定分类器类型中获取它们的对象(因为每个分类算法都将它们隐藏在不同的位置)。要创建这个分析器,我们只需让分类器进行如下操作:
正如你所看到的,这甚至适用于我们的元分类器。同样,此分析器是一个处理对象,在使用DataSet调用时返回所需的灵敏度。
为什么我们从分类器得到28幅敏感图?支持向量机为二值分类问题建立了一个模型。为了能够处理这个8类数据集,数据在内部被分割成所有可能的二进制问题(其中正好有28个)。提取了所有这些局部问题的灵敏度。
如果您对这个级别的细节不感兴趣,我们可以将映射组合成一个,就像我们以前对DataSet示例所做的那样。一种可行的算法可能是取任一地图中绝对灵敏度的每个特征最大值。得到的地图将表明特征对于某些部分分类的重要性。
>>> sens_comb = sens.get_mapped(maxofabs_sample())
你可能已经注意到了我们最近计算全脑敏感度图的方法中的一些不完善之处。我们从完整的数据集中派生出来,而不是从数据的交叉验证拆分中派生出来。用一种元测量来纠正这一点很容易。元度量类似于元分类器:采取基本度量、向其添加处理步骤并表现为度量本身的度量。我们想要使用的元度量是重复度量。
我们重新创建我们的基本灵敏度分析器,这一次自动应用后处理步骤,将所有部分分类的灵敏度映射组合在一起。最后,我们将其插入到元度量中,该度量使用nfoldparder生成的分区来提取数据集的每个折叠的训练部分。之后,我们可以运行分析器,并得到另一个数据集,这一次每个交叉验证分割都有一个灵敏度映射。我们可以用类似的方式组合这些地图,但是让我们来看看Anova特征选择的稳定性。
使用映射重叠助手,我们可以很容易地计算出在所有数据集分割中具有非零敏感性的特征的分数。
这可能是数据处理的结束。然而,通过使用元度量来计算灵敏度映射,我们失去了一种方便的方法来访问底层分类器的总体性能。为了再次获得它的访问权限,同时获得敏感度,我们可以稍微扭曲处理流水线。
我想这是值得解释的。我们用一个新的东西-分裂分类器来包装我们的特征选择分类器。这是另一个元分类器,它执行数据集的拆分,并在每个数据集上分别运行训练(和预测)。它可以在内部有效地执行交叉验证分析,我们要求它计算一个混淆矩阵。下一步是为这个元分类器获得一个灵敏度分析器(这一次没有后处理)。一旦我们得到了这一点,我们可以运行分析,并从所有内部培训的分类器获得灵敏度地图。此外,元灵敏度分析器还允许访问其内部元分类器,为我们提供混乱的统计信息。嗯!
虽然我们正在讨论,但值得注意的是,上述场景可以进一步扩展。我们可以向分类器中添加更多的选择或预处理步骤,例如将数据投影到pca组件上,并将分类器限制为每个拆分的前10个组件。pymvpa提供更复杂的元分类器(例如)在某些分析场景中可能非常有用。
4.7.4结束语
我们已经看到,灵敏度分析是一种有用的方法来定位信息,它比探照灯分析更少限制和要求。具体来说,我们可以使用它来发现分布在整个特征集(如全脑)的信号,但我们也可以用它进行基于ROI的分析。这是较少的计算要求,因为我们只训练分类器的一组特征,而不是数千,从而大大减少了所需的CPU时间。
然而,也有一些警告。在构造的模型中,灵敏度是表征特征重要性的一种更为直接的度量手段,但接近分类器的裸金属也存在问题。根据实际的分类算法和数据预处理,当不同的分类器进行比较时,敏感性可能意味着完全不同。例如,流行的svm算法通过识别最难建模的数据样本来解决分类问题。提取的灵敏度反映了这一特性。其他算法,如“高斯朴素贝叶斯”(GNB),对样本在每个类别中的分布进行了假设。GNB敏感性看起来可能完全不同,即使GNB和svm分类器都具有相当的精度。
还应注意的是,即使敏感性来自相同的算法并且只是在不同的数据集拆分上计算,也不能彼此直接比较。在分析中,必须先将它们正常化。例如,pymvpa提供了l1normed()和l2normed(),它们可以与fxmapper一起使用以将其用作后处理步骤。
在本教程中,我们还触及了另一个重要主题:特性选择。我们在分类之前进行了基于ANOVA的特征选择,以帮助svm获得可接受的性能。人们可能会想,这是否是一个聪明的想法,因为在多变量分析之前的单变量特征选择步骤与识别多变量信号的目标有些矛盾。只有功能将保留,以显示自己的一些信号。如果这是一个特定分析的问题,pymvpa为特性选择提供了许多多元的选择。有一个递归特性选择(RFE)的实现和一个帮助器来简化对一个典型用例(Splitrfe)的调用。所有的分类器灵敏度都可以用来选择f。
有了这些构建块,就可以运行相当复杂的分析。然而,解释结果可能并不总是直截了当的。在下一篇教程中,我们将着手消除以前执行的所有分析的另一个限制。我们将超越空间分析,探索时间维度。