Hi @FabiR
I found that you are not applying the haxby_dataset.mask_vt[0] mask to your fmri data in your sklearn approach. This masking is automatically done in nilearn approach when you provide mask_img in the Decoder class.
The lack of masking in sklearn approach is also apparent from its long run time (more voxels → more time to fit the decoder). It is still running for me: been >40 mins, nilearn one runs in a minute or 2, even when I parallelize it across the CV splits.
So, the large accuracy difference between nilearn vs. sklearn (0.6 vs 0.2) can simply be because of that.
When I applied the mask, I got similar accuracies with standardization for both approaches:
Nilearn, with standardization: 0.6636432350718064
Sklearn, with standardization: 0.6640211640211641
But without standardization – it is quite weird
Nilearn, without std.: 0.8784958427815571
Sklearn, without std.: 0.6031746031746033
So, the weird thing is why Nilearn decoder gives such high accuracy without the standardization. I do not have an answer yet. Still looking into it…
Here’s my code for your ref:
from sklearn.linear_model import LogisticRegression
from sklearn.preprocessing import StandardScaler
from tqdm import tqdm
from nilearn import datasets
from nilearn import plotting
from nilearn.image import mean_img
import numpy as np
import pandas as pd
from nilearn.image import index_img
from sklearn.model_selection import LeaveOneGroupOut
from nilearn.decoding import Decoder
from sklearn.model_selection import cross_validate
from sklearn.pipeline import Pipeline
from nilearn.maskers import NiftiMasker
haxby_dataset = datasets.fetch_haxby()
fmri_filename = haxby_dataset.func[0]
mask_filename = haxby_dataset.mask_vt[0]
###############################################################################
# Load behavioral information
behavioral = pd.read_csv(haxby_dataset.session_target[0], delimiter=" ")
conditions = behavioral["labels"]
runs = behavioral["chunks"]
condition_mask = conditions.isin(
["bottle", "cat", "chair", "face", "house", "scissors", "shoe"]
)
# Apply mask to chunks, runs and func data
fmri_niimgs = index_img(fmri_filename, condition_mask)
conditions = conditions[condition_mask]
runs = runs[condition_mask]
# print len and shape
print(len(conditions))
print(fmri_niimgs.shape[-1])
print(len(runs))
####### NILEARN DECODING ####################################
Xs_nilearn = []
cv = LeaveOneGroupOut()
for standardize_true_false in [True, False]:
if standardize_true_false:
standardize_true_false = "zscore_sample"
# explicitly avoid screening even though providing a small mask has the same effect
screening_percentile = 100
decoder_nilearn = Decoder(
estimator="logistic_l2",
mask=mask_filename,
cv=cv,
scoring="accuracy",
standardize=standardize_true_false,
n_jobs=10, # parallelizing
screening_percentile=screening_percentile,
)
decoder_nilearn.fit(fmri_niimgs, conditions, groups=runs)
print(f"Decoder results with standardize={standardize_true_false}")
print(
"Mean decoding accuracy:",
np.mean(list(decoder_nilearn.cv_scores_.values())),
)
print("Chance level:", 1.0 / len(np.unique(conditions)))
Xs_nilearn.append(
(
decoder_nilearn.masker_.fit_transform(fmri_niimgs),
standardize_true_false,
)
)
####### NILEARN MASKING WITHOUT STANDARDIZATION AND SKLEARN DECODING, STANDARD SCALING ####################################
Xs_sklearn = []
for standardize_true_false in [True, False]:
# mask data before decoding
masker = NiftiMasker(mask_img=mask_filename, standardize=False)
X = masker.fit_transform(fmri_niimgs)
y = np.array(conditions)
groups = np.array(runs)
logo = LeaveOneGroupOut()
clf = LogisticRegression(solver="liblinear", penalty="l2")
scaler = StandardScaler()
if standardize_true_false:
pipeline = Pipeline([("scaler", scaler), ("clf", clf)])
else:
pipeline = Pipeline([("clf", clf)])
decoder_sklearn = cross_validate(
pipeline,
X=X,
y=y,
groups=groups,
cv=logo,
n_jobs=8,
scoring="accuracy",
)
print(
f"Mean score in mask {mask_filename} with standardize={standardize_true_false}:",
np.mean(decoder_sklearn["test_score"]),
)
if standardize_true_false:
X_to_append = scaler.fit_transform(X)
else:
X_to_append = X.copy()
Xs_sklearn.append((X_to_append, standardize_true_false))
# check if the two standardization (or lack of them) gives similar results
for (X_nilearn, _), (X_sklearn, _) in zip(Xs_nilearn, Xs_sklearn):
print(np.allclose(X_nilearn, X_sklearn))
# returns False
# True