pymvpa_test

#Now I download file into pymvpa

from mvpa.suite import *

#filename of source fMRI timeseries image for each run

fmri_src1 = os.path.join(pymvpa_dataroot, '/home/cabiubuntu/temiadelore/s101CS_12runs/4D_s101
CS.nii.gz')

# load the samples attributes as usual and preserve the
# literal labels

attr = SampleAttributes(
                os.path.join(pymvpa_dataroot,
                             '/home/cabiubuntu/temiadelore/s101CS_12runs/attributes.txt'),
                        literallabels=True)

#setting up timeseries events
evs = [ev for ev in attr.toEvents()
            if ev['label'] in ['AI', 'VI', 'AC', 'VC', 'NULL']]
           

#Since we might want to look at the sensitivity profile ranging from just a little before to a little after each event,
#I will extract starting from ove volume prior to the event and one
#volume after the events i.e. two consecutive volumes are extracted a as sample for each event

verbose(1, "Segmenting timeseries into events")
fmri_src_ev = ERNiftiDataset(samples=fmri_src,
                    events=evs,
                    mask=None,
                    evconv = True,
                    enforce_dim=None,
                    labels_map={'AI': 1,
                                'VI': 2,
                                'AC': 3,
                                'VC': 4,
                                'NULL' : 0})

#Begin Preprocessing 

# detrend on full timeseries
detrend(fmri_src_ev, perchunk=True, model='linear')


#We perform normalization by Z-scoring the data 
# using rest/NULL as baseline

zscore(fmri_src_ev, perchunk=True)

#Using meta-classifier
#using k nearest neighbor as clasifier
from mvpa.clfs.knn import kNN
from mvpa.datasets.splitters import NFoldSplitter, OddEvenSplitter
from mvpa.clfs.meta import SplitClassifier, FeatureSelectionClassifier

#Performing a cross-validated sensitivity analysis.
#NFoldSplitter used to split the dataset prior to training
#cvtype = 1 is the leave one out cross validation split

clf = FeatureSelectionClassifier(
           kNN(k=5),
           SensitivityBasedFeatureSelection(OneWayAnova(transformer=N.abs),
                FractionTailSelector(0.05, mode='select', tail='upper')),
	  enable_states = ['feature_ids'],descr="kNN on 5%(ANOVA)")

sclf = SplitClassifier(clf, NFoldSplitter(cvtype = 1),
       enable_states=['confusion', 'training_confusion','feature_ids','trained_dataset'])

#Now train data
sclf.train(fmri_src_ev)

#The features chosen to predict test data
final_dataset = fmri_src_ev.selectFeatures(sclf.feature_ids)

 
#Checking the performance of the classifier
print sclf.confusion.asstring(description=True) \

#or plot the performance results
import pylab as P
sclf.confusion.plot() \
P.show() \

#Calculate the cross validated transfer error of the trained classifier
from mvpa.algorithms.cvtranserror import CrossValidatedTransferError
from mvpa.datasets.splitters import NFoldSplitter

clf = KNN(k=5) # you can change classifier type
terr = TransferError(clf)
cvterr = CrossValidatedTransferError(terr,
                                      NFoldSplitter(cvtype=1),
                                      enable_states=['confusion'])
error = cvterr(fmri_src_ev)

pymvpa_testupdate

#Now I download file into pymvpa

from mvpa.suite import *

#filename of source fMRI timeseries image for each run

fmri_src = os.path.join(pymvpa_dataroot, '/home/cabiubuntu/temiadelore/s101CS_12runs/4D_s101CS.nii.gz')
mask = os.path.join(pymvpa_dataroot, '/home/cabiubuntu/temiadelore/ROIS/AudVis_ROI.nii')
T1 = os.path.join(pymvpa_dataroot, '/home/cabiubuntu/temiadelore/s101CS_12runs/TI.nii.gz')

# load the samples (each scan)  attributes as usual and preserve the
# literal labels

attr = SampleAttributes(
                os.path.join(pymvpa_dataroot,
                             '/home/cabiubuntu/temiadelore/s101CS_12runs/attributes.txt'),
                        literallabels=True)

# actual labels are not important here, could be 'labels=1'
pre_fmri_src = NiftiDataset(samples=fmri_src, labels=attr.labels, chunks=attr.chunks)

# convert to floats
pre_fmri_src.setSamplesDType('float32')

#Begin Preprocessing

# detrend on full timeseries
detrend(pre_fmri_src, perchunk=True, model='linear')


# zscore dataset relative to baseline ('rest') mean
zscore(pre_fmri_src, perchunk=True, baselinelables='NULL', targetdtype='float32')

#load the reference function
from mvpa.misc.io.base import *
ev = SampleAttributes(os.path.join(pymvpa_dataroot,'/home/cabiubuntu/temiadelore/s101CS_12runs/events.txt'),header=['labels','onset','chunks'])
evs = ev.toEvents()

#Segmenting timeseries into events

fmri_src_ev = ERNiftiDataset(samples=pre_fmri_src.map2Nifti(),
                    events=evs,
                    mask=mask,
                    evconv = False,
                    enforce_dim=None,
                    labels_map={'AI': 1,
                                'VI': 2,
                                'AC': 3,
                                'VC': 4,
                                'NULL' : 0})

# remove baseline 'rest' samples from dataset for analysis
fmri_src_ev = fmri_src_ev.selectSamples(N.array([l != 0 for l in fmri_src_ev.labels],
                                       dtype='bool'))

#get dataset summary
fmri_src_ev.summary()

#Using meta-classifier to select features and choose classifier type
#using k nearest neighbor as classifier (k = 5 , number of neighboring voxels used)
from mvpa.clfs.knn import kNN
from mvpa.datasets.splitters import NFoldSplitter, OddEvenSplitter
from mvpa.clfs.meta import SplitClassifier, FeatureSelectionClassifier

#Performing a cross-validated sensitivity analysis.
#NFoldSplitter used to split the dataset prior to training
#cvtype = 1 is the leave one out cross validation split

clf = FeatureSelectionClassifier(
           kNN(k=5),
           SensitivityBasedFeatureSelection(OneWayAnova(transformer=N.abs),
                FractionTailSelector(0.05, mode='select', tail='upper')),
	  enable_states = ['feature_ids'],descr="kNN on 5%(ANOVA)")


sclf = SplitClassifier(clf, NFoldSplitter(cvtype = 1),
       enable_states=['confusion', 'training_confusion','feature_ids','trained_dataset'])

#need to register pickler/unpickler functions for ufuncs

import copy_reg

import numpy

def ufunc_pickler(ufunc):
    return ufunc.__name__

def ufunc_unpickler(name):
    import numpy
    return getattr(numpy, name)

copy_reg.pickle(numpy.ufunc, ufunc_pickler, ufunc_unpickler)

import cPickle

cPickle.dumps(numpy.add)
'cnumpy.core.umath\nadd\np1\n.'

cPickle.loads(cPickle.dumps(numpy.add))

#Now train data
sclf.train(fmri_src_ev)

#The features chosen to predict test data
final_dataset = fmri_src_ev.selectFeatures(sclf.feature_ids)

#convert to nfti image features selected
nftifeaures = fmri_src_ev.map2Nifti(final_dataset)

#visualizing/p;lotting features space on subject's brain
fig = P.figure(figsize=(12, 4), facecolor='white')

plotMRI(background=T1, background_mask=None, cmap_bg='gray', overlay=nftifeatures, overlay_mask=mask, cmap_overlay='autumn', vlim=(None,None), vlim_type=None, do_stretch_colors=True, add_info=True, add_hist=True, add_colorbar=True, fig=fig, interactive=None, nrows=None, ncolumns=None)

P.title('Feature Selection map')

#Checking the performance of the classifier
print sclf.confusion.asstring(description=True) \

#or plot the performance results
import pylab as P
sclf.confusion.plot() \
P.show() \

#Calculate the cross validated transfer error of the trained classifier
from mvpa.algorithms.cvtranserror import CrossValidatedTransferError
from mvpa.datasets.splitters import NFoldSplitter

clf = KNN(k=5) 
terr = TransferError(clf)
cvterr = CrossValidatedTransferError(terr,
                                      NFoldSplitter(cvtype=1),
                                      enable_states=['confusion'])
error = cvterr(fmri_src_ev)