%Jan 2017, Etemadi and Lu
% The RB of C
% Data format:
% % #run,# sampled edges,CW,W,CW in g,Psi_g,Omega_g,CW checked in G
%Example
%2,1021,6,33,0,6,13,1

clear all;clf('reset');

data_name1={'Wiki-Italian'};
data_name={'Ego-facebook';'CA-GrQc';'Wiki-Vote';'AstroPh';'CA-CondMat';'HepPh';'Enron-email';'Brightkite';'Facebook';'Epinions';'Slashdot-Zoo';'Prosper';'Livemocha';'Douban';'Gowalla';'Libimseti';'Digg';'Web-Stanford';'Dblp-Coau';...
    'Amazon';'Actor';'Web-NotreDame';'Citeseer';'Dogster';'Catster';'Web-Berkeley';'Web-Google';'Youtube';'Dblp';'Hyves';	'Wiki-Polish';'Trec-wt10g';'Wiki-Portuguese';'Wiki-Japanese';'Pokec';'As-skitter';...
    'Wiki-Italian';'Wiki-En';'Hudong';'Hollywood';'Baidu' ;'FLicker';'Flixster';'Wiki-Russian';'Wiki-Franch';'Orkut';'Wiki-German';'USpatent';'LiveJournal';...
    'Orkut2';'DBpedia';'Web-Arabic';'Gsh-2015';'Twitter';'MicrosoftAG';'Friendster'};
