RadialMazeTurns

RadialMazeTurns - Compute distribution of successive turn angles in radial maze.

    In order to detect stereotypic egocentric strategies that a rat could employ to
  solve a radial maze task, this function determines the distribution (over trials)
  of successive turn angles: if the rat does follow an egocentric strategy, then
  successive turn angles are expected to reproduce a stereotypical sequence.

  USAGE

    distribution = RadialMazeTurns(data,mode)

    data           list of [rat day trial arm configuration group] tuples
    mode           optional complementary analysis; possibilities include
                   'good' (include only trials with <=2 errors), 'preferred'
                   (include only trials starting with the preferred arm) and
                   'non-preferred' (include only trials starting with any
                   arm except the preferred arm) (default = 'preferred')

  NOTES

    The data is assumed to be ordered, i.e. sorted by rat, configuration, day, and trial.

    The 'preferred' arm is the most visited arm on the first visit.

0001 function turnDistribution = RadialMazeTurns(data,mode)
0002 
0003 %RadialMazeTurns - Compute distribution of successive turn angles in radial maze.
0004 %
0005 %    In order to detect stereotypic egocentric strategies that a rat could employ to
0006 %  solve a radial maze task, this function determines the distribution (over trials)
0007 %  of successive turn angles: if the rat does follow an egocentric strategy, then
0008 %  successive turn angles are expected to reproduce a stereotypical sequence.
0009 %
0010 %  USAGE
0011 %
0012 %    distribution = RadialMazeTurns(data,mode)
0013 %
0014 %    data           list of [rat day trial arm configuration group] tuples
0015 %    mode           optional complementary analysis; possibilities include
0016 %                   'good' (include only trials with <=2 errors), 'preferred'
0017 %                   (include only trials starting with the preferred arm) and
0018 %                   'non-preferred' (include only trials starting with any
0019 %                   arm except the preferred arm) (default = 'preferred')
0020 %
0021 %  NOTES
0022 %
0023 %    The data is assumed to be ordered, i.e. sorted by rat, configuration, day, and trial.
0024 %
0025 %    The 'preferred' arm is the most visited arm on the first visit.
0026 
0027 % Copyright (C) 2009-2011 by Michaël Zugaro
0028 %
0029 % This program is free software; you can redistribute it and/or modify
0030 % it under the terms of the GNU General Public License as published by
0031 % the Free Software Foundation; either version 3 of the License, or
0032 % (at your option) any later version.
0033 
0034 global RAT DAY TRIAL CONF GROUP ARM;
0035 
0036 if nargin < 2,
0037     mode = 'preferred';
0038 else
0039     mode = lower(mode);
0040 end
0041 
0042 nTrials = 3;
0043 nArms = 8;
0044 turnAngles = (-3:4)'*45;
0045 armAngles = mod(360-(2:9)'*45,360)-180; % arm 1 is at 12:00, arm 3 at 15:00, etc.
0046 
0047 % Column names for the 'data' matrix
0048 RAT = 1;
0049 DAY = 2;
0050 TRIAL = 3;
0051 CONF = 4;
0052 GROUP = 5;
0053 ARM = 6;
0054 arms = data(:,ARM);
0055 trials = data(:,TRIAL);
0056 rats = data(:,RAT);
0057 configurations = data(:,CONF);
0058 days = data(:,DAY);
0059 groups = data(:,GROUP);
0060 
0061 if nargin < 1,
0062   error('Incorrect number of parameters (type ''help <a href="matlab:help RadialMazeTurns">RadialMazeTurns</a>'' for details).');
0063 end
0064 
0065 
0066 % Turn angles are differences in successive arm angles (numbers)
0067 turns = diff(data(:,ARM));
0068 % Keep in [1,8]
0069 turns = mod(4-turns,8);
0070 turns(turns==0) = 8;
0071 turns = [0;turns];
0072 
0073 % Differences in succesive arm angles are meaningless when they correspond to different rats
0074 % (or configurations, or days, or trials): we will have to remove those cases, so we store this
0075 % information for later use
0076 change_rat = diff(data(:,RAT)) ~= 0;
0077 change_configuration = diff(data(:,CONF)) ~= 0;
0078 change_day = diff(data(:,DAY)) ~= 0;
0079 change_trial = diff(data(:,TRIAL)) ~= 0;
0080 change = change_rat|change_configuration|change_day|change_trial;
0081 change = [true;change];
0082 valid = ~change;
0083 
0084 % Trials longer than this limit will be discarded from the final analysis
0085 maxVisits = 5;
0086 
0087 % Group days 1-5, 6-10 and 11-15
0088 nBlocks = 3;
0089 blocks = Bin(days,nBlocks);
0090 
0091 % For each trial, label visits from 1 to N
0092 newTrial = [change;false];
0093 nVisits = diff([1;find(newTrial)]);
0094 visits = ones(size(change));
0095 visits(newTrial) = -nVisits+1;
0096 visits = cumsum(visits);
0097 
0098 P = [];
0099 nGroups = max(groups);
0100 for group = 1:nGroups,
0101 
0102     ratsInThisGroup = unique(rats(groups==group));
0103 
0104     nRats = length(ratsInThisGroup);
0105     for r = 1:nRats,
0106         thisRat = rats==ratsInThisGroup(r);
0107         fig = figure('position',[1285 52 1272 890]);
0108         c1 = uicontainer('position',[0 0 1 0.85]);
0109         c2 = uicontainer('position',[0 0.85 1 0.08]);
0110         c3 = uicontainer('position',[0 0.93 1 0.07]);
0111 
0112         nConfigurations = max(configurations);
0113         for configuration = 1:nConfigurations,
0114 
0115             thisConfiguration = configurations==configuration;
0116 
0117             for block = 1:nBlocks+1,  % iteration nBlocks+1 will provide more details about last block
0118 
0119                 thisBlock = blocks==min([block nBlocks]);
0120 
0121                 % Find most visited arm on first visit
0122                 nVisitsInEachArm = hist(arms(thisRat&thisConfiguration&visits==1),1:8);
0123                 [~,preferredArm] = max(nVisitsInEachArm);
0124 
0125                 % Compute distributions of turns (across trials and days) for successive visits
0126                 T = turns(valid&thisRat&thisConfiguration&thisBlock);
0127                 N = visits(valid&thisRat&thisConfiguration&thisBlock)-1; % turn number
0128                 if block == nBlocks+1,
0129                     if strcmp(mode,'preferred'),
0130                         % Discard trials starting with non-preferred arm
0131                         a = arms([valid(2:end,1);logical(1)]&thisRat&thisConfiguration&thisBlock);
0132                         i = length(N);
0133                         while i > 1,
0134                             n = N(i);
0135                             if a(i-n+1) ~= preferredArm,
0136                                 T(i-n+1:i)=[];
0137                                 N(i-n+1:i)=[];
0138                                 a(i-n+1:i)=[];
0139                             end
0140                             i = i - n;
0141                         end
0142                     elseif strcmp(mode,'non-preferred'),
0143                         % Discard trials starting with preferred arm
0144                         a = arms([valid(2:end,1);logical(1)]&thisRat&thisConfiguration&thisBlock);
0145                         i = length(N);
0146                         while i > 1,
0147                             n = N(i);
0148                             if a(i-n+1) == preferredArm,
0149                                 T(i-n+1:i)=[];
0150                                 N(i-n+1:i)=[];
0151                                 a(i-n+1:i)=[];
0152                             end
0153                             i = i - n;
0154                         end
0155                     else
0156                         % Discard long trials
0157                         i = length(N);
0158                         while i > 1,
0159                             n = N(i);
0160                             if n > maxVisits-1,
0161                                 T(i-n+1:i)=[];
0162                                 N(i-n+1:i)=[];
0163                             end
0164                             i = i - n;
0165                         end
0166                     end
0167                     if length(T) < 2, continue; end
0168                 end
0169                 turnDistribution = Accumulate([T N],1,[8 max(N)]);
0170                 nTurns = size(turnDistribution,2);
0171 
0172                 % Plot arm histograms
0173                 A = arms(thisRat&thisConfiguration&thisBlock);
0174                 V = visits(thisRat&thisConfiguration&thisBlock);
0175                 if block == nBlocks+1,
0176                     if strcmp(mode,'preferred'),
0177                         % Discard trials starting with non-preferred arm
0178                         i = length(V);
0179                         a = arms([valid(1,2:end);logical(1)]&thisRat&thisConfiguration&thisBlock);
0180                         while i > 1,
0181                             n = V(i);
0182                             if a(i-n+1) ~= preferredArm,
0183                                 V(i-n+1:i)=[];
0184                                 A(i-n+1:i)=[];
0185                                 a(i-n+1:i)=[];
0186                             end
0187                             i = i - n;
0188                         end
0189                     elseif strcmp(mode,'non-preferred'),
0190                         % Discard trials starting with preferred arm
0191                         a = arms([valid(2:end,1);logical(1)]&thisRat&thisConfiguration&thisBlock);
0192                         i = length(N);
0193                         while i > 1,
0194                             n = N(i);
0195                             if a(i-n+1) == preferredArm,
0196                                 V(i-n+1:i)=[];
0197                                 A(i-n+1:i)=[];
0198                                 a(i-n+1:i)=[];
0199                             end
0200                             i = i - n;
0201                         end
0202                     else
0203                         % Discard long trials
0204                         i = length(V);
0205                         while i > 1,
0206                             n = V(i);
0207                             if n > maxVisits,
0208                                 V(i-n+1:i)=[];
0209                                 A(i-n+1:i)=[];
0210                             end
0211                             i = i - n;
0212                         end
0213                     end
0214                 end
0215                 armDistribution = Accumulate([A V],1,[8 max(V)]);
0216                 nVisits = size(armDistribution,2);
0217                 for i = 1:nVisits,
0218                     subplot(1,nVisits*(nBlocks+1),i+(block-1)*nVisits,'parent',c3);hold on;
0219                     %  Setup axes
0220                     m = 1.5*max(armDistribution(:,i));
0221                     tt = 0:0.01:2*pi;
0222                     fill(m*cos(tt),m*sin(tt),'w');
0223                     %  Plot
0224                     h = rose(armAngles(A(V==i))*pi/180,(0:22.5:360)*pi/180);
0225                     set(h,'color','k');
0226                     xlim([-1 1]*m);
0227                     ylim([-1 1]*m);
0228                     axis square;axis off;
0229                 end
0230 
0231                 % Plot turn histograms
0232                 for i = 1:nTurns,
0233                     subplot(1,nTurns*(nBlocks+1),i+(block-1)*nTurns,'parent',c2);hold on;
0234                     % Setup axes
0235                     m = 1.5*max(turnDistribution(:,i));
0236                     tt = 0:0.01:2*pi;
0237                     fill(m*cos(tt),m*sin(tt),'w');
0238                     % Plot
0239                     h = rose(turnAngles(T(N==i))*pi/180,(0:22.5:360)*pi/180);
0240 %                      if sum(N==i) >= 1,
0241 %                          H = circ_raotest(turnAngles(T(N==i))*pi/180);
0242 %                          if H == 1, set(h,'color','r'); end
0243 %                      end
0244                     p = circ_rtest(turnAngles(T(N==i))*pi/180);
0245                     if p < 0.05, set(h,'color','r'); end
0246                     if block == nBlocks+1,
0247                         P = [P;p i group];
0248                     end
0249                     xlim([-1 1]*m);
0250                     ylim([-1 1]*m);
0251                     axis square;axis off;
0252                 end
0253 
0254 %  %                  subplot(2,4,block,'parent',c1);
0255 %  %                  colormap(jet(6));
0256 %  %                  td = [turnDistribution;turnDistribution(end,:)];td = [td td(:,end)];
0257 %  %                  pcolor(0:size(turnDistribution,2),(-3:5)*45,td);
0258 %  %                  caxis([0 6]);
0259 %  %                  set(gca,'ytick',(-3:4)*45);
0260 %  %  %                  set(gca,'xtick',[],'ytick',[]);
0261 %  %                  shading('flat');
0262 
0263                 % Plot distribution of turns
0264                 subplot(2,4,block,'parent',c1);
0265                 colormap(jet(6));
0266                 image(1:size(turnDistribution,2),(-3:4)*45,turnDistribution+1);
0267                 caxis([0 6]);set(gca,'ydir','normal','ytick',(-3:4)*45);
0268 
0269                 % Plot successive turns (across trials and days)
0270                 subplot(2,4,block,'parent',c1);
0271                 hold on;
0272                 start = [find(N==1);length(N)+1];
0273                 for i = 1:length(start)-1,
0274                     k = start(i):start(i+1)-1;
0275                     plot(N(k),turnAngles(T(k))+rand(length(k),1)*15,'w','LineWidth',2);
0276                 end
0277 
0278                 % Plot turn(i+1) = f(turn(i))
0279                 subplot(2,4,4+block,'parent',c1);
0280                 hold on;
0281                 start = [find(N==1);length(N)+1];
0282                 for i = 1:length(start)-1,
0283                     x = turnAngles(T(start(i):start(i+1)-2));
0284                     y = turnAngles(T(start(i)+1:start(i+1)-1));
0285                     plot(x+rand(length(x),1)*15,y+rand(length(x),1)*15,'k.-','LineWidth',1,'MarkerSize',12);
0286                 end
0287 
0288                 % Fix axes
0289                 subplot(2,4,4+block,'parent',c1);xlim([-155 200]);ylim([-155 200]);
0290                 set(gca,'xtick',turnAngles,'ytick',turnAngles);
0291             end
0292         end
0293         xlabel(['Group ' int2str(group) ' - Rat ' int2str(ratsInThisGroup(r)) ' [' mode ']']);
0294     end
0295 end
0296 
0297 figure;hold on;
0298 colors = 'rbkgmy';
0299 for group = 1:nGroups,
0300     p = P(P(:,3)==group,:);
0301     plot((p(:,3)-1)*10+p(:,2),p(:,1),colors(group),'marker','.','markersize',16,'linestyle','none');
0302 end
0303 xlim([0 nGroups*10]);
0304 PlotHVLines(0.05,'h',':');