dsSpecdiff

dsSpecDiff - Scans two DynaSim specs for differences

 Usage:
   [pop_diff,conn_diff] = dsSpecDiff(spec1,spec2)

 Inputs:
    - spec1, spec2        : DynaSim model specifications to be compared

 Outputs:
    - pop_diff            : Structure summarizing differences between populations
    - conn_diff           : Structure summarizing differences between connections
    - verbose             : Bool - show or hide text output

 Dependencies:
   Requires the nDDims class, which should be party of DynaSim. If not,
   get it here https://github.com/davestanley/nDDims

0001 function [pop_diff,conn_diff] = dsSpecDiff(spec1,spec2,verbose)
0002 %dsSpecDiff - Scans two DynaSim specs for differences
0003 %
0004 % Usage:
0005 %   [pop_diff,conn_diff] = dsSpecDiff(spec1,spec2)
0006 %
0007 % Inputs:
0008 %    - spec1, spec2        : DynaSim model specifications to be compared
0009 %
0010 % Outputs:
0011 %    - pop_diff            : Structure summarizing differences between populations
0012 %    - conn_diff           : Structure summarizing differences between connections
0013 %    - verbose             : Bool - show or hide text output
0014 %
0015 % Dependencies:
0016 %   Requires the nDDims class, which should be party of DynaSim. If not,
0017 %   get it here https://github.com/davestanley/nDDims
0018 
0019 debug_mode = 0;
0020 
0021 if nargin < 3
0022     verbose = true;
0023 end
0024 
0025 
0026 if ~exist('MDD','class')
0027     fprintf('This function requires the class MDD. Make sure this is in your MATLAB path.\n');
0028     fprintf('It can be downloaded from the following links: \n.');
0029     fprintf('1. https://github.com/davestanley/nDDims \n');
0030     fprintf('2. https://www.mathworks.com/matlabcentral/fileexchange/61656-multidimensional-dictionaries\n');
0031     error('Missing dependency');
0032 end
0033 
0034 % % % % % % % % % % % % % % % % % % % %
0035 % % % % % % DO POPULATIONS % % % % % %
0036 % % % % % % % % % % % % % % % % % % % %
0037 % Build tables summarizing all information for both models
0038 [val1,pop_name1,property1] = extract_pop_lists(spec1.populations);
0039 [val2,pop_name2,property2] = extract_pop_lists(spec2.populations);
0040 
0041 % For example, this is how the table for spec1 should look
0042 if debug_mode 
0043     % Arrange and view one of the tables
0044     table = horzcat(pop_name1,property1,val1);
0045     table
0046 end
0047 
0048 % Merge spec1 and spec2 tables together into one huge table. Create a new
0049 % column to identify spec1 vs spec2 entries
0050 val = vertcat(val1,val2);
0051 pop_names = vertcat(pop_name1,pop_name2);
0052 properties = vertcat(property1,property2);
0053 specID = vertcat(1*ones(length(val1),1), 2*ones(length(val2),1));       % Identifier saying whether entry is from spec1 or spec2.
0054 
0055 if debug_mode 
0056     % Arrange and view the huge table
0057     table = horzcat(pop_name1,property1,specID,val1);
0058     table
0059 end
0060 clear val1 val2 pop_name1 pop_name2 propert1 property2
0061 
0062 % Build specs structure.
0063 % population name x entry x spec1 or spec2
0064 nd = MDD;
0065 nd = nd.importDataTable(val,{pop_names,properties,specID});
0066 nd.axis(1).name = 'Populations';
0067 nd.axis(2).name = 'Properties';
0068 nd.axis(3).name = 'SpecID';
0069 
0070 if debug_mode
0071     nd.printAxisInfo
0072 end
0073 
0074 if verbose; run_comparison(nd); end
0075 if nargout > 0; pop_diff = return_comparison(nd); end
0076 
0077 if verbose; fprintf('\n\n'); end
0078 
0079 % % % % % % % % % % % % % % % % % % % %
0080 % % % % % % DO CONNECTIONS % % % % % %
0081 % % % % % % % % % % % % % % % % % % % %
0082 if isfield(spec1, 'connections')
0083     [val1,pop_name1,property1] = extract_conn_lists(spec1.connections);
0084 else
0085     val1 = {}; pop_name1 = {}; property1 = {};
0086 end
0087    
0088 if isfield(spec2, 'connections')
0089     [val2,pop_name2,property2] = extract_conn_lists(spec2.connections);
0090 else
0091     val2 = {}; pop_name2 = {}; property2 = {};
0092 end
0093 
0094 % Merge spec1 and spec2 tables together into one huge table. Create a new
0095 % column to identify spec1 vs spec2 entries
0096 val = vertcat(val1,val2);
0097 pop_names = vertcat(pop_name1,pop_name2);
0098 properties = vertcat(property1,property2);
0099 specID = vertcat(1*ones(length(val1),1), 2*ones(length(val2),1));       % Identifier saying whether entry is from spec1 or spec2.
0100 
0101 clear val1 val2 pop_name1 pop_name2 propert1 property2
0102 
0103 % Build specs structure.
0104 % Connection name x entry x spec1 or spec2
0105 if ~isempty(val)
0106     nd = MDD;
0107     nd = nd.importDataTable(val,{pop_names,properties,specID});
0108     nd.axis(1).name = 'Connections';
0109     nd.axis(2).name = 'Properties';
0110     nd.axis(3).name = 'SpecID';
0111     
0112     if debug_mode
0113         nd.printAxisInfo
0114     end
0115     
0116     if verbose; run_comparison(nd); end
0117     if nargout > 1; conn_diff = return_comparison(nd); end
0118     
0119 else
0120     conn_diff = struct(); 
0121 end
0122 
0123 end % naub fn
0124 
0125 
0126 %% Local Fn
0127 function [val,pop_name,properties] = extract_pop_lists(s)
0128     % This function builds a huge table describing the properties of all of
0129     % the mechanisms in each population
0130     
0131     % Initialize output variables to empty cells
0132     Nparams = arrayfun(@(x) length(x.parameters),s) / 2;
0133     Nproperties = Nparams + 3;      % Add 3 for: size, equations, mechanism_list
0134     val = cell(sum(Nproperties),1);
0135     pop_name = val;
0136     properties = val;
0137     
0138     N = length(s);
0139     
0140     k=0;
0141     for i = 1:N     % Loop over populations
0142         k=k+1; field = 'size'; pop_name{k} = s(i).name; properties{k} = ['spec_' field]; val{k} = s(i).(field);
0143         k=k+1; field = 'equations'; pop_name{k} = s(i).name; properties{k} = ['spec_' field]; val{k} = s(i).(field);
0144         k=k+1; field = 'mechanism_list'; pop_name{k} = s(i).name; properties{k} = ['spec_' field]; val{k} = s(i).(field);
0145         param_names = s(i).parameters(1:2:end);
0146         param_vals = s(i).parameters(2:2:end);
0147         for j = 1:length(param_names)
0148             k=k+1; pop_name{k} = s(i).name; properties{k} = param_names{j}; val{k} = param_vals{j};
0149         end
0150     end
0151 end
0152 
0153 
0154 function [val,pop_name,properties] = extract_conn_lists(s)
0155     % This function builds a huge table describing the properties of all of
0156     % the mechanisms in each connection
0157     
0158     % Initialize output variables to empty cells
0159     Nparams = arrayfun(@(x) length(x.parameters),s) / 2;
0160     Nproperties = Nparams + 1;      % Add 1 for: mechanism_list
0161     val = cell(sum(Nproperties),1);
0162     pop_name = val;
0163     properties = val;
0164     
0165     N = length(s);
0166     
0167     k=0;
0168     for i = 1:N     % Loop over connections
0169         k=k+1; field = 'mechanism_list'; pop_name{k} = s(i).direction; properties{k} = ['spec_' field]; val{k} = s(i).(field);
0170         param_names = s(i).parameters(1:2:end);
0171         param_vals = s(i).parameters(2:2:end);
0172         for j = 1:length(param_names)
0173             k=k+1; pop_name{k} = s(i).direction; properties{k} = param_names{j}; val{k} = param_vals{j};
0174         end
0175     end
0176 end
0177 
0178 
0179 
0180 function out = mycompare(x,y)
0181     if isempty(x) && isempty(y)
0182         out = -3;                   % Missing from both
0183     elseif isempty(x) && ~isempty(y)
0184         out = -1;                   % Missing from spec1
0185     elseif ~isempty(x) && isempty(y)
0186         out = -2;                   % Missing from spec2
0187     else
0188         % Parameter is present in both spec1 and spec2!!
0189         if isnumeric(x) && isnumeric(y)
0190             out = x == y;       % 1 codes for same, 0 codes for different
0191         elseif ischar(x) && ischar(y)
0192             out = strcmp(x,y);
0193         elseif iscell(x) && iscell(y)
0194             if length(x) == length(y)
0195                 temp = cellfun(@(x2,y2) mycompare(x2,y2),x,y);
0196                 out = any(temp == 1);
0197             else
0198                 out = -4;               % Both cells but different lengths!
0199             end
0200         else
0201             out = -5;                   % If we reach this, spec1 and spec2 both have an entry, but the variable types are different.
0202         end
0203     end
0204 
0205     out = double(out);
0206 end
0207 
0208 
0209 function fprintf_cells(string,cells2print)
0210     if ~isempty(cells2print)
0211         cells2print = cellfunu(@(x) [x ' '], cells2print);
0212         fprintf([string ' ' cells2print{:} '\n']);
0213     end
0214 
0215 end
0216 
0217 
0218 function [str_merged,ind1,ind2] = compare_string_cells(string1,string2)
0219     str_merged = vertcat(string1(:),string2(:));
0220     str_merged = unique(str_merged);
0221     
0222     Nstr = length(str_merged);
0223     
0224     ind1 = false(1,Nstr);
0225     ind2 = false(1,Nstr);
0226     for i = 1:length(str_merged)
0227         ind1(i) = any(strcmp(str_merged{i},string1));
0228         ind2(i) = any(strcmp(str_merged{i},string2));
0229     end
0230     
0231 end
0232 
0233 
0234 function run_comparison(nd)
0235     % Now, we start doing the actual comparison between the two specs. This
0236     % first section looks at the populations in both spec1 and spec2, and
0237     % identifies those that mutually present.
0238     allpops = nd.axis(1).values';
0239     Npops = length(allpops);
0240     spec1_pops = nd.axis(1).values(any(~cellfun(@isempty,nd.data(:,:,1)),2))';
0241     spec2_pops = nd.axis(1).values(any(~cellfun(@isempty,nd.data(:,:,2)),2))';
0242 %     spec1_pops = {spec1.populations.name};
0243 %     spec2_pops = {spec2.populations.name};
0244 
0245     ind1 = false(1,Npops);
0246     ind2 = false(1,Npops);
0247 
0248     for i = 1:length(allpops)
0249         ind1(i) = any(strcmp(allpops{i},spec1_pops));
0250         ind2(i) = any(strcmp(allpops{i},spec2_pops));
0251     end
0252 
0253     name = lower(nd.axis(1).name);
0254     fprintf('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n');
0255     fprintf(['~~~~~~~~~~~~~~~~~~~~~~~~Comparison of ' name ' ~~~~~~~~~~~~~~~~~~~~~~~~~\n']);
0256     fprintf('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n');
0257     fprintf_cells('Populations only in spec1:', allpops(ind1 & ~ind2));
0258     fprintf_cells('Populations only in spec2:', allpops(~ind1 & ind2));
0259     fprintf_cells('Populations in both:', allpops(ind1 & ind2));
0260     fprintf('\n');
0261 
0262     % Now that we know the populations mutually present in both models, we will
0263     % proceed to compare the parameters for these populations.
0264 
0265     inds_both = find(ind1 & ind2);        % Indices of populations present in both spec1 and spec2
0266 
0267     % First, build matrix summarizing differences between all populations in
0268     % both models
0269     params1 = nd.data(:,:,1);
0270     params2 = nd.data(:,:,2);
0271 
0272     inds = cellfunu(@(x,y) mycompare(x,y), params1,params2);
0273     inds2 = cell2mat(inds);
0274 
0275     % Loop through each population present in both models and compare
0276     % mechanisms
0277     for i = inds_both
0278         
0279         inds_curr = inds2(i,:);                     % Comparison indices for current population
0280         inds_different = inds_curr ~= 1;            % Tracks whether there are differences of some sort...
0281         inds_different(inds_curr == -3) = 0;        % If these differences are due to the mechanism being missing from both specs, ignore.
0282         
0283         if any(inds_different)          % If there are differences of some sort, show comparison
0284                                                                 %
0285             fprintf('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n');
0286             fprintf(['Comparing mechs in ' name ' ' nd.axis(1).values{i} ' \n']);
0287             
0288 
0289             % List mechanisms common to both populations, or those missing from one
0290             % but present in the other
0291             nd1 = nd.subset(i,'spec_mechanism_list',1);  % Mechanisms in spec1
0292             nd2 = nd.subset(i,'spec_mechanism_list',2);  % Mechanisms in spec2
0293             [str_merged,ind1,ind2] = compare_string_cells(nd1.data{:},nd2.data{:});
0294             fprintf_cells('Mechs only in spec1:', str_merged(ind1 & ~ind2));
0295             fprintf_cells('Mechs only in spec2:', str_merged(~ind1 & ind2));
0296             % fprintf_cells('Mechs in both:', str_merged(ind1 & ind2));
0297 
0298 
0299             % Print parameters common to both population, or those missing from one
0300             % but present in the other
0301             temp = inds_curr == -1; fprintf_cells(['Parameters present in spec1 but missing in spec2:' ], nd.axis(2).values(temp)); 
0302             temp = inds_curr == -2; fprintf_cells(['Parameters present in spec2 but missing in spec1:' ], nd.axis(2).values(temp)); 
0303             temp = inds_curr == -4; fprintf_cells(['Parameters are incomparable because they are cell arrays of different lengths:' ], nd.axis(2).values(temp)); 
0304             temp = inds_curr == -5; fprintf_cells(['Parameters are incomparable because they are different data types (e.g. one is string, one is numeric):' ], nd.axis(2).values(temp)); 
0305             %temp = inds_curr == 1; fprintf_cells(['Parameters are identical in spec1 and spec2:' ], nd.axis(2).values(temp));
0306             temp = inds_curr == 0; fprintf_cells(['Parameters differing between ' nd.axis(1).values{i} ':' ], nd.axis(2).values(temp)); 
0307 
0308             % For the mutual parameters, list ones that have different values
0309             ind_diff = find(inds_curr == 0);
0310             for j = ind_diff
0311                 val1 = nd.data{i,j,1};
0312                 val2 = nd.data{i,j,2};
0313                 if isnumeric(val1); val1 = num2str(val1); end
0314                 if isnumeric(val2); val2 = num2str(val2); end
0315                 if ischar(val1) && ischar(val2)
0316                     fprintf(['Mechanism ' nd.axis(2).values{j} ' is ' val1 ' for spec1 and ' val2 ' for spec2 \n']);
0317                 end
0318             end
0319             
0320             fprintf('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~\n');
0321         end
0322 
0323     end
0324 end
0325 
0326 function s = return_comparison(nd)
0327 
0328     sz = size(nd.data);
0329     Npop = sz(1);
0330     
0331     for i = 1:Npop
0332         s(i).name = nd.axis(1).values(i);
0333         temp = horzcat(nd.axis(2).values, nd.data(i,:,1)', nd.data(i,:,2)');
0334         inds = cellfun(@(x,y) isempty(x) && isempty(y),temp(:,2), temp(:,3));      % Parameters that are abscent from both pops
0335         temp = temp(~inds,:);                                                      % Remove these parameters
0336         s(i).parameter_names_values = temp;
0337         
0338         str = 'spec_size'; ind = strcmp(nd.axis(2).values,str); if any(ind); s(i).([str '1']) = nd.data{i,ind,1}; s(i).([str '2']) = nd.data{i,ind,2}; end
0339         str = 'spec_equations'; ind = strcmp(nd.axis(2).values,str); if any(ind); s(i).([str '1']) = nd.data{i,ind,1}; s(i).([str '2']) = nd.data{i,ind,2}; end
0340         str = 'spec_mechanism_list'; ind = strcmp(nd.axis(2).values,str); if any(ind); s(i).([str '1']) = nd.data{i,ind,1}; s(i).([str '2']) = nd.data{i,ind,2}; end
0341     end
0342     
0343 
0344 end