Chapter 6: evidence for zero-slope model

Here we will again vary the slope used to simulate data and compute evidence for the zero-slope model in each case.

clear plotFLAG=1;

COL=[.55*[1 1 1]; .3 .4 .85; .8 .5 .4];

N=6; thetas=linspace(-100,100,N); sig=30;

Nb=301; blist=linspace(-.25,1.25,Nb); %data slope

Ns=501; slist=linspace(0,50,Ns+1)'; slist=slist(2:end); %noise level

diffs=[diff(slist(1:2)) diff(blist(1:2))];

Nrot=22; rotlist=linspace(-pi/6,pi/4+pi/6,Nrot); %rotations in radians

%keep noise samples identical across simulations for clearer comparison

eorig=nrand(0,sig,[1 N]); eorig=eorig-mean(eorig);

dorig=[thetas(:) thetas(:)];

margL=nan(Nrot,2); sest=zeros(Nrot,2); best=zeros(Nrot,2);

for irot=1:Nrot, rotnow=rotlist(irot);

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

%%% DATA SIM using rotation ratio %%%

dnow=rot(dorig,rotnow-pi/4); Ov=ones(Nb,1);

slistmat=slist*ones(1,Nb); slistn2=slist.^(-2); p1blistmat2=pi*(1+ones(Ns,1)*blist.^2)/2;

tnow=dnow(:,1)'; dnow=dnow(:,2)'+eorig;

likeFN0=-(N/2)*slistn2*mean(dnow.^2,2) -(N+1)*log(slist); %[Nsx1] output

likeFN1=-(N/2)*slistn2*mean((Ov*dnow-blist'*tnow).^2,2)'-(N+1)*log(slistmat)-log(p1blistmat2); %[NsxNb] output

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

%%% param likelihood by condition %%%

Lall0=likeFN0;

Lall1=likeFN1;

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

%%% verify global sigma %%%

sest(irot,1)=logpeakval([Lall0 slist]);

sest(irot,2)=logpeakval([logsum(Lall1,2,diffs(2)) slist]);

best(irot,2)=logpeakval([logsum(Lall1,1,diffs(1)) blist']);

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

%%% Model likelihoods %%%

margL(irot,1)=logsum(Lall0,1,diffs(1));

margL(irot,2)=logsum(Lall1,[1 2],diffs); end

Elist=10*(margL(:,1)-margL(:,2));

figure(1); clf; hold on

plot(tan(rotlist),Elist,'ko-','MarkerFaceColor','k','MarkerSize',11);

xlabel('data slope','FontName','Arial','FontSize',15);

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

figure(30); clf; hold on; plot(tan(rotlist),margL(:,1),'ko:','MarkerFaceColor','k'); plot(tan(rotlist),margL(:,2),'ko:','MarkerFaceColor',COL(1,:))

xlabel('data slope','FontName','Arial','FontSize',15);

ylabel('mod likelihood','FontName','Arial','FontSize',15);

The evidence favoring a nonzero slope is simply the negative of these values.