%Election poll function main() close all %ATTENTION before using fopen remove the first row of comments in 'USA2012_3.dat' fid = fopen('USA2012_3.dat'); C = textscan(fid,'%d%s%f%f%f%f%s'); fclose(fid); %C{2} State name %C{3} percent votes candidate A %C{4} percent votes candidate B %C{5} Sample Size %C{6} Electoral votes %Bayes prior a1=1; a2=1; a3=1; %IDM s=2; %EXPLW c=2; %number of people changes vote Np=400000; for i=1:size(C{1},1) N=C{5}(i); n1=(C{3}(i)/100)*N; n2=(C{4}(i)/100)*N; %Dirichlet Kernel %dirkern=inline('x(1)^(a(1)-1) * x(2)^(a(2)-1) * (1-x(1)-x(2))^(a(3)-1)','x','a'); %Gamble g=inline('x(:,1)-x(:,2)>0'); %Bayesian Estimator expv=mc_int([n1+a1,n2+a2,N-n1-n2+a3],g,Np); %posterior probability of winning for A in each state P(i,1)=expv; % %IDM Estimator % sv=s; % expv1=mc_int([n1+sv,n2,N-n1-n2],g,Np); % %upper probability for A % UPIDM(i,1)=expv1; % expv2=mc_int([n1,n2+sv,N-n1-n2],g,Np); % %lower probability for A % LPIDM(i,1)=expv2; %EXPLW Estimator cv=c; expv3=mc_int([n1+cv,n2-cv,N-n1-n2],g,Np); %upper posterior probability of winning for Romney in each state UPEXPLW(i,1)=expv3; expv4=mc_int([n1-cv,n2+cv,N-n1-n2],g,Np); %lower posterior probability of winning for Romney in each state LPEXPLW(i,1)=expv4; end %Lower electoral map for Romney (best case for Obama) D=[]; for i=1:51 if LPEXPLW(i)<0.5 D{i}= 'Obama'; else D{i}= 'Romney' ; end end DL=struct('state',C{2}(:),'winner',D(:)); %Upper electoral map for Romney D=[]; for i=1:51 if UPEXPLW(i)<0.5 D{i}= 'Obama'; else D{i}= 'Romney' ; end end DU=struct('state',C{2}(:),'winner',D(:)); %Posterior distribution of the Electoral Votes for i=1:size(C{1},1) %Bayesian Estimator U(i,:)=(P(i,1)*ones(1,Np)>rand(1,Np))*(C{6}(i)); % % %IDM Estimator % UL(i,:)=(LPIDM(i,1)*ones(1,Np)>rand(1,Np))*(C{6}(i)); % UU(i,:)=(UPIDM(i,1)*ones(1,Np)>rand(1,Np))*(C{6}(i)); %EXPLW Estimator ULE(i,:)=(LPEXPLW(i,1)*ones(1,Np)>rand(1,Np))*(C{6}(i)); UUE(i,:)=(UPEXPLW(i,1)*ones(1,Np)>rand(1,Np))*(C{6}(i)); end figure ResEB=sum(U,1); hist(ResEB,20); median_BayesianEst=median(ResEB) mean_BayesianEst=mean(ResEB) std_BayesianEst=std(ResEB) % % % % figure % ResEL=sum(UL,1); % hist(ResEL,100); % median(ResEL) % mean(ResEL) % std(ResEL) % % figure % ResEU=sum(UU,1); % hist(ResEU,100); % median(ResEU) % mean(ResEU) % std(ResEU) % % figure ResELL=sum(ULE,1); hist(ResELL,20); lower_median=median(ResELL) lower_mean=mean(ResELL) lower_std=std(ResELL) % figure ResEUU=sum(UUE,1); hist(ResEUU,20); h = findobj(gca,'Type','patch'); set(h,'FaceColor','red') upper_median=median(ResEUU) upper_mean=mean(ResEUU) upper_std=std(ResEUU) i=find(ResELL>269); %lower probability of winning for Romney disp('lower probability of winning for Romney') LP=length(i)/length(ResELL) i=find(ResEUU>269); %upper probability of winning for Romney disp('upper probability of winning for Romney') UP=length(i)/length(ResEUU) % csvwrite('bayesian.dat',ResEB'); % csvwrite('LPIDM.dat',ResEL'); % csvwrite('UPIDM.dat',ResEU'); % csvwrite('LPEXPLW.dat',ResELL'); % csvwrite('UPEXPLW.dat',ResEUU'); end function expv=mc_int(a,g,N) %Samples from Dirichlet distribution % We use the method from p. 582 of: % A. Gelman, J. B. Carlin, H. S. Stern and D. B. Rubin: % "Bayesian Data Analysis" (2nd ed.), Chapman & Hall, 2004. % Which is: draw y1,...,yD from independent gamma distributions with % common scale and shape parameters a1,...,aD; for each i the % sample is xi = yi/sum(yj). %num of samples D = length(a); p = gamrnd(repmat(a,N,1),1,N,D); % Samples from gamma p = p ./ repmat(sum(p,2),1,D); % Normalisation expv=sum(g(p))/N; end