ROI-based decoding analysis in the Haxby dataset

In this script we reproduce the data analysis conducted by Haxby et al.[1].

Specifically, we look at decoding accuracy for different objects in three different masks:

• the full ventral stream (mask_vt),

• the house selective areas (mask_house)

• and the face selective areas (mask_face).

The masks were defined via a standard GLM-based analysis.

# We ignore some warnings that would otherwise
# be thrown when reading images
# or to tell us that some masks contain few voxels.
import warnings

warnings.filterwarnings(
    "ignore", message="The provided image has no sform in its header."
)
warnings.filterwarnings(
    "ignore", message="The decoding model will be trained only on"
)

Load and prepare the data

We fetch the data of a single subject for analysis. We also find the names of the different categories of stimuli, and identify in which run they were presented.

import pandas as pd

from nilearn import datasets

haxby_dataset = datasets.fetch_haxby()
func_filename = haxby_dataset.func[0]

labels = pd.read_csv(haxby_dataset.session_target[0], sep=" ")

stimuli = labels["labels"]

task_mask = stimuli != "rest"

categories = stimuli[task_mask].unique()
run_labels = labels["chunks"][task_mask]

[fetch_haxby] Dataset directory found: /home/runner/work/nilearn/nilearn/nilearn_data/haxby2001

We index the volumes that do NOT correspond to a rest condition.

from nilearn.image import index_img

task_data = index_img(func_filename, task_mask)

Decoding on the different masks

The classifier used here is a Support Vector Classifier (SVC). We use Decoder with a "svc_l1" estimator because it is intra subject. The mask of the region of interest is passed directly to the Decoder.

We will be doing a ‘leave one run out’ for cross validation by using sklearn.model_selection.LeaveOneGroupOut and passing the run labels at fit time.

We will use Decoder with a "dummy_classifier" to estimate a baseline.

import numpy as np
from sklearn.model_selection import LeaveOneGroupOut

from nilearn.decoding import Decoder

mask_names = ["mask_vt", "mask_face", "mask_house"]

mask_scores = {}
mask_chance_scores = {}

for mask_name in mask_names:
    print(f"\nWorking on {mask_name}")

    mask_filename = haxby_dataset[mask_name][0]

    mask_scores[mask_name] = {}
    mask_chance_scores[mask_name] = {}

    for category in categories:
        print(f"\tProcessing {category}")

        classification_target = stimuli[task_mask] == category

        decoder = Decoder(
            estimator="svc_l1",
            cv=LeaveOneGroupOut(),
            mask=mask_filename,
            scoring="roc_auc",
            screening_percentile=100,
            standardize="zscore_sample",
        )
        decoder.fit(task_data, classification_target, groups=run_labels)

        mask_scores[mask_name][category] = decoder.cv_scores_[1]

        mean = np.mean(mask_scores[mask_name][category])
        std = np.std(mask_scores[mask_name][category])
        print(f"  Scores: {mean:1.2f} +- {std:1.2f}")

        dummy_classifier = Decoder(
            estimator="dummy_classifier",
            cv=LeaveOneGroupOut(),
            mask=mask_filename,
            scoring="roc_auc",
            screening_percentile=100,
            standardize="zscore_sample",
        )
        dummy_classifier.fit(
            task_data, classification_target, groups=run_labels
        )

        mask_chance_scores[mask_name][category] = dummy_classifier.cv_scores_[
            1
        ]

Working on mask_vt
        Processing scissors
  Scores: 0.92 +- 0.05
        Processing face
  Scores: 0.98 +- 0.03
        Processing cat
  Scores: 0.96 +- 0.04
        Processing shoe
  Scores: 0.92 +- 0.07
        Processing house
  Scores: 1.00 +- 0.00
        Processing scrambledpix
  Scores: 0.99 +- 0.01
        Processing bottle
  Scores: 0.89 +- 0.08
        Processing chair
  Scores: 0.93 +- 0.04

Working on mask_face
        Processing scissors
  Scores: 0.70 +- 0.16
        Processing face
  Scores: 0.90 +- 0.06
        Processing cat
  Scores: 0.76 +- 0.12
        Processing shoe
  Scores: 0.75 +- 0.14
        Processing house
  Scores: 0.71 +- 0.15
        Processing scrambledpix
  Scores: 0.87 +- 0.09
        Processing bottle
  Scores: 0.70 +- 0.12
        Processing chair
  Scores: 0.65 +- 0.07

Working on mask_house
        Processing scissors
  Scores: 0.83 +- 0.08
        Processing face
  Scores: 0.90 +- 0.07
        Processing cat
  Scores: 0.86 +- 0.09
        Processing shoe
  Scores: 0.82 +- 0.12
        Processing house
  Scores: 1.00 +- 0.00
        Processing scrambledpix
  Scores: 0.96 +- 0.05
        Processing bottle
  Scores: 0.86 +- 0.10
        Processing chair
  Scores: 0.90 +- 0.10

We make a simple bar plot to summarize the results

import matplotlib.pyplot as plt

from nilearn.plotting import show

plt.figure(constrained_layout=True)

tick_position = np.arange(len(categories))
plt.xticks(tick_position, categories, rotation=45)

for color, mask_name in zip("rgb", mask_names, strict=False):
    score_means = [
        np.mean(mask_scores[mask_name][category]) for category in categories
    ]
    plt.bar(
        tick_position, score_means, label=mask_name, width=0.25, color=color
    )

    score_chance = [
        np.mean(mask_chance_scores[mask_name][category])
        for category in categories
    ]
    plt.bar(
        tick_position,
        score_chance,
        width=0.25,
        edgecolor="k",
        facecolor="none",
    )

    tick_position = tick_position + 0.2

plt.ylabel("Classification accuracy (AUC score)")
plt.xlabel("Visual stimuli category")
plt.ylim(0.3, 1)
plt.legend(loc="upper right")
plt.title("Category-specific classification accuracy for different masks")

show()

References

[1]

James V. Haxby, M. Ida Gobbini, Maura L. Furey, Alumit Ishai, Jennifer L. Schouten, and Pietro Pietrini. Distributed and overlapping representations of faces and objects in ventral temporal cortex. Science, 293(5539):2425–2430, 2001. doi:10.1126/science.1063736.

Estimated memory usage: 1077 MB