EV creation

Creation of EV Files

The purpose of this script is to isolate specific trial types in order to generate the contrasts in Level 1 analysis. EV creation scripts are typically written using Jupiter Notebook -- this allows quick creation and manipulation of EV file output.

• Navigate to the following folder: /home/data/madlab/scripts/tools/sample_scripts

• Open the sample script create_ev.ipynb

Always make sure to create a Markdown cell to explain the nature of the present experiment. Keep it simple, but make sure to provide sufficient notes for future researchers/programmers. And always comment your code!

Cell 1: Imports various libraries/tools

import os
from os.path import join
import numpy as np
import pandas as pd
from glob import glob

Cell 2: Creation of EV Files

• Create arrays containing subject IDs and set/run information:

#subject ids
subs = ['WMAZE_001', 'WMAZE_002', 'WMAZE_004', 'WMAZE_005', 'WMAZE_006',
        'WMAZE_007', 'WMAZE_008', 'WMAZE_009', 'WMAZE_010', 'WMAZE_012',
        'WMAZE_017', 'WMAZE_018', 'WMAZE_019', 'WMAZE_020', 'WMAZE_021',
        'WMAZE_022', 'WMAZE_023', 'WMAZE_024', 'WMAZE_026', 'WMAZE_027']

#three sets
stim_sets = ['set1', 'set2', 'set3']

• Create a loop to grab subject behavioral files:

#loop to grab correct 6 run text files for each subject
for sub in subs:
    # Array containing path to behavior files
    sub_dir = '/home/data/madlab/data/mri/wmaze/scanner_behav/{0}/'.format(sub)
    # Array containing current sub's 6 behavior file runs
    dir_file = glob(join(sub_dir, '{0}_wmazebl_2015*.txt'.format(sub)))    
    # Sort current sub's txt files in order of run
    dir_file.sort() 

• Create a loop to enumerate through set/run array. In this example, run1 and run2 will be processed in parallel for each set. Please note the removal of the last 3 trials from each run.

    #loop through each of the set types 
    for i, curr_set in enumerate(stim_sets):
        # Create dataframe for text files to extract EVS
        run1 = pd.read_table(dir_file[i * 2])
        run2 = pd.read_table(dir_file[i * 2 + 1])
        
        #removal of the last 3 trials to avoid scanner artifact
        run1 = run1[:-3]
        run2 = run2[:-3]

• In the planned analysis, the contrast of interest was fixed before correct vs. incorrect conditional trials, so the following code was required to address the "before" part:

        #convert dataframes into numpy arrays 
        run1_trialtype = run1['TrialType'].values
        run1_correct = run1['Correct'].values
        run1_resp = run1['Resp'].values 
        run2_trialtype = run2['TrialType'].values
        run2_correct = run2['Correct'].values
        run2_resp = run2['Resp'].values
        
        # Shift TrialType column by 1 & insert dummy element into first index
        run1_trial_shift = run1_trialtype[1:] 
        run1_trial_shift = np.insert(run1_trial_shift, -1, -1)
        run2_trial_shift = run2_trialtype[1:]
        run2_trial_shift = np.insert(run2_trial_shift, -1, -1)

        # Shift Correct column up by 1 & insert dummy element into first index
        run1_correct_shift = run1_correct[1:]
        run1_correct_shift = np.insert(run1_correct_shift, -1, -1)
        run2_correct_shift = run2_correct[1:]
        run2_correct_shift = np.insert(run2_correct_shift, -1, -1)

• In the above code, the newly created variables containing the values of the pandas data frame column are sliced to remove the first item, after which a "dummy element" is inserted in the last index of the now-converted numpy array. Depending on your unique analysis, this may not be necessary

• The next step is to isolate the precise trials that will comprise each EV file type. You can do this in a number of ways... here you will see the use of np.where(), which grabs the indices of the targeted trials.

        r1_all_b4_B_corr = np.where(((r1_trial_shift == 'B') 
                                    & (r1_correct_shift == 1)) 
                                    & (r1_trialtype != 'BL'))
        r1_all_b4_B_incorr = np.where(((r1_trial_shift == 'B') 
                                      & (r1_correct_shift == 0)) 
                                      & (r1_trialtype != 'BL'))
        r1_all_remain = np.where((r1_trial_shift != 'B') 
                                 & (r1_resp != 'NR') 
                                 | (r1_trialtype == 'BL'))
        r1_nonresponse = np.where((r1_resp == 'NR'))
         
        r2_all_b4_B_corr = np.where(((r2_trial_shift == 'B') 
                                    & (r2_correct_shift == 1)) 
                                    & (r2_trialtype != 'BL'))
        r2_all_b4_B_incorr = np.where(((r2_trial_shift == 'B') 
                                      & (r2_correct_shift == 0))
                                      & (r2_trialtype != 'BL'))
        r2_all_remain = np.where((r2_trial_shift != 'B') 
                                 & (r2_resp != 'NR') 
                                 | (r2_trialtype == 'BL'))
        r2_nonresponse = np.where((r2_resp == 'NR'))

