% Matlab Demo Program 
% "Flipping a coin in your head without monitoring outcomes?"
% Comments on predicting free choices and a demo program."
% by Martin Lages, Stephanie C Boyle, and Katarzyna Jaworska (2013) 
% Frontiers in Psychology, 4: 925. 

% tests predictability of binary choices using 1-back and 2-back SVM: 
% linear suppor vector machine on one or two preceding responses 

% try
clear all;

% Key mappings, unified across operating systems
% You may have to edit this depending on OS and platform
KbName('UnifyKeyNames');
leftkey = KbName('LeftArrow');
rightkey = KbName('RightArrow');
escape = KbName('ESCAPE');


pa=0; pa2=0;      % prediction accuracies
ab=zeros(1,100);  % response sequence



flag2=0; tc=0; nc=0; ttc=0; tnc=0; ntc=0; nnc=0; 

bias=0;
j=0;

%fprintf('Prediction is very difficult, especially if it''s about the future - Nils Bohr  \n');
fprintf('Can you beat the Predict-a-tron by generating a random response sequence? \n');
fprintf('Try to create a random sequence of left and right arrow key presses. \n');
fprintf('You win if response patterns predict less than 60 responses in 100 trials. \n');
fprintf('Press the left and right arrow key at your own pace now: \n');


% trial loop
while j<100 && bias==0 
    
       j=j+1;  % trial index
       
    
    while flag2==0
         
        [keyIsDown, secs, keyCode] = KbCheck;	
        WaitSecs(0.05);
        %ListenChar(2);
		  
        if keyCode(leftkey) 
            ab(j)=1; nc=nc+1;  
            break;  
  
        elseif keyCode(rightkey)  
            ab(j)=3;  tc=tc+1; 
            break;  
    
        elseif keyCode(escape)
            flag2=1;
            break;
            
        else
            % wrong key but
            % sounds depend on OS 
       
       end	
       
    end
   
    while KbCheck; end     

    fprintf('.');   % indicate trial
    if mod(j,10)==0 
         fprintf('\n');
    end  
     
        
    if j>1 % make sure there was at least one response
        % distinguish four possible response sequences
        if ab(j-1)==3 && ab(j)==3
            ttc=ttc+1;  
        elseif ab(j-1)==3 && ab(j)==1
            tnc=tnc+1; 
        elseif ab(j-1)==1 && ab(j)==3
            ntc=ntc+1;
        elseif ab(j-1)==1 && ab(j)==1
            nnc=nnc+1; 
        end
    end
    
    
end  % end of 100 trial loop 


if bias==0 & j==100
    
    
    % mean right response rate
    fprintf('\n Left key: %.2f \n', nc/j); 
    fprintf(' Right key: %.2f \n', tc/j); 
 
           
    % compute repetition and alternation
    fprintf('\n Alternated Left key p(L|R): %.2f \n', (ntc)/(j-1));
    fprintf(' Alternated Right key p(R|L): %.2f \n', (tnc)/(j-1));
    fprintf(' Sum of alternated keys: %.2f \n', (tnc+ntc)/(j-1));
    
    
    fprintf('\n Repeated Left key p(L|L): %.2f \n', nnc/(j-1));
    fprintf(' Repeated Right key p(R|R): %.2f \n', ttc/(j-1));
    fprintf(' Sum of repeated keys: %.2f \n', (ttc+nnc)/(j-1)); 
        
    
    % run svm 1-back
    cc=1; 
    for i=1:(j-1)  % 1-back 
          data(cc,:) = ab(i);   
          groups(cc,:) = ab(i+1);          
          cc=cc+1; 
    end        
    
        
    % Create a 10-fold cross-validation to compute classification error.
       indices = crossvalind('Kfold',groups,10);
       cp = classperf(groups);
       
       for i = 1:10
           test = (indices == i); 
           train = ~test;
           
           svmStruct = svmtrain(data(train, :), groups(train));
           class = svmclassify(svmStruct, data(test,:));
            
           classperf(cp,class,test);
       end
       classperf(cp,class,test);
       %cp.ErrorRate
       pa=cp.CorrectRate;
 
    
    if pa>=.6 % failed 1-back SVM
       fprintf('\n You lost against the Predict-a-tron (1-back) with %.3f \n', pa);   
       
    else 
    
        % run svm 2-back    
        cc2=1;   
        for i=1:(j-2)   % 2-back
             data2(cc2,:) = ab(i:(i+1));
            groups2(cc2,:) = ab(i+2);
            cc2=cc2+1; 
        end        
    
        
        % Create a 10-fold cross-validation to compute classification error.
       indices2 = crossvalind('Kfold',groups2,10);
       cp2 = classperf(groups2);
       
       for i = 1:10
           test2 = (indices2 == i); 
           train2 = ~test2;
           
           svmStruct2 = svmtrain(data2(train2, :), groups2(train2));
           class2 = svmclassify(svmStruct2, data2(test2,:));
           
           classperf(cp2,class2,test2);
       end
       
       classperf(cp2,class2,test2);
       %cp.ErrorRate
       pa2=cp2.CorrectRate;      
       
       if pa2>=.6 % failed 2-back SVM
          fprintf('\n You lost against the Predict-a-tron (2-back) with %.3f \n', pa2);
       elseif pa2<.6
          fprintf('\n You won against the Predict-a-tron with %.3f (1-back) and %.3f (2-back)\n', pa, pa2);
       end
    
    end % end of if clause
    
    % Walf-Wolfowitz (Runs) test (non-parametric
    fprintf('\n Results of Wald-Wolfowitz (Runs) test on all data');
    [h,p,stats]=runstest(ab, mean(ab))
    
    % Lilliefors test (unspecified standard normal)
    fprintf('\n Results of Lilliefors test on all data');
    [h,p,ksstat]=lillietest(ab) 
    
   
end   % of bias/j=100 clause       


% catch
%   psychrethrow(psychlasterror);
%   clear all;
%   close all;
% end