Chapter 6: evidence for signal equality

To examine the changes in evidence for a difference in Vernier bias with change in data mean, we simulate a series of comparison datasets whose mean changes with each simulation while keeping the no-delay dataset constant at zero mean.

clear

%%% Model likelihoods and evidence

sigprime=.3; N=10; Nvern=51; murange=[-2*sigprime 2*sigprime];

e0=randn([1 N])*sigprime; e0=e0-mean(e0); mu0=0;

e1=randn([1 N])*sigprime; e1=e1-mean(e1);

Lpri0=-log(diff(murange)); Lpri1=-2*log(diff(murange)); %normalized log-priors

Nmu=301; Nsig=201;

mulist=linspace(murange(1),murange(2),Nmu)';

slist=linspace(0,1,Nsig+1); slist=slist(2:end);

diffs0=[diff(mulist(1:2)) diff(slist(1:2))];

diffs1=[diff(mulist(1:2))*[1 1] diff(slist(1:2))];

Nx=31; xlist=linspace(.9*mulist(1),.9*mulist(end),Nx);

mumat=mulist*ones(1,N); Ltmp0=nan(Nx,1); Ltmp1=nan(Nx,1); EV=nan(Nx,1);

for ix=1:Nx, Lall0=zeros(Nmu,1,Nsig); Lall1=zeros(Nmu,Nmu,Nsig);

mutrue=xlist(ix);

%%%%%%%%R%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%% DATA SIM using rotation ratio %%%

d1=ones(Nmu,1)*(e0+mu0);

d2=ones(Nmu,1)*(e1+mutrue);

%output is Nb x 1

likeFN0=@(snow) [-(N/2)*((mean((d1-mumat).^2,2)/snow^2)+(mean((d2-mumat).^2,2)/snow^2) )-(2*N+1)*log(snow)];

%output is Nb x Nb

likeFN1=@(snow) [-(N/2)*((mean((d1-mumat).^2,2)/snow^2)*ones(1,Nmu)+ones(Nmu,1)*(mean((d2-mumat).^2,2)/snow^2)')-(2*N+1)*log(snow)];

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

%%% param likelihood by condition %%%

for is=1:Nsig, snow=slist(is);

Lall0(:,1,is)=squeeze(Lall0(:,1,is))+likeFN0(snow);

Lall1(:,:,is)=squeeze(Lall1(:,:,is))+likeFN1(snow); end

%%% compute likelihoods

Ltmp0(ix)=logsum(squeeze(Lall0),[1 2],diffs0); %unif prior range cancels

Ltmp1(ix)=logsum(Lall1,[1 2 3],diffs1)-10*log10(diff(murange));

%%% compute evidence

EV(ix)=10*(Ltmp0(ix)-Ltmp1(ix)); end

%find EV zero-crossings

Ltmp0mat=[Ltmp0(1:end-1) Ltmp0(2:end)]; Ltmp1mat=[Ltmp1(1:end-1) Ltmp1(2:end)];

inds=[find(and(Ltmp0mat(:,1)<Ltmp1mat(:,1),Ltmp0mat(:,2)>Ltmp1mat(:,2))) find(and(Ltmp0mat(:,1)>Ltmp1mat(:,1),Ltmp0mat(:,2)<Ltmp1mat(:,2)))];

xcrit=[interp1(Ltmp0(inds(1)+[0:1]),xlist(inds(1)+[0:1]),mean(Ltmp1(inds(1)+[0:1])))...

interp1(Ltmp0(inds(2)+[0:1]),xlist(inds(2)+[0:1]),mean(Ltmp1(inds(2)+[0:1])))];

%%% plotting

COL={.4*[1 1 1],.8*[1 1 1]};

xlist=xlist/sigprime;

figure(4); clf; hold on

plot(xlist,Ltmp0,'ko','MarkerFaceColor',COL{1},'MarkerSize',10)

plot(xlist,Ltmp1,'ko','MarkerFaceColor',COL{2},'MarkerSize',10)

xlabel(['data mean (\sigma\prime units)'],'FontSize',15);

ylabel('model likelihoods','FontName','Arial','FontSize',15);

figure(5); clf; hold on

EVrange=[-40 20];

plot(xcrit(1)*[1 1],EVrange,'k:')

plot(xcrit(2)*[1 1],EVrange,'k:')

plot(xlist(EV>0),EV(EV>0),'ko','MarkerFaceColor',COL{1},'MarkerSize',12,'LineWidth',1);

plot(xlist(EV<0),EV(EV<0),'ko','MarkerFaceColor',COL{2},'MarkerSize',12,'LineWidth',1);

plot(xlist([1 end]),[0 0],'k--','LineWidth',2);

axis([murange/sigprime EVrange]);

xlabel(['data mean (\sigma\prime units)'],'FontSize',15);

ylabel('evidence (db)','FontName','Arial','FontSize',15);