• The EV files are a 2-D matrix of (n,3) -- three columns of (1) onset time, (2) duration, and (3) amplitude

        r1_onsets = run1['StimOnset']
        r2_onsets = run2['StimOnset']
       
        #obtain the onset time with indices from np.where() 
        #don't forget the [0] because np.where() returns a tuple
        r1_b4_B_corr_onsets = r1_onsets.values[r1_all_b4_B_corr[0]] 
        r1_b4_B_incorr_onsets = r1_onsets.values[r1_all_b4_B_incorr[0]]
        r1_remain_onsets = r1_onsets.values[r1_all_remain[0]]
        r1_nonresponse_onsets = r1_onsets.values[r1_nonresponse[0]]
        
        r2_b4_B_corr_onsets = r2_onsets.values[r2_all_b4_B_corr[0]] 
        r2_b4_B_incorr_onsets = r2_onsets.values[r2_all_b4_B_incorr[0]]
        r2_remain_onsets = r2_onsets.values[r2_all_remain[0]]
        r2_nonresponse_onsets = r2_onsets.values[r2_nonresponse[0]]
             
        #use v-stack to create matrix containing *ALL* onsets, 
        # durations, and amplitudes in vertical columns 
        r1_mtrx = np.vstack((r1_onsets,
                             # Numpy array filled with 3's
                             np.ones(len(r1_onsets)) * 3.0,
                             # Numpy array filled with 1's
                             np.ones(len(r1_onsets)))).T

        r1_b4_B_corr_mtrx = np.vstack((r1_b4_B_corr_onsets, 
                                     np.ones(len(r1_b4_B_corr_onsets))*3.0,
                                     np.ones(len(r1_b4_B_corr_onsets)))).T  
        r1_b4_B_incorr_mtrx = np.vstack((run1_b4_B_incorr_onsets,
                                       np.ones(len(r1_b4_B_incorr_onsets))*3.0,
                                       np.ones(len(r1_b4_B_incorr_onsets)))).T
        r1_all_remain_mtrx = np.vstack((r1_all_remain_onsets,
                                      np.ones(len(r1_all_remain_onsets))*3.0,
                                      np.ones(len(r1_all_remain_onsets)))).T
        r1_nonresponse_mtrx = np.vstack((r1_nonresponse_onsets,
                                       np.ones(len(r1_nonresponse_onsets))*3.0,
                                       np.ones(len(r1_nonresponse_onsets)))).T
                     
        r2_mtrx = np.vstack((r2_onsets,
                           np.ones(len(r2_onsets))*3.0,
                           np.ones(len(r2_onsets)))).T 
        r2_b4_B_corr_mtrx = np.vstack((r2_before_B_corr_onsets, 
                                     np.ones(len(r2_b4_B_corr_onsets))*3.0,
                                     np.ones(len(r2_b4_B_corr_onsets)))).T  
        r2_b4_B_incorr_mtrx = np.vstack((r2_b4_B_incorr_onsets,
                                       np.ones(len(r2_b4_B_incorr_onsets))*3.0,
                                       np.ones(len(r2_b4_B_incorr_onsets)))).T
        r2_all_remain_mtrx = np.vstack((r2_all_remain_onsets,
                                      np.ones(len(r2_all_remain_onsets))*3.0,
                                      np.ones(len(r2_all_remain_onsets)))).T
        r2_nonresponse_mtrx = np.vstack((r2_nonresponse_onsets,
                                       np.ones(len(r2_nonresponse_onsets))*3.0,
                                       np.ones(len(r2_nonresponse_onsets)))).T

• Now you only need to create a new directory and save the EV files you have created:

        # If the output directory does not exist
        if not os.path.exists(join(sub_dir, 'MRthesis/', 'model3_drop3/', 'EVs/')):
            # Create it
            os.makedirs(join(sub_dir, 'MRthesis/', 'model3_drop3/', 'EVs/')) 

        # Run 1 matrix    
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}.txt'.format(i * 2 + 1), 
                   r1_mtrx, delimiter='\t', fmt='%.4f')

        # A before CORRECT B ONLY
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_before_B_corr.txt'.format(i * 2 + 1), 
                   r1_before_B_corr_mtrx, delimiter='\t', fmt='%.4f')

        # A before INCORRECT B ONLY
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_before_B_incorr.txt'.format(i * 2 + 1), 
                   r1_before_B_incorr_mtrx, delimiter='\t', fmt='%.4f')
        
        # All-remaining ONLY
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_all_remain.txt'.format(i * 2 + 1), 
                   r1_all_remain_mtrx, delimiter='\t', fmt='%.4f')
                   
        # Nonresponse ONLY
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_nonresponse.txt'.format(i * 2 + 1), 
                   r1_nonresponse_mtrx, delimiter='\t', fmt='%.4f')   
        
        # Run 2 files
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}.txt'.format(i * 2 + 2), 
                   r2_mtrx, delimiter='\t', fmt='%.4f')
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_before_B_corr.txt'.format(i * 2 + 2), 
                   r2_before_B_corr_mtrx, delimiter='\t', fmt='%.4f')
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_before_B_incorr.txt'.format(i * 2 + 2), 
                   r2_before_B_incorr_mtrx, delimiter='\t', fmt='%.4f')
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_all_remain.txt'.format(i * 2 + 2), 
                   r2_all_remain_mtrx, delimiter='\t', fmt='%.4f')
        np.savetxt(sub_dir + '/MRthesis/' + 'model3_drop3/' + 'EVs/' + 
                   'r{0}_nonresponse.txt'.format(i * 2 + 2), 
                   r2_nonresponse_mtrx, delimiter='\t', fmt='%.4f')

• Once the script is complete, hit "Run" in Jupiter Notebook. Make sure to save when finished!

Review the output EV files

• Go the the output directory and open one of each type of EV file you have created.
• Compare the EV files' onset times to the original behavioral files to ensure that you grabbed the intended trials. This is the most likely point of error in your analyses, so take the extra time to double-check

Final Notes

• Each script to produce EV files will need to be customized to the precise analysis you wish to conduct.
• The EV files each represent a separate regressor for the general linear model in the Level 1 script.