size_subplot=size(data_name,1);
for i=1: size_subplot
    str1=data_name(i);
    path=char(strcat('Results_EG_Bias_p/',str1,'-raw-data-plot-REG.txt'));
    properties=dlmread(path,',',[5 0 5 9]);
    C=properties(1); N=properties(2); M=properties(3); triangles=properties(4); K=properties(5); CW=properties(6); W=properties(7);KW=properties(9);KWT=properties(10);
    Omega=2*KWT+6*triangles;
    Psi=KW;
    CW=3*triangles;
    
    iteration=dlmread(path,',',[6 0 6 0]);
    data=dlmread(path,',',9);
    point_number=1.0*size(data,1)/(iteration+1);
    
    for j=1:point_number
        % Compute ob_RB
        start=(j-1)*(iteration+1)+1;
        p(j)=data(start,1);
        p1=p(j);
        
        start=start+1;
        ob_W=data(start:start+iteration-1,4);
        ob_Lambda_g=mean(ob_W);
        ob_CW=data(start:start+iteration-1,3);
        ob_Delta_g=mean(ob_CW);
        cw=mean(ob_CW);
        obs_c=(ob_CW)./(ob_W);
        ob_RB(j)=(mean(obs_c)-C)/C;
        
        ob_Psi1=data(start:start+iteration-1,7);
        ob_Psi_g=mean(ob_Psi1);
        
        
        ob_Omega1=data(start:start+iteration-1,8);
        ob_Omega_g=mean(ob_Omega1);
        
        
        % estimate RB using Delta_g, Lambda_g, Omega_g, Psi_g in sample graph
        RB_aprox_g(j)=2*ob_Psi_g/(ob_Lambda_g.^2)-ob_Omega_g/(ob_Delta_g*ob_Lambda_g);
        
        
        % estimate RB
        RB_aprox_general(j)=(2.0/p1)*(Psi/(W^2))-(1.0/p1)*(Omega/(CW*W));
        
        obs_c_correction=(ob_CW)./(ob_W.*(1+RB_aprox_g(j)));
        ob_RB_correction(j)=(mean(obs_c_correction)-C)/C;
        
    end
    
    h1=subplot(8,7,i);
    start=1;
    end1=point_number;
    
    line1=plot(p(start:end1),ob_RB(start:end1),'r-o','MarkerSize',5,'LineWidth',1.0);
    hold on
    line2=plot(p(start:end1),RB_aprox_g(start:end1),'k','LineWidth',1.0);
    %hold on
    %line3=plot(p(start:end1),ob_RB_correction(start:end1),'k-^','LineWidth',1.0);
    %h=refline([0 0]);
    max_point=max(max(ob_RB),max(RB_aprox_g));
    min_point=min(min(ob_RB),min(RB_aprox_g));
    if min(ob_RB)>0
        ylim([min_point-min_point/10 max_point+max_point/21]);
    elseif i==2
        ylim([min_point+min_point/50 max_point-max_point/30]);
    else
        ylim([min_point+min_point/50 max_point+max_point/10]);
    end
    xlim([p(1)-p(1)/12 p(floor(point_number))+p(floor(point_number))/60]);
    box off;
    p = get(h1, 'position');
    if mod (i,7)==1
        p(1)=p(1)-0.05;
    elseif mod (i,7)==2
        p(1)=p(1)-0.09;
    elseif mod (i,7)==3
        p(1)=p(1)-0.135;
    elseif mod (i,7)==4
        p(1)=p(1)-0.18;
    elseif mod (i,7)==5
        p(1)=p(1)-0.225;
    elseif mod (i,7)==6
        p(1)=p(1)-0.27;
    elseif mod (i,7)==0
        p(1)=p(1)-0.315;
        %      elseif mod (i,8)==0
        %         p(1)=p(1)-0.19;
    end
    if i<=7
        p(2)=p(2)+0.042;
        
    elseif i<=14
        p(2)=p(2)+0.027;
        
    elseif i<=21
        p(2)=p(2)+0.018;
        
    elseif i<=28
        p(2)=p(2)+0.006;
        
    elseif i<=35
        p(2)=p(2)-0.009;
    elseif i<=42
        p(2)=p(2)-0.018;
    elseif i<=49
        p(2)=p(2)-0.029;
    else
        p(2)=p(2)-0.05;
    end
    p(3)=p(3)-0.03;  %+0.036;
    p(4)=p(4)+0.005;
    set(h1, 'position', p);
    %h=refline([0 0]);
    j = findobj(gca,'Type','text');
    set(j,'Position', [0 0 500 500]);
    %name=strcat('(',num2str(i),')',{' '},str1);
    t=title(str1,'FontSize',8,'Interpreter','Latex');
    p1=get(t,'Position');
    set(t,'Position',[p1(1) p1(2)  p1(3)]);
    
    if i==1||i==8 || i==15 || i==22 ||i==29||i==36||i==43||i==50
        ylabel('RB','FontSize',8,'Interpreter','Latex');
    end
    if i>=50
        xlabel('$p$','Interpreter','Latex','FontSize',8);
    end
    set(gca,'FontSize',3,'XTickMode','auto');
    set(gca,'XTickMode','manual','YTickMode','manual');
    set(gca,'FontSize',6);
    NumTicks = 4;
    L = get(gca,'XLim');
    set(gca,'XTick',linspace(L(1),L(2),NumTicks));
    n=get(gca,'Xtick');
    b=size(n,2);
    z=n(1,1:b-1);
    format shortG;
    set(gca,'Xtick',[z'],'Xticklabel',num2str(get(gca,'XTick')','%1.0e'));
    
end
%'Est. RB based on $\frac{var(\Lambda_g)}{\Lambda^2p^2}-\frac{cov(\Delta_g,\Lambda_g)}{\Delta\Lambda p^2}$',...
hL = legend([line1 line2 ],{'Obs. ','Est. ','Obs. corrected $\widehat{\mathcal{C}}^+$'},'FontSize',11, 'LineWidth', 4,'interpreter','latex','Orientation','horizontal');
newPosition = [0.32 0 0.005 0.04];
newUnits = 'normalized';
set(hL,'Position', newPosition,'Units', newUnits);
legend('boxoff');
set(gcf,'color','w');
export_fig RB_G_LargeBias_p56 -painters -eps

