Steps to use the model: 1. In the file loaddata.m let variable nt equal the number of trials (currently nt=150) 2. In the file loaddata.m define the directory for the data file (currently ''c:\data_file.txt') 3. The data file should be designed as follows : serial number (1 to number of trials) choice (1,2,3, or 4) a redundent field gain loss for example: 1 2 2 100 0 2 2 2 100 0 3 2 2 100 0 4 2 2 100 0 5 2 2 100 0 6 2 2 100 0 7 2 2 100 0 8 2 2 100 0 9 2 2 100 1250 10 1 1 100 0 . . . same for subject 2 etc... 4. After the model is run (it could take several hours, depending on the number of subjects, computing power etc.) the results show a a table with: a. subject number b. model fit c. paramater a (recency) d. parameter b1 (attention to losses) e. paramater c. %%%%************************************************************************** %This function loads the data into a global variable and specifies global information about the experiment %NEEDS TO BE CHANGED WHEN DATA CHANGES function loaddata %CONSTANT global data nt ns %Experiment specific data load 'c:\data_file.txt'; data = data_file; %loading and declaration of data; type in the file your data is stored in and declare it data nt = 150; %number of trials ns = size(data,1)/nt; %number of subjects; might have a different formula, depending on the data %%%%******************************************* clear clear global %MANUAL: %All data and subject specific information is dealt with in the three files: 1) loaddata (loading of constants); %2)subjectdata(loading of subject specific data 3)subjectinfo (additional information to be read from data) %CONSTANT: NEED TO BE SPECIFIED BY USER ONCE; AFTER THAT THEY HAVE TO REMAIN CONSTANT global model np sr %Model specific: Model specifier, number of parameters global lba uba lbb1 ubb1 lbc ubc %Bounds for rescaling of raster=Bounds for parameters of the model global data nt ns %Data specific: data, columns, number of trials, number of subjects %VARIABLE DATA: global rp %current point in raster from which the search starts global rpc1 rpc2 rpc4 %raster point coordinates; just used for passing on the % global deck %Initialization of some variables needed for recursions; storage for results obtained by loope1 parmvc=[]; %initialize collection of parameter matrices from all raster points chivc=[]; %initialize collection of chi vectors from all raster points sm=[]; %initialize collection of subject matrices from all raster points maxchivc=[]; %initialize maxparmvc=[]; %initialize %-----------------------DEFINE CONSTANT GLOBAL VARIABLE--------------------- %Model specific constants model = 7; %model specification np = 3; %number of parameters %Data specific: loaddata file NEEDS CHANGE IF DATA CHANGES loaddata %Bounds for Rescaling: CONSTANT %Lower and upper bounds of parameters a, b1, b2, and c for the rescaling process lba=0; uba=1;% learning rate or recency lbb1=0; ubb1=1;% relative attention wgt to loss lbc=-5; ubc=5; % sensitivity %----------------------------------------RASTER--------------------------------- %Creation of raster p1a = [0.15 0.5 0.85]; % learning rate and recency p1b1 = [0.15 0.5 0.85]; % relative attention wgt to losses p1c = [-1 0 1]; % sensitivity sr=size(p1a,2); %size of raster; raster has sr^np points %Rescale raster to be used by minsearch for pi = 1:sr p0a(pi) = scale(uba,lba,p1a(pi)); p0b1(pi) = scale(ubb1,lbb1,p1b1(pi)); p0c(pi) = scale(ubc,lbc,p1c(pi)); end %---------------------------------------Model---------------------------------- 'Begin model with np parameters' model np %For each point in the raster the optimization process is done here if np==3 %error check for rpc1 = 1:sr for rpc2 = 1:sr for rpc4 = 1:sr rp=[p0a(rpc1) p0b1(rpc2) p0c(rpc4)]' %raster point in rescaled coordinates %MAIN PART [parmv1, chiv1]=loope1; %Result of the optimization process starting at this raster point %chiv1 is a column vector with ns entries %parmv1 is a ns x np matrix; one row per particant; one column per parameter parmvc=[parmvc parmv1]; %Collects all parameter-matrices via horizontal concatination from each of the raster points chivc=[chivc chiv1]; %Collects all chi-column-vectors via horizontal concatination from each of the raster points save gbt1mod7p4cnp chivc parmvc; %Storing; in case the program crashes end end end else 'Error: number of parameters does not fit model' end %------------------------------------Loading of best parameters-------------------------- mm=size(chivc,1); %should be equal to the number of subjects for m = 1:mm %Find best chi for subject m [c,i] = max (chivc(m,:)); %look for the best chi for subject m in the chi vectors for all raster points maxchivc = [maxchivc; c]; %add the best chi to the maxchi column vector (one entry per subject) %Find the parameters that produced the best chi j = ((i-1)*np+1):((i-1)*np+np); %vector containing the indices where to find the paramters for the best chi maxparmvc=[maxparmvc; parmvc(m,j)]; %Collects parameters with best chi %maxparmvc will have one row per subject and one column per parameter end %Rescaling the returned parameters of minsearch to the model parameters; same as in fitr1 for m = 1:mm tmaxparmvc(m,1) = scalei(uba,lba,maxparmvc(m,1)); tmaxparmvc(m,2) = scalei(ubb1,lbb1,maxparmvc(m,2)); tmaxparmvc(m,3) = scalei(ubc,lbc,maxparmvc(m,3)); end %Loads additionals information about the subjects for n=1:ns sv=subjectinfo(n); %returns subjectinfo vector (row) sm=[sm; sv]; %collects all column subject vectors into a matrix with one column per subject %sm will have one row per subject and one column per additional information about the subject end res = [sm maxchivc tmaxparmvc] %res will have one row per participant % for unconstrain save results1 sm maxchivc tmaxparmvc; save results2.txt res -ascii 'End model with np parameters' model np %%%%%%%%%********************************************************************* %Call main function for each of the "ns" subjects %No need for change, even if data or model changes function [parmv, chiv]=loope1; global nt np ns % for desperation meter global rpc1 rpc2 rpc4 sr %global deck chiv=[]; parmv=[]; %Init %looping through all subjects for i = 1:ns %Printout for status check 'Subject:' i %MAIN PART MAIN PARt [parm, chi] = main1(i); %optimization process is called for the i-th subject %Printout for status check [i rpc4+(rpc2-1)*sr+(rpc1-1)*(sr^2) sr^np ] %print out subject number and the associated chi %Collection of the data for all subjects parmv = [parmv ; parm']; %collect all parameter vectors (parmv) to a concatinated matrix chiv = [chiv; chi]; %collect all chi of all subjects to a chi column vector end %%%%%%%%%%********************************************************************* %Estimates parameters and computes chi sq for a single subject with subject number/id=i function [parm, chi] = main1(i); %CONSTANT global nt %VARIABLE global rp %Point in raster from which the search starts global x y %Subject specific data; needed for training (introduced first here) %size of deck %global deck %if i < 51 % deck = 60; % else % deck = 40; % end %Load data into global variables [x,y]=subjectdata(i); %Specifies the subject specific data %Main Part-Call optimization which uses the loaded data [parm, chi1]=fminsearch('fitr1',rp); %--------------------------Comparison with other Model---------------- mp = sum(y)'/nt; chi2=bin(mp) ; % 'Bin Chi Sq minus Model Chi Sq ' chi = chi2-chi1; %%%%*************************************************************************** % searches for best chi square for one particular participant % x y represent the data (see main1) function chi = fitr1(parm); %CONSTANT global nt np %number of trials and parameters (no need for ns) global lba uba lbb1 ubb1 lbc ubc %Needed for rescaling of raster %VARIABLE global x y %subject data; specified in main1 and need for global passing % ************************************************************************ %Rescaling: Needs to be changed only if np changes if np==3 aa=parm(1,1); bb1 = parm(2,1); cc = parm(3,1); %dummy variables %Rescaling back to model parameters (we are now in the inner part of the optimization process) a = scalei(uba,lba,aa); b1 = scalei(ubb1,lbb1,bb1); c = scalei(ubc,lbc,cc); else 'Error number of parameters different from implementation' end % ************************************************************************ chi= 0; %initialize v = [0;0;0;0]; %initialize close_gate = 0; %-----------------------------------MODEL----------------------- %begin of training for tt = 2:nt t = tt-1; yt = y(t,:)'; % current choice just made this choice ytt= y(tt,:)'; % the next choice one step ahead we wish to predict xt = x(t,:); % current payoff just experienced; one row in x, for trial t win = abs(xt(2))/100; % Eldad loss = abs(xt(3))/100; % determine loss and rescale it v = v + a * yt.* ((1-b1)*win - b1*loss - v); gam = (t/10)^c; % sensitivity s = exp(gam*v)+ 0.0000000001; % strength p = s./sum(s); % prob of choice for each deck pp = .0001 + .9998*p; chi = chi + log(pp)'*ytt; end %end of training %-----------------------------------------MODEL---------------------------------- chi = -2*chi;% -2 times log likelihood %%%%%%%%%%%%**************************************************888888888888888888888 function chi=bin(parm); % chi = bin(parm) % computes chi square for binomial model %CONSTANT global model nt %VARIABLE global x y a=parm(1,1); b = parm(2,1); c=parm(3,1); d = 1-(a+b+c); chi=0; for t = 2:nt yt = y(t,:)'; p = [a b c d]'; pp = .0001 + .9998*p; chi = chi + log(pp)'*yt; end % end of training chi = -2*chi; %%%%%%%%%%%************************************************************ function [xx,yy]=subjectdata(sid) %reads out the subject specific data of subject number "subid", which is needed for the training loop in fitr1 global data nt ns %TASK: Write functions that extract x and y from the data %xx (nt x 2 matrix) contains the choice in one column and the loss in the second column, row each trial (row) %yy is a nt x 4 matrix, that contains a 4d vector at each column which represents the choice at that trial in the form (0;1;0;0) %---------------------------------------Determine the subject specific data needed for the training--------- %------------------------------ONE ALTERNATIVE--------------------------------------- xx = data(((sid-1)*nt+1):((sid)*nt),[2 4 5]); %choice and loss columns of all 100 trials + version (1 or 2, disparity) yy = (ones(nt,1)*[1 2 3 4]) == (xx(:,1)*[1 1 1 1]); %decision matrix; one 4d row-vector for each trial; a one denotes the choice %== is an elementwise comparison that returns true or false % Eldad : Changed the size of XX matrix to 3 columns: Choice (put into yy), wins and losses %%%%%%%%%%%*************************************************************** function sv=subjectinfo(sid) %returns subject information of subject number "sid" that is not needed for the training %but for further analysis of the data global data nt %Add the data you want to read out here and build a ROW vector sv %If you want another structure than a ROW vector you need to change the last part in start.m ser_num=data(nt*(sid-1)+1,2); %version = data(nt*(sid-1)+1,3); %sv=[sid ser_num version]; sv=[sid ser_num]; %%%%%%%%%%%*********************************************************************** function back=scale(ub,lb,value) %used to scale the interval (lb,ub) to the reals back = -log((ub-lb)/(value-lb) - 1); %% **************************** function back=scalei(ub,lb,value) %used to scale the reals to the interval (lb,ub) back=(ub-lb)/(1+exp(-value))+lb;