function [SWN,SWP,global_CC,global_rand_CC,global_latt_CC,CC,L_real,L_rand,L_latt,L,degree,deltaC,deltaL,degree_rand,degree_latt] = small_world_sporns(signf_edges_0,dim,neuron)
% small_world_analysis: small world analysis using Sporns real-based
% functions. 
%
% INPUTS:
%   signf_edges_0: significant connections between neurons
%   dim: number of neurons
%   neuron: array that contains the exact neurons we are working with
%
% OUTPUTS:
%   SWN: small-worldness coefficient
%   SWP: small-world propensity coefficient
%   global_CC: clustering coefficient (transitivity) of the real network
%   global_rand_CC: clustering coefficient (transitivity) of the random network
%   global_latt_CC: clustering coefficient (transitivity) of the lattice network
%   CC: parameter needed to compute the small-worldness coefficient
%   L_real: shortest path length of the real network
%   L_rand: shortest path length of the random network
%   L_latt: shortest path length of the lattice network
%   L: parameter needed to compute the small-worldness coefficient
%   degree: degree of connectivity of the real network
%   deltaC: parameter needed to compute the small-world propensity coefficient
%   deltaL: parameter needed to compute the small-world propensity coefficient
%   degree_rand: degree of connectivity of the random network
%   degree_latt: degree of connectivity of the lattice network

addpath('BCT\2017_01_15_BCT');

%create the connection matrix
for i = 1:size(signf_edges_0,1)
    conn_matrix(signf_edges_0(i,1),signf_edges_0(i,2)) = 1;    
end

%make the connection matrix unidirectional 
for i = 1:dim
    for j = 1:dim
        if conn_matrix(neuron(i),neuron(j))~=conn_matrix(neuron(j),neuron(i))
            conn_matrix(neuron(i),neuron(j)) = 1;
            conn_matrix(neuron(j),neuron(i)) = 1;
        end
    end
end

%unidirectional edges
uni_edges  = [];
uni_conn_matrix = triu(conn_matrix);
aux = 1;
for i = 1:size(uni_conn_matrix,1)
    for j = 1:size(uni_conn_matrix,2)
        if uni_conn_matrix(i,j)==1
            uni_edges(aux,1)=i;
            uni_edges(aux,2)=j;
            aux = aux+1;
        end
    end
end

%unidirectional global CC
global_CC = clustering_coefficient_v2(uni_edges,dim,neuron)

%DoC
[degree,mean_degree,hubs] = DoC(conn_matrix);

%statistics of the DoC
[cdf_degree,x] = ecdf(degree);
max_degree     = max(degree);
min_degree     = min(degree(degree~=0));
median_degree  = median(degree(degree~=0));

% shortest path length
DG = sparse(signf_edges_0(:,1),signf_edges_0(:, 2),1); %unweight

[short_path_edges ] = graphallshortestpaths(DG,'directed',false);
L_real              = sum(short_path_edges(short_path_edges~=Inf))/(dim*(dim-1));
median_pathlen      = median(median(short_path_edges));
max_pathlen         = max(max(short_path_edges));
min_pathlen         = min(min(short_path_edges));


%% random network

%create a random network using Sporns real-based function
[G1,~]  = randmio_und(conn_matrix,1);
und_G1  = triu(G1);
conn_comp_rand = [];

%make the random network unidirectional
for i = 1:size(und_G1,1)
    for j = 1:size(und_G1,2)
        if und_G1(i,j) == 1
            conn_comp_rand(i,j) = 1;
            conn_comp_rand(j,i) = 1;
        end
    end
end

%compute the connected components of the random network. 
%Connected component: network that has at least 3 nodes connected
comp_rand = connected_components_v1(conn_comp_rand, size(conn_comp_rand,1));
delete_comp = find(cell2mat(comp_rand(:,2))<=3);
comp_rand(delete_comp,:) = [];
dim_rand = comp_rand{1,2};
nodes_rand = comp_rand{1,3};

