Chapter 6: multisource bias

To understand how the evidence for a bias follows the size of the true bias, we will simulate a dataset and vary its mean from left to right of no bias, computing and plotting the evidence in each case.

%%% Biased Gaussian

%%% [models compare unbiased vs. unknown bias]

clear %#ok<*UNRCH>

N=10; S=5;

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

sig=.1;

Ns=201; slist=linspace(0,10,Ns)'; slist=slist(2:end); sn2=slist.^(-2);

err=nrand([S,N])*sig;

for s=1:S, err(s,:)=err(s,:)-mean(err(s,:)); end

brange=[-2*sig 2*sig]; %don't expect any large biases

Nb=101; blist=linspace(brange(1),brange(2),Nb);

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

Nx=31; xlist=linspace(.9*brange(1),.9*brange(2),Nx);

ML1=zeros(Nx,1); ML2=zeros(Nx,1);

for ix=1:Nx, xnow=xlist(ix);

Lall1=zeros(1); Lall2=zeros(Nb,1);

for s=1:S,

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

%%% DATA SIM %%%

dnow=err(s,:)+xnow; Ov=ones(Nb,1); ON=ones(1,N);

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

%%% param likelihood by condition %%%

likeFN0=-(N/2)*sn2*mean(( dnow ).^2,2)'-(N+1)*log(slist); %output is Ns x 1

likeFNb=-(N/2)*sn2*mean((Ov*dnow-blist'*ON).^2,2)'-(N+1)*log(slist); %output is Ns x Nb

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

%%% marginalize source-specific sig %%%

Lall1=Lall1+logsum(likeFN0,1,diffs(1));

Lall2=Lall2+logsum(likeFNb,1,diffs(1)); end

ML1(ix)=Lall1; ML2(ix)=logsum(Lall2,1,diffs(2)); end

%%% compute evidence

EV=10*(ML2-ML1);

%%% plotting

figure(4); clf; hold on

ind=findnearestN(xlist,.05,1);

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

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

xlabel(['data mean'],'FontSize',15);

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

figure(5); clf; hold on

plot(.05*[1 1],[-10 EV(ind)],'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([brange -50 400]);

xlabel(['data mean'],'FontSize',15);

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