Hi Alle, I think that what has changed between 2014 and today is that the coefficients coef are now a 2D array (number of hyperplanes x number of features). In your case, the first direction is of length one, so you could just do:
coef = clf.coef_[0] and your script should work. The code of the Abraham paper cannot be updated, because papers are not living objects. However, we are maintaining a package that encodes all these patterns in higher-level construct: http://nilearn.github.io/ It might be a good idea to use this package, as it is maintained and has quality assurance. Best, Gaël On Thu, May 04, 2017 at 05:52:27PM +0200, Alle Meije Wink wrote: > I have a script to classify MRI perfusion maps from healthy subjects and > patients. For the file IO and the classifier I have started with the example > code in Abraham et al 2014 [https://arxiv.org/pdf/1412.3919.pdf]. > I use the same classifier as in the paper to produce a back-projected map of > classification weights, which I then want to 'unmask' like in the paper: > coef=clf.coef_ > coef=featureselection.inverse_transform(coef) > and > map_name='weights_check.nii.gz' > wmap=np.zeros(mask.shape, dtype=X.dtype) > wmap[mask]=coef > img=nb.Nifti1Image(wmap,np.eye(4)) > img.to_filename(map_name) > But the line "wmap[mask]=coef" throws an error "ValueError: boolean index > array > should have 1 dimension". I tried the example code from the paper and that > works. > Is the 'coef' array of back-projected SVM weights in some way different than > the masked input image? Or am I doing something else wrong? The error suggests > that the mask array is the problem. > The complete script is attached. > Many thanks for your help! > Alle Meije > # -*- coding: utf-8 -*- > """ > classify MR images from controls (CON) and patients (MCI) > """ > import os > import numpy as np > import nibabel as nb > import sklearn as sl > from sklearn.feature_selection import f_classif > from sklearn import svm > def subj_lists( rootdir='/data/a.wink/MR', starts=['CON','MCI'] ): > # make and empty list for each patient group > slists=[] > for i in range(0,len(starts)): > slists.append([]); > # add subjects to the group based on > for root, dirs, files in os.walk(rootdir): > for fname in files: > for i in range(0,len(starts)): > if fname.startswith(starts[i]): > slists[i].append(os.path.join(root,fname)) > print('finished building subject lists') > return slists > def mk_mask ( subj_lists=[] ): > # check whether mask can be loaded > mname='mask_check.nii.gz' > if os.path.isfile(mname): > msk=nb.load(mname).get_data() > # or must be built (values >20% in the subject images and >20% grey > matter in the template) > else: > num_im=0; > for i in range(0,len(subj_lists)): > for fname in subj_lists[i]: > imdata=nb.load(fname).get_data() > if 'msk' not in locals(): > msk=np.absolute(imdata) > else: > msk=msk+np.absolute(imdata) > num_im+=1 > msk/=num_im > msk/=np.amax(msk) > msk[msk<.2]=0 > msk[msk>0]=1 > grey=nb.load('/usr/local/spm8/apriori/grey.nii').get_data() > grey[grey<.2]=0 > grey[grey>0]=1 > msk=msk*grey > img=nb.Nifti1Image(msk,np.eye(4)) > img.to_filename(mname) > msk=msk.astype(bool) > print('finished building mask') > return msk, mname > def load_images( subj_lists=[], mask=[] ): > # load all the images, build a matrix X or subjects (rows) * > inmask-voxels (columns) > num_im=0; > for i in range(0,len(subj_lists)): > for fname in subj_lists[i]: > imdata=nb.load(fname).get_data() > num_im+=1 > print '\r' + str(num_im) +'\r', > if 'the_matrix' not in locals(): > the_matrix=imdata[mask].T > else: > the_matrix=np.vstack((the_matrix,imdata[mask].T)) > print('finished building matrix X') > return the_matrix > def main(): > # get the input data > subjlists = subj_lists() > y=np.concatenate( ([-1 for _ in range(len(subjlists[0]))], > [ 1 for _ in range(len(subjlists[1]))]) ) > mask, mas_name = mk_mask(subjlists) > X = load_images(subjlists,mask) > print "number of subjects, voxels: %d, %d" % X.shape > # select features > featureselection=sl.feature_selection.SelectKBest(f_classif, k=8000) > X_reduced=featureselection.fit_transform(X,y) > # classify > clf=svm.SVC(kernel='linear') > clf.fit(X_reduced,y) > # make discrimination map > coef=clf.coef_ > coef=featureselection.inverse_transform(coef) > map_name='weights_check.nii.gz' > wmap=np.zeros(mask.shape, dtype=X.dtype) > wmap[mask]=coef > img=nb.Nifti1Image(wmap,np.eye(4)) > img.to_filename(map_name) > if __name__ == "__main__": > main() > _______________________________________________ > scikit-learn mailing list > scikit-learn@python.org > https://mail.python.org/mailman/listinfo/scikit-learn -- Gael Varoquaux Researcher, INRIA Parietal NeuroSpin/CEA Saclay , Bat 145, 91191 Gif-sur-Yvette France Phone: ++ 33-1-69-08-79-68 http://gael-varoquaux.info http://twitter.com/GaelVaroquaux _______________________________________________ scikit-learn mailing list scikit-learn@python.org https://mail.python.org/mailman/listinfo/scikit-learn
