clear

% s = 973554513;

% randn('state',s)

X = [];

th = 1;

h = .01;

x = 0;

dd = 1;std=1;

t=0;

while abs(x) < th

t=t+h;

x = x + dd.*h + std*sqrt(h).*randn(1);

X = [X ; [t x]];

end % a single trial

n = size(X,1);

Y = (n+1):(n+100);

Y = [ h*Y' zeros(100,1)];

X = [X ; Y];

N = size(X,1);

plot(X(:,1),X(:,2),'.',X(:,1),zeros(N,1),'-',...

X(:,1),th*ones(N,1),'-',X(:,1),-th*ones(N,1),'-');

xlabel('Time')

ylabel('Preference State')

axis( [0 t+100*h -3*th 3*th] )

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% generate data for random walk simulation

clear

n = 300; % trials

ns = 20; % subjects

ng=3; % groups

d = [-.75 ; .75]; % noise signal

th =1; z = .2; stdd=.5; sth=.1;

% d = [-.5 ; .5];

% th = 1; z = .2; stdd=0; sth=0; mt = 0;

type = 2;

if type ==1

datam = zeros(n,4,ns,ng);

% n x 4 :

% n trials

% 4 cond's: noise-up noise-down sig-up sig-down

% nstim-Sresp nstim-Nresp sstim-Sresp sstim-Nresp

for age = 1:ng

age

mt = (age/ng);

for sub = 1:ns

thr = th + sth*age*randn(1);

choice = simulate(n,d,thr,z,stdd,mt);

datam(:,:,sub,age) = choice;

end

end

else

load simdat;

end

% run data_anal;

run summary;

% run test_diff

% save simdat datam Ym Mdat

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function choice = simulate(n,d,th,z,stdd,mt)

% [choice, Rt] = simulate(n,d,th,z,stdd,mt,sth)

% produces a sample of n simulations from a diffusion process

h = .0001;

choice = zeros(n,2,2);

for c = 1:2 % pick a d value

for k = 1:n % run trials

k

t = mt;

x = z;

dd = d(c) + stdd*randn(1);

while abs(x) < th

t=t+h;

x = x + dd.*h + sqrt(h).*randn(1);

end % a single trial

choice(k,:,c) = t.*[(x>0) (x<0)];

% [up = resp sig; down = resp noise]

end % all trials

end % both d's

choice = [choice(:,:,1) choice(:,:,2)];

% trial, cond, d

% cond = noise-up noise-down sig-up sig-down

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%summary statistics

% Ym(age,subj, 4 conditions, 7 measures)

% 4 conditions: (noise up, noise down, sig up, sig down)

% 7 measures: (4 n's per quartile, 3 quartiles)

Z=squeeze(sum(Ym(:,:,:,1:4),4)); % sum over four n's

C1 = Z(:,:,1)./(Z(:,:,1)+Z(:,:,2)); % prob up for noise

C2 = Z(:,:,3)./(Z(:,:,3)+Z(:,:,4)); % prob up for sig

mC1 = mean(C1,2); % mean over subj for each group , noise

mC2 = mean(C2,2); % mean over subj for each group , sig

' for each group, prob up given noise, prob up given sig'

Pc = [mC1 mC2]

' mrt per group for ( up noise, down noise, up sig, down sig)'

mrt = squeeze(mean(Ym(:,:,:,6),2)) % mean over subj of medians

Nq = (n/4)*[Pc 1-Pc];

Nq = cat(3,Nq,Nq,Nq,Nq);

Mq = squeeze(mean(Ym(:,:,:,5:7),2));

Mdat = cat(3,Nq,Mq);

% (mrt(:,1:2)+mrt(:,3:4))/2 % group by resp

% m1 = mrt(:,1).*Pc(:,1) + mrt(:,2).*(1-Pc(:,1));

% m2 = mrt(:,3).*Pc(:,2) + mrt(:,4).*(1-Pc(:,2));

% (m1+m2)/2