Test Neurolearn Tools

This notebook tests the prediction and apply mask functions of the neurolearn tools
https://github.com/ljchang/neurolearn

Written by Luke Chang


In [1]:
import os
import pandas as pd
import matplotlib.pyplot as plt
import time
import nibabel as nb
import numpy as np
import sklearn
from nltools.analysis import Predict, apply_mask
%matplotlib inline

# Github Project Root
basedir = os.path.abspath('..')

# Set Paths
dest_dir = "/Users/lukechang/Dropbox/NEPA/Test_Analysis"

# Separate image types and sort
combined_df = pd.read_csv('%s/metadata.csv' % dest_dir, encoding='utf8')
neg_list = combined_df.ix[combined_df.name_image.str.contains('Neg'),].sort(columns='name_image')
neu_list = combined_df.ix[combined_df.name_image.str.contains('Neu'),].sort(columns='name_image')

Test Predict Function


In [2]:
tic = time.time() #Start Timer

# Load data using nibabel
neg_file_list = [dest_dir + '/resampled/00' + str(x) + '.nii.gz' for x in neg_list.image_id]
neu_file_list = [dest_dir + '/resampled/00' + str(x) + '.nii.gz' for x in neu_list.image_id]
dat = nb.funcs.concat_images(neg_file_list+neu_file_list)
Y = np.array([1]*len(neg_file_list)+[0]*len(neu_file_list))
sublist = np.array([x[:-7].split('/')[-1] for x in neg_file_list] *2)

# Test SVM with kfold xVal
negvneu = Predict(dat,Y,algorithm='svm',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5}, **{'kernel':"linear"})
negvneu.predict()

# # Test SVM with LOSO xVal
# negvneu = Predict(dat,Y,algorithm='svm',subject_id = sublist, output_dir=dest_dir, cv_dict = {'loso':sublist}, **{'kernel':"linear"})
# negvneu.predict()

# # Test SVM with kfold xVal & Platt scaling
# negvneu = Predict(dat,Y,algorithm='svm',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5}, **{'kernel':'linear','probability':True})
# negvneu.predict()