%shortest path length
G1_sparse   = sparse(G1);
[dist_rand] = graphallshortestpaths(G1_sparse,'directed',false);
L_rand      = (sum(dist_rand(dist_rand~=Inf)))/(dim_rand*(dim_rand-1))

%find the connected edges. 
%Create 2 columns matrix with these connections. 
%Column 1: neuron A. Column 2: neuron B.
conn_edges_rand = [];
aux = 1;
for i = 1:size(und_G1,1)
    for j = 1:size(und_G1,2)
        if und_G1(i,j)==1
            conn_edges_rand(aux,1)=i;
            conn_edges_rand(aux,2)=j;
            aux = aux+1;
        end
    end
end

%clustering coefficient
global_rand_CC = clustering_coefficient_v2(conn_edges_rand,dim_rand,nodes_rand)

%degree of connectivity
[degree_rand,mean_degree_rand,~] = DoC(G1);
%statistics of the degree of connectivity
max_degree_rand    = max(degree_rand);
min_degree_rand    = min(degree_rand(degree_rand~=0));
median_degree_rand = median(degree_rand(degree_rand~=0));

%% lattice network

%create a lattice network using Sporns real-based function
[G0,~,~,~] = latmio_und(conn_matrix,1);
und_G0 = triu(G0);
conn_comp_latt = [];

%make the lattice network unidirectional
for i = 1:size(und_G0,1)
    for j = 1:size(und_G0,2)
        if und_G0(i,j) == 1
            conn_comp_latt(i,j) = 1;
            conn_comp_latt(j,i) = 1;
        end
    end
end

%compute the connected components of the lattice network. 
%Connected component: network that has at least 3 nodes connected
comp_latt = connected_components_v1(conn_comp_latt, size(conn_comp_latt,1));
delete_comp = find(cell2mat(comp_latt(:,2))<=3);
comp_latt(delete_comp,:) = [];
dim_latt = comp_latt{1,2};
nodes_latt = comp_latt{1,3};

%shortest path length
G0_sparse   = sparse(G0);
[dist_latt] = graphallshortestpaths(G0_sparse,'directed',false);
L_latt      = (sum(dist_latt(dist_latt~=Inf)))/(dim_latt*(dim_latt-1))

%find the connected edges. 
%Create 2 columns matrix with these connections. 
%Column 1: neuron A. Column 2: neuron B.
conn_edges_latt  = [];
aux = 1;
for i = 1:size(und_G0,1)
    for j = 1:size(und_G0,2)
        if und_G0(i,j)==1
            conn_edges_latt(aux,1)=i;
            conn_edges_latt(aux,2)=j;
            aux = aux+1;
        end
    end
end

%clustering coefficient
global_latt_CC = clustering_coefficient_v2(conn_edges_latt,dim_latt,nodes_latt)

%degree of connectivity
[degree_latt,mean_degree_latt,~] = DoC(G0);
%statistics of the degree of connectivity
max_degree_latt    = max(degree_latt);
min_degree_latt    = min(degree_latt(degree_latt~=0));
median_degree_latt = median(degree_latt(degree_latt~=0));

%% Small World Ness

CC  = global_CC/global_rand_CC;
L   = L_real/L_rand;

%small-worldness coefficient
SWN = CC/L

%% Small World Propensity

deltaC = (global_latt_CC - global_CC)/(global_latt_CC-global_rand_CC);

%limit the coefficient deltaC between [0,1]
if deltaC > 1
    deltaC = 1;
elseif deltaC < 0
    deltaC = 0;
end

deltaL = (L_real - L_rand)/(L_latt - L_rand);
%limit the coefficient deltaL between [0,1]
if deltaL > 1
    deltaL = 1;
elseif deltaL < 0
    deltaL = 0;
end

%small-world propensity coefficient
SWP = 1 - sqrt((power(deltaC,2) + power(deltaL,2))/2)
end