Mind Reading

The goal of so-called "mind reading" techniques is to use machine learning algorithms to decode brain states (as indexed by global patterns of cortical/subcortical activation) into their associated experimental conditions. Where this technique uses supervised learning, two things are required: First, a series of input patterns must be generated from fMRI activations. Second, a target pattern must be generated for each of the input patterns. This target pattern will indicate the experimental condition to be associated with each of the input patterns (e.g., REST (0) vs TASK (1)).

Determine Inputs

Source data for input patterns are obtained using the same methods used for connectivity analyses:

1. Load surface-space time series data (loadFSTS)
2. Remove any non-linear trends in the data (NLDetrendFSTS)
   • Note which, if any, rows you choose to drop from the data (e.g., it is customary to drop the first ~8-10 seconds of volumes from fMRI time series)
3. Normalize the time series data (normalizeMatrix)
4. Binarize and scale the normalized time series data (binarizeMatrix)

Determine Targets

Classifying Task vs Baseline Blocks

If the goal is simply to distinguish task block from rest periods, use findBlockBoundaries as follows:

%As mentioned above, keep in mind any samples dropped from the inputs:
droprows=[1 2 3 4]; %this vector was used in the call to NLDetrendFSTS;

sample_rate=2.047; %an fMRI study, with a TR=2.047 seconds
b=findBlockBoundaries([], sample_rate);
targets=zeros(1,b(end)); %1 zero for each volume -- default=baseline
for block=1:size(b,1)
 targets(b(block,1):b(block,2))=1; %block volumes get a '1'
end

targets(droprows)=[]; %remove from the schedule the dropped time points
%the first 4 × 2.047 seconds have been dropped from the schedule and from the input data
%as if this period never happened

It is possible that the number of samples in the input data may exceed the number of targets. For example, the experiment script may terminate at the 6:00 mark, but the neuroimaging may continue for an additional few seconds, ending at, say the 6:04 mark. Any input samples extending beyond the timestamp associated with the final target value obtained using the above method would be undefined, though they could be construed as being associated with non-task events. To avoid ambiguity, however, it might be advisable to drop input samples that were obtained after the last target.

Classifying Task Blocks

If blocks are associated with different tasks or conditions to be classified, the schedule vector, schedule can be created similarly, but with some modification. Here's some examples.

sample_rate=2.047; %each sample spans 2.047 seconds in this fMRI example
expinfo=load('LDT_Sub_1004_Run_11_18-Apr-2016.mat');
t=cell2mat({expinfo.data.timestamp});
vols=floor(t/sample_rate)+1; %convert the timestamps into volume numbers
cond=double(cell2mat({expinfo.data.conditon})); %what condition is each trial? 
%p.s., note the typo on "conditon"
b=double(cell2mat({expinfo.data.block})); %what block is each trial? 
bnums=unique(b); %what are the different blocks?
lookup=[0,1;0,2]; %condition codes that make up each block condition
bcodes=[1,2];%code assigned to each block condition
schedule=zeros(1,vols(end));%blank schedule preallocated for each volume

%%This loop will iterate through all the numbered blocks and use the lookup
%%table to determine which block code to assign each block, depending on
%%the individual trial conditions present in that block.
%%Then, all individual volumes that belong to that block will get assigned
%%that condition code. Anything not assigned a block code will remain '0'
%%(or 'rest'/'baseline')
for i=1:length(bnums)
   idx=find(b==bnums(i));%get indices of current block
   codes=unique(cond(idx));%what conditions are represented in this block?
   blockcondition=bcodes(ismember(lookup, codes, 'rows'));%lookup the blockcode (in bcodes) that matches the conditions in this block
   firstvol=vols(idx(1));
   lastvol=vols(idx(end));
   schedule(firstvol:lastvol)=blockcondition;
end

Match Targets to Inputs

Input values will come from time series data, imported and scaled as described here. If any volumes have been discarded from the input time series, the same volumes will need to be deleted from the schedule of targets. Assuming we have a schedule and a set of inputs for a single run, we need to match each input vector to the corresponding condition in the schedule vector and produce a MikeNet example file. This will use a Matlab function, tentatively called mindReadingXFiles.

mindReadingXFiles('inputs', SCALED, 'targets', SCHEDULE, 'window', 'all', ...
'clamp', 2, 'tdelay', 2, 'prefix', 'wholething');

This function will generate one or more .ex files with the specified filename prefix. A sliding window is used to generate the examples in each file. A window size of w (4 is the default) will present w consecutive input patterns, on successive time steps (e.g., when the window is set to 4, the first example will present the first, second, third and fourth activation patterns; the second example will present the second, third, fourth and fifth activation patterns, etc.). The target for each input pattern is the corresponding row from the TARGETS vector. The clamp parameter sets the number of time steps that each input pattern needs to be clamped for, and the tdelay parameter sets the time offset from the start of an input pattern presentation to check the output activation. The tdelay should be selected to allow for inputs in a multilayer network to have an opportunity to influence the output layer.

Measuring Classification Accuracy

A testing module takes as input the continuous time-series for an entire run (i.e., with a window set to all) and generates an output value for each unit in the classification output layer for every time point in the example file. Note that the number of time points will likely exceed the number of volumes, effectively oversampling the time series. There is a mathematical relationship between the oversampling and the clamp/target parameters, but I'm not positive what it is just now.

Graphing Performance

Load the network output into a matrix, A. The first two columns contain information about the trial and the time point. The information you will want to display will be found in A(:,[2:end]).