# Test Ridge Classifier with kfold xVal 
negvneu = Predict(dat,Y,algorithm='ridgeClassifier',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
negvneu.predict()

# # Test Ridge Classifier with LOSO xVal
# negvneu = Predict(dat,Y,algorithm='ridgeClassifier',subject_id = sublist, output_dir=dest_dir, cv_dict = {'loso':sublist})
# negvneu.predict()

# # Test Ridge Classifier with kfold xVal + grid search for regularization
# negvneu = Predict(dat,Y,algorithm='ridgeClassifierCV',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5},**{'alphas':np.linspace(.1, 10, 5)})
# negvneu.predict()

# # Test RandomForest Classifier with kfold xVal - FAILS: NO COEF
# negvneu = Predict(dat,Y,algorithm='randomforestClassifier',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
# negvneu.predict()

# Test Logistic Regression with kfold xVal
negvneu = Predict(dat,Y,algorithm='logistic',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
negvneu.predict()

# # Test LASSO with kfold xVal: - FAILS : empty mask
# negvneu = Predict(dat,Y,algorithm='lasso',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
# negvneu.predict()

# # Test LASSO with kfold xVal + grid search for regularization
# negvneu = Predict(dat,Y,algorithm='lassoCV',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5},**{'alphas':np.linspace(.1, 10, 5)})
# negvneu.predict()

# Test Ridge with kfold xVal
negvneu = Predict(dat,Y,algorithm='ridge',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
negvneu.predict()

# # Test Ridge with kfold xVal + grid search for regularization
# negvneu = Predict(dat,Y,algorithm='ridgeCV',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5},**{'alphas':np.linspace(.1, 10, 5)})
# negvneu.predict()

# # Test Linear Regression with kfold xVal
# negvneu = Predict(dat,Y,algorithm='linear',output_dir=dest_dir)
# negvneu.predict()

# Test LASSO-PCR with kfold xVal
negvneu = Predict(dat,Y,algorithm='lassopcr',subject_id = sublist, output_dir=dest_dir, cv_dict = {'kfolds':5})
negvneu.predict()

print 'Elapsed: %.2f seconds' % (time.time() - tic) #Stop timer


overall CV accuracy: 0.78
overall CV accuracy: 0.80
overall CV accuracy: 0.82
overall Root Mean Squared Error: 0.54
overall Correlation: 0.61
overall Root Mean Squared Error: 0.46
overall Correlation: 0.63
Elapsed: 68.37 seconds

Test Apply Mask Function


In [9]:
tic = time.time() #Start Timer

# Load data using nibabel
neg_file_list = [dest_dir + '/resampled/00' + str(x) + '.nii.gz' for x in neg_list.image_id]
neu_file_list = [dest_dir + '/resampled/00' + str(x) + '.nii.gz' for x in neu_list.image_id]
dat = nb.funcs.concat_images(neg_file_list+neu_file_list)
pines = nb.load(os.path.abspath(dest_dir + '/resampled/00' + str(combined_df.image_id[0]) + '.nii.gz'))
Y = np.array([1]*len(neg_file_list)+[0]*len(neu_file_list))
sublist = np.array([x[:-7].split('/')[-1] for x in neg_file_list] *2)

pexpd = apply_mask(data=dat, weight_map=pines, output_dir=dest_dir, method='dot_product', save_output=True)
pexpc = apply_mask(data=dat, weight_map=pines, output_dir=dest_dir, method='correlation', save_output=True)

plt.subplot(2, 1, 1)
plt.plot(pexpd)
plt.title('Pattern Expression')
plt.ylabel('Dot Product')

plt.subplot(2, 1, 2)
plt.plot(pexpc)
plt.xlabel('Subject')
plt.ylabel('Correlation')

print 'Elapsed: %.2f seconds' % (time.time() - tic) #Stop timer


Elapsed: 2.93 seconds

Compare weight maps: Similarity


In [39]:
# Create LASSO-PCR
tic = time.time() #Start Timer

import glob
from nilearn.input_data import NiftiMasker     
from nilearn.image import mean_img
from sklearn.metrics import pairwise_distances
import seaborn as sns
from nilearn.plotting import plot_stat_map, plot_img, plot_epi
resource_dir = os.path.abspath('/Users/lukechang/Github/neurolearn/nltools/resources')

# Load data using nibabel
file_list = glob.glob(dest_dir + '/*nii.gz')
dat = nb.funcs.concat_images(file_list)
mask = nb.load(os.path.join(resource_dir,'MNI152_T1_2mm_brain_mask_dil.nii.gz'))
nifti_masker = NiftiMasker(mask_img=mask)
data = nifti_masker.fit_transform(dat)

# Plot Similarity Matrix
alg = [os.path.basename(x).split('_')[0] for x in file_list]
sim = 1 - pairwise_distances(data, metric='correlation')
sns.heatmap(sim, linewidths=0, square=True, vmin=0, vmax=1,annot=True,yticklabels=alg,xticklabels=alg, cmap='RdBu_r')

# plot Mean
mn = mean_img(dat)
overlay_img = nb.load(os.path.join(resource_dir,'MNI152_T1_2mm_brain.nii.gz'))

fig1 = plot_stat_map(mn, overlay_img, title='Mean', cut_coords=range(-40, 40, 10), display_mode='z')
# fig2 = plot_img(mn, title='Mean', cut_coords=range(-40, 40, 10), display_mode='z')
# fig2 = plot_epi(mn, title='Mean', cut_coords=range(-40, 40, 10), display_mode='z')


print 'Elapsed: %.2f seconds' % (time.time() - tic) #Stop timer


Elapsed: 1.43 seconds

In [26]:


In [ ]: