initial commit

Author: stanleyhchan

Figure1.m

@@ -0,0 +1,54 @@
+clear all
+close all
+clc
+
+G1 = [1 0; 0 1];
+G2 = [0 1; 1 0];
+G3 = 0.5;
+
+G1 = imresize(G1,50,'nearest');
+G2 = imresize(G2,50,'nearest');
+G3 = imresize(G3,100,'nearest');
+
+Fig = figure(1);
+G1(:,1) = 0;
+G1(1,:) = 0;
+G1(:,end) = 0;
+G1(end,:) = 0;
+imshow(G1);
+
+text(22,25,'0','color',[0 0 0], 'fontsize',20);
+text(72,25,'1','color',[0.99 0.99 0.99],'backgroundcolor','k', 'fontsize',20);
+text(22,75,'1','color',[0.99 0.99 0.99],'backgroundcolor','k', 'fontsize',20);
+text(72,75,'0','color',[0 0 0], 'fontsize',20);
+
+Fig_prop = get(gcf);
+boundbox = Fig_prop.PaperPosition;
+width    = boundbox(3);
+height   = boundbox(4);
+set(gcf, 'PaperPosition', [0, 0, width, height]);
+set(gcf, 'PaperSize', [width height]);
+% saveas(Fig, 'twin_eg_w.eps','epsc');
+
+
+Fig = figure(2);
+G2(:,1) = 0;
+G2(1,:) = 0;
+G2(:,end) = 0;
+G2(end,:) = 0;
+
+imshow(G2);
+text(22,25,'1','color',[0.99 0.99 0.99],'backgroundcolor','k', 'fontsize',20);
+text(72,25,'0','color',[0 0 0], 'fontsize',20);
+text(22,75,'0','color',[0 0 0], 'fontsize',20);
+text(72,75,'1','color',[0.99 0.99 0.99],'backgroundcolor','k', 'fontsize',20);
+
+Fig_prop = get(gcf);
+boundbox = Fig_prop.PaperPosition;
+width    = boundbox(3);
+height   = boundbox(4);
+set(gcf, 'PaperPosition', [0, 0, width, height]);
+set(gcf, 'PaperSize', [width height]);
+% saveas(Fig, 'twin_eg_w1.eps','epsc');
+
+

Figure2.m

@@ -0,0 +1,79 @@
+clear all
+close all
+clc
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This routine reproduces the results of Figure 2
+%
+% Stanley Chan
+% Harvard University
+% Dec 28, 2013
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%========= Load Package ===========%
+addpath(genpath('./deconvtv_v1/'));
+
+
+
+%========= Experiment Setup =======%
+n = 500;                                    % Number of nodes
+caseid = 4;                                 % Graphon ID
+wtrue = gen_graphon(n,caseid);              % Generate a graphon
+G = construct_a_graph_from_P(wtrue,n,1);    % Sample a random graph
+pidx = randperm(n);                         % Randomly permute col/row
+G = G(pidx,pidx);
+
+
+%========= Sort and Smooth ========%
+[M0 N0] = size(G);
+n = M0;
+h = round(log(n));
+
+% Empirical Degree Sorting
+d = mean(G);
+[~, pos] = sort(d,'ascend');
+A = G(pos,pos);
+
+% Histogram
+H = imfilter(A, ones(h)/h^2, 'symmetric');
+H = H(1:h:end, 1:h:end);
+
+% TV
+opts.beta    = [1 1 0];
+opts.print   = false;
+opts.method  = 'l2';
+[M N] = size(H);
+M2   = round(M/2);
+N2   = round(N/2);
+Hpad = padarray(H, [M2 N2], 'symmetric');
+
+out  = deconvtv(Hpad,1,7.5,opts);
+what = out.f;
+wtv  = what(M2+1:end-M2, N2+1:end-N2);
+
+west = imresize(wtv,[M0,N0],'nearest');
+
+
+
+
+
+%========= Display Results ========%
+Fig = figure(1);
+imagesc(G); axis image; colormap(gray);
+export_fig Illustration_G.eps -transparent
+
+Fig = figure(2);
+imagesc(A); axis image; colormap(gray);
+export_fig Illustration_A.eps -transparent
+
+
+Fig = figure(3);
+imagesc(H); axis image;
+export_fig Illustration_H.eps -transparent
+
+Fig = figure(4);
+imagesc(wtv); axis image;
+export_fig Illustration_wtv.eps -transparent
+
+
+

Figure3.m

@@ -0,0 +1,66 @@
+clear all
+close all
+clc
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This routine reproduces the results of Figure 3
+%
+% Stanley Chan
+% Harvard University
+% Dec 28, 2013
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+%========= Load Package ===========%
+addpath(genpath('./deconvtv_v1/'));
+
+
+%========= Experiment Setup =======%
+n = 1000;                                    % Number of nodes
+caseid = 1;                                 % Graphon ID
+wtrue = gen_graphon(n,caseid);              % Generate a graphon
+G = construct_a_graph_from_P(wtrue,n,1);    % Sample a random graph
+pidx = randperm(n);                         % Randomly permute col/row
+G = G(pidx,pidx);
+
+
+
+
+
+%========= Sort and Smooth (Proposed) ====%
+t0    = tic;
+wsas  = sort_and_smooth(G);
+t_sas = toc(t0);
+
+%========= USVT (Chatterjee 2012) ========%
+t1     = tic;
+wusvt  = usvt(G);
+t_usvt = toc(t1);
+
+%========= SBA (Airoldi et al. 2013) =====%
+t2    = tic;
+h_opt = oracle_h(G,wtrue);
+wsba  = stochastic_block(G,h_opt);
+t_sba = toc(t2);
+
+
+
+
+
+%========= Display Result =====%
+Fig = figure(1);
+imagesc(wtrue); axis image; axis off;
+export_fig wtrue_1.eps -transparent
+
+Fig = figure(2);
+imagesc(wsas); axis image; axis off;
+export_fig wsas_1.eps -transparent
+
+Fig = figure(3);
+imagesc(wusvt); axis image; axis off;
+export_fig wusvt_1.eps -transparent
+
+Fig = figure(4);
+imagesc(wsba); axis image; axis off;
+export_fig wsba_1.eps -transparent
+
+

Figure4.m

@@ -0,0 +1,64 @@
+clear all
+close all
+clc
+
+addpath(genpath('./deconvtv_v1/'));
+n_set = round(logspace(log10(100),log10(2500),50));
+
+mse_sas  = zeros(10,length(n_set));
+mse_usvt = zeros(10,length(n_set));
+mse_sba  = zeros(10,length(n_set));
+
+t_sas  = zeros(10,length(n_set));
+t_usvt = zeros(10,length(n_set));
+t_sba  = zeros(10,length(n_set));
+
+
+% Main Loop
+for cidx = 1:1
+    caseid = cidx;
+    fprintf('cidx = %3g, \n', cidx);
+    
+    for nidx=1:length(n_set)
+        fprintf('nidx = %3g, \n', nidx);
+        
+        % Generate a graphon
+        n     = n_set(nidx);
+        wtrue = gen_graphon(n,caseid);
+        
+        % Sample a random graph from the graphon
+        G = construct_a_graph_from_P(wtrue,n,1);
+        pidx = randperm(n);
+        G = G(pidx,pidx);
+        
+        % Apply Sort and Smooth
+        t0 = tic;
+        wsas = sort_and_smooth(G);
+        t_sas(cidx,nidx) = toc(t0);
+        
+        % Apply USVT
+        t1 = tic;
+        wusvt = usvt(G);
+        t_usvt(cidx,nidx) = toc(t1);
+        
+        % Apply Stochastic Blockmodel Approximation
+        t2 = tic;
+        wsba = stochastic_block(G,round(log(n)));
+        t_sba(cidx,nidx) = toc(t2);
+        
+        % Record MSE
+        mse_sas(cidx,nidx)  = mean((wsas(:) -wtrue(:)).^2);
+        mse_usvt(cidx,nidx) = mean((wusvt(:)-wtrue(:)).^2);
+        mse_sba(cidx,nidx)  = mean((wsba(:) -wtrue(:)).^2);
+    end
+end
+
+Fig = figure(1);
+plot(n_set(1:2:end), mean(t_sas(:,1:2:end)),  '-x',  'LineWidth',2); hold on; 
+plot(n_set(1:2:end), mean(t_usvt(:,1:2:end)), 'r-o', 'LineWidth',2); 
+plot(n_set(1:2:end), mean(t_sba(:,1:2:end)),  'k-^',  'LineWidth',2);
+grid on;
+xlabel('n');
+ylabel('Runtime (Seconds)');
+legend('SAS (proposed)','USVT','SBA','Location','NW');
+% export_fig runtime.eps -transparent

Figure5.m

@@ -0,0 +1,62 @@
+clear all
+close all
+clc
+
+%%%%%%%%%%%%%%%%%%%%%%%%%%%
+% This program tests the Stanford Network Data
+%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+addpath(genpath('./deconvtv_v1/'));
+
+load('./data/ca-AstroPh.mat');
+G = Problem.A;
+clear Problem;
+
+% Reduce data size and symmetrize data
+% G = G(1:10:end,1:10:end);
+G = G - diag(diag(G)); % remove diagonal
+G = (G + G')/2>0;
+
+
+% tic
+% Empirical Degree Sorting and rearrangement
+d = mean(G);
+[~, pos] = sort(d,'descend');
+A = G(pos,pos);
+d = mean(A);
+
+
+% Parameter
+n = size(A,1);
+h = 20;                 % h = 20 for ca-Astroph, h = 30 for soc-Epinion
+k = round(n/h); % number of bins
+
+binwidth = h*ones(1,k);
+cum = cumsum(binwidth);                           % cumulative sum of bin widths
+binwidth(cum>n) = [];                             % remove oversized bins
+binwidth = [binwidth n-sum(binwidth)];            % compensate the last bin
+binstart = cumsum(binwidth);                      % record the starting bin ID
+
+% Histogram
+H = zeros(k,k);
+for i=1:length(binstart)-1
+    fprintf('i = %3g / %3g \n', i,k);
+    for j=1:length(binstart)-1
+        H(i,j) = (nnz(A(binstart(i):binstart(i+1), binstart(j):binstart(j+1)))>0);
+    end
+end
+
+% TV minimization
+opts.beta    = [1 1 0];
+opts.print   = true;
+opts.method  = 'l2';
+[M N] = size(H);
+M2   = 10;
+N2   = 10;
+Hpad = padarray(H, [M2 N2], 'symmetric');
+out  = deconvtv(Hpad,1,2.5,opts);
+wtv  = out.f;
+wtv  = wtv(M2+1:end-M2, N2+1:end-N2);
+
+imagesc(wtv); axis image; axis off;
+% export_fig Graphon_soc_Epinions1.eps -transparent

Readme.txt

@@ -0,0 +1,61 @@
+================================================
+Sort and Smooth (SAS) Algorithm for Consistent Graphon Estimation
+
+
+================================================
+This MATLAB package is a supplement to the paper 
+
+S. H. Chan and E. M. Airoldi, "A Consistent Histogram Estimator for Exchangeable Graph Models", in Proceedings of International Conference on Machine Learning, 2014.
+
+
+================================================
+Content:
+
+1. Construct Graphs from a Graphon
+	Method 1: [G P u] = construct_a_graph(w,n,T)
+		Input:  w - a Graphon 
+			n - number of nodes
+			T - number of observations
+		Output: G - graph (size nxnxT)
+			P - probability of each node
+			u - label indices
+
+	Method 2: G = construct_a_graph_from_P(P,n,T)
+		Input:  P - probability of each node
+			n - number of nodes
+			T - number of observations
+		Output: G - graph (size nxnxT)
+
+2. Sorting and Smoothing (sort_and_smooth.m)
+	Input:  G    - a graph
+	Output: west - estimated graphon
+	
+	Algorithm dependency: ./deconvtv_v1/
+
+
+3. Results reported in the paper
+	Figure1.m - plot a twin graphon
+	Figure2.m - display an example of the SAS algorithm
+	Figure3.m - results of SAS, USVT and SBA for graphons no. 5 and no. 10
+	Figure4.m - runtime plot
+	Figure5.m - graphon estimation of soc-Epinion1 and ca-astroph network
+	Table2.m  - mean squared error (average and standard deviation) of SAS, USVT and SBA.
+
+4. Compared Methods
+	(i)   stochastic_block.m (Stochastic Blockmodel Approximation, Airoldi et al. 2013)
+	(ii)  usvt.m             (Universal Singular Value Thresholding, Chatterjee 2012)
+
+
+	
+References
+[1] E. M. Airoldi, T. B. Costa, and S. H. Chan, "Stochastic blockmodel approximation of a graphon: Theory and consistent estimation", Advances in Neural Information Processing Systems. ArXiv: 1311.1731. 2013.
+
+[2] S. Chatterjee. Matrix estimation by universal singular value thresholding. ArXiv:1212.1247. 2012.
+
+
+
+
+================================================
+Please report bugs to Stanley Chan [email protected]
+
+Last update: January 10, 2014