510 lines
26 KiB
Matlab
510 lines
26 KiB
Matlab
function [ide_dynamic, ide_reducedform, ide_moments, ide_spectrum, ide_minimal] = identification_checks_via_subsets(ide_dynamic, ide_reducedform, ide_moments, ide_spectrum, ide_minimal, totparam_nbr, modparam_nbr, options_ident,error_indicator)
|
|
%[ide_dynamic, ide_reducedform, ide_moments, ide_spectrum, ide_minimal] = identification_checks_via_subsets(ide_dynamic, ide_reducedform, ide_moments, ide_spectrum, ide_minimal, totparam_nbr, modparam_nbr, options_ident,error_indicator)
|
|
% -------------------------------------------------------------------------
|
|
% Finds problematic sets of paramters via checking the necessary rank condition
|
|
% of the Jacobians for all possible combinations of parameters. The rank is
|
|
% computed via an inbuild function based on the SVD, similar to matlab's
|
|
% rank. The idea is that once we have the Jacobian for all parameters, we
|
|
% can easily set up Jacobians containing all combinations of parameters by
|
|
% picking the relevant columns/elements of the full Jacobian. Then the rank
|
|
% of these smaller Jacobians indicates whether this paramter combination is
|
|
% identified or not. To speed up computations:
|
|
% (1) single parameters are removed from possible higher-order sets,
|
|
% (2) for parameters that are collinear, i.e. rank failure for 2 element sets,
|
|
% we replace the second parameter by the first one, and then compute
|
|
% higher-order combinations [uncommented]
|
|
% (3) all lower-order problematic sets are removed from higher-order sets
|
|
% by an inbuild function
|
|
% (4) we could replace nchoosek by a mex version, e.g. VChooseK
|
|
% (https://de.mathworks.com/matlabcentral/fileexchange/26190-vchoosek) as
|
|
% nchoosek could be the bottleneck in terms of speed (and memory for models
|
|
% with totparam_nbr > 150)
|
|
% =========================================================================
|
|
% INPUTS
|
|
% ide_reducedform: [structure] Containing results from identification
|
|
% analysis based on the reduced-form solution (Ratto
|
|
% and Iskrev, 2011). If either options_ident.no_identification_reducedform
|
|
% or error_indicator.identification_reducedform are 1 then the search for problematic parameter sets will be skipped
|
|
% ide_moments: [structure] Containing results from identification
|
|
% analysis based on moments (Iskrev, 2010). If either
|
|
% options_ident.no_identification_moments or error_indicator.identification_moments are 1 then the search for
|
|
% problematic parameter sets will be skipped
|
|
% ide_spectrum: [structure] Containing results from identification
|
|
% analysis based on the spectrum (Qu and Tkachenko, 2012).
|
|
% If either options_ident.no_identification_spectrum or error_indicator.identification_spectrum are 1 then the search for
|
|
% problematic parameter sets will be skipped
|
|
% ide_minimal: [structure] Containing results from identification
|
|
% analysis based on the minimal state space system
|
|
% (Komunjer and Ng, 2011). If either options_ident.no_identification_minimal or
|
|
% error_indicator.identification_minimal are 1 then the search for problematic parameter sets will be skipped
|
|
% totparam_nbr: [integer] number of estimated stderr, corr and model parameters
|
|
% numzerotolrank: [double] tolerance level for rank compuations
|
|
% error_indicator [structure] indicators whether objects could be computed
|
|
% -------------------------------------------------------------------------
|
|
% OUTPUTS
|
|
% ide_reducedform, ide_moments, ide_spectrum, ide_minimal are augmented by the
|
|
% following fields:
|
|
% * problpars: [1 by totparam_nbr] cell with the following structure for j=1:totparam_nbr
|
|
% problpars{j}: [nonidentified_j_set_parameters_nbr by j]
|
|
% matrix with j collinear parameters in each row
|
|
% * rank: [integer] rank of Jacobian
|
|
% -------------------------------------------------------------------------
|
|
% This function is called by
|
|
% * identification_analysis.m
|
|
% =========================================================================
|
|
% Copyright (C) 2019-2021 Dynare Team
|
|
%
|
|
% This file is part of Dynare.
|
|
%
|
|
% Dynare is free software: you can redistribute it and/or modify
|
|
% it under the terms of the GNU General Public License as published by
|
|
% the Free Software Foundation, either version 3 of the License, or
|
|
% (at your option) any later version.
|
|
%
|
|
% Dynare is distributed in the hope that it will be useful,
|
|
% but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
% GNU General Public License for more details.
|
|
%
|
|
% You should have received a copy of the GNU General Public License
|
|
% along with Dynare. If not, see <https://www.gnu.org/licenses/>.
|
|
% =========================================================================
|
|
|
|
%% initialize output objects and get options
|
|
no_identification_dynamic = 0; %always compute dynamic
|
|
no_identification_reducedform = options_ident.no_identification_reducedform;
|
|
no_identification_moments = options_ident.no_identification_moments;
|
|
no_identification_spectrum = options_ident.no_identification_spectrum;
|
|
no_identification_minimal = options_ident.no_identification_minimal;
|
|
tol_rank = options_ident.tol_rank;
|
|
max_dim_subsets_groups = options_ident.max_dim_subsets_groups;
|
|
dynamic_problpars = cell(1,max_dim_subsets_groups);
|
|
reducedform_problpars = cell(1,max_dim_subsets_groups);
|
|
moments_problpars = cell(1,max_dim_subsets_groups);
|
|
spectrum_problpars = cell(1,max_dim_subsets_groups);
|
|
minimal_problpars = cell(1,max_dim_subsets_groups);
|
|
|
|
indtotparam = 1:totparam_nbr; %initialize index of parameters
|
|
|
|
%% Prepare Jacobians and check rank
|
|
% initialize linear rational expectations model
|
|
if ~no_identification_dynamic
|
|
dDYNAMIC = ide_dynamic.dDYNAMIC;
|
|
dDYNAMIC(ide_dynamic.ind_dDYNAMIC,:) = dDYNAMIC(ide_dynamic.ind_dDYNAMIC,:)./ide_dynamic.norm_dDYNAMIC; %normalize
|
|
if totparam_nbr > modparam_nbr
|
|
dDYNAMIC = [zeros(size(ide_dynamic.dDYNAMIC,1),totparam_nbr-modparam_nbr) dDYNAMIC]; %add derivatives wrt stderr and corr parameters
|
|
end
|
|
if strcmp(tol_rank,'robust')
|
|
rank_dDYNAMIC = rank(full(dDYNAMIC)); %compute rank with imposed tolerance level
|
|
else
|
|
rank_dDYNAMIC = rank(full(dDYNAMIC),tol_rank); %compute rank with imposed tolerance level
|
|
end
|
|
ide_dynamic.rank = rank_dDYNAMIC;
|
|
% check rank criteria for full Jacobian
|
|
if rank_dDYNAMIC == totparam_nbr
|
|
% all parameters are identifiable
|
|
no_identification_dynamic = 1; %skip in the following
|
|
indparam_dDYNAMIC = [];
|
|
else
|
|
% there is lack of identification
|
|
indparam_dDYNAMIC = indtotparam; %initialize for nchoosek
|
|
end
|
|
else
|
|
indparam_dDYNAMIC = []; %empty for nchoosek
|
|
end
|
|
|
|
% initialize for reduced form solution criteria
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
dREDUCEDFORM = ide_reducedform.dREDUCEDFORM;
|
|
dREDUCEDFORM(ide_reducedform.ind_dREDUCEDFORM,:) = dREDUCEDFORM(ide_reducedform.ind_dREDUCEDFORM,:)./ide_reducedform.norm_dREDUCEDFORM; %normalize
|
|
if strcmp(tol_rank,'robust')
|
|
rank_dREDUCEDFORM = rank(full(dREDUCEDFORM)); %compute rank with imposed tolerance level
|
|
else
|
|
rank_dREDUCEDFORM = rank(full(dREDUCEDFORM),tol_rank); %compute rank with imposed tolerance level
|
|
end
|
|
ide_reducedform.rank = rank_dREDUCEDFORM;
|
|
% check rank criteria for full Jacobian
|
|
if rank_dREDUCEDFORM == totparam_nbr
|
|
% all parameters are identifiable
|
|
no_identification_reducedform = 1; %skip in the following
|
|
indparam_dREDUCEDFORM = [];
|
|
else
|
|
% there is lack of identification
|
|
indparam_dREDUCEDFORM = indtotparam; %initialize for nchoosek
|
|
end
|
|
else
|
|
indparam_dREDUCEDFORM = []; %empty for nchoosek
|
|
end
|
|
|
|
% initialize for moments criteria
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
dMOMENTS = ide_moments.dMOMENTS;
|
|
dMOMENTS(ide_moments.ind_dMOMENTS,:) = dMOMENTS(ide_moments.ind_dMOMENTS,:)./ide_moments.norm_dMOMENTS; %normalize
|
|
if strcmp(tol_rank,'robust')
|
|
rank_dMOMENTS = rank(full(dMOMENTS)); %compute rank with imposed tolerance level
|
|
else
|
|
rank_dMOMENTS = rank(full(dMOMENTS),tol_rank); %compute rank with imposed tolerance level
|
|
end
|
|
ide_moments.rank = rank_dMOMENTS;
|
|
% check rank criteria for full Jacobian
|
|
if rank_dMOMENTS == totparam_nbr
|
|
% all parameters are identifiable
|
|
no_identification_moments = 1; %skip in the following
|
|
indparam_dMOMENTS = [];
|
|
else
|
|
% there is lack of identification
|
|
indparam_dMOMENTS = indtotparam; %initialize for nchoosek
|
|
end
|
|
else
|
|
indparam_dMOMENTS = []; %empty for nchoosek
|
|
end
|
|
|
|
% initialize for spectrum criteria
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
dSPECTRUM = ide_spectrum.tilda_dSPECTRUM; %tilda dSPECTRUM is normalized dSPECTRUM matrix in identification_analysis.m
|
|
%alternative normalization
|
|
%dSPECTRUM = ide_spectrum.dSPECTRUM;
|
|
%dSPECTRUM(ide_spectrum.ind_dSPECTRUM,:) = dSPECTRUM(ide_spectrum.ind_dSPECTRUM,:)./ide_spectrum.norm_dSPECTRUM; %normalize
|
|
if strcmp(tol_rank,'robust')
|
|
rank_dSPECTRUM = rank(full(dSPECTRUM)); %compute rank with imposed tolerance level
|
|
else
|
|
rank_dSPECTRUM = rank(full(dSPECTRUM),tol_rank); %compute rank with imposed tolerance level
|
|
end
|
|
ide_spectrum.rank = rank_dSPECTRUM;
|
|
% check rank criteria for full Jacobian
|
|
if rank_dSPECTRUM == totparam_nbr
|
|
% all parameters are identifiable
|
|
no_identification_spectrum = 1; %skip in the following
|
|
indparam_dSPECTRUM = [];
|
|
else
|
|
% lack of identification
|
|
indparam_dSPECTRUM = indtotparam; %initialize for nchoosek
|
|
end
|
|
else
|
|
indparam_dSPECTRUM = []; %empty for nchoosek
|
|
end
|
|
|
|
% initialize for minimal system criteria
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
dMINIMAL = ide_minimal.dMINIMAL;
|
|
dMINIMAL(ide_minimal.ind_dMINIMAL,:) = dMINIMAL(ide_minimal.ind_dMINIMAL,:)./ide_minimal.norm_dMINIMAL; %normalize
|
|
dMINIMAL_par = dMINIMAL(:,1:totparam_nbr); %part of dMINIMAL that is dependent on parameters
|
|
dMINIMAL_rest = dMINIMAL(:,(totparam_nbr+1):end); %part of dMINIMAL that is independent of parameters
|
|
if strcmp(tol_rank,'robust')
|
|
rank_dMINIMAL = rank(full(dMINIMAL)); %compute rank via SVD see function below
|
|
else
|
|
rank_dMINIMAL = rank(full(dMINIMAL),tol_rank); %compute rank via SVD see function below
|
|
end
|
|
ide_minimal.rank = rank_dMINIMAL;
|
|
dMINIMAL_fixed_rank_nbr = size(dMINIMAL_rest,2);
|
|
% check rank criteria for full Jacobian
|
|
if rank_dMINIMAL == totparam_nbr + dMINIMAL_fixed_rank_nbr
|
|
% all parameters are identifiable
|
|
no_identification_minimal = 1; %skip in the following
|
|
indparam_dMINIMAL = [];
|
|
else
|
|
% lack of identification
|
|
indparam_dMINIMAL = indtotparam; %initialize for nchoosek
|
|
end
|
|
else
|
|
indparam_dMINIMAL = []; %empty for nchoosek
|
|
end
|
|
|
|
%% Check single parameters
|
|
for j=1:totparam_nbr
|
|
if ~no_identification_dynamic
|
|
% Columns correspond to single parameters, i.e. full rank would be equal to 1
|
|
if strcmp(tol_rank,'robust')
|
|
if rank(full(dDYNAMIC(:,j))) == 0
|
|
dynamic_problpars{1} = [dynamic_problpars{1};j];
|
|
end
|
|
else
|
|
if rank(full(dDYNAMIC(:,j)),tol_rank) == 0
|
|
dynamic_problpars{1} = [dynamic_problpars{1};j];
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
% Columns correspond to single parameters, i.e. full rank would be equal to 1
|
|
if strcmp(tol_rank,'robust')
|
|
if rank(full(dREDUCEDFORM(:,j))) == 0
|
|
reducedform_problpars{1} = [reducedform_problpars{1};j];
|
|
end
|
|
else
|
|
if rank(full(dREDUCEDFORM(:,j)),tol_rank) == 0
|
|
reducedform_problpars{1} = [reducedform_problpars{1};j];
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
% Columns correspond to single parameters, i.e. full rank would be equal to 1
|
|
if strcmp(tol_rank,'robust')
|
|
if rank(full(dMOMENTS(:,j))) == 0
|
|
moments_problpars{1} = [moments_problpars{1};j];
|
|
end
|
|
else
|
|
if rank(full(dMOMENTS(:,j)),tol_rank) == 0
|
|
moments_problpars{1} = [moments_problpars{1};j];
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
% Diagonal values correspond to single parameters, absolute value needs to be greater than tolerance level
|
|
if abs(dSPECTRUM(j,j)) < tol_rank
|
|
spectrum_problpars{1} = [spectrum_problpars{1};j];
|
|
end
|
|
end
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
% Columns of dMINIMAL_par correspond to single parameters, needs to be augmented with dMINIMAL_rest (part that is independent of parameters),
|
|
% full rank would be equal to 1+dMINIMAL_fixed_rank_nbr
|
|
if strcmp(tol_rank,'robust')
|
|
if rank(full([dMINIMAL_par(:,j) dMINIMAL_rest])) == dMINIMAL_fixed_rank_nbr
|
|
minimal_problpars{1} = [minimal_problpars{1};j];
|
|
end
|
|
else
|
|
if rank(full([dMINIMAL_par(:,j) dMINIMAL_rest]),tol_rank) == dMINIMAL_fixed_rank_nbr
|
|
minimal_problpars{1} = [minimal_problpars{1};j];
|
|
end
|
|
end
|
|
end
|
|
end
|
|
|
|
% Check whether lack of identification is only due to single parameters
|
|
if ~no_identification_dynamic
|
|
if size(dynamic_problpars{1},1) == (totparam_nbr - rank_dDYNAMIC)
|
|
%found all nonidentified parameter sets
|
|
no_identification_dynamic = 1; %skip in the following
|
|
else
|
|
%still parameter sets that are nonidentified
|
|
indparam_dDYNAMIC(dynamic_problpars{1}) = []; %remove single unidentified parameters from higher-order sets of indparam
|
|
end
|
|
end
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
if size(reducedform_problpars{1},1) == (totparam_nbr - rank_dREDUCEDFORM)
|
|
%found all nonidentified parameter sets
|
|
no_identification_reducedform = 1; %skip in the following
|
|
else
|
|
%still parameter sets that are nonidentified
|
|
indparam_dREDUCEDFORM(reducedform_problpars{1}) = []; %remove single unidentified parameters from higher-order sets of indparam
|
|
end
|
|
end
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
if size(moments_problpars{1},1) == (totparam_nbr - rank_dMOMENTS)
|
|
%found all nonidentified parameter sets
|
|
no_identification_moments = 1; %skip in the following
|
|
else
|
|
%still parameter sets that are nonidentified
|
|
indparam_dMOMENTS(moments_problpars{1}) = []; %remove single unidentified parameters from higher-order sets of indparam
|
|
end
|
|
end
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
if size(spectrum_problpars{1},1) == (totparam_nbr - rank_dSPECTRUM)
|
|
%found all nonidentified parameter sets
|
|
no_identification_spectrum = 1; %skip in the following
|
|
else
|
|
%still parameter sets that are nonidentified
|
|
indparam_dSPECTRUM(spectrum_problpars{1}) = []; %remove single unidentified parameters from higher-order sets of indparam
|
|
end
|
|
end
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
if size(minimal_problpars{1},1) == (totparam_nbr + dMINIMAL_fixed_rank_nbr - rank_dMINIMAL)
|
|
%found all nonidentified parameter sets
|
|
no_identification_minimal = 1; %skip in the following
|
|
else
|
|
%still parameter sets that are nonidentified
|
|
indparam_dMINIMAL(minimal_problpars{1}) = []; %remove single unidentified parameters from higher-order sets of indparam
|
|
end
|
|
end
|
|
|
|
|
|
%% check higher order (j>1) parameter sets
|
|
%get common parameter indices from which to sample higher-order sets using nchoosek (we do not want to run nchoosek three times), most of the times indparamdMOMENTS, indparamdSPECTRUM, and indparamdMINIMAL are equal anyways
|
|
indtotparam = unique([indparam_dDYNAMIC indparam_dREDUCEDFORM indparam_dMOMENTS indparam_dSPECTRUM indparam_dMINIMAL]);
|
|
|
|
for j=2:min(length(indtotparam),max_dim_subsets_groups) % Check j-element subsets
|
|
h = dyn_waitbar(0,['Brute force collinearity for ' int2str(j) ' parameters.']);
|
|
%Step1: get all possible unique subsets of j elements
|
|
if ~no_identification_dynamic ...
|
|
|| (~no_identification_reducedform && ~error_indicator.identification_reducedform)...
|
|
|| (~no_identification_moments && ~error_indicator.identification_moments)...
|
|
|| (~no_identification_spectrum && ~error_indicator.identification_spectrum)...
|
|
|| (~no_identification_minimal && ~error_indicator.identification_minimal)
|
|
indexj=nchoosek(int16(indtotparam),j); % int16 speeds up nchoosek
|
|
% One could also use a mex version of nchoosek to speed this up, e.g.VChooseK from https://de.mathworks.com/matlabcentral/fileexchange/26190-vchoosek
|
|
end
|
|
|
|
%Step 2: remove already problematic sets and initialize rank vector
|
|
if ~no_identification_dynamic
|
|
indexj_dDYNAMIC = RemoveProblematicParameterSets(indexj,dynamic_problpars);
|
|
rankj_dDYNAMIC = zeros(size(indexj_dDYNAMIC,1),1);
|
|
else
|
|
indexj_dDYNAMIC = [];
|
|
end
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
indexj_dREDUCEDFORM = RemoveProblematicParameterSets(indexj,reducedform_problpars);
|
|
rankj_dREDUCEDFORM = zeros(size(indexj_dREDUCEDFORM,1),1);
|
|
else
|
|
indexj_dREDUCEDFORM = [];
|
|
end
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
indexj_dMOMENTS = RemoveProblematicParameterSets(indexj,moments_problpars);
|
|
rankj_dMOMENTS = zeros(size(indexj_dMOMENTS,1),1);
|
|
else
|
|
indexj_dMOMENTS = [];
|
|
end
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
indexj_dSPECTRUM = RemoveProblematicParameterSets(indexj,spectrum_problpars);
|
|
rankj_dSPECTRUM = zeros(size(indexj_dSPECTRUM,1),1);
|
|
else
|
|
indexj_dSPECTRUM = [];
|
|
end
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
indexj_dMINIMAL = RemoveProblematicParameterSets(indexj,minimal_problpars);
|
|
rankj_dMINIMAL = zeros(size(indexj_dMINIMAL,1),1);
|
|
else
|
|
indexj_dMINIMAL = [];
|
|
end
|
|
|
|
%Step3: Check rank criteria on submatrices
|
|
k_dDYNAMIC=0; k_dREDUCEDFORM=0; k_dMOMENTS=0; k_dSPECTRUM=0; k_dMINIMAL=0; %initialize counters
|
|
maxk = max([size(indexj_dDYNAMIC,1), size(indexj_dREDUCEDFORM,1), size(indexj_dMOMENTS,1), size(indexj_dMINIMAL,1), size(indexj_dSPECTRUM,1)]);
|
|
for k=1:maxk
|
|
if ~no_identification_dynamic
|
|
k_dDYNAMIC = k_dDYNAMIC+1;
|
|
if k_dDYNAMIC <= size(indexj_dDYNAMIC,1)
|
|
dDYNAMIC_j = dDYNAMIC(:,indexj_dDYNAMIC(k_dDYNAMIC,:)); % pick columns that correspond to parameter subset
|
|
if strcmp(tol_rank,'robust')
|
|
rankj_dDYNAMIC(k_dDYNAMIC,1) = rank(full(dDYNAMIC_j)); %compute rank with imposed tolerance
|
|
else
|
|
rankj_dDYNAMIC(k_dDYNAMIC,1) = rank(full(dDYNAMIC_j),tol_rank); %compute rank with imposed tolerance
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
k_dREDUCEDFORM = k_dREDUCEDFORM+1;
|
|
if k_dREDUCEDFORM <= size(indexj_dREDUCEDFORM,1)
|
|
dREDUCEDFORM_j = dREDUCEDFORM(:,indexj_dREDUCEDFORM(k_dREDUCEDFORM,:)); % pick columns that correspond to parameter subset
|
|
if strcmp(tol_rank,'robust')
|
|
rankj_dREDUCEDFORM(k_dREDUCEDFORM,1) = rank(full(dREDUCEDFORM_j)); %compute rank with imposed tolerance
|
|
else
|
|
rankj_dREDUCEDFORM(k_dREDUCEDFORM,1) = rank(full(dREDUCEDFORM_j),tol_rank); %compute rank with imposed tolerance
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
k_dMOMENTS = k_dMOMENTS+1;
|
|
if k_dMOMENTS <= size(indexj_dMOMENTS,1)
|
|
dMOMENTS_j = dMOMENTS(:,indexj_dMOMENTS(k_dMOMENTS,:)); % pick columns that correspond to parameter subset
|
|
if strcmp(tol_rank,'robust')
|
|
rankj_dMOMENTS(k_dMOMENTS,1) = rank(full(dMOMENTS_j)); %compute rank with imposed tolerance
|
|
else
|
|
rankj_dMOMENTS(k_dMOMENTS,1) = rank(full(dMOMENTS_j),tol_rank); %compute rank with imposed tolerance
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
k_dMINIMAL = k_dMINIMAL+1;
|
|
if k_dMINIMAL <= size(indexj_dMINIMAL,1)
|
|
dMINIMAL_j = [dMINIMAL_par(:,indexj_dMINIMAL(k_dMINIMAL,:)) dMINIMAL_rest]; % pick columns in dMINIMAL_par that correspond to parameter subset and augment with parameter-indepdendet part dMINIMAL_rest
|
|
if strcmp(tol_rank,'robust')
|
|
rankj_dMINIMAL(k_dMINIMAL,1) = rank(full(dMINIMAL_j)); %compute rank with imposed tolerance
|
|
else
|
|
rankj_dMINIMAL(k_dMINIMAL,1) = rank(full(dMINIMAL_j),tol_rank); %compute rank with imposed tolerance
|
|
end
|
|
end
|
|
end
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
k_dSPECTRUM = k_dSPECTRUM+1;
|
|
if k_dSPECTRUM <= size(indexj_dSPECTRUM,1)
|
|
dSPECTRUM_j = dSPECTRUM(indexj_dSPECTRUM(k_dSPECTRUM,:),indexj_dSPECTRUM(k_dSPECTRUM,:)); % pick rows and columns that correspond to parameter subset
|
|
if strcmp(tol_rank,'robust')
|
|
rankj_dSPECTRUM(k_dSPECTRUM,1) = rank(full(dSPECTRUM_j)); % Compute rank with imposed tol
|
|
else
|
|
rankj_dSPECTRUM(k_dSPECTRUM,1) = rank(full(dSPECTRUM_j),tol_rank); % Compute rank with imposed tol
|
|
end
|
|
end
|
|
end
|
|
dyn_waitbar(k/maxk,h)
|
|
end
|
|
|
|
%Step 4: Compare rank conditions for all possible subsets. If rank condition is violated, then the corresponding numbers of the parameters are stored
|
|
if ~no_identification_dynamic
|
|
dynamic_problpars{j} = indexj_dDYNAMIC(rankj_dDYNAMIC < j,:);
|
|
end
|
|
if ~no_identification_reducedform && ~error_indicator.identification_reducedform
|
|
reducedform_problpars{j} = indexj_dREDUCEDFORM(rankj_dREDUCEDFORM < j,:);
|
|
end
|
|
if ~no_identification_moments && ~error_indicator.identification_moments
|
|
moments_problpars{j} = indexj_dMOMENTS(rankj_dMOMENTS < j,:);
|
|
end
|
|
if ~no_identification_minimal && ~error_indicator.identification_minimal
|
|
minimal_problpars{j} = indexj_dMINIMAL(rankj_dMINIMAL < (j+dMINIMAL_fixed_rank_nbr),:);
|
|
end
|
|
if ~no_identification_spectrum && ~error_indicator.identification_spectrum
|
|
spectrum_problpars{j} = indexj_dSPECTRUM(rankj_dSPECTRUM < j,:);
|
|
end
|
|
% % Optional Step 5: % remove redundant 2-sets, eg. if the problematic sets are [(p1,p2);(p1,p3);(p2,p3)], then the unique problematic parameter sets are actually only [(p1,p2),(p1,p3)]
|
|
% if j == 2
|
|
% for jj=1:max([size(dynamic_problpars{2},1), size(reducedform_problpars{2},1), size(moments_problpars{2},1), size(spectrum_problpars{2},1), size(minimal_problpars{2},1)])
|
|
% if jj <= size(dynamic_problpars{2},1)
|
|
% dynamic_problpars{2}(dynamic_problpars{2}(jj,2)==dynamic_problpars{2}(:,1)) = dynamic_problpars{2}(jj,1);
|
|
% end
|
|
% if jj <= size(reducedform_problpars{2},1)
|
|
% reducedform_problpars{2}(reducedform_problpars{2}(jj,2)==reducedform_problpars{2}(:,1)) = reducedform_problpars{2}(jj,1);
|
|
% end
|
|
% if jj <= size(moments_problpars{2},1)
|
|
% moments_problpars{2}(moments_problpars{2}(jj,2)==moments_problpars{2}(:,1)) = moments_problpars{2}(jj,1);
|
|
% end
|
|
% if jj <= size(spectrum_problpars{2},1)
|
|
% spectrum_problpars{2}(spectrum_problpars{2}(jj,2)==spectrum_problpars{2}(:,1)) = spectrum_problpars{2}(jj,1);
|
|
% end
|
|
% if jj <= size(minimal_problpars{2},1)
|
|
% minimal_problpars{2}(minimal_problpars{2}(jj,2)==minimal_problpars{2}(:,1)) = minimal_problpars{2}(jj,1);
|
|
% end
|
|
% end
|
|
% dynamic_problpars{2} = unique(dynamic_problpars{2},'rows');
|
|
% reducedform_problpars{2} = unique(reducedform_problpars{2},'rows');
|
|
% moments_problpars{2} = unique(moments_problpars{2},'rows');
|
|
% spectrum_problpars{2} = unique(spectrum_problpars{2},'rows');
|
|
% minimal_problpars{2} = unique(minimal_problpars{2},'rows');
|
|
% % in indparam we replace the second parameter of problematic 2-sets by the observational equivalent first parameter to speed up nchoosek
|
|
% idx2 = unique([dynamic_problpars{2}; reducedform_problpars{2}; moments_problpars{2}; spectrum_problpars{2}; minimal_problpars{2}],'rows');
|
|
% if ~isempty(idx2)
|
|
% indtotparam(ismember(indtotparam,idx2(:,2))) = [];
|
|
% end
|
|
% end
|
|
dyn_waitbar_close(h);
|
|
end
|
|
|
|
%% Save output variables
|
|
if ~isempty(dynamic_problpars{1})
|
|
dynamic_problpars{1}(ismember(dynamic_problpars{1},1:(totparam_nbr-modparam_nbr))) = []; % get rid of stderr and corr variables for dynamic
|
|
end
|
|
ide_dynamic.problpars = dynamic_problpars;
|
|
ide_reducedform.problpars = reducedform_problpars;
|
|
ide_moments.problpars = moments_problpars;
|
|
ide_spectrum.problpars = spectrum_problpars;
|
|
ide_minimal.problpars = minimal_problpars;
|
|
|
|
%% Auxiliary functions
|
|
function idx = RemoveProblematicParameterSets(idx,problparset)
|
|
% Remove already problematic parameters
|
|
% INPUTS:
|
|
% * idx: complete index of possible combinations
|
|
% * problparset [cell] of all lower order combinations that are problematic
|
|
% * iset [integer] number of elements in set to consider
|
|
% OUTPUTS:
|
|
% * idx: index of possible combinations without already
|
|
% problematic lower order sets
|
|
iset = size(idx,2);
|
|
for iii=1:(iset-1)
|
|
if ~isempty(problparset{iii})
|
|
for kkk=1:size(problparset{iii},1)
|
|
idx((sum(ismember(idx,problparset{iii}(kkk,:)),2)==iii),:) = [];
|
|
end
|
|
end
|
|
end
|
|
end%RemoveProblematicParameterSets ed
|
|
|
|
|
|
end %main function end |