diff --git a/matlab/+bvar/density.m b/matlab/+bvar/density.m
index ba69bf9bbc1d23f5d4c8f5da99ef2801ebb967fb..9e1a219db796c4346d07e8052cbf89f0b70d4bc2 100644
--- a/matlab/+bvar/density.m
+++ b/matlab/+bvar/density.m
@@ -34,7 +34,7 @@ global oo_
 oo_.bvar.log_marginal_data_density=NaN(maxnlags,1);
 
 for nlags = 1:maxnlags
-    [ny, nx, posterior, prior] = bvar.toolbox(nlags);
+    [ny, ~, posterior, prior] = bvar.toolbox(nlags);
     oo_.bvar.posterior{nlags}=posterior;
     oo_.bvar.prior{nlags}=prior;
 
@@ -75,8 +75,8 @@ function w = matrictint(S, df, XXi)
 
 k=size(XXi,1);
 ny=size(S,1);
-[cx,p]=chol(XXi);
-[cs,q]=chol(S);
+cx = chol(XXi);
+cs = chol(S);
 
 if any(diag(cx)<100*eps)
     error('singular XXi')
diff --git a/matlab/+bvar/forecast.m b/matlab/+bvar/forecast.m
index 6fb0d2ca9d4d26ffc54e8b3e1075564bb5c04903..1e44f163f0a16e809c171b162fe669d2000f22f1 100644
--- a/matlab/+bvar/forecast.m
+++ b/matlab/+bvar/forecast.m
@@ -33,7 +33,7 @@ global options_ oo_ M_
 if options_.forecast == 0
     error('bvar.forecast: you must specify "forecast" option')
 end
-[ny, nx, posterior, prior, forecast_data] = bvar.toolbox(nlags);
+[ny, nx, posterior, ~, forecast_data] = bvar.toolbox(nlags);
 
 sims_no_shock = NaN(options_.forecast, ny, options_.bvar_replic);
 sims_with_shocks = NaN(options_.forecast, ny, options_.bvar_replic);
diff --git a/matlab/+bvar/irf.m b/matlab/+bvar/irf.m
index 69074e2cceabc39d9a0e09b601584ce602b80168..f2a216940395dec52ba808847d78449805826e80 100644
--- a/matlab/+bvar/irf.m
+++ b/matlab/+bvar/irf.m
@@ -35,7 +35,7 @@ if nargin==1
     identification = 'Cholesky';
 end
 
-[ny, nx, posterior, prior] = bvar.toolbox(nlags);
+[ny, nx, posterior] = bvar.toolbox(nlags);
 
 S_inv_upper_chol = chol(inv(posterior.S));
 
diff --git a/matlab/+estimate/nls.m b/matlab/+estimate/nls.m
index a7673004aec884b61058760d0b70738442a89e40..39422a75496de484053a55214ceaabfaaa74a55b 100644
--- a/matlab/+estimate/nls.m
+++ b/matlab/+estimate/nls.m
@@ -292,23 +292,23 @@ end
 %
 
 if is_gauss_newton
-    [params1, SSR, exitflag] = gauss_newton(resfun, params0);
+    [params1, SSR] = gauss_newton(resfun, params0);
 elseif is_lsqnonlin
     if ismember('levenberg-marquardt', varargin)
         % Levenberg Marquardt does not handle boundary constraints.
-        [params1, SSR, ~, exitflag] = lsqnonlin(resfun, params0, [], [], optimset(varargin{:}));
+        [params1, SSR] = lsqnonlin(resfun, params0, [], [], optimset(varargin{:}));
     else
-        [params1, SSR, ~, exitflag] = lsqnonlin(resfun, params0, bounds(:,1), bounds(:,2), optimset(varargin{:}));
+        [params1, SSR] = lsqnonlin(resfun, params0, bounds(:,1), bounds(:,2), optimset(varargin{:}));
     end
 else
     % Estimate the parameters by minimizing the sum of squared residuals.
-    [params1, SSR, exitflag] = dynare_minimize_objective(ssrfun, params0, ...
-                                                         minalgo, ...
-                                                         options_, ...
-                                                         bounds, ...
-                                                         parameter_names, ...
-                                                         [], ...
-                                                         []);
+    [params1, SSR] = dynare_minimize_objective(ssrfun, params0, ...
+                                               minalgo, ...
+                                               options_, ...
+                                               bounds, ...
+                                               parameter_names, ...
+                                               [], ...
+                                               []);
 end
 
 % Revert local modifications to the options.
diff --git a/matlab/+gsa/stability_mapping.m b/matlab/+gsa/stability_mapping.m
index 0da8eefb20ba53b54d9492b7385bdf96087a49dd..8badb3bc7d23f45c1a520b5883942b27e122c422 100644
--- a/matlab/+gsa/stability_mapping.m
+++ b/matlab/+gsa/stability_mapping.m
@@ -89,7 +89,7 @@ lpmat0=zeros(Nsam,0);
 xparam1=[];
 
 %% prepare prior bounds
-[~,~,~,lb,ub,~] = set_prior(estim_params_,M_,options_); %Prepare bounds
+[~,~,~,lb,ub] = set_prior(estim_params_,M_,options_); %Prepare bounds
 if ~isempty(bayestopt_) && any(bayestopt_.pshape > 0)
     % Set prior bounds
     bounds = prior_bounds(bayestopt_, options_.prior_trunc);
@@ -600,4 +600,4 @@ end
 skipline(1);
 
 xparam1=x0;
-save([OutputDirectoryName filesep 'prior_ok.mat'],'xparam1');
\ No newline at end of file
+save([OutputDirectoryName filesep 'prior_ok.mat'],'xparam1');
diff --git a/matlab/+identification/analysis.m b/matlab/+identification/analysis.m
index c1a387e2a5e0e612c4b555f717d8a1e9ae2e19b0..84059027c16c6c8a386362f849d0883d4e851f58 100644
--- a/matlab/+identification/analysis.m
+++ b/matlab/+identification/analysis.m
@@ -206,7 +206,7 @@ if info(1) == 0 %no errors in solution
                 options_ident_local.no_identification_spectrum = 1; %do not recompute dSPECTRUM
                 options_ident_local.ar = nlags;     %store new lag number
                 options_.ar      = nlags;           %store new lag number
-                [~, ~, ~, ~, ~, ~, MOMENTS, dMOMENTS, ~, ~, ~, ~] = identification.get_jacobians(estim_params_, M_, options_, options_ident_local, indpmodel, indpstderr, indpcorr, indvobs, oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+                [~, ~, ~, ~, ~, ~, MOMENTS, dMOMENTS] = identification.get_jacobians(estim_params_, M_, options_, options_ident_local, indpmodel, indpstderr, indpcorr, indvobs, oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
 
                 ind_dMOMENTS = (find(max(abs(dMOMENTS'),[],1) > tol_deriv)); %new index with non-zero rows
             end
diff --git a/matlab/+identification/checks.m b/matlab/+identification/checks.m
index 03895b6ad13bb9674835e079a8f6c59cbf1c056a..015575a8bf05d5843a7cc5f2926edda796c25b03 100644
--- a/matlab/+identification/checks.m
+++ b/matlab/+identification/checks.m
@@ -86,7 +86,7 @@ else
     Xparnonzero = Xpar(:,ind1); % focus on non-zero columns
 end
 
-[eu, ee2, ee1] = svd( [Xparnonzero Xrest], 0 );
+[~, ~, ee1] = svd( [Xparnonzero Xrest], 0 );
 condX = cond([Xparnonzero Xrest]);
 rankX = rank(X, tol_rank);
 icheck = 0; %initialize flag which is equal to 0 if we already found all single parameters that are not identified
@@ -118,7 +118,7 @@ if icheck
     else
         Xparnonzero = Xpar(:,ind1); % focus on non-zero columns
     end
-    [eu, ee2, ee1] = svd( [Xparnonzero Xrest], 0 );
+    [~, ~, ee1] = svd( [Xparnonzero Xrest], 0 );
     condX = cond([Xparnonzero Xrest]);
     rankX = rank(X,tol_rank);
 end
diff --git a/matlab/+identification/get_minimal_state_representation.m b/matlab/+identification/get_minimal_state_representation.m
index deb6d43bbca09eab2c394753ce6cc0a284433a59..24d18bca49536b79c98ec9d19d6fb352a662c891 100644
--- a/matlab/+identification/get_minimal_state_representation.m
+++ b/matlab/+identification/get_minimal_state_representation.m
@@ -209,7 +209,7 @@ function [mSYS,U] = minrealold(SYS,tol)
     a = SYS.A;
     b = SYS.B;
     c = SYS.C;
-    [ns,nu] = size(b);
+    [ns,~] = size(b);
     [am,bm,cm,U,k] = ControllabilityStaircaseRosenbrock(a,b,c,tol);
     kk = sum(k);
     nu = ns - kk;
@@ -219,7 +219,7 @@ function [mSYS,U] = minrealold(SYS,tol)
     cm = cm(:,nu+1:ns);
     ns = ns - nu;
     if ns
-        [am,bm,cm,t,k] = ObservabilityStaircaseRosenbrock(am,bm,cm,tol);
+        [am,bm,cm,~,k] = ObservabilityStaircaseRosenbrock(am,bm,cm,tol);
         kk = sum(k);
         nu = ns - kk;
         nn = nn + nu;
@@ -242,8 +242,8 @@ end
 
 function [abar,bbar,cbar,t,k] = ControllabilityStaircaseRosenbrock(a, b, c, tol)
     % Controllability staircase algorithm of Rosenbrock, 1968
-    [ra,ca] = size(a);
-    [rb,cb] = size(b);
+    [ra,~] = size(a);
+    [~,cb] = size(b);
     ptjn1 = eye(ra);
     ajn1 = a;
     bjn1 = b;
@@ -255,8 +255,8 @@ function [abar,bbar,cbar,t,k] = ControllabilityStaircaseRosenbrock(a, b, c, tol)
         tol = ra*norm(a,1)*eps;
     end
     for jj = 1 : ra
-        [uj,sj,vj] = svd(bjn1);
-        [rsj,csj] = size(sj);
+        [uj,sj] = svd(bjn1);
+        [rsj,~] = size(sj);
         %p = flip(eye(rsj),2);
         p = eye(rsj);
         p = p(:,end:-1:1);
@@ -264,7 +264,7 @@ function [abar,bbar,cbar,t,k] = ControllabilityStaircaseRosenbrock(a, b, c, tol)
         uj = uj*p;
         bb = uj'*bjn1;
         roj = rank(bb,tol);
-        [rbb,cbb] = size(bb);
+        [rbb,~] = size(bb);
         sigmaj = rbb - roj;
         sigmajn1 = sigmaj;
         k(jj) = roj;
diff --git a/matlab/+identification/numerical_objective.m b/matlab/+identification/numerical_objective.m
index 84c1695b7f2b2c2b6f21c5888397b36cddcf83eb..7237ec3a9c6b2e98879ba827cd74d3109b8628aa 100644
--- a/matlab/+identification/numerical_objective.m
+++ b/matlab/+identification/numerical_objective.m
@@ -78,7 +78,7 @@ else
 end
 
 %% compute Kalman transition matrices and steady state with updated parameters
-[dr,info,M_.params] = compute_decision_rules(M_,options_,dr, steady_state, exo_steady_state, exo_det_steady_state);
+[dr,~,M_.params] = compute_decision_rules(M_,options_,dr, steady_state, exo_steady_state, exo_det_steady_state);
 options_ = rmfield(options_,'options_ident');
 pruned = pruned_SS.pruned_state_space_system(M_, options_, dr, indvar, nlags, useautocorr, 0);
 
diff --git a/matlab/+identification/run.m b/matlab/+identification/run.m
index fac0c5d22080c5555eee2d793694d6a953b4b5aa..a65204dd67f1c21457fef675b25cf0e1f726efe8 100644
--- a/matlab/+identification/run.m
+++ b/matlab/+identification/run.m
@@ -302,7 +302,7 @@ options_.mode_compute = 0;
 options_.plot_priors = 0;
 options_.smoother = 1;
 options_.options_ident = [];
-[dataset_, dataset_info, xparam1, hh, M_, options_, oo_, estim_params_, bayestopt_, bounds] = dynare_estimation_init(M_.endo_names, fname, 1, M_, options_, oo_, estim_params_, bayestopt_);
+[~, dataset_info, ~, ~, M_, options_, oo_, estim_params_, bayestopt_] = dynare_estimation_init(M_.endo_names, fname, 1, M_, options_, oo_, estim_params_, bayestopt_);
 
 % set method to compute identification Jacobians (kronflag). Default:0
 options_ident = set_default_option(options_ident,'analytic_derivation_mode', options_.analytic_derivation_mode); % if not set by user, inherit default global one
@@ -470,7 +470,7 @@ if iload <=0
     end
     options_ident.tittxt = parameters; %title text for graphs and figures
     % perform identification analysis for single point
-    [ide_moments_point, ide_spectrum_point, ide_minimal_point, ide_hess_point, ide_reducedform_point, ide_dynamic_point, derivatives_info_point, info, error_indicator_point] = ...
+    [ide_moments_point, ide_spectrum_point, ide_minimal_point, ide_hess_point, ide_reducedform_point, ide_dynamic_point, ~, info, error_indicator_point] = ...
         identification.analysis(M_,options_,oo_,bayestopt_,estim_params_,params, indpmodel, indpstderr, indpcorr, options_ident, dataset_info, prior_exist, 1); %the 1 at the end implies initialization of persistent variables
     if info(1)~=0
         % there are errors in the solution algorithm
@@ -487,7 +487,7 @@ if iload <=0
                 params = Prior.draw();
                 options_ident.tittxt = 'Random_prior_params'; %title text for graphs and figures
                 % perform identification analysis
-                [ide_moments_point, ide_spectrum_point, ide_minimal_point, ide_hess_point, ide_reducedform_point, ide_dynamic_point, derivatives_info_point, info, error_indicator_point] = ...
+                [ide_moments_point, ide_spectrum_point, ide_minimal_point, ide_hess_point, ide_reducedform_point, ide_dynamic_point, ~, info, error_indicator_point] = ...
                     identification.analysis(M_,options_,oo_,bayestopt_,estim_params_,params, indpmodel, indpstderr, indpcorr, options_ident, dataset_info, prior_exist, 1);
             end
         end
@@ -542,7 +542,7 @@ if iload <=0
         end
         options_ident.tittxt = []; % clear title text for graphs and figures
         % run identification analysis
-        [ide_moments, ide_spectrum, ide_minimal, ide_hess, ide_reducedform, ide_dynamic, ide_derivatives_info, info, error_indicator] = ...
+        [ide_moments, ide_spectrum, ide_minimal, ~, ide_reducedform, ide_dynamic, ~, info, error_indicator] = ...
             identification.analysis(M_,options_,oo_,bayestopt_,estim_params_,params, indpmodel, indpstderr, indpcorr, options_MC, dataset_info, prior_exist, 0); % the 0 implies that we do not initialize persistent variables anymore
 
         if iteration==0 && info(1)==0 % preallocate storage in the first admissable run
@@ -911,7 +911,7 @@ if SampleSize > 1
                 fprintf('\nTesting %s.\n',tittxt);
                 if ~iload
                     options_ident.tittxt = tittxt; %title text for graphs and figures
-                    [ide_moments_max, ide_spectrum_max, ide_minimal_max, ide_hess_max, ide_reducedform_max, ide_dynamic_max, derivatives_info_max, info_max, error_indicator_max] = ...
+                    [ide_moments_max, ide_spectrum_max, ide_minimal_max, ide_hess_max, ide_reducedform_max, ide_dynamic_max, ~, ~, error_indicator_max] = ...
                         identification.analysis(M_,options_,oo_,bayestopt_,estim_params_,pdraws(jmax,:), indpmodel, indpstderr, indpcorr, options_ident, dataset_info, prior_exist, 1); %the 1 at the end initializes some persistent variables
                     save([IdentifDirectoryName '/' fname '_identif.mat'], 'ide_hess_max', 'ide_moments_max', 'ide_spectrum_max', 'ide_minimal_max','ide_reducedform_max', 'ide_dynamic_max', 'jmax', '-append');
                 end
@@ -978,4 +978,4 @@ end
 %reset warning state
 warning_config;
 
-fprintf('\n==== Identification analysis completed ====\n\n')
\ No newline at end of file
+fprintf('\n==== Identification analysis completed ====\n\n')
diff --git a/matlab/+mom/mode_compute_gmm_smm.m b/matlab/+mom/mode_compute_gmm_smm.m
index 4b934bc512c1f6595d3e72b7e490638c8c1de002..2c1fc5ef6772b646a175701162e0190a3aefaebf 100644
--- a/matlab/+mom/mode_compute_gmm_smm.m
+++ b/matlab/+mom/mode_compute_gmm_smm.m
@@ -112,8 +112,8 @@ for stage_iter = 1:size(options_mom_.mom.weighting_matrix,1)
             else
                 options_mom_.mom.compute_derivs = false;
             end
-            [xparam1, fval, exitflag] = dynare_minimize_objective(objective_function, xparam0, options_mom_.optimizer_vec{optim_iter}, options_mom_, [Bounds.lb Bounds.ub], bayestopt_.name, bayestopt_, [],...
-                                                                  Bounds, oo_, estim_params_, M_, options_mom_);
+            [xparam1, fval] = dynare_minimize_objective(objective_function, xparam0, options_mom_.optimizer_vec{optim_iter}, options_mom_, [Bounds.lb Bounds.ub], bayestopt_.name, bayestopt_, [],...
+                                                        Bounds, oo_, estim_params_, M_, options_mom_);
             if options_mom_.mom.vector_output
                 fval = fval'*fval;
             end
@@ -126,4 +126,4 @@ for stage_iter = 1:size(options_mom_.mom.weighting_matrix,1)
     end
     options_mom_.vector_output = false;
     [~, ~, ~,~,~, oo_] = feval(objective_function, xparam1, Bounds, oo_, estim_params_, M_, options_mom_); % get oo_.mom.model_moments for iterated GMM/SMM to compute optimal weighting matrix
-end
\ No newline at end of file
+end
diff --git a/matlab/+mom/run.m b/matlab/+mom/run.m
index 3994ffbe91832fec7782800c64d54dae4e335e88..b5579aeb49ce6bb483a50334113b0c7e86ff9a5b 100644
--- a/matlab/+mom/run.m
+++ b/matlab/+mom/run.m
@@ -403,7 +403,7 @@ end
 % Build dataset
 if strcmp(options_mom_.mom.mom_method,'GMM') || strcmp(options_mom_.mom.mom_method,'SMM')
     % Check if datafile has same name as mod file
-    [~,name,~] = fileparts(options_mom_.datafile);
+    [~,name] = fileparts(options_mom_.datafile);
     if strcmp(name,M_.fname)
         error('method_of_moments: ''datafile'' and mod file are not allowed to have the same name; change the name of the ''datafile''!')
     end
diff --git a/matlab/+occbin/DSGE_smoother.m b/matlab/+occbin/DSGE_smoother.m
index 5195ec5f55b6340847cc0363eae4eff7d475c83b..d2bcdc6a42b8dbc26f1f36f27b499d13587bfe9f 100644
--- a/matlab/+occbin/DSGE_smoother.m
+++ b/matlab/+occbin/DSGE_smoother.m
@@ -118,7 +118,7 @@ occbin_options.first_period_occbin_update = options_.occbin.smoother.first_perio
 occbin_options.opts_regime = opts_simul; % this builds the opts_simul options field needed by occbin.solver
 occbin_options.opts_regime.binding_indicator = options_.occbin.likelihood.init_binding_indicator;
 occbin_options.opts_regime.regime_history=options_.occbin.likelihood.init_regime_history;
-[alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T0,R0,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0, diffuse_steps] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,occbin_options);%     T1=TT;
+[alphahat,etahat,epsilonhat,ahat,SteadyState,trend_coeff,aK,T0,R0,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_,alphahat0,state_uncertainty0] = DsgeSmoother(xparam1,gend,Y,data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,occbin_options);%     T1=TT;
 
 oo_.occbin.smoother.realtime_regime_history = oo_.occbin.smoother.regime_history;
 regime_history = oo_.occbin.smoother.regime_history;
diff --git a/matlab/+occbin/IVF_posterior.m b/matlab/+occbin/IVF_posterior.m
index 3a8b018c98ea9b38e5a2f406f453d22c5725a51b..69c6c110ddb34815d31c8df8701ae1e24c5fe29a 100644
--- a/matlab/+occbin/IVF_posterior.m
+++ b/matlab/+occbin/IVF_posterior.m
@@ -71,7 +71,7 @@ end
 
 if ~isempty(xparam1)
     M_ = set_all_parameters(xparam1,estim_params_,M_);
-    [fval,info,exit_flag,Q,H]=check_bounds_and_definiteness_estimation(xparam1, M_, estim_params_, BoundsInfo);
+    [fval,info,exit_flag]=check_bounds_and_definiteness_estimation(xparam1, M_, estim_params_, BoundsInfo);
     if info(1)
         return
     end
@@ -81,7 +81,7 @@ err_index=options_.occbin.likelihood.IVF_shock_observable_mapping; % err_index=
 COVMAT1 = M_.Sigma_e(err_index,err_index);
 
 % Linearize the model around the deterministic steady state and extract the matrices of the state equation (T and R).
-[T,R,SteadyState,info,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state,'restrict');
+[~,~,SteadyState,info,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state,'restrict');
 
 % Return, with endogenous penalty when possible, if dynare_resolve issues an error code (defined in resol).
 if info(1)
@@ -195,4 +195,4 @@ end
 
 % remember that the likelihood has already been multiplied by -1
 % hence, posterior is -1 times the log of the prior
-fval = like+prior;
\ No newline at end of file
+fval = like+prior;
diff --git a/matlab/+occbin/forecast.m b/matlab/+occbin/forecast.m
index 56cc1999fd5ed0db7793ed49678dfd669a1b0b29..07c8a94f139f8a6b3c6d7a6202a1b240475cd512 100644
--- a/matlab/+occbin/forecast.m
+++ b/matlab/+occbin/forecast.m
@@ -28,7 +28,7 @@ if opts.replic
     effective_exo_nbr=  length(ishock);
     effective_exo_names = M_.exo_names(ishock);
     effective_Sigma_e = M_.Sigma_e(ishock,ishock);
-    [U,S,V] = svd(effective_Sigma_e);
+    [U,S] = svd(effective_Sigma_e);
     if opts.qmc
         opts.replic =2^(round(log2(opts.replic+1)))-1;
         SHOCKS_ant =   qmc_sequence(forecast*effective_exo_nbr, int64(1), 1, opts.replic)';
diff --git a/matlab/+occbin/kalman_update_algo_1.m b/matlab/+occbin/kalman_update_algo_1.m
index e8d4a487d519de49780123d7a23641a571dda72a..fe4a71220af6fcad370862f50d6f81a28a1ed62d 100644
--- a/matlab/+occbin/kalman_update_algo_1.m
+++ b/matlab/+occbin/kalman_update_algo_1.m
@@ -292,8 +292,7 @@ error_flag = out.error_flag;
 if ~error_flag && niter>options_.occbin.likelihood.max_number_of_iterations && ~isequal(regimes_(1),regimes0(1))
     error_flag = 1;
     if M_.occbin.constraint_nbr==1 % try some other regime
-        [ll, il]=sort(lik_hist);
-        [ll, il]=sort(regime_end);
+        [~, il]=sort(regime_end);
         rr=regime_hist(il(2:3));
         newstart=1;
         if length(rr{1}.regimestart)>1
diff --git a/matlab/+occbin/kalman_update_algo_3.m b/matlab/+occbin/kalman_update_algo_3.m
index b07a6d2276740e1cffedcb1bc369c2e271084e77..2ed4c603ba59a3789e89545df01968be2e94499f 100644
--- a/matlab/+occbin/kalman_update_algo_3.m
+++ b/matlab/+occbin/kalman_update_algo_3.m
@@ -287,7 +287,7 @@ if error_flag==0 && niter>options_.occbin.likelihood.max_number_of_iterations &&
     
   if M_.occbin.constraint_nbr==1
     % try some other regime before giving up
-    [ll, il]=sort(regime_end);
+    [~, il]=sort(regime_end);
     rr=regime_hist(il(2:3));
     newstart=1;
     if length(rr{1}(1).regimestart)>1
diff --git a/matlab/+pac/+estimate/init.m b/matlab/+pac/+estimate/init.m
index 891770b5ab1afedce359306b8f33091fa6994e21..b81b31d460862d36a45372bf8c2bfa0fc6756d55 100644
--- a/matlab/+pac/+estimate/init.m
+++ b/matlab/+pac/+estimate/init.m
@@ -18,7 +18,7 @@ function [pacmodl, lhs, rhs, pnames, enames, xnames, rname, pid, eid, xid, pname
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 % Get the original equation to be estimated
-[LHS, RHS] = get_lhs_and_rhs(eqname, M_, true);
+[~, RHS] = get_lhs_and_rhs(eqname, M_, true);
 
 % Check that the equation has a PAC expectation term.
 if ~contains(RHS, 'pac_expectation', 'IgnoreCase', true)
diff --git a/matlab/+pac/+estimate/iterative_ols.m b/matlab/+pac/+estimate/iterative_ols.m
index cbf5655fe42cbe696c2d976b2028b5bef03b37a5..3808648c929fb7ac3bf1eccb3b5e1eab402e6cac 100644
--- a/matlab/+pac/+estimate/iterative_ols.m
+++ b/matlab/+pac/+estimate/iterative_ols.m
@@ -41,7 +41,7 @@ function iterative_ols(eqname, params, data, range)
 
 global M_ oo_ options_
 
-[pacmodl, ~, rhs, ~, ~, ~, rname, ~, ~, ~, ~, ipnames_, params, data, ~] = ...
+[pacmodl, ~, rhs, ~, ~, ~, rname, ~, ~, ~, ~, ipnames_, params, data] = ...
     pac.estimate.init(M_, oo_, eqname, params, data, range);
 
 % Set initial condition.
diff --git a/matlab/+pac/+estimate/nls.m b/matlab/+pac/+estimate/nls.m
index 5fa8a626ecf13b0fc5949d5c804edf3979b00b90..873a14783d50c6bad288fb924d0e1ec346a27afc 100644
--- a/matlab/+pac/+estimate/nls.m
+++ b/matlab/+pac/+estimate/nls.m
@@ -200,23 +200,23 @@ if isnan(ssr0) || isinf(ssr0) || ~isreal(ssr0)
 end
 
 if is_gauss_newton
-    [params1, SSR, exitflag] = gauss_newton(resfun, params0);
+    [params1, SSR] = gauss_newton(resfun, params0);
 elseif is_lsqnonlin
     if ismember('levenberg-marquardt', varargin)
         % Levenberg Marquardt does not handle boundary constraints.
-        [params1, SSR, ~, exitflag] = lsqnonlin(resfun, params0, [], [], optimset(varargin{:}));
+        [params1, SSR] = lsqnonlin(resfun, params0, [], [], optimset(varargin{:}));
     else
-        [params1, SSR, ~, exitflag] = lsqnonlin(resfun, params0, bounds(:,1), bounds(:,2), optimset(varargin{:}));
+        [params1, SSR] = lsqnonlin(resfun, params0, bounds(:,1), bounds(:,2), optimset(varargin{:}));
     end
 else
     % Estimate the parameters by minimizing the sum of squared residuals.
-    [params1, SSR, exitflag] = dynare_minimize_objective(ssrfun, params0, ...
-                                                      minalgo, ...
-                                                      options_, ...
-                                                      bounds, ...
-                                                      parameter_names, ...
-                                                      [], ...
-                                                      []);
+    [params1, SSR] = dynare_minimize_objective(ssrfun, params0, ...
+                                               minalgo, ...
+                                               options_, ...
+                                               bounds, ...
+                                               parameter_names, ...
+                                               [], ...
+                                               []);
 end
 
 % Revert local modifications to the options.
diff --git a/matlab/+pac/check.m b/matlab/+pac/check.m
index 1bf282abb8bfb4e25eb4d6cb06fad4e2150a92f8..dc6b49a6728268ddbd7216019a96d00dbbf94a3e 100644
--- a/matlab/+pac/check.m
+++ b/matlab/+pac/check.m
@@ -40,7 +40,7 @@ end
 errorcode = 0;
 
 % Get the original equation to be estimated
-[LHS, RHS] = get_lhs_and_rhs(eqname, M_, true);
+[~, RHS] = get_lhs_and_rhs(eqname, M_, true);
 
 % Check that the equation has a PAC expectation term.
 if ~contains(RHS, 'pac_expectation', 'IgnoreCase', true)
diff --git a/matlab/AIM/dynAIMsolver1.m b/matlab/AIM/dynAIMsolver1.m
index ca1b44f06a621d70f871911ca7b898334fef5e67..78d770cd37890fd847ae37dec7ed1502ca788d70 100644
--- a/matlab/AIM/dynAIMsolver1.m
+++ b/matlab/AIM/dynAIMsolver1.m
@@ -89,7 +89,7 @@ if leads ==0
 end
 %disp('gensysToAMA:running ama');
 try % try to run AIM
-    [bb,rts,ia,nexact,nnumeric,lgroots,aimcode] =...
+    [bb,rts,~,~,~,~,aimcode] =...
         SPAmalg(theAIM_H,neq, lags,leads,condn,uprbnd);
 catch
     err = lasterror;
diff --git a/matlab/aggregate.m b/matlab/aggregate.m
index 979a9d93067131c8b96b147b4e9177d92b662ca1..b8ae195b04f2703bd52f722f7a79f1098d40c1d7 100644
--- a/matlab/aggregate.m
+++ b/matlab/aggregate.m
@@ -84,7 +84,7 @@ for i=1:length(varargin)
             end
             % Ensure that the equation tag name matches the LHS variable.
             eqtagname = regexp(model{j}, 'name=''(\w*)''', 'match');
-            [lhs, ~] = getequation(model{j+1});
+            lhs = getequation(model{j+1});
             endovar = getendovar(lhs);
             eqtagname_ = strcat('name=''', endovar{1}, '''');
             if ~isempty(eqtagname)
diff --git a/matlab/backward/backward_model_forecast.m b/matlab/backward/backward_model_forecast.m
index 836737402ef6bf8c1b8707c95709b63923144149..505ee9596f72d194ee7d78cd06c0f79345907693 100644
--- a/matlab/backward/backward_model_forecast.m
+++ b/matlab/backward/backward_model_forecast.m
@@ -63,7 +63,7 @@ for i=1:M_.exo_nbr
 end
 
 % Set up initial conditions
-[initialcondition, periods, innovations, options_local, M_local, oo_local, endonames, exonames, dynamic_resid, dynamic_g1, y] = ...
+[initialcondition, periods, innovations, options_local, M_local, oo_local, endonames, ~, dynamic_resid, dynamic_g1] = ...
     simul_backward_model_init(initialcondition, periods, options_, M_, oo_, zeros(periods, M_.exo_nbr));
 
 % Get vector of indices for the selected endogenous variables.
@@ -110,9 +110,9 @@ if withuncertainty
     for i=1:B
         innovations = transpose(sigma*randn(M_.exo_nbr, periods));
         if options_.linear
-            [ysim__, xsim__, errorflag] = simul_backward_linear_model_(initialcondition, periods, options_local, M_local, oo_local, innovations, dynamic_resid, dynamic_g1);
+            [ysim__, ~, errorflag] = simul_backward_linear_model_(initialcondition, periods, options_local, M_local, oo_local, innovations, dynamic_resid, dynamic_g1);
         else
-            [ysim__, xsim__, errorflag] = simul_backward_nonlinear_model_(initialcondition, periods, options_local, M_local, oo_local, innovations, dynamic_resid, dynamic_g1);
+            [ysim__, ~, errorflag] = simul_backward_nonlinear_model_(initialcondition, periods, options_local, M_local, oo_local, innovations, dynamic_resid, dynamic_g1);
         end
         if errorflag
             error('Simulation failed.')
diff --git a/matlab/backward/backward_model_irf.m b/matlab/backward/backward_model_irf.m
index eadf6001676d33ec93eab080ebccb16fa274e59c..fd60176bb639448bec5e84319ce5d25521e6bd9c 100644
--- a/matlab/backward/backward_model_irf.m
+++ b/matlab/backward/backward_model_irf.m
@@ -139,7 +139,7 @@ if ~isempty(innovationbaseline)
 end
 
 % Set up initial conditions
-[initialcondition, periods, Innovations, options_local, M_local, oo_local, endonames, exonames, dynamic_resid, dynamic_g1, y] = ...
+[initialcondition, periods, Innovations, options_local, M_local, oo_local, endonames, exonames, dynamic_resid, dynamic_g1] = ...
     simul_backward_model_init(initialcondition, periods, options_, M_, oo_, Innovations);
 
 % Get the covariance matrix of the shocks.
diff --git a/matlab/backward/simul_backward_model_init.m b/matlab/backward/simul_backward_model_init.m
index f7f8d4a14765e2378bfbc8658abe7fa450ead07a..c77176a1a571d9524c53ff9eac09979cb410a116 100644
--- a/matlab/backward/simul_backward_model_init.m
+++ b/matlab/backward/simul_backward_model_init.m
@@ -57,7 +57,7 @@ if isempty(initialconditions)
                                 vertcat(M_.endo_names(1:M_.orig_endo_nbr), M_.exo_names));
 end
 
-[initialconditions, info] = checkdatabase(initialconditions, M_, false, true);
+initialconditions = checkdatabase(initialconditions, M_, false, true);
 
 % Test if the first argument contains all the lagged endogenous variables
 endonames = M_.endo_names;
diff --git a/matlab/backward/simul_backward_nonlinear_model_.m b/matlab/backward/simul_backward_nonlinear_model_.m
index 3aa4d180ded5b63520ce5c570f020bf71fc2dbd5..51ef7eedfe573685baf4f56ad9d5f008b0f7169e 100644
--- a/matlab/backward/simul_backward_nonlinear_model_.m
+++ b/matlab/backward/simul_backward_nonlinear_model_.m
@@ -143,7 +143,7 @@ for it = initialconditions.nobs+(1:samplesize)
         %
         % Evaluate and check the residuals
         %
-        [r, J] = dynamic_backward_model_for_simulation(ytm, dynamic_resid, dynamic_g1, ytm, x, M_.params, oo_.steady_state, M_.dynamic_g1_sparse_rowval, M_.dynamic_g1_sparse_colval, M_.dynamic_g1_sparse_colptr);
+        r = dynamic_backward_model_for_simulation(ytm, dynamic_resid, dynamic_g1, ytm, x, M_.params, oo_.steady_state, M_.dynamic_g1_sparse_rowval, M_.dynamic_g1_sparse_colval, M_.dynamic_g1_sparse_colptr);
         residuals_evaluating_to_nan = isnan(r);
         residuals_evaluating_to_inf = isinf(r);
         residuals_evaluating_to_complex = ~isreal(r);
diff --git a/matlab/clear_persistent_variables.m b/matlab/clear_persistent_variables.m
index 32223192bd90c50541718c4a61b7951cfb58d705..09775d39c3c8236dd9eda575c8f92a5cc007852d 100644
--- a/matlab/clear_persistent_variables.m
+++ b/matlab/clear_persistent_variables.m
@@ -52,7 +52,7 @@ if writelistofroutinestobecleared
                 end
             end
         end
-        [~, list_of_functions, ~] = cellfun(@fileparts, list_of_files, 'UniformOutput',false);
+        [~, list_of_functions] = cellfun(@fileparts, list_of_files, 'UniformOutput',false);
         cellofchar2mfile(sprintf('%slist_of_functions_to_be_cleared.m', DYNARE_FOLDER), list_of_functions)
     end
     return
diff --git a/matlab/cli/+cli/+evaluate/smoother.m b/matlab/cli/+cli/+evaluate/smoother.m
index 94a2bb8f778dcd39cd327558b2a189874da28fa2..8fb020cb6fc55fd4f787e960279e0c24357853cd 100644
--- a/matlab/cli/+cli/+evaluate/smoother.m
+++ b/matlab/cli/+cli/+evaluate/smoother.m
@@ -35,4 +35,4 @@ end
 
 parameters = strrep(parameters, ' ', '_');
 
-[oo_, M_, options_, bayestopt_, Smoothed_variables_declaration_order_deviation_form] = evaluate_smoother(parameters, varlist, M_, oo_, options_, bayestopt_, estim_params_);
\ No newline at end of file
+[oo_, M_, options_, bayestopt_] = evaluate_smoother(parameters, varlist, M_, oo_, options_, bayestopt_, estim_params_);
diff --git a/matlab/cli/prior.m b/matlab/cli/prior.m
index 9ed74fde1907293527bfc785a386808bf55ca531..6d4825b39153a1b59649269d9c99bd1cd64110c0 100644
--- a/matlab/cli/prior.m
+++ b/matlab/cli/prior.m
@@ -61,7 +61,7 @@ if (size(estim_params_.var_endo,1) || size(estim_params_.corrn,1))
 end
 
 % Fill or update bayestopt_ structure
-[xparam1, estim_params_, BayesOptions, lb, ub, M_local] = set_prior(estim_params_, M_, options_);
+[~, estim_params_, ~, lb, ub, M_local] = set_prior(estim_params_, M_, options_);
 % Set restricted state space
 options_plot_priors_old=options_.plot_priors;
 options_.plot_priors=0;
diff --git a/matlab/convergence_diagnostics/mcmc_diagnostics.m b/matlab/convergence_diagnostics/mcmc_diagnostics.m
index 1c1b1abb72f91cc1c6d8f9f873cbe3b26b564e95..6872074fbecf75cd1a2df41bde77c1e278760b2e 100644
--- a/matlab/convergence_diagnostics/mcmc_diagnostics.m
+++ b/matlab/convergence_diagnostics/mcmc_diagnostics.m
@@ -272,7 +272,7 @@ else
     end
     NamFileInput={[M_.dname '/metropolis/'],[ModelName '_mh*_blck*.mat']};
 
-    [fout, nBlockPerCPU, totCPU] = masterParallel(options_.parallel, 1, npar,NamFileInput,'mcmc_diagnostics_core', localVars, [], options_.parallel_info);
+    [fout, ~, totCPU] = masterParallel(options_.parallel, 1, npar,NamFileInput,'mcmc_diagnostics_core', localVars, [], options_.parallel_info);
     UDIAG = fout(1).UDIAG;
     for j=2:totCPU
         UDIAG = cat(3,UDIAG ,fout(j).UDIAG);
diff --git a/matlab/convergence_diagnostics/raftery_lewis.m b/matlab/convergence_diagnostics/raftery_lewis.m
index 48c6b4941aa330d4a81aa1c118d3af1ac1a468cc..324a524e07688e087c3bea38721cddcae7b6864f 100644
--- a/matlab/convergence_diagnostics/raftery_lewis.m
+++ b/matlab/convergence_diagnostics/raftery_lewis.m
@@ -93,7 +93,7 @@ for nv = 1:n_vars % big loop over variables
     % Find thinning factor for which first-order Markov Chain is preferred to second-order one
     while(bic > 0)
         thinned_chain=work(1:k_thin_current_var:n_runs,1);
-        [g2, bic] = first_vs_second_order_MC_test(thinned_chain);
+        [~, bic] = first_vs_second_order_MC_test(thinned_chain);
         k_thin_current_var = k_thin_current_var+1;
     end
 
@@ -108,11 +108,11 @@ for nv = 1:n_vars % big loop over variables
     beta = transition_matrix(2,1)/(transition_matrix(2,1)+transition_matrix(2,2));  %prob of going from 2 to 1
 
     kmind=k_thin_current_var;
-    [g2, bic]=independence_chain_test(thinned_chain);
+    [~, bic]=independence_chain_test(thinned_chain);
 
     while(bic > 0)
         thinned_chain=work(1:kmind:n_runs,1);
-        [g2, bic] = independence_chain_test(thinned_chain);
+        [~, bic] = independence_chain_test(thinned_chain);
         kmind = kmind+1;
     end
 
diff --git a/matlab/discretionary_policy/discretionary_policy_1.m b/matlab/discretionary_policy/discretionary_policy_1.m
index f8c60afadc8459dabb9c5bc1c7bbdf8e5fcaa2f0..30f7ca0480c634c2134cbf561b76f1b20b8872a9 100644
--- a/matlab/discretionary_policy/discretionary_policy_1.m
+++ b/matlab/discretionary_policy/discretionary_policy_1.m
@@ -49,7 +49,7 @@ else
 end
 params=M_.params;
 
-[U,Uy,W] = feval([M_.fname,'.objective.static'],zeros(M_.endo_nbr,1),[], M_.params);
+[~,Uy,W] = feval([M_.fname,'.objective.static'],zeros(M_.endo_nbr,1),[], M_.params);
 if any(any(isnan(Uy)))
     info = 64 ; %the derivatives of the objective function contain NaN
     return;
diff --git a/matlab/distributions/mode_and_variance_to_mean.m b/matlab/distributions/mode_and_variance_to_mean.m
index aabfb4ae6d140011782907a168b77a75cb0991ea..e4a892829bb2d50eb8d7270a3e22c3d837fafcdf 100644
--- a/matlab/distributions/mode_and_variance_to_mean.m
+++ b/matlab/distributions/mode_and_variance_to_mean.m
@@ -121,7 +121,7 @@ if (distribution==3)% Inverted Gamma 1 distribution
         nu = 2;
         s  = 3*(m*m);
     else
-        [mu, parameters] = mode_and_variance_to_mean(m,s2,2);
+        [~, parameters] = mode_and_variance_to_mean(m,s2,2);
         nu = sqrt(parameters(1));
         nu2 = 2*nu;
         nu1 = 2;
diff --git a/matlab/distributions/rand_inverse_wishart.m b/matlab/distributions/rand_inverse_wishart.m
index 59a5845e38b96e3d3c9829adfb5d559074a25687..13def973208ab4f7f7e19b94e8cfcd7b4e435ced 100644
--- a/matlab/distributions/rand_inverse_wishart.m
+++ b/matlab/distributions/rand_inverse_wishart.m
@@ -48,6 +48,6 @@ X = randn(v, m) * H_inv_upper_chol;
 % G = inv(X'*X);
 
 % Rather compute inv(X'*X) using the SVD
-[U,S,V] = svd(X, 0);
+[~,S,V] = svd(X, 0);
 SSi = 1 ./ (diag(S) .^ 2);
 G = (V .* repmat(SSi', m, 1)) * V';
\ No newline at end of file
diff --git a/matlab/distributions/weibull_specification.m b/matlab/distributions/weibull_specification.m
index 0c5a9b128402bca9a796660a6545d362a3a661be..8330f8815d0d00c94232ca791fb0983b26bcf2a7 100644
--- a/matlab/distributions/weibull_specification.m
+++ b/matlab/distributions/weibull_specification.m
@@ -56,7 +56,7 @@ eqn = @(k) gammaln(1+2./k) - 2*gammaln(1+1./k) - log(1+sigma2/mu2);
 eqn2 = @(k) eqn(k).*eqn(k);
 
 kstar = fminbnd(eqn2, 1e-9, 100);
-[shape, fval, exitflag]  = fzero(eqn, kstar);
+[shape, ~, exitflag]  = fzero(eqn, kstar);
 
 if exitflag<1
     shape = NaN;
diff --git a/matlab/dseries b/matlab/dseries
index 65ea12c7a5c0ac29ba42428e0ec5d991b6dd7f5e..4fceb85e3a1d10e8040b60bcb5cf79a3c533956c 160000
--- a/matlab/dseries
+++ b/matlab/dseries
@@ -1 +1 @@
-Subproject commit 65ea12c7a5c0ac29ba42428e0ec5d991b6dd7f5e
+Subproject commit 4fceb85e3a1d10e8040b60bcb5cf79a3c533956c
diff --git a/matlab/dynare.m b/matlab/dynare.m
index 84ff83a087b9a01d915477fb148472c5cab9e47e..7f13270a29ba5874807bf8e3802923461b5b084e 100644
--- a/matlab/dynare.m
+++ b/matlab/dynare.m
@@ -240,7 +240,7 @@ end
 % https://forum.dynare.org/t/issue-with-dynare-preprocessor-4-6-1/15448/1
 if ~fast
     if ispc && ~isoctave && exist(['+',fname(1:end-4)],'dir')
-        [~,~]=rmdir(['+', fname(1:end-4)],'s');
+        rmdir(['+', fname(1:end-4)],'s');
     end
 end
 
diff --git a/matlab/endogenous_prior_restrictions.m b/matlab/endogenous_prior_restrictions.m
index 8d045cf43a2db6ade3e80c6ec1f7720b7d3a8546..ef1987ff928a1974b87a7bbf8e527359446c38dd 100644
--- a/matlab/endogenous_prior_restrictions.m
+++ b/matlab/endogenous_prior_restrictions.m
@@ -51,7 +51,7 @@ if ~isempty(endo_prior_restrictions.irf)
     end
     varlist = M_.endo_names(dr.order_var);
     if isempty(T)
-        [T,R,SteadyState,infox,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
+        [T,R,~,~,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
     else % check if T and R are given in the restricted form!!!
         if size(T,1)<length(varlist)
             varlist = varlist(dr.restrict_var_list);
@@ -65,7 +65,7 @@ if ~isempty(endo_prior_restrictions.irf)
         end
         if ~varlistok
             varlist = M_.endo_names(dr.order_var);
-            [T,R,SteadyState,infox,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
+            [T,R,~,~,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
         end
     end
     NT=1;
@@ -138,7 +138,7 @@ if ~isempty(endo_prior_restrictions.moment)
         end
     end
     options_.ar = max(abs(NTmin),NTmax);
-    [gamma_y,stationary_vars] = th_autocovariances(dr, ivar, M_, options_,1);
+    gamma_y = th_autocovariances(dr, ivar, M_, options_,1);
     for t=NTmin:NTmax
         RR = gamma_y{abs(t)+1};
         if t==0
@@ -171,4 +171,4 @@ if ~isempty(endo_prior_restrictions.moment)
     if any(info_moment)
         info=[49, info(2) + sum(info_moment(:,2))];
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/ep/ep_accuracy_check.m b/matlab/ep/ep_accuracy_check.m
index f973c2cd8e04c083f8a9588ab92646bbced695c7..2feb42e1e07f91bd866ca5946d4127ddf082c063 100644
--- a/matlab/ep/ep_accuracy_check.m
+++ b/matlab/ep/ep_accuracy_check.m
@@ -20,14 +20,14 @@ function e = ep_accuracy_check(M_,options_,oo_)
 
 endo_simul = oo_.endo_simul;
 n = size(endo_simul,2);
-[initialconditions, innovations, pfm, ep, verbosity, options_, oo_] = ...
+[~, innovations, pfm, ~, ~, options_, oo_] = ...
     extended_path_initialization([], n-1, [], options_, M_, oo_);
 
 options_.ep.accuracy.stochastic.order = options_.ep.stochastic.order;
-[nodes,weights,nnodes] = setup_integration_nodes(options_.ep.accuracy,pfm);
+[nodes,weights] = setup_integration_nodes(options_.ep.accuracy,pfm);
 
 e = zeros(M_.endo_nbr,n);
 for i=1:n
     e(:,i) = euler_equation_error(endo_simul(:,i),oo_.exo_simul, ...
                                   innovations, M_, options_,oo_,pfm,nodes,weights);
-end
\ No newline at end of file
+end
diff --git a/matlab/ep/euler_equation_error.m b/matlab/ep/euler_equation_error.m
index bd735522b6996b3d3e334aad0de70aeb9bf51410..3bb304a5697dfe7a67b054bb6ed844da3960deb0 100644
--- a/matlab/ep/euler_equation_error.m
+++ b/matlab/ep/euler_equation_error.m
@@ -20,31 +20,29 @@ function e = euler_equation_error(y0,x,innovations,M_,options_,oo_,pfm,nodes,wei
 
 dynamic_model = str2func([M_.fname '.dynamic']);
 ep = options_.ep;
-[y1, info_convergence, endogenousvariablespaths] = extended_path_core(ep.periods, ...
-                                                  M_.endo_nbr, M_.exo_nbr, ...
-                                                  innovations.positive_var_indx, ...
-                                                  x, ep.init, y0, oo_.steady_state, ...
-                                                  0, ...
-                                                  ep.stochastic.order, M_, ...
-                                                  pfm, ep.stochastic.algo, ...
-                                                  ep.solve_algo, ...
-                                                  ep.stack_solve_algo, ...
-                                                  options_.lmmcp, options_, oo_, ...
-                                                  []);
+y1 = extended_path_core(ep.periods, ...
+                        M_.endo_nbr, M_.exo_nbr, ...
+                        innovations.positive_var_indx, ...
+                        x, ep.init, y0, oo_.steady_state, ...
+                        0, ...
+                        ep.stochastic.order, M_, ...
+                        pfm, ep.stochastic.algo, ...
+                        ep.solve_algo, ...
+                        ep.stack_solve_algo, ...
+                        options_.lmmcp, options_, oo_, ...
+                        []);
 i_pred = find(M_.lead_lag_incidence(1,:));
 i_fwrd = find(M_.lead_lag_incidence(3,:));
 x1 = [x(2:end,:); zeros(1,M_.exo_nbr)];
 for i=1:length(nodes)
     x2 = x1;
     x2(2,:) = x2(2,:) + nodes(i,:);
-    [y2, info_convergence, endogenousvariablespaths] = ...
-        extended_path_core(ep.periods, M_.endo_nbr, M_.exo_nbr, ...
-                           innovations.positive_var_indx, x2, ep.init, ...
-                           y1, oo_.steady_state, 0, ...
-                           ep.stochastic.order, M_, pfm, ep.stochastic.algo, ...
-                           ep.solve_algo, ep.stack_solve_algo, options_.lmmcp, ...
-                           options_, oo_, []);
-
+    y2 = extended_path_core(ep.periods, M_.endo_nbr, M_.exo_nbr, ...
+                            innovations.positive_var_indx, x2, ep.init, ...
+                            y1, oo_.steady_state, 0, ...
+                            ep.stochastic.order, M_, pfm, ep.stochastic.algo, ...
+                            ep.solve_algo, ep.stack_solve_algo, options_.lmmcp, ...
+                            options_, oo_, []);
     z = [y0(i_pred); y1; y2(i_fwrd)];
     res(:,i) = dynamic_model(z,x,M_.params,oo_.steady_state,2);
 end
diff --git a/matlab/ep/extended_path_initialization.m b/matlab/ep/extended_path_initialization.m
index 54daed1186a96b52316e25db71140a42f68a7880..ea4f235c82a0deb45bda03983933b480825214cc 100644
--- a/matlab/ep/extended_path_initialization.m
+++ b/matlab/ep/extended_path_initialization.m
@@ -78,7 +78,7 @@ dr = struct();
 if ep.init
     options_.order = 1;
     oo_.dr=set_state_space(dr,M_);
-    [oo_.dr,Info,M_.params] = resol(0,M_,options_,oo_.dr,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+    [oo_.dr,~,M_.params] = resol(0,M_,options_,oo_.dr,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
 end
 
 % Do not use a minimal number of perdiods for the perfect foresight solver (with bytecode and blocks)
diff --git a/matlab/ep/extended_path_mc.m b/matlab/ep/extended_path_mc.m
index 43eebd994873da28a649c6f1d11ab9da1be2beb9..fcfe284b54c0270eea6db5b74443ca5db2a81a73 100644
--- a/matlab/ep/extended_path_mc.m
+++ b/matlab/ep/extended_path_mc.m
@@ -36,7 +36,7 @@ function Simulations = extended_path_mc(initialconditions, samplesize, replic, e
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[initialconditions, innovations, pfm, ep, verbosity, options_, oo_] = ...
+[initialconditions, innovations, pfm, ep, ~, options_, oo_] = ...
     extended_path_initialization(initialconditions, samplesize, exogenousvariables, options_, M_, oo_);
 
 % Check the dimension of the first input argument
diff --git a/matlab/estimation/PosteriorIRF.m b/matlab/estimation/PosteriorIRF.m
index c07fc92313968cc3fe6f81d6fb2c5b42b714ed23..0bd16e606525a8dfeb4fd6787b88b9ae8737154c 100644
--- a/matlab/estimation/PosteriorIRF.m
+++ b/matlab/estimation/PosteriorIRF.m
@@ -200,7 +200,7 @@ if isnumeric(options_.parallel)
     nosaddle = fout.nosaddle;
 else
     % Parallel execution!
-    [nCPU, totCPU, nBlockPerCPU] = distributeJobs(options_.parallel, 1, B);
+    [~, totCPU, nBlockPerCPU] = distributeJobs(options_.parallel, 1, B);
     for j=1:totCPU-1
         nfiles = ceil(nBlockPerCPU(j)/MAX_nirfs_dsge);
         NumberOfIRFfiles_dsge(j+1) =NumberOfIRFfiles_dsge(j)+nfiles;
diff --git a/matlab/estimation/dsge_likelihood.m b/matlab/estimation/dsge_likelihood.m
index 50c7e8585cc35c19429489142c6ccc8f14e1a64a..114edec8a9f28acfb06429aee6374f3c7ea97292 100644
--- a/matlab/estimation/dsge_likelihood.m
+++ b/matlab/estimation/dsge_likelihood.m
@@ -464,7 +464,7 @@ if analytic_derivation
     AHess = [];
     iv = dr.restrict_var_list;
     if nargin<13 || isempty(derivatives_info)
-        [A,B,nou,nou,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
+        [~,~,~,~,dr, M_.params] = dynare_resolve(M_,options_,dr, endo_steady_state, exo_steady_state, exo_det_steady_state);
         if ~isempty(estim_params_.var_exo)
             indexo=estim_params_.var_exo(:,1);
         else
diff --git a/matlab/estimation/dsge_simulated_theoretical_conditional_variance_decomposition.m b/matlab/estimation/dsge_simulated_theoretical_conditional_variance_decomposition.m
index e70663ac4abc92f5876362d552c6fc09ad99f30d..6e052b6ec6c506af9d502c62a11188088b793920 100644
--- a/matlab/estimation/dsge_simulated_theoretical_conditional_variance_decomposition.m
+++ b/matlab/estimation/dsge_simulated_theoretical_conditional_variance_decomposition.m
@@ -81,7 +81,7 @@ options_.ar = 0;
 NumberOfSavedElementsPerSimulation = nvar*M_.exo_nbr*length(Steps);
 MaXNumberOfConditionalDecompLines = ceil(options_.MaxNumberOfBytes/NumberOfSavedElementsPerSimulation/8);
 
-[ME_present,observable_pos_requested_vars,index_subset,index_observables]=check_measurement_error_requested_vars(M_,options_,ivar);
+[ME_present,observable_pos_requested_vars] = check_measurement_error_requested_vars(M_,options_,ivar);
 
 if ME_present && ~isempty(observable_pos_requested_vars)
     nobs_ME=length(observable_pos_requested_vars);
@@ -132,7 +132,7 @@ for file = 1:NumberOfDrawsFiles
             dr = temp.pdraws{linee,2};
         else
             M_=set_parameters_locally(M_,temp.pdraws{linee,1});
-			[dr,info,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+			[dr,~,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
         end
 
         M_ = set_measurement_errors(temp.pdraws{linee,1},temp.estim_params_,M_);
diff --git a/matlab/estimation/dsge_simulated_theoretical_correlation.m b/matlab/estimation/dsge_simulated_theoretical_correlation.m
index 7db0ae3646528f0e36e662e8abf0cf8db095d923..9dc90c44e8f2077c2a5e9ce84621278bcb343ea5 100644
--- a/matlab/estimation/dsge_simulated_theoretical_correlation.m
+++ b/matlab/estimation/dsge_simulated_theoretical_correlation.m
@@ -36,7 +36,7 @@ function [nvar,vartan,CorrFileNumber] = dsge_simulated_theoretical_correlation(S
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 nvar = length(ivar);
-[ivar,vartan, options_] = get_variables_list(options_, M_);
+[~,vartan] = get_variables_list(options_, M_);
 
 % Get informations about the _posterior_draws files.
 if strcmpi(type,'posterior')
diff --git a/matlab/estimation/dsge_simulated_theoretical_covariance.m b/matlab/estimation/dsge_simulated_theoretical_covariance.m
index 93d9e1fef045992b9c929bae3098199ee546e5c8..62cbcabd495687e8b2b7955427d882c84e3d27c8 100644
--- a/matlab/estimation/dsge_simulated_theoretical_covariance.m
+++ b/matlab/estimation/dsge_simulated_theoretical_covariance.m
@@ -127,7 +127,7 @@ for file = 1:NumberOfDrawsFiles
             dr = temp.pdraws{linee,2};
         else
             M_=set_parameters_locally(M_,temp.pdraws{linee,1});
-            [dr,info,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+            [dr,~,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
         end
         if ~options_.pruning
             tmp = th_autocovariances(dr,ivar,M_,options_,nodecomposition);
diff --git a/matlab/estimation/dsge_simulated_theoretical_variance_decomposition.m b/matlab/estimation/dsge_simulated_theoretical_variance_decomposition.m
index 38863c0ac038aad84e3588d1d6e01cd80ed978c2..57e6508cdb4faa1b90a5f8238effcab843f0d707 100644
--- a/matlab/estimation/dsge_simulated_theoretical_variance_decomposition.m
+++ b/matlab/estimation/dsge_simulated_theoretical_variance_decomposition.m
@@ -137,10 +137,10 @@ for file = 1:NumberOfDrawsFiles
             dr = temp.pdraws{linee,2};
         else
             M_=set_parameters_locally(M_,temp.pdraws{linee,1});
-            [dr,info,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+            [dr,~,M_.params] = compute_decision_rules(M_,options_,oo_.dr, oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
         end
         if file==1 && linee==1
-            [tmp, stationary_vars] = th_autocovariances(dr,ivar,M_,options_,nodecomposition);
+            [~, stationary_vars] = th_autocovariances(dr,ivar,M_,options_,nodecomposition);
             if isempty(stationary_vars)
                 fprintf('\ndsge_simulated_theoretical_variance_decomposition:: All requested endogenous variables have a unit root and thus infinite variance.\n')
                 fprintf('dsge_simulated_theoretical_variance_decomposition:: No decomposition is performed.\n')
diff --git a/matlab/estimation/dynare_estimation_1.m b/matlab/estimation/dynare_estimation_1.m
index 05e8f949556e759510c692b611eb5b53d258ef71..56904f5be5b3e9549f03cebd2133014ce8dc3bb6 100644
--- a/matlab/estimation/dynare_estimation_1.m
+++ b/matlab/estimation/dynare_estimation_1.m
@@ -201,14 +201,14 @@ if isequal(options_.mode_compute,0) && isempty(options_.mode_file) && ~options_.
                 end
             else
                 if options_.occbin.smoother.status
-                    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = occbin.DSGE_smoother(xparam1,gend,transpose(data),data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,dataset_,dataset_info);
+                    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = occbin.DSGE_smoother(xparam1,gend,transpose(data),data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_,dataset_,dataset_info);
                     if oo_.occbin.smoother.error_flag(1)==0
                         [oo_]=store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
                     else
                         fprintf('\nOccbin: smoother did not succeed. No results will be written to oo_.\n')
                     end
                 else
-                    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = DsgeSmoother(xparam1,gend,transpose(data),data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_);
+                    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = DsgeSmoother(xparam1,gend,transpose(data),data_index,missing_value,M_,oo_,options_,bayestopt_,estim_params_);
                     [oo_]=store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
                 end
             end
@@ -437,7 +437,7 @@ if issmc(options_) || (any(bayestopt_.pshape>0) && options_.mh_replic) ||  (any(
         CutSample(M_, options_, dispString);
     end
     if options_.mh_posterior_mode_estimation || (issmc(options_) && options_.smc_posterior_mode_estimation)
-        [~, covariance, posterior_mode, ~] = compute_posterior_covariance_matrix(bayestopt_.name, M_.fname, M_.dname, options_);
+        [~, covariance, posterior_mode] = compute_posterior_covariance_matrix(bayestopt_.name, M_.fname, M_.dname, options_);
         oo_ = fill_mh_mode(posterior_mode, sqrt(diag(covariance)), M_, options_, estim_params_, oo_, 'posterior');
         %reset qz_criterium
         options_.qz_criterium = qz_criterium_old;
diff --git a/matlab/estimation/dynare_estimation_init.m b/matlab/estimation/dynare_estimation_init.m
index 8238f7a9a4e2cedc713f0b49fc167ce0f67f622f..789af16fd1b6cd45f660e0bf2c1bdf83c7b2ab85 100644
--- a/matlab/estimation/dynare_estimation_init.m
+++ b/matlab/estimation/dynare_estimation_init.m
@@ -345,7 +345,7 @@ if options_.analytic_derivation
         else
             steadystate_check_flag = 1;
         end
-        [tmp1, params] = evaluate_steady_state(oo_.steady_state,[oo_.exo_steady_state; oo_.exo_det_steady_state],M_local,options_,steadystate_check_flag);
+        [~, params] = evaluate_steady_state(oo_.steady_state,[oo_.exo_steady_state; oo_.exo_det_steady_state],M_local,options_,steadystate_check_flag);
         change_flag=any(find(params-M_local.params));
         if change_flag
             skipline()
@@ -374,7 +374,7 @@ bayestopt_.jscale(k) = options_.mh_jscale;
 
 % Build the dataset
 if ~isempty(options_.datafile)
-    [pathstr,name,ext] = fileparts(options_.datafile);
+    [~,name] = fileparts(options_.datafile);
     if strcmp(name,M_.fname)
         error('Data-file and mod-file are not allowed to have the same name. Please change the name of the data file.')
     end
@@ -383,7 +383,7 @@ end
 if isnan(options_.first_obs)
     options_.first_obs=1;
 end
-[dataset_, dataset_info, newdatainterfaceflag] = makedataset(options_, options_.dsge_var*options_.dsge_varlag, gsa_flag);
+[dataset_, dataset_info] = makedataset(options_, options_.dsge_var*options_.dsge_varlag, gsa_flag);
 
 %set options for old interface from the ones for new interface
 if ~isempty(dataset_)
diff --git a/matlab/estimation/execute_prior_posterior_function.m b/matlab/estimation/execute_prior_posterior_function.m
index 520f02c211c2a8b22adc6836629784cd08dfc833..f23d4664babdb2bb206837bc6671e1a6abeb5253 100644
--- a/matlab/estimation/execute_prior_posterior_function.m
+++ b/matlab/estimation/execute_prior_posterior_function.m
@@ -34,7 +34,7 @@ function oo_=execute_prior_posterior_function(posterior_function_name,M_,options
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[directory,basename,extension] = fileparts(posterior_function_name);
+[~,basename,extension] = fileparts(posterior_function_name);
 if isempty(extension)
     extension = '.m';
 end
@@ -64,7 +64,7 @@ elseif strcmpi(type,'prior')
     % Get informations about the prior distribution.
     if isempty(bayestopt_)
         if ~isempty(estim_params_) && ~(isfield(estim_params_,'nvx') && (size(estim_params_.var_exo,1)+size(estim_params_.var_endo,1)+size(estim_params_.corrx,1)+size(estim_params_.corrn,1)+size(estim_params_.param_vals,1))==0)
-            [xparam1,estim_params_,bayestopt_,lb,ub,M_] = set_prior(estim_params_,M_,options_);
+            [~,estim_params_,bayestopt_,~,~,M_] = set_prior(estim_params_,M_,options_);
         else
             error('The prior distributions are not properly set up.')
         end
diff --git a/matlab/estimation/initial_estimation_checks.m b/matlab/estimation/initial_estimation_checks.m
index fd03649ee619b4c76bf12ec59974ce40ce7b25e2..88e6014fc105e24b0c08a7b7c85e1e6a92110cb2 100644
--- a/matlab/estimation/initial_estimation_checks.m
+++ b/matlab/estimation/initial_estimation_checks.m
@@ -172,7 +172,7 @@ end
 % display warning if some parameters are still NaN
 test_for_deep_parameters_calibration(M_);
 
-[lnprior,~,~,info]= priordens(xparam1,bayestopt_.pshape,bayestopt_.p6,bayestopt_.p7,bayestopt_.p3,bayestopt_.p4);
+[~,~,~,info]= priordens(xparam1,bayestopt_.pshape,bayestopt_.p6,bayestopt_.p7,bayestopt_.p3,bayestopt_.p4);
 if any(info)
     fprintf('The prior density evaluated at the initial values is Inf for the following parameters: %s\n',bayestopt_.name{info,1})
     error('The initial value of the prior is -Inf')
diff --git a/matlab/estimation/maximize_prior_density.m b/matlab/estimation/maximize_prior_density.m
index 5c47e1a3c600e59383ea523498999c43ba6bde1e..ef6babca55e969b9326685b5b4100bdcb61da41c 100644
--- a/matlab/estimation/maximize_prior_density.m
+++ b/matlab/estimation/maximize_prior_density.m
@@ -46,7 +46,7 @@ function [xparams, lpd, hessian_mat] = ...
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[xparams, lpd, exitflag, hessian_mat] = dynare_minimize_objective('minus_logged_prior_density', ...
+[xparams, lpd, ~, hessian_mat] = dynare_minimize_objective('minus_logged_prior_density', ...
                                                                   iparams, ...
                                                                   options_.mode_compute, ...
                                                                   options_, ...
diff --git a/matlab/estimation/model_comparison.m b/matlab/estimation/model_comparison.m
index ac27031d3ba8759eb8b6d21fc483f4393ef7f05d..e39ec2b4bf4d8759ef0bbd32ef2db3e0a4982227 100644
--- a/matlab/estimation/model_comparison.m
+++ b/matlab/estimation/model_comparison.m
@@ -110,7 +110,7 @@ end
 % In order to avoid overflow, we divide the numerator and the denominator
 % of the Posterior Odds Ratio by the largest Marginal Posterior Density
 lmpd = log(ModelPriors)+MarginalLogDensity;
-[maxval,k] = max(lmpd);
+maxval = max(lmpd);
 elmpd = exp(lmpd-maxval);
 
 % Now I display the posterior probabilities.
diff --git a/matlab/estimation/optimize_prior.m b/matlab/estimation/optimize_prior.m
index ace2283e6d1947e7e44e0a6e7b7ed85720b6c879..9e57741b4020fccac3b4f0cde553f149f859bfc3 100644
--- a/matlab/estimation/optimize_prior.m
+++ b/matlab/estimation/optimize_prior.m
@@ -37,7 +37,7 @@ while look_for_admissible_initial_condition
     xinit = xparam1+scale*randn(size(xparam1));
     if all(xinit>Prior.p3) && all(xinit<Prior.p4)
         M_ = set_all_parameters(xinit, estim_params_, M_);
-        [dr, INFO, M_, oo_] = resol(0, M_, options_, oo_);
+        [~, INFO, M_, oo_] = resol(0, M_, options_, oo_);
         if ~INFO(1)
             look_for_admissible_initial_condition = false;
         end
@@ -52,8 +52,7 @@ while look_for_admissible_initial_condition
 end
 
 % Maximization of the prior density
-[xparams, lpd, hessian_mat] = ...
-    maximize_prior_density(xinit, pnames, options_, M_, Prior, estim_params_, oo_);
+xparams = maximize_prior_density(xinit, pnames, options_, M_, Prior, estim_params_, oo_);
 
 % Display results.
 skipline(2)
diff --git a/matlab/estimation/pm3.m b/matlab/estimation/pm3.m
index d08b068304bc008e0b55d1867adc2eccf0e77faa..165fdb1ad133e6798adacea0ceaa32c1dc6d42aa 100644
--- a/matlab/estimation/pm3.m
+++ b/matlab/estimation/pm3.m
@@ -336,7 +336,7 @@ if ~options_.nograph && ~options_.no_graph.posterior
                 fout = pm3_core(localVars,1,nvar,0);
             else
                 globalVars = [];
-                [fout, nvar0, totCPU] = masterParallel(options_.parallel, 1, nvar, [],'pm3_core', localVars,globalVars, options_.parallel_info);
+                [~, nvar0] = masterParallel(options_.parallel, 1, nvar, [],'pm3_core', localVars,globalVars, options_.parallel_info);
             end
         end
     else
diff --git a/matlab/estimation/posterior_analysis.m b/matlab/estimation/posterior_analysis.m
index 2bc14e27c8f1d5ad0243bc38fa0a53c4396b0f0b..a0487d21c9c252bb5d7ad4572bb30a46731987d8 100644
--- a/matlab/estimation/posterior_analysis.m
+++ b/matlab/estimation/posterior_analysis.m
@@ -52,14 +52,14 @@ end
 switch type
   case 'variance'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [nvar,vartan] = ...
             dsge_simulated_theoretical_covariance(SampleSize,arg3,M_,options_,oo_,'posterior');
     end
     oo_ = covariance_mc_analysis(SampleSize,'posterior',M_.dname,M_.fname,...
                                  vartan,nvar,arg1,arg2,options_.mh_conf_sig,oo_,options_);
   case 'decomposition'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [~,vartan] = ...
             dsge_simulated_theoretical_variance_decomposition(SampleSize,M_,options_,oo_,'posterior');
     end
     oo_ = variance_decomposition_mc_analysis(SampleSize,'posterior',M_.dname,M_.fname,...
@@ -77,14 +77,14 @@ switch type
     end
   case 'correlation'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [nvar,vartan] = ...
             dsge_simulated_theoretical_correlation(SampleSize,arg3,M_,options_,oo_,'posterior');
     end
     oo_ = correlation_mc_analysis(SampleSize,'posterior',M_.dname,M_.fname,...
                                   vartan,nvar,arg1,arg2,arg3,options_.mh_conf_sig,oo_,M_,options_);
   case 'conditional decomposition'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [~,vartan] = ...
             dsge_simulated_theoretical_conditional_variance_decomposition(SampleSize,arg3,M_,options_,oo_,'posterior');
     end
     oo_ = conditional_variance_decomposition_mc_analysis(SampleSize,'posterior',M_.dname,M_.fname,...
diff --git a/matlab/estimation/posterior_sampler_iteration.m b/matlab/estimation/posterior_sampler_iteration.m
index 7a9dc61750a097fb3a8eb59ad9be0f1ea0fdcc7f..efb6e147054788e837ac2b383f5ae49abeb14eae 100644
--- a/matlab/estimation/posterior_sampler_iteration.m
+++ b/matlab/estimation/posterior_sampler_iteration.m
@@ -91,9 +91,9 @@ switch posterior_sampling_method
         blocked_draws_counter=blocked_draws_counter+1;
         nxopt=length(indices(blocks==block_iter,1)); %get size of current block
         par_start_current_block=current_draw(indices(blocks==block_iter,1));
-        [xopt_current_block, fval, exitflag, hess_mat_optimizer, options_, Scale] = dynare_minimize_objective(@TaRB_optimizer_wrapper,par_start_current_block,sampler_options.mode_compute,options_,[mh_bounds.lb(indices(blocks==block_iter,1),1) mh_bounds.ub(indices(blocks==block_iter,1),1)],bayestopt_.name,bayestopt_,[],...
-                                                          current_draw,indices(blocks==block_iter,1),TargetFun,...% inputs for wrapper
-                                                          varargin{:}); %inputs for objective
+        [xopt_current_block, ~, ~, ~, options_] = dynare_minimize_objective(@TaRB_optimizer_wrapper,par_start_current_block,sampler_options.mode_compute,options_,[mh_bounds.lb(indices(blocks==block_iter,1),1) mh_bounds.ub(indices(blocks==block_iter,1),1)],bayestopt_.name,bayestopt_,[],...
+                                                                            current_draw,indices(blocks==block_iter,1),TargetFun,...% inputs for wrapper
+                                                                            varargin{:}); %inputs for objective
         %% covariance for proposal density
         hessian_mat = reshape(hessian('TaRB_optimizer_wrapper',xopt_current_block, ...
                                       options_.gstep,...
@@ -186,4 +186,4 @@ switch posterior_sampling_method
         par = last_draw;
         logpost = last_posterior;
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/estimation/print_table_prior.m b/matlab/estimation/print_table_prior.m
index b9253af4c7e9a2facb5d1bd6bbfd62dd078f9ae9..fa7220a1a27022c99d6d672c835447db3832e348 100644
--- a/matlab/estimation/print_table_prior.m
+++ b/matlab/estimation/print_table_prior.m
@@ -52,7 +52,7 @@ options_.prior_trunc = prior_trunc_backup ;
 RESIZE = false;
 
 for i=1:size(bayestopt_.name,1)
-    [Name,~] = get_the_name(i,1,M_,estim_params_,options_.varobs);
+    Name = get_the_name(i,1,M_,estim_params_,options_.varobs);
     if length(Name)>size(T1,2)
         resize = true;
     else
@@ -163,4 +163,4 @@ if ~isnumeric(UpperBound)
 else
     format_string = [ format_string , ' %6.4f \t'];
 end
-format_string = [ format_string , ' %6.4f \t %6.4f'];
\ No newline at end of file
+format_string = [ format_string , ' %6.4f \t %6.4f'];
diff --git a/matlab/estimation/prior_analysis.m b/matlab/estimation/prior_analysis.m
index b41f4c53a11ae5f7f96ad053a3100d5336b8dfcd..51c6805d4f6727e150a68f92384c217220a4b894 100644
--- a/matlab/estimation/prior_analysis.m
+++ b/matlab/estimation/prior_analysis.m
@@ -52,28 +52,28 @@ end
 switch type
   case 'variance'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [nvar,vartan] = ...
             dsge_simulated_theoretical_covariance(SampleSize,M_,options_,oo_,'prior');
     end
     oo_ = covariance_mc_analysis(SampleSize,'prior',M_.dname,M_.fname,...
                                  vartan,nvar,arg1,arg2,options_.mh_conf_sig,oo_,options_);
   case 'decomposition'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [~,vartan] = ...
             dsge_simulated_theoretical_variance_decomposition(SampleSize,M_,options_,oo_,'prior');
     end
     oo_ = variance_decomposition_mc_analysis(SampleSize,'prior',M_.dname,M_.fname,...
                                              M_.exo_names,arg2,vartan,arg1,options_.mh_conf_sig,oo_,options_);
   case 'correlation'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [nvar,vartan] = ...
             dsge_simulated_theoretical_correlation(SampleSize,arg3,M_,options_,oo_,'prior');
     end
     oo_ = correlation_mc_analysis(SampleSize,'prior',M_.dname,M_.fname,...
                                   vartan,nvar,arg1,arg2,arg3,options_.mh_conf_sig,oo_,M_,options_);
   case 'conditional decomposition'
     if nargin==narg1
-        [nvar,vartan,NumberOfFiles] = ...
+        [~,vartan] = ...
             dsge_simulated_theoretical_conditional_variance_decomposition(SampleSize,arg3,M_,options_,oo_,'prior');
     end
     oo_ = conditional_variance_decomposition_mc_analysis(SampleSize,'prior',M_.dname,M_.fname,...
@@ -86,4 +86,4 @@ switch type
     end
   otherwise
     disp('Not yet implemented')
-end
\ No newline at end of file
+end
diff --git a/matlab/estimation/prior_posterior_statistics_core.m b/matlab/estimation/prior_posterior_statistics_core.m
index b98d1e950d6b6f7ca1d6cc8e63385af7def29ea4..2f112308addb5dca70385fb35b4bf1ade6388333 100644
--- a/matlab/estimation/prior_posterior_statistics_core.m
+++ b/matlab/estimation/prior_posterior_statistics_core.m
@@ -368,7 +368,7 @@ for b=fpar:B
             stock_smoothed_uncert(dr.order_var,dr.order_var,:,irun(13)) = state_uncertainty;
         end
     else
-        [~,~,SteadyState,info] = dynare_resolve(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
+        [~,~,SteadyState] = dynare_resolve(M_,options_,oo_.dr,oo_.steady_state,oo_.exo_steady_state,oo_.exo_det_steady_state);
     end
     stock_param(irun(5),:) = deep;
     stock_logpo(irun(5),1) = logpo;
@@ -594,4 +594,4 @@ for iter=1:horizon
 end
 
 yf(dr.order_var,:,:) = yf;
-yf=permute(yf,[2 1 3]);
\ No newline at end of file
+yf=permute(yf,[2 1 3]);
diff --git a/matlab/estimation/priordens.m b/matlab/estimation/priordens.m
index f68376794b55ef4e49dcad1e24fa2c43ba59bbcc..f285606510ed9d0ffccc5ca8326360feac98a134 100644
--- a/matlab/estimation/priordens.m
+++ b/matlab/estimation/priordens.m
@@ -95,9 +95,9 @@ if tt1
         return
     end
     if nargout == 2
-        [tmp, dlprior(id1)]=lpdfgbeta(x(id1),p6(id1),p7(id1),p3(id1),p4(id1));
+        [~, dlprior(id1)]=lpdfgbeta(x(id1),p6(id1),p7(id1),p3(id1),p4(id1));
     elseif nargout == 3
-        [tmp, dlprior(id1), d2lprior(id1)]=lpdfgbeta(x(id1),p6(id1),p7(id1),p3(id1),p4(id1));
+        [~, dlprior(id1), d2lprior(id1)]=lpdfgbeta(x(id1),p6(id1),p7(id1),p3(id1),p4(id1));
     end
 end
 
@@ -110,18 +110,18 @@ if tt2
         return
     end
     if nargout == 2
-        [tmp, dlprior(id2)]=lpdfgam(x(id2)-p3(id2),p6(id2),p7(id2));
+        [~, dlprior(id2)]=lpdfgam(x(id2)-p3(id2),p6(id2),p7(id2));
     elseif nargout == 3
-        [tmp, dlprior(id2), d2lprior(id2)]=lpdfgam(x(id2)-p3(id2),p6(id2),p7(id2));
+        [~, dlprior(id2), d2lprior(id2)]=lpdfgam(x(id2)-p3(id2),p6(id2),p7(id2));
     end
 end
 
 if tt3
     logged_prior_density = logged_prior_density + sum(lpdfnorm(x(id3),p6(id3),p7(id3))) ;
     if nargout == 2
-        [tmp, dlprior(id3)]=lpdfnorm(x(id3),p6(id3),p7(id3));
+        [~, dlprior(id3)]=lpdfnorm(x(id3),p6(id3),p7(id3));
     elseif nargout == 3
-        [tmp, dlprior(id3), d2lprior(id3)]=lpdfnorm(x(id3),p6(id3),p7(id3));
+        [~, dlprior(id3), d2lprior(id3)]=lpdfnorm(x(id3),p6(id3),p7(id3));
     end
 end
 
@@ -134,9 +134,9 @@ if tt4
         return
     end
     if nargout == 2
-        [tmp, dlprior(id4)]=lpdfig1(x(id4)-p3(id4),p6(id4),p7(id4));
+        [~, dlprior(id4)]=lpdfig1(x(id4)-p3(id4),p6(id4),p7(id4));
     elseif nargout == 3
-        [tmp, dlprior(id4), d2lprior(id4)]=lpdfig1(x(id4)-p3(id4),p6(id4),p7(id4));
+        [~, dlprior(id4), d2lprior(id4)]=lpdfig1(x(id4)-p3(id4),p6(id4),p7(id4));
     end
 end
 
@@ -166,9 +166,9 @@ if tt6
         return
     end
     if nargout == 2
-        [tmp, dlprior(id6)]=lpdfig2(x(id6)-p3(id6),p6(id6),p7(id6));
+        [~, dlprior(id6)]=lpdfig2(x(id6)-p3(id6),p6(id6),p7(id6));
     elseif nargout == 3
-        [tmp, dlprior(id6), d2lprior(id6)]=lpdfig2(x(id6)-p3(id6),p6(id6),p7(id6));
+        [~, dlprior(id6), d2lprior(id6)]=lpdfig2(x(id6)-p3(id6),p6(id6),p7(id6));
     end
 end
 
@@ -181,9 +181,9 @@ if tt8
         return
     end
     if nargout==2
-        [tmp, dlprior(id8)] = lpdfgweibull(x(id8),p6(id8),p7(id8));
+        [~, dlprior(id8)] = lpdfgweibull(x(id8),p6(id8),p7(id8));
     elseif nargout==3
-        [tmp, dlprior(id8), d2lprior(id8)] = lpdfgweibull(x(id8),p6(id8),p7(id8));
+        [~, dlprior(id8), d2lprior(id8)] = lpdfgweibull(x(id8),p6(id8),p7(id8));
     end
 end
 
diff --git a/matlab/estimation/recursive_moments.m b/matlab/estimation/recursive_moments.m
index 93509e51d88d36b5eec344dec5372a67f9254cb0..cb1badcef570d67a9d46e89764cd9ec4bcd27dd7 100644
--- a/matlab/estimation/recursive_moments.m
+++ b/matlab/estimation/recursive_moments.m
@@ -36,7 +36,7 @@ function [mu,sigma,offset] = recursive_moments(m0,s0,data,offset)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[T,n] = size(data);
+[T,~] = size(data);
 
 for t = 1:T
     tt = t+offset;
diff --git a/matlab/estimation/smc/dsmh.m b/matlab/estimation/smc/dsmh.m
index 5382639d598176c03aba2bbf2dbef2c216bdd4ec..dea9cbdd5eae83f31c5ceda38086734a79a4ddf1 100644
--- a/matlab/estimation/smc/dsmh.m
+++ b/matlab/estimation/smc/dsmh.m
@@ -96,7 +96,7 @@ TeX = options_.TeX;
 
 str = sprintf(' Param. \t Lower Bound (95%%) \t Mean \t Upper Bound (95%%)');
 for l=1:npar
-    [name,~] = get_the_name(l,TeX,M_,estim_params_,options_.varobs);
+    name = get_the_name(l,TeX,M_,estim_params_,options_.varobs);
     str = sprintf('%s\n %s \t\t %5.4f \t\t %7.5f \t\t %5.4f', str, name, lb95_xparam(l), mean_xparam(l), ub95_xparam(l));
 end
 disp([str])
@@ -104,8 +104,6 @@ disp('')
 
 %% Plot parameters densities
 
-[nbplt,nr,nc,lr,lc,nstar] = pltorg(npar);
-
 if TeX
     fidTeX = fopen([M_.fname '_param_density.tex'],'w');
     fprintf(fidTeX,'%% TeX eps-loader file generated by DSMH.m (Dynare).\n');
diff --git a/matlab/get_ar_ec_matrices.m b/matlab/get_ar_ec_matrices.m
index fef4dc52a744c340da690bd4556314ec882b04dd..38ddaa049c579afa29582cbb45562c1e63dcb0ea 100644
--- a/matlab/get_ar_ec_matrices.m
+++ b/matlab/get_ar_ec_matrices.m
@@ -68,7 +68,7 @@ else
 end
 
 %% Call Dynamic Function
-[junk, g1] = feval([M_.fname '.dynamic'], ...
+[~, g1] = feval([M_.fname '.dynamic'], ...
     ones(max(max(M_.lead_lag_incidence)), 1), ...
     ones(1, M_.exo_nbr), ...
     M_.params, ...
@@ -159,12 +159,12 @@ for i = 1:length(lhs)
         if g1col ~= 0 && any(g1(:, g1col))
             if rhsvars{i}.arRhsIdxs(j) > 0
                 % Fill AR
-                [lag, ndiffs] = findLagForVar(var, -rhsvars{i}.lags(j), 0, lhs);
+                lag = findLagForVar(var, -rhsvars{i}.lags(j), 0, lhs);
                 oo_.(model_type).(model_name).ar(i, rhsvars{i}.arRhsIdxs(j), lag) = ...
                     oo_.(model_type).(model_name).ar(i, rhsvars{i}.arRhsIdxs(j), lag) + g1(i, g1col);
             elseif rhsvars{i}.ecRhsIdxs(j) > 0
                 % Fill EC
-                [lag, ndiffs] = findLagForVar(var, -rhsvars{i}.lags(j), 0, ecRhsVars);
+                lag = findLagForVar(var, -rhsvars{i}.lags(j), 0, ecRhsVars);
                 if lag==1
                     if size(oo_.(model_type).(model_name).ec, 3) < lag
                         oo_.(model_type).(model_name).ec(i, rhsvars{i}.ecRhsIdxs(j), lag) = 0;
diff --git a/matlab/get_file_extension.m b/matlab/get_file_extension.m
index cff2eaf0e34e70b76a0441ef755e5a836e644495..92d41fbcc4d93b8d258d4bb331a10e5fa50e16bf 100644
--- a/matlab/get_file_extension.m
+++ b/matlab/get_file_extension.m
@@ -28,7 +28,7 @@ function ext = get_file_extension(file)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[dir, fname, ext] = fileparts(file);
+[~, ~, ext] = fileparts(file);
 
 if ~isempty(ext)
     % Removes the leading dot.
diff --git a/matlab/histvalf_initvalf.m b/matlab/histvalf_initvalf.m
index 48d62d67454677073a66f419280e8af109fb13ec..7885ac4e82c30aa3bd1e095a03a022fe1bceb5f7 100644
--- a/matlab/histvalf_initvalf.m
+++ b/matlab/histvalf_initvalf.m
@@ -58,7 +58,7 @@ if isfield(options, 'datafile')
 end
 
 if datafile
-    [directory,basename,extension] = fileparts(datafile);
+    [~,basename,extension] = fileparts(datafile);
     % Auto-detect extension if not provided
     if isempty(extension)
         if exist([basename '.m'],'file')
diff --git a/matlab/kalman/evaluate_smoother.m b/matlab/kalman/evaluate_smoother.m
index be08936b8f51746fddd3ed570f14b3d490eb0005..f093701df8e428693db326b38b37710e69bc5026 100644
--- a/matlab/kalman/evaluate_smoother.m
+++ b/matlab/kalman/evaluate_smoother.m
@@ -118,11 +118,11 @@ if options_.occbin.smoother.status
             bayestopt_.mf = bayestopt_.smoother_var_list(bayestopt_.smoother_mf);
         end
     else
-        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = ...
+        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = ...
             occbin.DSGE_smoother(parameters,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_,dataset_,dataset_info);
     end
 else
-    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = ...
+    [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = ...
         DsgeSmoother(parameters,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_);
 end
 if ~(options_.occbin.smoother.status && options_.occbin.smoother.inversion_filter)
diff --git a/matlab/kalman/likelihood/kalman_filter.m b/matlab/kalman/likelihood/kalman_filter.m
index b85d3edada24e724618ce41793e2e148317f5b0b..369c086013495de1c3277cb56f2186223acba355 100644
--- a/matlab/kalman/likelihood/kalman_filter.m
+++ b/matlab/kalman/likelihood/kalman_filter.m
@@ -273,7 +273,7 @@ if t <= last
             Hess = Hess + tmp{3};
         end
     else
-        [tmp, likk(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log_dF, Z, pp, Zflag);
+        [~, likk(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log_dF, Z, pp, Zflag);
     end
 end
 
diff --git a/matlab/kalman/likelihood/kalman_filter_fast.m b/matlab/kalman/likelihood/kalman_filter_fast.m
index e7815f52ee33de9e5ddb0492e9e7f1068713d790..f457f42360da92bf972429d1de566efb38747567 100644
--- a/matlab/kalman/likelihood/kalman_filter_fast.m
+++ b/matlab/kalman/likelihood/kalman_filter_fast.m
@@ -218,7 +218,7 @@ if t <= last
             Hess = Hess + tmp{3};
         end
     else
-        [tmp, likk(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log(dF), Z, pp, Zflag);
+        [~, likk(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log(dF), Z, pp, Zflag);
     end
 end
 
diff --git a/matlab/kalman/likelihood/missing_observations_kalman_filter.m b/matlab/kalman/likelihood/missing_observations_kalman_filter.m
index bf4407e617d8b9802595b27a49a18c9cf9a78433..47a707b5511f7603dcb1676c1e1f1fd710ef95c0 100644
--- a/matlab/kalman/likelihood/missing_observations_kalman_filter.m
+++ b/matlab/kalman/likelihood/missing_observations_kalman_filter.m
@@ -294,7 +294,7 @@ lik(1:s) = .5*lik(1:s);
 
 % Call steady state Kalman filter if needed.
 if t<=last
-    [tmp, lik(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log_dF, Z, pp, Zflag);
+    [~, lik(s+1:end)] = kalman_filter_ss(Y, t, last, a, T, K, iF, log_dF, Z, pp, Zflag);
 end
 
 % Compute minus the log-likelihood.
diff --git a/matlab/kalman/likelihood/univariate_kalman_filter.m b/matlab/kalman/likelihood/univariate_kalman_filter.m
index 89de1da5d317812e0ae64ad95f661395756f806f..50cdc6020ffc533028c4ef5c37ad7c5434ce5939 100644
--- a/matlab/kalman/likelihood/univariate_kalman_filter.m
+++ b/matlab/kalman/likelihood/univariate_kalman_filter.m
@@ -260,7 +260,7 @@ if t <= last
             Hess = Hess + tmp{3};
         end
     else
-        [tmp, lik(s+1:end,:)] = univariate_kalman_filter_ss(Y,t,last,a,P,kalman_tol,T,H,Z,pp,Zflag);
+        [~, lik(s+1:end,:)] = univariate_kalman_filter_ss(Y,t,last,a,P,kalman_tol,T,H,Z,pp,Zflag);
     end
 end
 
diff --git a/matlab/kalman/missing_DiffuseKalmanSmootherH3_Z.m b/matlab/kalman/missing_DiffuseKalmanSmootherH3_Z.m
index 23d57a8c43ea8a54ae15a1ce9eba9133c4213ad5..d3acfa72ac5e33d5d234c7b5d2ccea41179b3fc4 100644
--- a/matlab/kalman/missing_DiffuseKalmanSmootherH3_Z.m
+++ b/matlab/kalman/missing_DiffuseKalmanSmootherH3_Z.m
@@ -509,7 +509,7 @@ while notsteady && t<smpl
             opts_simul.init_regime = []; %regimes_(t);
             opts_simul.waitbar=0;
             options_.occbin.simul=opts_simul;
-            [~, out, ss] = occbin.solver(M_,options_,dr,steady_state,exo_steady_state,exo_det_steady_state);
+            [~, out] = occbin.solver(M_,options_,dr,steady_state,exo_steady_state,exo_det_steady_state);
         end
         for jnk=1:nk
             if filter_covariance_flag
diff --git a/matlab/kalman/save_display_classical_smoother_results.m b/matlab/kalman/save_display_classical_smoother_results.m
index a0becb01f57d0b3e5f7bff7d11c3ed9ba0e48e06..48dcad0b5b31a0e7eae91bbe043acf70331dac18 100644
--- a/matlab/kalman/save_display_classical_smoother_results.m
+++ b/matlab/kalman/save_display_classical_smoother_results.m
@@ -52,20 +52,20 @@ if options_.occbin.smoother.status && options_.occbin.smoother.inversion_filter
     end
 else
     if options_.occbin.smoother.status
-        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = occbin.DSGE_smoother(xparam1,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_,dataset_,dataset_info);
+        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = occbin.DSGE_smoother(xparam1,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_,dataset_,dataset_info);
         if oo_.occbin.smoother.error_flag(1)==0
-            [oo_,yf]=store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
+            oo_ = store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
         else
             fprintf('\nOccbin: smoother did not succeed. No results will be written to oo_.\n')
         end
     else
-        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,T,R,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = DsgeSmoother(xparam1,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_);
-        [oo_,yf]=store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
+        [atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,~,~,P,PK,decomp,Trend,state_uncertainty,oo_,bayestopt_] = DsgeSmoother(xparam1,dataset_.nobs,transpose(dataset_.data),dataset_info.missing.aindex,dataset_info.missing.state,M_,oo_,options_,bayestopt_,estim_params_);
+        oo_ = store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
     end
     [oo_,yf]=store_smoother_results(M_,oo_,options_,bayestopt_,dataset_,dataset_info,atT,innov,measurement_error,updated_variables,ys,trend_coeff,aK,P,PK,decomp,Trend,state_uncertainty);
 end
 if ~options_.nograph
-    [nbplt,nr,nc,lr,lc,nstar] = pltorg(M_.exo_nbr);
+    [nbplt,nr,nc,~,~,nstar] = pltorg(M_.exo_nbr);
     if ~exist([M_.dname '/graphs'],'dir')
         mkdir(M_.dname,'graphs');
     end
@@ -133,7 +133,7 @@ if estim_params_.nvn
         end
     end
     if ~options_.nograph
-        [nbplt,nr,nc,lr,lc,nstar] = pltorg(number_of_plots_to_draw);
+        [nbplt,nr,nc,~,~,nstar] = pltorg(number_of_plots_to_draw);
         if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
             fidTeX = fopen([M_.dname, '/graphs/' M_.fname '_SmoothedObservationErrors.tex'],'w');
             fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation_1.m (Dynare).\n');
@@ -194,7 +194,7 @@ end
 %%  Historical and smoothed variabes
 %%
 if ~options_.nograph
-    [nbplt,nr,nc,lr,lc,nstar] = pltorg(n_varobs);
+    [nbplt,nr,nc,~,~,nstar] = pltorg(n_varobs);
     if options_.TeX && any(strcmp('eps',cellstr(options_.graph_format)))
         fidTeX = fopen([M_.dname, '/graphs/' M_.fname '_HistoricalAndSmoothedVariables.tex'],'w');
         fprintf(fidTeX,'%% TeX eps-loader file generated by dynare_estimation_1.m (Dynare).\n');
diff --git a/matlab/latex/collect_latex_files.m b/matlab/latex/collect_latex_files.m
index 536a3b808894559bcf0bf1a02e391b40bffca2a2..e5d776febb0780919e4dad8fab51bfd9aebab7ff 100644
--- a/matlab/latex/collect_latex_files.m
+++ b/matlab/latex/collect_latex_files.m
@@ -43,7 +43,7 @@ fprintf(fid,'%s \n','\begin{document}');
 %% Include LaTeX files from <fname>/latex/ directory, except the standalone ones
 TeX_Files=dir([M_.dname filesep 'latex' filesep '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(f_name, 'dynamic') && ...
             ~strcmp(f_name, 'static') && ...
             ~strcmp(f_name, 'original') && ...
@@ -55,7 +55,7 @@ end
 %% Output directory
 TeX_Files=dir([M_.dname filesep 'Output' filesep  M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/Output' '/',f_name,'}']);
     end
@@ -64,7 +64,7 @@ end
 %% graphs directory
 TeX_Files=dir([M_.dname filesep 'graphs' filesep  M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/graphs' '/',f_name,'}']);
     end
@@ -73,7 +73,7 @@ end
 %% Identification directory
 TeX_Files=dir([M_.dname filesep 'identification' filesep  M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/identification' '/',f_name,'}']);
     end
@@ -83,7 +83,7 @@ end
 %% Identification/Output directory
 TeX_Files=dir([M_.dname filesep 'identification' filesep 'Output' filesep M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/identification/Output' '/',f_name,'}']);
     end
@@ -92,7 +92,7 @@ end
 %% GSA directory
 TeX_Files=dir([M_.dname filesep 'gsa' filesep  M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/gsa' '/',f_name,'}']);
     end
@@ -101,7 +101,7 @@ end
 %% GSA/Output directory
 TeX_Files=dir([M_.dname filesep 'gsa' filesep 'Output' filesep  M_.fname '*.tex']);
 for ii=1:length(TeX_Files)
-    [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+    [~,f_name] = fileparts(TeX_Files(ii).name);
     if ~strcmp(TeX_Files(ii).name,f_name_binder)
         fprintf(fid,'%s \n',['\include{', M_.dname '/gsa/Output' '/',f_name,'}']);
     end
@@ -122,14 +122,14 @@ for level1_iter = 1:numsubdir_level1
     for level2_iter = 1:numsubdir_level2
         TeX_Files=dir([M_.dname filesep 'gsa' filesep dirinfo_parent(level1_iter).name filesep  dirinfo_subfolder(level2_iter).name filesep M_.fname '*.tex']);
         for ii=1:length(TeX_Files)
-            [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+            [~,f_name] = fileparts(TeX_Files(ii).name);
             if ~strcmp(TeX_Files(ii).name,f_name_binder)
                 fprintf(fid,'%s \n',['\include{', M_.dname '/gsa/',dirinfo_parent(level1_iter).name '/'  dirinfo_subfolder(level2_iter).name ,'/',f_name,'}']);
             end
         end
         TeX_Files=dir([M_.dname filesep 'gsa' filesep dirinfo_parent(level1_iter).name filesep  dirinfo_subfolder(level2_iter).name filesep 'Output' filesep  M_.fname '*.tex']);
         for ii=1:length(TeX_Files)
-            [pathstr,f_name,ext] = fileparts(TeX_Files(ii).name);
+            [~,f_name] = fileparts(TeX_Files(ii).name);
             if ~strcmp(TeX_Files(ii).name,f_name_binder)
                 fprintf(fid,'%s \n',['\include{', M_.dname '/gsa/', dirinfo_parent(level1_iter).name '/'  dirinfo_subfolder(level2_iter).name, '/Output' '/',f_name,'}']);
             end
diff --git a/matlab/lmmcp/dyn_lmmcp.m b/matlab/lmmcp/dyn_lmmcp.m
index aa63fa789f15fd145ccf53b6b065bc041db4ebab..87526e9d3b6cfd760ddf957078094fb43830e6e7 100644
--- a/matlab/lmmcp/dyn_lmmcp.m
+++ b/matlab/lmmcp/dyn_lmmcp.m
@@ -57,7 +57,7 @@ x = endo_simul(:);
 
 model_dynamic = str2func([M_.fname,'.dynamic']);
 z = x(find(lead_lag_incidence'));
-[res,A] = model_dynamic(z, exo_simul, params, steady_state,2);
+[~,A] = model_dynamic(z, exo_simul, params, steady_state,2);
 nnzA = nnz(A);
 
 LB = repmat(lb,periods,1);
diff --git a/matlab/load_m_file_data_legacy.m b/matlab/load_m_file_data_legacy.m
index 8b2afe4faa2a552416c1e4d8729c9b5ba49a0acd..782dba58ae37d4020c149dcd3cf0dea6bf28a2fb 100644
--- a/matlab/load_m_file_data_legacy.m
+++ b/matlab/load_m_file_data_legacy.m
@@ -18,7 +18,7 @@ function o2WysrOISH  = load_m_file_data_legacy(datafile, U7ORsJ0vy3)
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
 cXDHdrXnqo5KwwVpTRuc6OprAW = datafile(1:end-2);
-[pathtocXDHdrXnqo5KwwVpTRuc6OprAW,cXDHdrXnqo5KwwVpTRuc6OprAW,~] = fileparts(cXDHdrXnqo5KwwVpTRuc6OprAW);
+[pathtocXDHdrXnqo5KwwVpTRuc6OprAW,cXDHdrXnqo5KwwVpTRuc6OprAW] = fileparts(cXDHdrXnqo5KwwVpTRuc6OprAW);
 
 if ~isempty(pathtocXDHdrXnqo5KwwVpTRuc6OprAW)
     % We need to change directory, first we keep the current directory in memory...
diff --git a/matlab/long_run_variance.m b/matlab/long_run_variance.m
index 9cee0a0c5a89a6dc2c47146918b4c1d0360762c8..a3c02e406ea2448b4f632d32f94b7eb950916898 100644
--- a/matlab/long_run_variance.m
+++ b/matlab/long_run_variance.m
@@ -31,7 +31,7 @@ function sigma = long_run_variance(data,band)
 verbose = 1;
 
 if nargin<2
-    [T,m] = size(data);
+    [T,~] = size(data);
     band = ceil(T^(1/4));
     if verbose
         disp(['Bandwidth parameter is equal to ' num2str(band) '.'])
diff --git a/matlab/missing/mex/mjdgges/mjdgges.m b/matlab/missing/mex/mjdgges/mjdgges.m
index 4b1ef0648c87790382e96d6449e680fbdb1723b4..67ce0dc81288ffe978620b7f9052b36b7bfa493d 100644
--- a/matlab/missing/mex/mjdgges/mjdgges.m
+++ b/matlab/missing/mex/mjdgges/mjdgges.m
@@ -60,7 +60,7 @@ try
     eigval = ordeig(ss, tt);
     select = abs(eigval) < qz_criterium;
     sdim = sum(select);
-    [ss, tt, qq, zz] = ordqz(ss, tt, qq, zz, select);
+    [ss, tt, ~, zz] = ordqz(ss, tt, qq, zz, select);
     eigval = ordeig(ss, tt);
 catch
     info = 1; % Not as precise as lapack's info!
diff --git a/matlab/model_diagnostics.m b/matlab/model_diagnostics.m
index c17ac45780ca3510bde6e8c9233dcbcfcfd58459..70953988398697efc0c75d48cabdf93ad2e4a9e3 100644
--- a/matlab/model_diagnostics.m
+++ b/matlab/model_diagnostics.m
@@ -153,10 +153,10 @@ exo = [oo_.exo_steady_state; oo_.exo_det_steady_state];
 for b=1:nb
     if options_.bytecode
         if nb == 1
-            [res, jacob] = bytecode(M_, options_, dr.ys, exo, M_.params, dr.ys, 1, exo, ...
+            [~, jacob] = bytecode(M_, options_, dr.ys, exo, M_.params, dr.ys, 1, exo, ...
                                     'evaluate', 'static');
         else
-            [res, jacob] = bytecode(M_, options_, dr.ys, exo, M_.params, dr.ys, 1, exo, ...
+            [~, jacob] = bytecode(M_, options_, dr.ys, exo, M_.params, dr.ys, 1, exo, ...
                                     'evaluate', 'static', 'block_decomposed', ['block=' ...
                                 int2str(b)]);
         end
@@ -169,7 +169,7 @@ for b=1:nb
                 M_.block_structure_stat.block(b).g1_sparse_colptr, T);
             n_vars_jacob=size(jacob,2);
         else
-            [res, T_order, T] = feval([M_.fname '.sparse.static_resid'], dr.ys, exo, M_.params);
+            [~, T_order, T] = feval([M_.fname '.sparse.static_resid'], dr.ys, exo, M_.params);
             jacob = feval([M_.fname '.sparse.static_g1'], dr.ys, exo, M_.params, M_.static_g1_sparse_rowval, M_.static_g1_sparse_colval, M_.static_g1_sparse_colptr, T_order, T);
             n_vars_jacob=M_.endo_nbr;
         end
diff --git a/matlab/moments.m b/matlab/moments.m
index 86d5d8e8e4e4ae05894c42c760fec1b7bf788c35..2f5c3c3f0c57ca466701415521f78095c1b09e65 100644
--- a/matlab/moments.m
+++ b/matlab/moments.m
@@ -40,7 +40,7 @@ switch order
     if round(order)-order
         error('The second input argument (order) has to be an integer!')
     end
-    [T,n] = size(X);
+    [~,n] = size(X);
     c = mean(X);
     m = zeros(n,1);
     for i=1:n
diff --git a/matlab/moments/UnivariateSpectralDensity.m b/matlab/moments/UnivariateSpectralDensity.m
index 862dde5e612b3b2f67bac80688c330d26b2828f8..dcab15e93ae061d08e1a7bd01c3ea0fc28bd9e97 100644
--- a/matlab/moments/UnivariateSpectralDensity.m
+++ b/matlab/moments/UnivariateSpectralDensity.m
@@ -93,7 +93,7 @@ for i=M_.maximum_lag:-1:2
 end
 
 [A,B] = kalman_transition_matrix(oo_.dr,ikx',1:nx);
-[vx, u] =  lyapunov_symm(A,B*M_.Sigma_e*B',options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold,[],options_.debug);
+[~, u] =  lyapunov_symm(A,B*M_.Sigma_e*B',options_.lyapunov_fixed_point_tol,options_.qz_criterium,options_.lyapunov_complex_threshold,[],options_.debug);
 iky = iv(ivar);
 if ~isempty(u)
     iky = iky(find(any(abs(ghx(iky,:)*u) < options_.schur_vec_tol,2)));
diff --git a/matlab/moments/correlation_mc_analysis.m b/matlab/moments/correlation_mc_analysis.m
index 548860c9f5fdc49d10fbc4334029e4a4d4ef43bf..ecf58c4a458b6f179f8ca7351c57f28f99793529 100644
--- a/matlab/moments/correlation_mc_analysis.m
+++ b/matlab/moments/correlation_mc_analysis.m
@@ -61,8 +61,6 @@ if isfield(oo_,[TYPE 'TheoreticalMoments'])
                         % INITIALIZATION:
                         oo_ = initialize_output_structure(var1,var2,nar,type,oo_);
                         delete([PATH fname '_' TYPE 'Correlations*'])
-                        [nvar,vartan,NumberOfFiles] = ...
-                            dsge_simulated_theoretical_correlation(SampleSize,nar,M_,options_,oo_,type);
                     else
                         if ~isnan(temporary_structure_2(nar))
                             %Nothing to do.
@@ -143,4 +141,4 @@ switch moment
     oo_.([type, 'TheoreticalMoments']).dsge.correlation.(moment).(var1).(var2)(lag,1) = {result};
   otherwise
     disp('fill_output_structure:: Unknown field!')
-end
\ No newline at end of file
+end
diff --git a/matlab/moments/disp_moments.m b/matlab/moments/disp_moments.m
index 7cc3dad8a91607e8bde7e0c1ea6d580dc23aaccd..9ef54069e44688ba8fa6c2d60f352a38f805c26b 100644
--- a/matlab/moments/disp_moments.m
+++ b/matlab/moments/disp_moments.m
@@ -244,9 +244,9 @@ end
 function y = get_filtered_time_series(y, m, options_)
 
 if options_.hp_filter && ~options_.one_sided_hp_filter  && ~options_.bandpass.indicator
-    [hptrend,y] = sample_hp_filter(y,options_.hp_filter);
+    [~,y] = sample_hp_filter(y,options_.hp_filter);
 elseif ~options_.hp_filter && options_.one_sided_hp_filter && ~options_.bandpass.indicator
-    [hptrend,y] = one_sided_hp_filter(y,options_.one_sided_hp_filter);
+    [~,y] = one_sided_hp_filter(y,options_.one_sided_hp_filter);
 elseif ~options_.hp_filter && ~options_.one_sided_hp_filter && options_.bandpass.indicator
     data_temp=dseries(y,'0q1');
     data_temp=baxter_king_filter(data_temp,options_.bandpass.passband(1),options_.bandpass.passband(2),options_.bandpass.K);
diff --git a/matlab/ms-sbvar/ms_mardd.m b/matlab/ms-sbvar/ms_mardd.m
index aa1d1985331df3e100ea590f2ea5061ea40f2432..851c6556d3a1d01892373b48c983eef360e8ef12 100644
--- a/matlab/ms-sbvar/ms_mardd.m
+++ b/matlab/ms-sbvar/ms_mardd.m
@@ -89,7 +89,7 @@ vlog_Y_a = -0.5*nvar*fss*log(2*pi) + fss*log(abs(det(A0xhat))) + Yexpt
 logMarLHres = 0;   % Initialize log of the marginal likelihood (restricted or constant parameters).
 for ki=1:fss   %ndobs+1:fss     % Forward recursion to get the marginal likelihood.  See F on p.19 and pp. 48-49.
                %----  Restricted log marginal likelihood function (constant parameters).
-    [A0l,A0u] = lu(A0xhat);
+    [~,A0u] = lu(A0xhat);
     ada = sum(log(abs(diag(A0u))));   % log|A0|
     termexp = y(ki,:)*A0xhat - phi(ki,:)*Apxhat;   % 1-by-nvar
     logMarLHres = logMarLHres - (0.5*nvar)*log(2*pi) + ada - 0.5*termexp*termexp';  % log value
@@ -148,7 +148,7 @@ for k=1:nvar
             %------ Computing p(a0_k|Y,a_others) at some point such as the peak along the dimensions of indx_ks.
             Vk = Tinv{k}\Wcell{k};  %V_k on p.71 of Forecast (II).
             gbeta = Vk\bk;  % inv(V_k)*b_k on p.71 of Forecast (II) where alpha_k = b_k in our notation.
-            [Vtq,Vtr]=qr(Vk',0);  %To get inv(V_k)'*inv(V_k) in (*) on p.71 of Forecast (II).
+            [~,Vtr]=qr(Vk',0);  %To get inv(V_k)'*inv(V_k) in (*) on p.71 of Forecast (II).
                                   %
             vlog(draws) = 0.5*(fss+nk)*log(fss)-log(abs(det(Vk)))-0.5*(nk-1)*log(2*pi)-...
                 0.5*(fss+1)*log(2)-gammaln(0.5*(fss+1))+fss*log(abs(gbeta(1)))-...
@@ -164,12 +164,12 @@ for k=1:nvar
     else
         skipline()
         disp('The last(6th) column or equation in A0 with no Gibbs draws')
-        [A0gbs1, Wcell] = fn_gibbsrvar(A0gbs0,UT,nvar,fss,n0,indx_ks)
+        [~, Wcell] = fn_gibbsrvar(A0gbs0,UT,nvar,fss,n0,indx_ks)
         %------ See p.71, Forecast (II).
         %------ Computing p(a0_k|Y,a_others) at some point such as the peak along the dimensions of indx_ks.
         Vk = Tinv{k}\Wcell{k};  %V_k on p.71 of Forecast (II).
         gbeta = Vk\bk;  % inv(V_k)*b_k on p.71 of Forecast (II) where alpha_k = b_k in our notation.
-        [Vtq,Vtr]=qr(Vk',0);  %To get inv(V_k)'*inv(V_k) in (*) on p.71 of Forecast (II).
+        [~,Vtr]=qr(Vk',0);  %To get inv(V_k)'*inv(V_k) in (*) on p.71 of Forecast (II).
                               %
         vloglast = 0.5*(fss+nk)*log(fss)-log(abs(det(Vk)))-0.5*(nk-1)*log(2*pi)-...
             0.5*(fss+1)*log(2)-gammaln(0.5*(fss+1))+fss*log(abs(gbeta(1)))-...
diff --git a/matlab/ms-sbvar/ms_sbvar_setup.m b/matlab/ms-sbvar/ms_sbvar_setup.m
index 9eb6db3872dbb87e9b38ae80d0300ddc3e36313c..17e0f2299016febc8a58c69cdeacbed8715d0418 100644
--- a/matlab/ms-sbvar/ms_sbvar_setup.m
+++ b/matlab/ms-sbvar/ms_sbvar_setup.m
@@ -209,7 +209,7 @@ nexo=1;
 % Arranged data information, WITHOUT dummy obs when 0 after mu is used.
 % See fn_rnrprior_covres_dobs.m for using the dummy observations as part of
 % an explicit prior.
-[xtx,xty,yty,fss,phi,y,ncoef,xr,Bh] = fn_dataxy(nvar,options_.ms.nlags,xdgel,mu,0,nexo);
+[xtx,xty,yty,fss,phi,y,ncoef] = fn_dataxy(nvar,options_.ms.nlags,xdgel,mu,0,nexo);
 
 
 %======================================================================
@@ -239,7 +239,7 @@ if indxPrior
     %*** Obtains asymmetric prior (with no linear restrictions) with dummy observations as part of an explicit prior (i.e,
     %      reflected in Hpmulti and Hpinvmulti).  See Forecast II, pp.69a-69b for details.
     if 1  % Liquidity effect prior on both MS and MD equations.
-        [Pi,H0multi,Hpmulti,H0invmulti,Hpinvmulti] = fn_rnrprior_covres_dobs(nvar,q_m,options_.ms.nlags,xdgel,mu,indxDummy,hpmsmd,indxmsmdeqn);
+        [Pi,H0multi,Hpmulti] = fn_rnrprior_covres_dobs(nvar,q_m,options_.ms.nlags,xdgel,mu,indxDummy,hpmsmd,indxmsmdeqn);
     else
         [Pi,H0multi,Hpmulti,H0invmulti,Hpinvmulti] = fn_rnrprior(nvar,q_m,options_.ms.nlags,xdgel,mu);
     end
@@ -268,12 +268,12 @@ else
     crit = 1.0e-9;
     nit = 10000;
     %
-    [fhat,xhat,grad,Hhat,itct,fcount,retcodehat] = csminwel('fn_a0freefun',x,H0,'fn_a0freegrad',crit,nit,Ui,nvar,n0,fss,H0inv);
+    [~,xhat] = csminwel('fn_a0freefun',x,H0,'fn_a0freegrad',crit,nit,Ui,nvar,n0,fss,H0inv);
 
     A0hat = fn_tran_b2a(xhat,Ui,nvar,n0);
 
     xhat = fn_tran_a2b(A0hat,Ui,nvar,n0);
-    [Aphat,ghat] = fn_gfmean(xhat,Pmat,Vi,nvar,ncoef,n0,np);
+    Aphat = fn_gfmean(xhat,Pmat,Vi,nvar,ncoef,n0,np);
     if indxC0Pres
         Fhatur0P = Fhat;  % ur: unrestriced across A0 and A+
         for ki = 1:size(ixmC0Pres,1)   % loop through the number of equations in which
diff --git a/matlab/ms-sbvar/plot_ms_probabilities.m b/matlab/ms-sbvar/plot_ms_probabilities.m
index ff1b08422e7f7c1058701eea2918370d70e5912b..4a5d81c6e7284a790f2afdd60240976df29494fa 100644
--- a/matlab/ms-sbvar/plot_ms_probabilities.m
+++ b/matlab/ms-sbvar/plot_ms_probabilities.m
@@ -28,7 +28,7 @@ function plot_ms_probabilities(computed_probabilities, options_)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[T,num_grand_regimes] = size(computed_probabilities);
+[T,~] = size(computed_probabilities);
 num_chains = length(options_.ms.ms_chain);
 for i=1:num_chains
     chains(i).num_regimes = length(options_.ms.ms_chain(i).regime);
diff --git a/matlab/ms-sbvar/svar_global_identification_check.m b/matlab/ms-sbvar/svar_global_identification_check.m
index 5baf019fa8c753fed6f33c5e8d56d0904203bc00..92751789acb1e2a8cde6105a1ade0344c41eaa3b 100644
--- a/matlab/ms-sbvar/svar_global_identification_check.m
+++ b/matlab/ms-sbvar/svar_global_identification_check.m
@@ -48,7 +48,7 @@ end
 nvar = length(options_.varobs);   % number of endogenous variables
 nexo = 1;
 
-[Uiconst,Viconst,n0,np,ixmC0Pres,Qi,Ri] = exclusions(nvar,nexo,options_.ms );
+[~,~,~,~,~,Qi,Ri] = exclusions(nvar,nexo,options_.ms );
 
 % order column constraints by rank
 QQranks = zeros(nvar,2);
diff --git a/matlab/nonlinear-filters/auxiliary_particle_filter.m b/matlab/nonlinear-filters/auxiliary_particle_filter.m
index 7071f68f8e741d8e933380b004a00b6f27caf364..713864308fc673f0592e9e9f48fcda88a9ce796f 100644
--- a/matlab/nonlinear-filters/auxiliary_particle_filter.m
+++ b/matlab/nonlinear-filters/auxiliary_particle_filter.m
@@ -111,9 +111,9 @@ for t=1:sample_size
     if pruning
         yhat_ = bsxfun(@minus,StateVectors_,state_variables_steady_state_);
         if order == 2
-            [tmp, tmp_] = local_state_space_iteration_2(yhat,zeros(number_of_structural_innovations,number_of_particles),ghx,ghu,constant,ghxx,ghuu,ghxu,yhat_,steadystate,ThreadsOptions.local_state_space_iteration_2);
+            tmp = local_state_space_iteration_2(yhat,zeros(number_of_structural_innovations,number_of_particles),ghx,ghu,constant,ghxx,ghuu,ghxu,yhat_,steadystate,ThreadsOptions.local_state_space_iteration_2);
         elseif order == 3
-            [tmp, tmp_] = local_state_space_iteration_3(yhat_, zeros(number_of_structural_innovations,number_of_particles), ghx, ghu, ghxx, ghuu, ghxu, ghs2, ghxxx, ghuuu, ghxxu, ghxuu, ghxss, ghuss, steadystate, ThreadsOptions.local_state_space_iteration_3, pruning);
+            tmp = local_state_space_iteration_3(yhat_, zeros(number_of_structural_innovations,number_of_particles), ghx, ghu, ghxx, ghuu, ghxu, ghs2, ghxxx, ghuuu, ghxxu, ghxuu, ghxss, ghuss, steadystate, ThreadsOptions.local_state_space_iteration_3, pruning);
         else
             error('Pruning is not available for orders > 3');
         end
diff --git a/matlab/nonlinear-filters/fit_gaussian_mixture.m b/matlab/nonlinear-filters/fit_gaussian_mixture.m
index 9836fa59e6ffe33cf49fdf77388e5891d6f14c63..4bf25f25b00785ca980f64fd28dd72c40abe476b 100644
--- a/matlab/nonlinear-filters/fit_gaussian_mixture.m
+++ b/matlab/nonlinear-filters/fit_gaussian_mixture.m
@@ -17,7 +17,7 @@ function [StateMu,StateSqrtP,StateWeights] = fit_gaussian_mixture(X,X_weights,St
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[dim,Ndata] = size(X);
+[~,Ndata] = size(X);
 M = size(StateMu,2) ;
 if check                        % Ensure that covariances don't collapse
     MIN_COVAR_SQRT = sqrt(eps);
@@ -26,7 +26,7 @@ end
 eold = -Inf;
 for n=1:niters
     % Calculate posteriors based on old parameters
-    [prior,likelihood,marginal,posterior] = probability3(StateMu,StateSqrtP,StateWeights,X,X_weights);
+    [~,~,marginal,posterior] = probability3(StateMu,StateSqrtP,StateWeights,X,X_weights);
     e = sum(log(marginal));
     if (n > 1 && abs((e - eold)/eold) < crit)
         return;
@@ -40,7 +40,7 @@ for n=1:niters
         diffs = bsxfun(@minus,X,StateMu(:,j));
         tpost = (1/sqrt(new_pr(j)))*sqrt(posterior(j,:));
         diffs = bsxfun(@times,diffs,tpost);
-        [foo,tcov] = qr2(diffs',0);
+        [~,tcov] = qr2(diffs',0);
         StateSqrtP(:,:,j) = tcov';
         if check
             if min(abs(diag(StateSqrtP(:,:,j)))) < MIN_COVAR_SQRT
diff --git a/matlab/nonlinear-filters/mykmeans.m b/matlab/nonlinear-filters/mykmeans.m
index e7c424b00646f8ccf3dbcee70395b218c632cd54..60383ad95a98494f23d12a9f8146bec743cf9e23 100644
--- a/matlab/nonlinear-filters/mykmeans.m
+++ b/matlab/nonlinear-filters/mykmeans.m
@@ -32,7 +32,7 @@ for iter=1:300
     for i=1:g
         dist(i,:) = sum(bsxfun(@power,bsxfun(@minus,x,c(:,i)),2));
     end
-    [rien,ind] = min(dist) ;
+    [~,ind] = min(dist) ;
     if isequal(ind,indold)
         break ;
     end
diff --git a/matlab/nonlinear-filters/online_auxiliary_filter.m b/matlab/nonlinear-filters/online_auxiliary_filter.m
index 0405ff27351ad73b5133ff49ae72269543484872..7495d4bc66c5d5d1602419d193c4f2e8ac0f0e23 100644
--- a/matlab/nonlinear-filters/online_auxiliary_filter.m
+++ b/matlab/nonlinear-filters/online_auxiliary_filter.m
@@ -172,11 +172,11 @@ for t=1:sample_size
                 if pruning
                     yhat_ = bsxfun(@minus,StateVectors_(:,i),state_variables_steady_state_);
                     if order == 2
-                        [tmp, ~] = local_state_space_iteration_2(yhat, zeros(number_of_structural_innovations, 1), ghx, ghu, constant, ghxx, ghuu, ghxu, yhat_, steadystate, options_.threads.local_state_space_iteration_2);
+                        tmp = local_state_space_iteration_2(yhat, zeros(number_of_structural_innovations, 1), ghx, ghu, constant, ghxx, ghuu, ghxu, yhat_, steadystate, options_.threads.local_state_space_iteration_2);
                     elseif order == 3
-                        [tmp, tmp_] = local_state_space_iteration_3(yhat_, zeros(number_of_structural_innovations, 1), ghx, ghu, ghxx, ghuu, ghxu, ghs2, ghxxx, ghuuu, ghxxu, ghxuu, ghxss, ghuss, steadystate, options_.threads.local_state_space_iteration_3, pruning);
+                        tmp = local_state_space_iteration_3(yhat_, zeros(number_of_structural_innovations, 1), ghx, ghu, ghxx, ghuu, ghxu, ghs2, ghxxx, ghuuu, ghxxu, ghxuu, ghxss, ghuss, steadystate, options_.threads.local_state_space_iteration_3, pruning);
                     else
-                        error('Pruning is not available for orders > 3');
+                    error('Pruning is not available for orders > 3');
                     end
                 else
                     if order == 2
diff --git a/matlab/nonlinear-filters/reduced_rank_cholesky.m b/matlab/nonlinear-filters/reduced_rank_cholesky.m
index 4f65157fcbdcf2fc8f79a5c65a7f655ed64ad9b8..9afd8e72b5295f9fc60cf433d0af50dc9384b0c2 100644
--- a/matlab/nonlinear-filters/reduced_rank_cholesky.m
+++ b/matlab/nonlinear-filters/reduced_rank_cholesky.m
@@ -58,7 +58,7 @@ function T = reduced_rank_cholesky(X)
 if X_is_not_positive_definite
     n = length(X);
     [U,D] = eig(X);
-    [tmp,max_elements_indices] = max(abs(U),[],1);
+    [~,max_elements_indices] = max(abs(U),[],1);
     negloc = (U(max_elements_indices+(0:n:(n-1)*n))<0);
     U(:,negloc) = -U(:,negloc);
     D = diag(D);
diff --git a/matlab/optimal_policy/dyn_ramsey_static.m b/matlab/optimal_policy/dyn_ramsey_static.m
index ac4b4d23ca97347c95b3ba4f3a6c320a4db99805..8e99a14eba3e76078994280b1f36f662d3e45dc5 100644
--- a/matlab/optimal_policy/dyn_ramsey_static.m
+++ b/matlab/optimal_policy/dyn_ramsey_static.m
@@ -78,7 +78,7 @@ elseif options_.steadystate_flag
         end
     end
     ys_init(k_inst) = inst_val;
-    [xx,params] = evaluate_steady_state_file(ys_init,exo_ss,M_,options_,~options_.steadystate.nocheck); %run steady state file again to update parameters
+    [~,params] = evaluate_steady_state_file(ys_init,exo_ss,M_,options_,~options_.steadystate.nocheck); %run steady state file again to update parameters
     [~,~,steady_state] = nl_func(inst_val); %compute and return steady state
 else
     xx = ys_init(1:M_.orig_endo_nbr);
@@ -149,7 +149,7 @@ resids1 = y+A*mult;
 if inst_nbr == 1
     r1 = sqrt(resids1'*resids1);
 else
-    [q,r,e] = qr([A y]');
+    [~,r] = qr([A y]');
     k = size(A,1)+(1-inst_nbr:0);
     r1 = r(end,k)';
 end
@@ -167,7 +167,7 @@ end
 function result = check_static_model(ys,exo_ss,M_,options_)
 result = false;
 if (options_.bytecode)
-    [res, ~] = bytecode('static', M_, options, ys, exo_ss, M_.params, 'evaluate');
+    res = bytecode('static', M_, options, ys, exo_ss, M_.params, 'evaluate');
 else
     res = feval([M_.fname '.sparse.static_resid'], ys, exo_ss, M_.params);
 end
diff --git a/matlab/optimal_policy/mult_elimination.m b/matlab/optimal_policy/mult_elimination.m
index bc5d4d52503224cfdf5eec9bf39162a7c6cab4bf..8961c99c5f0181a2021ec78837efea9640772162 100644
--- a/matlab/optimal_policy/mult_elimination.m
+++ b/matlab/optimal_policy/mult_elimination.m
@@ -53,7 +53,7 @@ A22 = A2(il,:);
 B1 = B(nil,:);
 B2 = B(il,:);
 
-[Q1,R1,E1] = qr([A12; A22]);
+[Q1,R1] = qr([A12; A22]);
 n1 = sum(abs(diag(R1)) > 1e-8);
 
 Q1_12 = Q1(1:nm_nbr,n1+1:end);
diff --git a/matlab/optimization/analytic_gradient_wrapper.m b/matlab/optimization/analytic_gradient_wrapper.m
index 52800ad717e27d75e4f8bbac9b8588132bbe9083..8c0061bdab0aaecb8d7a7379fab1a601182a7811 100644
--- a/matlab/optimization/analytic_gradient_wrapper.m
+++ b/matlab/optimization/analytic_gradient_wrapper.m
@@ -31,7 +31,7 @@ function [fval, grad, hess, exit_flag]=analytic_gradient_wrapper(x, fcn, varargi
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[fval, info, exit_flag, grad, hess] = fcn(x, varargin{:});
+[fval, ~, exit_flag, grad, hess] = fcn(x, varargin{:});
 if size(grad,2)==1
     grad=grad'; %should be row vector for Matlab; exception lsqnonlin where Jacobian is required
 end
\ No newline at end of file
diff --git a/matlab/optimization/cmaes.m b/matlab/optimization/cmaes.m
index 71783eaf3900ea0d97301e4cf3f7791c6db74ff4..c08bf008fb2451ac3dffd57b1349c78ef6a540dd 100644
--- a/matlab/optimization/cmaes.m
+++ b/matlab/optimization/cmaes.m
@@ -684,7 +684,7 @@ while irun <= myeval(opts.Restarts) % for-loop does not work with resume
                 for namecell = filenames(:)'
                     name = namecell{:};
 
-                    [fid, err] = fopen(['./' filenameprefix name '.dat'], 'w');
+                    fid = fopen(['./' filenameprefix name '.dat'], 'w');
                     if fid < 1 % err ~= 0
                         warning(['could not open ' filenameprefix name '.dat']);
                         filenames(find(strcmp(filenames,name))) = [];
@@ -1151,7 +1151,7 @@ while irun <= myeval(opts.Restarts) % for-loop does not work with resume
                     % TODO: this is not in compliance with the paper Hansen&Ros2010,
                     %       where simply arnorms = arnorms(end:-1:1) ?
                     [arnorms, idxnorms] = sort(sqrt(sum(arzneg.^2, 1)));
-                    [ignore, idxnorms] = sort(idxnorms);  % inverse index
+                    [~, idxnorms] = sort(idxnorms);  % inverse index
                     arnormfacs = arnorms(end:-1:1) ./ arnorms;
                     % arnormfacs = arnorms(randperm(neg.mu)) ./ arnorms;
                     arnorms = arnorms(end:-1:1); % for the record
@@ -1596,7 +1596,7 @@ while irun <= myeval(opts.Restarts) % for-loop does not work with resume
             for namecell = filenames(:)'
                 name = namecell{:};
 
-                [fid, err] = fopen(['./' filenameprefix name '.dat'], 'a');
+                fid = fopen(['./' filenameprefix name '.dat'], 'a');
                 if fid < 1 % err ~= 0
                     warning(['could not open ' filenameprefix name '.dat']);
                 else
@@ -1981,7 +1981,7 @@ end
 
 % sort inar
 if nargin < 3 || isempty(idx)
-    [sar, idx] = sort(inar);
+    sar = sort(inar);
 else
     sar = inar(idx);
 end
@@ -2070,8 +2070,8 @@ end
 
 %%% compute rank changes into rankDelta
 % compute ranks
-[ignore, idx] = sort([arf1 arf2]);
-[ignore, ranks] = sort(idx);
+[~, idx] = sort([arf1 arf2]);
+[~, ranks] = sort(idx);
 ranks = reshape(ranks, lam, 2)';
 
 rankDelta = ranks(1,:) - ranks(2,:) - sign(ranks(1,:) - ranks(2,:));
@@ -2199,7 +2199,7 @@ end
 % plot fitness etc
 foffset = 1e-99;
 dfit = d.f(:,6)-min(d.f(:,6));
-[ignore, idxbest] = min(dfit);
+[~, idxbest] = min(dfit);
 dfit(dfit<1e-98) = NaN;
 subplot(2,2,1); hold off;
 dd = abs(d.f(:,7:8)) + foffset;
@@ -2256,7 +2256,7 @@ ax(2) = max(minxend, ax(2));
 axis(ax);
 
 % add some annotation lines
-[ignore, idx] = sort(d.x(end,6:end));
+[~, idx] = sort(d.x(end,6:end));
 % choose no more than 25 indices
 idxs = round(linspace(1, size(d.x,2)-5, min(size(d.x,2)-5, 25)));
 yy = repmat(NaN, 2, size(d.x,2)-5);
@@ -2300,7 +2300,7 @@ ax = axis;
 ax(2) = max(minxend, ax(2));
 axis(ax);
 % add some annotation lines
-[ignore, idx] = sort(d.std(end,6:end));
+[~, idx] = sort(d.std(end,6:end));
 % choose no more than 25 indices
 idxs = round(linspace(1, size(d.x,2)-5, min(size(d.x,2)-5, 25)));
 yy = repmat(NaN, 2, size(d.std,2)-5);
@@ -2380,15 +2380,14 @@ function f=fmixranks(x)
 N = size(x,1);
 f=(10.^(0*(0:(N-1))/(N-1))*x.^2).^0.5;
 if size(x, 2) > 1 % compute ranks, if it is a population
-    [ignore, idx] = sort(f);
-    [ignore, ranks] = sort(idx);
+    [~, idx] = sort(f);
     k = 9; % number of solutions randomly permuted, lambda/2-1
            % works still quite well (two time slower)
     for i = k+1:k-0:size(x,2)
         idx(i-k+(1:k)) = idx(i-k+randperm(k));
     end
     %disp([ranks' f'])
-    [ignore, ranks] = sort(idx);
+    [~, ranks] = sort(idx);
     %disp([ranks' f'])
     %pause
     f = ranks+1e-9*randn(1,1);
diff --git a/matlab/optimization/dynare_minimize_objective.m b/matlab/optimization/dynare_minimize_objective.m
index 18115ded7dca894dd39d88762c67b047eac66a54..ea80ac09139d1c8459f93f57f6bb612ae7fd9236 100644
--- a/matlab/optimization/dynare_minimize_objective.m
+++ b/matlab/optimization/dynare_minimize_objective.m
@@ -107,23 +107,23 @@ switch minimizer_algorithm
     if options_.analytic_derivation || (isfield(options_,'mom') && options_.mom.analytic_jacobian==1) %use wrapper
         func = @(x) analytic_gradient_wrapper(x,objective_function,varargin{:});
         if ~isoctave
-            [opt_par_values,fval,exitflag,output,lamdba,grad,hessian_mat] = ...
+            [opt_par_values,fval,exitflag,~,~,~,hessian_mat] = ...
                 fmincon(func,start_par_value,[],[],[],[],bounds(:,1),bounds(:,2),[],optim_options);
         else
             % Under Octave, use a wrapper, since fmincon() does not have an 11th
             % arg. Also, only the first 4 output arguments are available.
-            [opt_par_values,fval,exitflag,output] = ...
+            [opt_par_values,fval,exitflag] = ...
                 fmincon(func,start_par_value,[],[],[],[],bounds(:,1),bounds(:,2),[],optim_options);
         end
     else
         if ~isoctave
-            [opt_par_values,fval,exitflag,output,lamdba,grad,hessian_mat] = ...
+            [opt_par_values,fval,exitflag,~,~,~,hessian_mat] = ...
                 fmincon(objective_function,start_par_value,[],[],[],[],bounds(:,1),bounds(:,2),[],optim_options,varargin{:});
         else
             % Under Octave, use a wrapper, since fmincon() does not have an 11th
             % arg. Also, only the first 4 output arguments are available.
             func = @(x) objective_function(x,varargin{:});
-            [opt_par_values,fval,exitflag,output] = ...
+            [opt_par_values,fval,exitflag] = ...
                 fmincon(func,start_par_value,[],[],[],[],bounds(:,1),bounds(:,2),[],optim_options);
         end    
     end
@@ -177,7 +177,7 @@ switch minimizer_algorithm
     end
     sa_options.initial_step_length= sa_options.initial_step_length*ones(npar,1); %bring step length to correct vector size
     sa_options.step_length_c= sa_options.step_length_c*ones(npar,1); %bring step_length_c to correct vector size
-    [opt_par_values, fval,exitflag, n_accepted_draws, n_total_draws, n_out_of_bounds_draws, t, vm] =...
+    [opt_par_values, fval,exitflag] =...
         simulated_annealing(objective_function,start_par_value,sa_options,LB,UB,varargin{:});
   case 3
     if isoctave && ~user_has_octave_forge_package('optim')
@@ -269,7 +269,7 @@ switch minimizer_algorithm
         analytic_grad=[];
     end
     % Call csminwell.
-    [fval,opt_par_values,grad,inverse_hessian_mat,itct,fcount,exitflag] = ...
+    [fval,opt_par_values,~,inverse_hessian_mat,~,~,exitflag] = ...
         csminwel1(objective_function, start_par_value, H0, analytic_grad, crit, nit, numgrad, epsilon, Verbose, Save_files, varargin{:});
     hessian_mat=inv(inverse_hessian_mat);
   case 5
@@ -337,7 +337,7 @@ switch minimizer_algorithm
     hess_info.robust=robust;
     % here we force 7th input argument (flagg) to be 0, since outer product
     % gradient Hessian is handled in dynare_estimation_1
-    [opt_par_values,hessian_mat,gg,fval,invhess,new_rat_hess_info] = newrat(objective_function,start_par_value,bounds,analytic_grad,crit,nit,0,Verbose,Save_files,hess_info,prior_information.p2,options_.gradient_epsilon,parameter_names,varargin{:});    %hessian_mat is the plain outer product gradient Hessian
+    [opt_par_values,hessian_mat,~,fval,~,new_rat_hess_info] = newrat(objective_function,start_par_value,bounds,analytic_grad,crit,nit,0,Verbose,Save_files,hess_info,prior_information.p2,options_.gradient_epsilon,parameter_names,varargin{:});    %hessian_mat is the plain outer product gradient Hessian
     new_rat_hess_info.new_rat_hess_info = new_rat_hess_info;
     new_rat_hess_info.newratflag = newratflag;
     if options_.analytic_derivation 
@@ -453,7 +453,7 @@ switch minimizer_algorithm
     end
     warning('off','CMAES:NonfinitenessRange');
     warning('off','CMAES:InitialSigma');
-    [x, fval, COUNTEVAL, STOPFLAG, OUT, BESTEVER] = cmaes(func2str(objective_function),start_par_value,H0,cmaesOptions,varargin{:});
+    [~, ~, ~, ~, ~, BESTEVER] = cmaes(func2str(objective_function),start_par_value,H0,cmaesOptions,varargin{:});
     opt_par_values=BESTEVER.x;
     fval=BESTEVER.f;
   case 10
@@ -498,7 +498,7 @@ switch minimizer_algorithm
   case 11
     options_.cova_compute = 0;
     subvarargin = [varargin(1), varargin(3:6), varargin(8)];
-    [opt_par_values, stdh, lb_95, ub_95, med_param] = online_auxiliary_filter(start_par_value, subvarargin{:});
+    opt_par_values = online_auxiliary_filter(start_par_value, subvarargin{:});
   case 12
     if isoctave
         error('Option mode_compute=12 is not available under Octave')
@@ -512,10 +512,10 @@ switch minimizer_algorithm
     if ~isempty(options_.optim_opt)
         options_list = read_key_value_string(options_.optim_opt);
         SupportedListOfOptions = {'CreationFcn', 'Display', 'DisplayInterval', 'FunctionTolerance', ...
-                            'FunValCheck', 'HybridFcn', 'InertiaRange', 'InitialSwarmMatrix', 'InitialSwarmSpan', ...
-                            'MaxIterations', 'MaxStallIterations', 'MaxStallTime', 'MaxTime', ...
-                            'MinNeighborsFraction', 'ObjectiveLimit', 'OutputFcn', 'PlotFcn', 'SelfAdjustmentWeight', ...
-                            'SocialAdjustmentWeight', 'SwarmSize', 'UseParallel', 'UseVectorized'};
+                                  'FunValCheck', 'HybridFcn', 'InertiaRange', 'InitialSwarmMatrix', 'InitialSwarmSpan', ...
+                                  'MaxIterations', 'MaxStallIterations', 'MaxStallTime', 'MaxTime', ...
+                                  'MinNeighborsFraction', 'ObjectiveLimit', 'OutputFcn', 'PlotFcn', 'SelfAdjustmentWeight', ...
+                                  'SocialAdjustmentWeight', 'SwarmSize', 'UseParallel', 'UseVectorized'};
         for i=1:rows(options_list)
             if ismember(options_list{i,1}, SupportedListOfOptions)
                 particleswarmOptions = optimoptions(particleswarmOptions, options_list{i,1}, options_list{i,2});
@@ -537,7 +537,7 @@ switch minimizer_algorithm
         FVALS = zeros(particleswarmOptions.SwarmSize, 1);
         while p<=particleswarmOptions.SwarmSize
             candidate = rand(numberofvariables, 1).*(UB-LB)+LB;
-            [fval, info, exit_flag] = objfun(candidate);
+            [fval, ~, exit_flag] = objfun(candidate);
             if exit_flag
                 particleswarmOptions.InitialSwarmMatrix(p,:) = transpose(candidate);
                 FVALS(p) = fval;
@@ -552,7 +552,7 @@ switch minimizer_algorithm
     % Define penalized objective.
     objfun = @(x) penalty_objective_function(x, objective_function, penalty, varargin{:});
     % Minimize the penalized objective (note that the penalty is not updated).
-    [opt_par_values, fval, exitflag, output] = particleswarm(objfun, length(start_par_value), LB, UB, particleswarmOptions);
+    [opt_par_values, fval, exitflag] = particleswarm(objfun, length(start_par_value), LB, UB, particleswarmOptions);
     opt_par_values = opt_par_values(:);
   case 13
     % Matlab's lsqnonlin (Optimization toolbox needed).
@@ -586,15 +586,15 @@ switch minimizer_algorithm
             optim_options.SpecifyObjectiveGradient = true;
         end
         func = @(x) analytic_gradient_wrapper(x,objective_function,varargin{:});
-        [opt_par_values,Resnorm,fval,exitflag,OUTPUT,LAMBDA,JACOB] = ...
+        [opt_par_values,~,fval,exitflag] = ...
             lsqnonlin(func,start_par_value,bounds(:,1),bounds(:,2),optim_options);
     else
         if ~isoctave
-            [opt_par_values,Resnorm,fval,exitflag,OUTPUT,LAMBDA,JACOB] = lsqnonlin(objective_function,start_par_value,bounds(:,1),bounds(:,2),optim_options,varargin{:});
+            [opt_par_values,~,fval,exitflag] = lsqnonlin(objective_function,start_par_value,bounds(:,1),bounds(:,2),optim_options,varargin{:});
         else
             % Under Octave, use a wrapper, since lsqnonlin() does not have a 6th arg
             func = @(x)objective_function(x,varargin{:});
-            [opt_par_values,Resnorm,fval,exitflag,OUTPUT,LAMBDA,JACOB] = lsqnonlin(func,start_par_value,bounds(:,1),bounds(:,2),optim_options);
+            [opt_par_values,~,fval,exitflag] = lsqnonlin(func,start_par_value,bounds(:,1),bounds(:,2),optim_options);
         end
     end    
   case 101
@@ -656,7 +656,7 @@ switch minimizer_algorithm
         end
     end
     func = @(x)objective_function(x,varargin{:});
-    [opt_par_values,fval,exitflag,output] = simulannealbnd(func,start_par_value,bounds(:,1),bounds(:,2),optim_options);
+    [opt_par_values,fval,exitflag] = simulannealbnd(func,start_par_value,bounds(:,1),bounds(:,2),optim_options);
   otherwise
     if ischar(minimizer_algorithm)
         if exist(minimizer_algorithm)
diff --git a/matlab/optimization/dynare_solve.m b/matlab/optimization/dynare_solve.m
index 871737faab0dfdeba3c2a6a87c06c6cb7ebea23b..b41f82e73ed4820d8031ccd09ce34ee428410c4e 100644
--- a/matlab/optimization/dynare_solve.m
+++ b/matlab/optimization/dynare_solve.m
@@ -318,7 +318,7 @@ elseif options.solve_algo == 11
     mcp_data.func = f;
     mcp_data.args = varargin;
     try
-        [x, fval, jac, mu] = pathmcp(x,omcppath.lb,omcppath.ub,'mcp_func',omcppath.A,omcppath.b,omcppath.t,omcppath.mu0);
+        x = pathmcp(x,omcppath.lb,omcppath.ub,'mcp_func',omcppath.A,omcppath.b,omcppath.t,omcppath.mu0);
     catch
         errorflag = true;
     end
diff --git a/matlab/optimization/mr_gstep.m b/matlab/optimization/mr_gstep.m
index fd83ea50bfaf85ca268af9ff109df9cd8251166d..8063d0cf019ecc54a7bd4d8fa2ea313ea6810290 100644
--- a/matlab/optimization/mr_gstep.m
+++ b/matlab/optimization/mr_gstep.m
@@ -88,7 +88,7 @@ while i<n
                 b=[ones(7,1) xa xa.*xa./2]\fa;
                 gg(i)=x(i)*b(3)+b(2);
                 hh(i)=1/b(3);
-                [ff2, xx2, ~, ~] = csminit1(func0,x,penalty,f0,gg,0,diag(hh),Verbose,varargin{:});
+                [ff2, xx2] = csminit1(func0,x,penalty,f0,gg,0,diag(hh),Verbose,varargin{:});
                 if ff2<ff
                     ff=ff2;
                     xx=xx2;
diff --git a/matlab/optimization/numgrad5_.m b/matlab/optimization/numgrad5_.m
index 102bef77f7a0016509078cad98cd7aafdf229efb..551e75df547d39134d83503a46ab1717eedcb224 100644
--- a/matlab/optimization/numgrad5_.m
+++ b/matlab/optimization/numgrad5_.m
@@ -48,13 +48,13 @@ for i=1:n
     xiold = x(i);
     h = step_length_correction(xiold,scale,i)*delta;
     x(i) = xiold+h;
-    [f1,~,cost_flag1,] = penalty_objective_function(x, fcn, penalty, varargin{:});
+    [f1,~,~,] = penalty_objective_function(x, fcn, penalty, varargin{:});
     x(i) = xiold-h;
-    [f2,~,cost_flag2] = penalty_objective_function(x, fcn, penalty, varargin{:});
+    f2 = penalty_objective_function(x, fcn, penalty, varargin{:});
     x(i) = xiold+2*h;
-    [f3,~,cost_flag3] = penalty_objective_function(x, fcn, penalty, varargin{:});
+    f3 = penalty_objective_function(x, fcn, penalty, varargin{:});
     x(i) = xiold-2*h;
-    [f4,~,cost_flag4] = penalty_objective_function(x, fcn, penalty, varargin{:});
+    f4 = penalty_objective_function(x, fcn, penalty, varargin{:});
     if f0<f1 && f1<f3 && f0<f2 && f2<f4
         g0 = 0;
     else
@@ -68,4 +68,4 @@ for i=1:n
         g(i) = 0;
         badg = 1;
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/optimization/solve_one_boundary.m b/matlab/optimization/solve_one_boundary.m
index 876876d8c4fbc4515dedea851ba31447f94f3ae3..0ff3f1167bd977f1245bdd76c110f290812efe4e 100644
--- a/matlab/optimization/solve_one_boundary.m
+++ b/matlab/optimization/solve_one_boundary.m
@@ -182,8 +182,8 @@ for it_=start:incr:finish
                 g = (r'*g1)';
                 f = 0.5*r'*r;
                 p = -g1\r ;
-                [ya,f,r,check]=lnsrch1(ya,f,g,p,stpmax, ...
-                                       @lnsrch1_wrapper_one_boundary,nn, ...
+                ya = lnsrch1(ya,f,g,p,stpmax, ...
+                             @lnsrch1_wrapper_one_boundary,nn, ...
                                        nn, options_.solve_tolx, y_index_eq, fh, Block_Num, y, x, params, steady_state, T(:, it_), it_, M_);
                 dx = ya - y(y_index_eq, it_);
                 y(y_index_eq, it_) = ya;
diff --git a/matlab/parallel/AnalyseComputationalEnvironment.m b/matlab/parallel/AnalyseComputationalEnvironment.m
index 95f835b4270fd54ff1b86d6c6fb24bacee87bfd6..6722e45435ae21c78d10f47b5f0e265cce11e368 100644
--- a/matlab/parallel/AnalyseComputationalEnvironment.m
+++ b/matlab/parallel/AnalyseComputationalEnvironment.m
@@ -32,7 +32,7 @@ function [ErrorCode] = AnalyseComputationalEnvironment(DataInput, DataInputAdd)
 %   computer) is used.
 
 if ispc
-    [tempo, MasterName]=system('hostname');
+    [~, MasterName]=system('hostname');
     MasterName=deblank(MasterName);
 end
 
@@ -175,12 +175,12 @@ for Node=1:length(DataInput) % To obtain a recursive function remove the 'for'
 
         if Environment
             if OScallerWindows
-                [si1, de1]=system(['ping ', DataInput(Node).ComputerName]);
+                [si1]=system(['ping ', DataInput(Node).ComputerName]);
             else
-                [si1, de1]=system(['ping ', DataInput(Node).ComputerName, ' -c 4']);
+                [si1]=system(['ping ', DataInput(Node).ComputerName, ' -c 4']);
             end
         else
-            [si1, de1]=system(['ping ', DataInput(Node).ComputerName]);
+            [si1]=system(['ping ', DataInput(Node).ComputerName]);
         end
 
         if (si1)
@@ -387,15 +387,15 @@ for Node=1:length(DataInput) % To obtain a recursive function remove the 'for'
         else
             if ~strcmp(DataInput(Node).ComputerName,MasterName) % run on remote machine
                 if  strfind([DataInput(Node).MatlabOctavePath], 'octave') % Hybrid computing Matlab(Master)->Octave(Slaves) and Vice Versa!
-                    [NonServeS, NenServeD]=system(['start /B psexec \\',DataInput(Node).ComputerName,' -e -u ',DataInput(Node).UserName,' -p ',DataInput(Node).Password,' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' Tracing.m']);
+                    system(['start /B psexec \\',DataInput(Node).ComputerName,' -e -u ',DataInput(Node).UserName,' -p ',DataInput(Node).Password,' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' Tracing.m']);
                 else
-                    [NonServeS, NenServeD]=system(['start /B psexec \\',DataInput(Node).ComputerName,' -e -u ',DataInput(Node).UserName,' -p ',DataInput(Node).Password,' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' -nosplash -nodesktop -minimize -r Tracing']);
+                    system(['start /B psexec \\',DataInput(Node).ComputerName,' -e -u ',DataInput(Node).UserName,' -p ',DataInput(Node).Password,' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' -nosplash -nodesktop -minimize -r Tracing']);
                 end
             else % run on local machine via the network: user and passwd cannot be used!
                 if  strfind([DataInput(Node).MatlabOctavePath], 'octave') % Hybrid computing Matlab(Master)->Octave(Slaves) and Vice Versa!
-                    [NonServeS, NenServeD]=system(['start /B psexec \\',DataInput(Node).ComputerName,' -e ',' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' Tracing.m']);
+                    system(['start /B psexec \\',DataInput(Node).ComputerName,' -e ',' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' Tracing.m']);
                 else
-                    [NonServeS, NenServeD]=system(['start /B psexec \\',DataInput(Node).ComputerName,' -e ',' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' -nosplash -nodesktop -minimize -r Tracing']);
+                    system(['start /B psexec \\',DataInput(Node).ComputerName,' -e ',' -W ',DataInput(Node).RemoteDrive,':\',DataInput(Node).RemoteDirectory,'\',RemoteTmpFolder ' -low   ',DataInput(Node).MatlabOctavePath,' -nosplash -nodesktop -minimize -r Tracing']);
                 end
             end
 
@@ -570,7 +570,7 @@ for Node=1:length(DataInput) % To obtain a recursive function remove the 'for'
 
     % Questo controllo penso che si possa MIGLIORARE!!!!!
     if  isempty (RealCPUnbr) && Environment1==0
-        [si0, de0]=system(['psinfo \\',DataInput(Node).ComputerName]);
+        [~, de0]=system(['psinfo \\',DataInput(Node).ComputerName]);
     end
     RealCPUnbr=GiveCPUnumber(de0,Environment1);
 
@@ -629,4 +629,4 @@ for Node=1:length(DataInput) % To obtain a recursive function remove the 'for'
 
     disp(['Test for Cluster computation, computer ',DataInput(Node).ComputerName, ' ..... Passed!'])
     skipline(2)
-end
\ No newline at end of file
+end
diff --git a/matlab/parallel/InitializeComputationalEnvironment.m b/matlab/parallel/InitializeComputationalEnvironment.m
index 9d52cf015daedfcf143a141ddabd17469a966c10..ca29beb838093b92f44bbf38b689adb6f30b6b44 100644
--- a/matlab/parallel/InitializeComputationalEnvironment.m
+++ b/matlab/parallel/InitializeComputationalEnvironment.m
@@ -94,7 +94,7 @@ CPUWeightTemp=ones(1,lP)*(-1);
 CPUWeightTemp=CPUWeight;
 
 for i=1:lP
-    [NoNServes, mP]=max(CPUWeightTemp);
+    [~, mP]=max(CPUWeightTemp);
     NewPosition(i)=mP;
     CPUWeightTemp(mP)=-1;
 end
diff --git a/matlab/parallel/distributeJobs.m b/matlab/parallel/distributeJobs.m
index 336a898dbbfee084b7373e020054fb4f1c1d2e5c..9c010503662bfe82cbec597757de3c7807e63d4e 100644
--- a/matlab/parallel/distributeJobs.m
+++ b/matlab/parallel/distributeJobs.m
@@ -129,7 +129,7 @@ if SumOfJobs~=NumbersOfJobs
         % Many choices are possible:
         % - ... (see above).
 
-        [NonServe, VeryFast]= min(CPUWeight);
+        [~, VeryFast]= min(CPUWeight);
 
         while SumOfJobs<NumbersOfJobs
             JobsForNode(VeryFast)=JobsForNode(VeryFast)+1;
diff --git a/matlab/parallel/dynareParallelDelete.m b/matlab/parallel/dynareParallelDelete.m
index 21c14d3d52621b618f933bd4eee3ec595129ec79..2d4703c3ed2b31636c2493c77e326434276e748c 100644
--- a/matlab/parallel/dynareParallelDelete.m
+++ b/matlab/parallel/dynareParallelDelete.m
@@ -52,7 +52,7 @@ for indPC=1:length(Parallel)
         if ~isempty(directory)
             directory = [directory '/'];
         end
-        [~, ~] = system(['ssh ',ssh_token,username,Parallel(indPC).ComputerName,' ''/bin/bash --norc -c "rm -f ',directory,pname,fname,'"''']);
+        system(['ssh ',ssh_token,username,Parallel(indPC).ComputerName,' ''/bin/bash --norc -c "rm -f ',directory,pname,fname,'"''']);
     else
         fname_temp=['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',pname,fname];
         if exist(fname_temp,'file')
diff --git a/matlab/parallel/dynareParallelGetFiles.m b/matlab/parallel/dynareParallelGetFiles.m
index 349308c2c7159a65a22dcfbaaedbcc0e49eaf4be..5fc995d24d3bf7b362ea4a05c1caec2a51d5361b 100644
--- a/matlab/parallel/dynareParallelGetFiles.m
+++ b/matlab/parallel/dynareParallelGetFiles.m
@@ -61,13 +61,13 @@ for indPC=1:length(Parallel)
 
                     if isempty (FindAst)
 
-                        [NonServeL, NonServeR]= system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',NamFileInput{jfil,1}]);
+                        system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',NamFileInput{jfil,1}]);
 
                     else
 
                         filenameTemp=NamFileInput{jfil,2};
 
-                        [NotUsed, FlI]=system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' ls ',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',filenameTemp, ' 2> OctaveStandardOutputMessage.txt']);
+                        [~, FlI]=system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' ls ',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',filenameTemp, ' 2> OctaveStandardOutputMessage.txt']);
 
                         if isempty (FlI)
                             return
@@ -81,13 +81,13 @@ for indPC=1:length(Parallel)
                         for i=1: NumFileToCopy
                             Ni=num2str(i);
                             filenameTemp(1,AstPos)=Ni;
-                            [NonServeL, NonServeR]= system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},filenameTemp,' ',NamFileInput{jfil,1}]);
+                            system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},filenameTemp,' ',NamFileInput{jfil,1}]);
                         end
                     end
 
                 else
 
-                    [NonServeL, NonServeR]= system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',NamFileInput{jfil,1}]);
+                    system(['scp ',scp_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',NamFileInput{jfil,1}]);
                 end
 
             end
diff --git a/matlab/parallel/dynareParallelListAllFiles.m b/matlab/parallel/dynareParallelListAllFiles.m
index e91b131d7c81801291de113f525cb53c942357b5..bda2ae8a08c7f0fee278676c5da29e4fc6d5e953 100644
--- a/matlab/parallel/dynareParallelListAllFiles.m
+++ b/matlab/parallel/dynareParallelListAllFiles.m
@@ -44,7 +44,7 @@ if (~ispc || strcmpi('unix',Parallel.OperatingSystem))
     end
 
     % Get the data for the current remote directory.
-    [Flag, fL]=system(['ssh ',ssh_token,' ',' ',Parallel.UserName,'@',Parallel.ComputerName,' ls ',Parallel.RemoteDirectory,'/',PRCDir, ' -R -p -1']);
+    [~, fL]=system(['ssh ',ssh_token,' ',' ',Parallel.UserName,'@',Parallel.ComputerName,' ls ',Parallel.RemoteDirectory,'/',PRCDir, ' -R -p -1']);
 
     % Format the return value fL.
 
diff --git a/matlab/parallel/dynareParallelMkDir.m b/matlab/parallel/dynareParallelMkDir.m
index 3d908d81b8bdbf38503a77af972fb892abcdde40..a3f31e7aff71fc03fc60bf8ca8c9925909140bcb 100644
--- a/matlab/parallel/dynareParallelMkDir.m
+++ b/matlab/parallel/dynareParallelMkDir.m
@@ -42,9 +42,9 @@ for indPC=1:length(Parallel)
             else
                 ssh_token = '';
             end
-            [NonServeS, NonServeD]=system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' mkdir -p ',Parallel(indPC).RemoteDirectory,'/',PRCDir]);
+            system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' mkdir -p ',Parallel(indPC).RemoteDirectory,'/',PRCDir]);
         else
-            [NonServeS, NonServeD]=mkdir(['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir]);
+            mkdir(['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir]);
         end
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/parallel/dynareParallelRmDir.m b/matlab/parallel/dynareParallelRmDir.m
index cc5db203161a06535a4af13f7258dea46972d755..3c57742ec69900bc5cd230e9770a459e8d7bf5b4 100644
--- a/matlab/parallel/dynareParallelRmDir.m
+++ b/matlab/parallel/dynareParallelRmDir.m
@@ -65,10 +65,10 @@ for indPC=1:length(Parallel)
             else
                 ssh_token = '';
             end
-            [stat, NonServe] = system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' rm -fr ',Parallel(indPC).RemoteDirectory,'/',PRCDir,]);
+            system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' rm -fr ',Parallel(indPC).RemoteDirectory,'/',PRCDir,]);
             break
         else
-            [stat, mess, id] = rmdir(['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir],'s');
+            stat = rmdir(['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir],'s');
 
             if stat==1
                 break
@@ -81,4 +81,4 @@ for indPC=1:length(Parallel)
             end
         end
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/parallel/dynareParallelSendFiles.m b/matlab/parallel/dynareParallelSendFiles.m
index c0c49764a9d7a7df1a3ca9f10bb695e0b225de0b..3602493b68937dda129e1423ea7d469a21f3fca7 100644
--- a/matlab/parallel/dynareParallelSendFiles.m
+++ b/matlab/parallel/dynareParallelSendFiles.m
@@ -54,9 +54,9 @@ for indPC=1:length(Parallel)
             end
             for jfil=1:size(NamFileInput,1)
                 if ~isempty(NamFileInput{jfil,1})
-                    [NonServeL, NonServeR]=system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' mkdir -p ',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1}]);
+                    system(['ssh ',ssh_token,' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,' mkdir -p ',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1}]);
                 end
-                [NonServeL, NonServeR]=system(['scp ',scp_token,' ',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1}]);
+                system(['scp ',scp_token,' ',NamFileInput{jfil,1},NamFileInput{jfil,2},' ',Parallel(indPC).UserName,'@',Parallel(indPC).ComputerName,':',Parallel(indPC).RemoteDirectory,'/',PRCDir,'/',NamFileInput{jfil,1}]);
             end
         else
             for jfil=1:size(NamFileInput,1)
@@ -81,7 +81,7 @@ for indPC=1:length(Parallel)
                                 end
                             end
 
-                            [NonServeL, NonServeR]=system(['mkdir \\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInputTemp]);
+                            system(['mkdir \\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInputTemp]);
 
                         else
                             mkdir(['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInput{jfil,1}]);
@@ -107,7 +107,7 @@ for indPC=1:length(Parallel)
                         end
                     end
 
-                    [NonServeS, NonServeD]=system(['copy ',NamFileInputTemp, NamFileInput{jfil,2},' \\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInputTemp]);
+                    system(['copy ',NamFileInputTemp, NamFileInput{jfil,2},' \\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInputTemp]);
 
                 else
                     copyfile([NamFileInput{jfil,1},NamFileInput{jfil,2}],['\\',Parallel(indPC).ComputerName,'\',Parallel(indPC).RemoteDrive,'$\',Parallel(indPC).RemoteDirectory,'\',PRCDir,'\',NamFileInput{jfil,1}]);
diff --git a/matlab/parallel/masterParallel.m b/matlab/parallel/masterParallel.m
index ca0ec02bd5b564adb8e5f51c8019002c18b0ed6c..b569551b91535b0bfc9d1221aa525e07d876ae55 100644
--- a/matlab/parallel/masterParallel.m
+++ b/matlab/parallel/masterParallel.m
@@ -169,7 +169,7 @@ if parallel_recover ==0
     DyMo=pwd;
     % fInputVar.DyMo=DyMo;
     if ispc
-        [tempo, MasterName]=system('hostname');
+        [~, MasterName]=system('hostname');
         MasterName=deblank(MasterName);
     end
     % fInputVar.MasterName = MasterName;
@@ -893,7 +893,7 @@ switch Strategy
         end
 
         if isempty(dir('dynareParallelLogFiles'))
-            [A, B, C]=rmdir('dynareParallelLogFiles');
+            rmdir('dynareParallelLogFiles');
             mkdir('dynareParallelLogFiles');
         end
         try
@@ -911,7 +911,7 @@ switch Strategy
     delete(['temp_input.mat'])
     if newInstance
         if isempty(dir('dynareParallelLogFiles'))
-            [A, B, C]=rmdir('dynareParallelLogFiles');
+            rmdir('dynareParallelLogFiles');
             mkdir('dynareParallelLogFiles');
         end
     end
diff --git a/matlab/partial_information/PI_gensys.m b/matlab/partial_information/PI_gensys.m
index ec5b67a7975f226d0f5a6497c16a1efa399e4231..20a860ff45dd4767a626ba7eb180888911bd4ae4 100644
--- a/matlab/partial_information/PI_gensys.m
+++ b/matlab/partial_information/PI_gensys.m
@@ -76,7 +76,7 @@ try
         warning('off','MATLAB:nearlySingularMatrix');
         warning('off','MATLAB:singularMatrix');
         UAVinv=inv(C2); % i.e. inv(U02'*a1*V02)
-        [LastWarningTxt, LastWarningID]=lastwarn;
+        [~, LastWarningID]=lastwarn;
         if any(any(isinf(UAVinv)))
             singular=1;
         end
@@ -84,7 +84,7 @@ try
     if singular == 1 || strcmp('MATLAB:nearlySingularMatrix',LastWarningID) == 1 || ...
                  strcmp('MATLAB:illConditionedMatrix',LastWarningID)==1 || ...
                  strcmp('MATLAB:singularMatrix',LastWarningID)==1
-        [C1,C2,C3,C4, C5, F1, F2, F3, F4, F5, M1, M2, UAVinv, FL_RANK, V01, V02] = PI_gensys_singularC(C1,C2,C3,C4, C5, F1, F2, F3, F4, F5, V01, V02, 0);
+        [C1,~,C3,C4, C5, F1, F2, F3, F4, F5, ~, ~, UAVinv, FL_RANK, V01, V02] = PI_gensys_singularC(C1,C2,C3,C4, C5, F1, F2, F3, F4, F5, V01, V02, 0);
     end
     warning('on','MATLAB:singularMatrix');
     warning('on','MATLAB:nearlySingularMatrix');
@@ -239,7 +239,7 @@ for i=1:nn
 end
 div ;
 if ~zxz
-    [a, b, q, z]=qzdiv(div,a,b,q,z);
+    [a, b, ~, z]=qzdiv(div,a,b,q,z);
 end
 
 gev=[diag(a) diag(b)];
diff --git a/matlab/partial_information/PI_gensys_singularC.m b/matlab/partial_information/PI_gensys_singularC.m
index fbc77f3c5df7e34c9cac07e3b95b0ade053c1ce2..15a2a6f75ccee561fcbb43b4f96529a305bb9e8a 100644
--- a/matlab/partial_information/PI_gensys_singularC.m
+++ b/matlab/partial_information/PI_gensys_singularC.m
@@ -34,7 +34,7 @@ M1=[];M2=[]; UAVinv=[];
 % Find SVD of a0, and create partitions of U, S and V
 %
 
-[J0,K0,L0] = svd(C2in);
+[J0,K0] = svd(C2in);
 n=size(C2in,1);
 K_RANK=rank(K0);
 J2=J0(1:n,K_RANK+1:n);
@@ -71,7 +71,7 @@ try
         singular=1;
     else
         UAVinv=inv(C2);
-        [LastWarningTxt, LastWarningID]=lastwarn;
+        [~, LastWarningID]=lastwarn;
         if any(any(isinf(UAVinv)))
             singular=1;
         end
@@ -83,7 +83,7 @@ try
         [C1,C2,C3,C4, C5, F1, F2, F3, F4, F5, M1, M2, UAVinv, FL_RANK, V01, V02] = PI_gensys_singularC(C1,C2,C3,C4, C5, F1, F2, F3, F4, F5, V01, V02, level);
     end
 catch
-    [errmsg, errcode]=lasterr;
+    errmsg = lasterr;
     warning(['error callig PI_gensys_singularC: ' errmsg ],'errcode');
     error('errcode',['error callig PI_gensys_singularC: ' errmsg ]);
 end
diff --git a/matlab/partial_information/disc_riccati_fast.m b/matlab/partial_information/disc_riccati_fast.m
index 3a9283c394ae02b127aa5a1ed9ae0e62877e7838..dabc2a6531129ba4f643579dd3f18e0c1c44cc0f 100644
--- a/matlab/partial_information/disc_riccati_fast.m
+++ b/matlab/partial_information/disc_riccati_fast.m
@@ -84,7 +84,7 @@ clear X0 X1 Y0 Y1 P1 I INVPY;
 
 % Check that X is positive definite
 if flag_ch==1
-    [C,p]=chol(Z);
+    [~,p]=chol(Z);
     if p ~= 0
         error('Z is not positive definite')
     end
diff --git a/matlab/partial_information/dr1_PI.m b/matlab/partial_information/dr1_PI.m
index dc62187cfbddc7203a48476f133a051b4b729f2e..262de852845020504ff3d3aee958c3d0a4ca1857 100644
--- a/matlab/partial_information/dr1_PI.m
+++ b/matlab/partial_information/dr1_PI.m
@@ -124,7 +124,7 @@ else
         lq_instruments.tct_ruleids=tct_ruleids;
         %if isfield(lq_instruments,'xsopt_SS') %% changed by BY
         [~, lq_instruments.xsopt_SS,lq_instruments.lmopt_SS,s2,check] = opt_steady_get;%% changed by BY
-        [qc, DYN_Q] = QPsolve(lq_instruments, s2, check); %% added by BY
+        [~, DYN_Q] = QPsolve(lq_instruments, s2, check); %% added by BY
         z = repmat(lq_instruments.xsopt_SS,1,klen);
     else
         z = repmat(dr.ys,1,klen);
diff --git a/matlab/partial_information/qzdiv.m b/matlab/partial_information/qzdiv.m
index d36433b8943300e2a318972b50477d59a7f32a1a..4c802eb8e50b87942a3aa05e5a31ea35afb9782c 100644
--- a/matlab/partial_information/qzdiv.m
+++ b/matlab/partial_information/qzdiv.m
@@ -27,7 +27,7 @@ function [A,B,Q,Z] = qzdiv(stake,A,B,Q,Z)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[n, jnk] = size(A);
+[n, ~] = size(A);
 root = abs([diag(A) diag(B)]);
 root(:,1) = root(:,1)-(root(:,1)<1.e-13).*(root(:,1)+root(:,2));
 root(:,2) = root(:,2)./root(:,1);
diff --git a/matlab/perfect-foresight-models/basic_plan.m b/matlab/perfect-foresight-models/basic_plan.m
index 22b929862d3ed6e1c4a5499612b0704fd17f2718..0605a060f54e5c507af3ed57b374ed1ddd8d17b8 100644
--- a/matlab/perfect-foresight-models/basic_plan.m
+++ b/matlab/perfect-foresight-models/basic_plan.m
@@ -57,7 +57,7 @@ if ~isempty(plan.options_cond_fcst_.controlled_varexo)
     if ~isempty(common_var)
         common_date = intersect(date, plan.constrained_date_{common_var});
         if ~isempty(common_date)
-            [date_, i_date] = setdiff(date, common_date);
+            [~, i_date] = setdiff(date, common_date);
             value = value(i_date);
             if common_date.length > 1
                 the_dates = [cell2mat(strings(common_date(1))) ':' cell2mat(strings(common_date(end)))];
diff --git a/matlab/perfect-foresight-models/det_cond_forecast.m b/matlab/perfect-foresight-models/det_cond_forecast.m
index 2fa7b9b1af57009f7191a545d7901e7ddc9fb612..89b723e05fc56c40595255b7a4327fa1d9415f38 100644
--- a/matlab/perfect-foresight-models/det_cond_forecast.m
+++ b/matlab/perfect-foresight-models/det_cond_forecast.m
@@ -47,7 +47,7 @@ if ~isfield(oo_,'dr') || ~isfield(oo_.dr,'ghx')
     dr = struct();
     oo_.dr=set_state_space(dr,M_);
     options_.order = 1;
-    [oo_.dr,Info,M_.params] = resol(0,M_,options_,oo_.dr ,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
+    [oo_.dr,~,M_.params] = resol(0,M_,options_,oo_.dr ,oo_.steady_state, oo_.exo_steady_state, oo_.exo_det_steady_state);
     fprintf('done\n');
 end
 b_surprise = 0;
diff --git a/matlab/perfect-foresight-models/solve_two_boundaries_stacked.m b/matlab/perfect-foresight-models/solve_two_boundaries_stacked.m
index 829cd19a85d324efbfab1aaff277a61a264c5b90..d34ee78422313ee3e988376f8c85ef2cfe7e36dd 100644
--- a/matlab/perfect-foresight-models/solve_two_boundaries_stacked.m
+++ b/matlab/perfect-foresight-models/solve_two_boundaries_stacked.m
@@ -262,7 +262,7 @@ while ~(cvg || iter > options_.simul.maxit)
             g = (ra'*g1a)';
             f = 0.5*ra'*ra;
             p = -g1a\ra;
-            [yn,f,ra,check]=lnsrch1(ya,f,g,p,stpmax,@lnsrch1_wrapper_two_boundaries,nn,nn, options_.solve_tolx, fh, Block_Num, y, y_index,x, M_.params, steady_state, T, periods, Blck_size, M_);
+            yn = lnsrch1(ya,f,g,p,stpmax,@lnsrch1_wrapper_two_boundaries,nn,nn, options_.solve_tolx, fh, Block_Num, y, y_index,x, M_.params, steady_state, T, periods, Blck_size, M_);
             dx = ya - yn;
             y(y_index, y_kmin+(1:periods))=reshape(yn',length(y_index),periods);
         end
diff --git a/matlab/read_variables.m b/matlab/read_variables.m
index 8cf3d941c4e1109645a74e0f511e13c1f87cbf12..9cc92f8631d1b1b881e20e61869df65ca9cf35a9 100644
--- a/matlab/read_variables.m
+++ b/matlab/read_variables.m
@@ -89,7 +89,7 @@ switch (extension)
         dyn_data_01(:,dyn_i_01) = dyn_tmp_01;
     end
   case { '.xls', '.xlsx' }
-    [freq,init,data,varlist] = load_xls_file_data(fullname,xls_sheet,xls_range);
+    [~,~,data,varlist] = load_xls_file_data(fullname,xls_sheet,xls_range);
     for dyn_i_01=1:var_size_01
         iv = strmatch(strtrim(var_names_01{dyn_i_01}),varlist,'exact');
         if ~isempty(iv)
@@ -105,7 +105,7 @@ switch (extension)
         end
     end
   case '.csv'
-    [freq,init,data,varlist] = load_csv_file_data(fullname);
+    [~,~,data,varlist] = load_csv_file_data(fullname);
     for dyn_i_01=1:var_size_01
         iv = strmatch(var_names_01{dyn_i_01},varlist,'exact');
         if ~isempty(iv)
diff --git a/matlab/sample_hp_filter.m b/matlab/sample_hp_filter.m
index 93ed0566aa09c8dbf0974afcd299f3bf108536e9..f917aa56da4a05ba73793d99c99463bdd1a41e3c 100644
--- a/matlab/sample_hp_filter.m
+++ b/matlab/sample_hp_filter.m
@@ -29,7 +29,7 @@ function [hptrend,hpcycle] = sample_hp_filter(y,s)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[T,n] = size(y);
+[T,~] = size(y);
 
 if nargin<2 || isempty(s)
     s = 1600;
diff --git a/matlab/shock_decomposition/annualized_shock_decomposition.m b/matlab/shock_decomposition/annualized_shock_decomposition.m
index 4654adcc67f2db552a50eda01b6091e78ae47539..d9180cc31ef8c347f83d75f7c05408c09a4ed14b 100644
--- a/matlab/shock_decomposition/annualized_shock_decomposition.m
+++ b/matlab/shock_decomposition/annualized_shock_decomposition.m
@@ -129,7 +129,7 @@ if realtime_==0
         myopts=options_;
         myopts.plot_shock_decomp.type='qoq';
         myopts.plot_shock_decomp.realtime=0;
-        [z, ~] = plot_shock_decomposition(M_,oo_,myopts,[]);
+        z = plot_shock_decomposition(M_,oo_,myopts,[]);
     else
         z = oo_;
     end
@@ -158,7 +158,7 @@ if realtime_ && isstruct(oo_) && isfield(oo_, 'realtime_shock_decomposition')
         myopts.plot_shock_decomp.realtime=1;
         myopts.plot_shock_decomp.vintage=i;
         % retrieve quarterly shock decomp
-        [z, ~] = plot_shock_decomposition(M_,oo_,myopts,[]);
+        z = plot_shock_decomposition(M_,oo_,myopts,[]);
         zdim = size(z);
         z = z(i_var,:,:);
         if isstruct(aux)
@@ -185,7 +185,7 @@ if realtime_ && isstruct(oo_) && isfield(oo_, 'realtime_shock_decomposition')
             if qvintage_>i-4 && qvintage_<i
                 myopts.plot_shock_decomp.vintage=qvintage_;
                 % retrieve quarterly shock decomp
-                [z, ~] = plot_shock_decomposition(M_,oo_,myopts,[]);
+                z = plot_shock_decomposition(M_,oo_,myopts,[]);
                 z(:,:,end+1:zdim(3))=nan; % fill with nan's remaining time points to reach Q4
                 z = z(i_var,:,:);
                 if isstruct(aux)
diff --git a/matlab/shock_decomposition/expand_group.m b/matlab/shock_decomposition/expand_group.m
index ba953230a8fd670cd45ec311477d601e1da0054e..332fc32efd44c522874fea2d179ed4228b9bf105 100644
--- a/matlab/shock_decomposition/expand_group.m
+++ b/matlab/shock_decomposition/expand_group.m
@@ -31,7 +31,7 @@ if nargin<4
     no_graph=0;
 end
 filename = get(gcf,'filename');
-[filepath, name, ext]=fileparts(filename);
+filepath = fileparts(filename);
 M_ = evalin('base','M_');
 oo_ = evalin('base','oo_');
 options_ = evalin('base','options_');
diff --git a/matlab/shock_decomposition/initial_condition_decomposition.m b/matlab/shock_decomposition/initial_condition_decomposition.m
index c130869f937f777457f02991dfffb7d49130ae61..1004645c1c9ba2a3bda999e53ea3265ef787eb00 100644
--- a/matlab/shock_decomposition/initial_condition_decomposition.m
+++ b/matlab/shock_decomposition/initial_condition_decomposition.m
@@ -64,7 +64,7 @@ if isempty(varlist)
 end
 
 if ~isequal(varlist,0)
-    [i_var, nvar, index_uniques] = varlist_indices(varlist, M_.endo_names);
+    [~, ~, index_uniques] = varlist_indices(varlist, M_.endo_names);
     varlist = varlist(index_uniques);
 end
 
@@ -99,7 +99,7 @@ if ~isfield(oo_,'initval_decomposition') || isequal(varlist,0)
     with_epilogue = options_.initial_condition_decomp.with_epilogue;
     options_.selected_variables_only = 0; %make sure all variables are stored
     options_.plot_priors=0;
-    [oo_local,M,~,~,Smoothed_Variables_deviation_from_mean] = evaluate_smoother(parameter_set,varlist,M_,oo_,options_,bayestopt_,estim_params_);
+    [oo_local,~,~,~,Smoothed_Variables_deviation_from_mean] = evaluate_smoother(parameter_set,varlist,M_,oo_,options_,bayestopt_,estim_params_);
 
     % reduced form
     dr = oo_local.dr;
diff --git a/matlab/shock_decomposition/plot_shock_decomposition.m b/matlab/shock_decomposition/plot_shock_decomposition.m
index 4be4cdf8700b62fdfd0aea771980103085907806..9dcc361346bab45542ba075c50fb92f23f706e31 100644
--- a/matlab/shock_decomposition/plot_shock_decomposition.m
+++ b/matlab/shock_decomposition/plot_shock_decomposition.m
@@ -471,7 +471,7 @@ switch type
                 end
                 i_var0 = i_var;
                 steady_state_0 = steady_state;
-                [za, endo_names, endo_names_tex, steady_state, i_var, oo_] = ...
+                [za, endo_names, endo_names_tex, steady_state, i_var] = ...
                     annualized_shock_decomposition(z,M_, options_, i_var, t0, options_.nobs, realtime_, vintage_, steady_state,q2a);
                 if options_.plot_shock_decomp.interactive && ~isempty(options_.plot_shock_decomp.use_shock_groups)
                     mygroup = options_.plot_shock_decomp.use_shock_groups;
diff --git a/matlab/stochastic_solver/stochastic_solvers.m b/matlab/stochastic_solver/stochastic_solvers.m
index 4e34396d5e7d8179f8f38c6d23af5a6a9d5cbaad..6c6f985828a9638a8d5968ec4ea0a7c805d74cd5 100644
--- a/matlab/stochastic_solver/stochastic_solvers.m
+++ b/matlab/stochastic_solver/stochastic_solvers.m
@@ -120,7 +120,7 @@ elseif local_order == 2
     [~,jacobia_,hessian1] = feval([M_.fname '.dynamic'],z(iyr0),...
                                   exo_simul, ...
                                   M_.params, dr.ys, it_);
-    [infrow,infcol]=find(isinf(hessian1));
+    [infrow, ~] = find(isinf(hessian1));
     if options_.debug
         if ~isempty(infrow)
             fprintf('\nSTOCHASTIC_SOLVER: The Hessian of the dynamic model contains Inf.\n')
@@ -132,7 +132,7 @@ elseif local_order == 2
         info(1)=11;
         return
     end
-    [nanrow,nancol]=find(isnan(hessian1));
+    [nanrow, ~] = find(isnan(hessian1));
     if options_.debug
         if ~isempty(nanrow)
             fprintf('\nSTOCHASTIC_SOLVER: The Hessian of the dynamic model contains NaN.\n')
@@ -343,8 +343,8 @@ end
 
 if options_.loglinear
     % this needs to be extended for order=2,3
-    [il,il1,ik,k1] = indices_lagged_leaded_exogenous_variables(dr.order_var,M_);
-    [illag,illag1,iklag,klag1] = indices_lagged_leaded_exogenous_variables(dr.order_var(M_.nstatic+(1:M_.nspred)),M_);
+    [il,~,ik,k1] = indices_lagged_leaded_exogenous_variables(dr.order_var,M_);
+    [illag,~,iklag,klag1] = indices_lagged_leaded_exogenous_variables(dr.order_var(M_.nstatic+(1:M_.nspred)),M_);
     if ~isempty(ik)
         if M_.nspred > 0
             dr.ghx(ik,iklag) = repmat(1./dr.ys(k1),1,length(klag1)).*dr.ghx(ik,iklag).* ...
diff --git a/matlab/user_has_matlab_license.m b/matlab/user_has_matlab_license.m
index b2eaf20c422fff4183f3839b6e1da2d9ccb59c6b..354b03344af8e1b85681f615e1295c32b8f357f8 100644
--- a/matlab/user_has_matlab_license.m
+++ b/matlab/user_has_matlab_license.m
@@ -28,7 +28,7 @@ function [hasLicense] = user_has_matlab_license(toolbox)
 % You should have received a copy of the GNU General Public License
 % along with Dynare.  If not, see <https://www.gnu.org/licenses/>.
 
-[hasLicense, ~] = license('checkout',toolbox);
+hasLicense = license('checkout',toolbox);
 
 if ~hasLicense
     return
diff --git a/matlab/utilities/dataset/makedataset.m b/matlab/utilities/dataset/makedataset.m
index f48e5837f4b1a5cd3ffed339df1374fd9c4050b5..4303b3df8b25a9db60b56e0823fbc269f379bd2f 100644
--- a/matlab/utilities/dataset/makedataset.m
+++ b/matlab/utilities/dataset/makedataset.m
@@ -91,7 +91,7 @@ if ~isempty(datafile)
     datafile_extension = get_file_extension(datafile);
     if isempty(datafile_extension)
         available_extensions = {}; j = 1;
-        [datafilepath,datafilename,datafileext] = fileparts(datafile);
+        [datafilepath, datafilename] = fileparts(datafile);
         if isempty(datafilepath)
             datafilepath = '.';
         end
diff --git a/matlab/utilities/doc/dynInfo.m b/matlab/utilities/doc/dynInfo.m
index 434cfbfc6eb9ad4ccb907872cee167ef87c4bd3e..c4c105200ab7235b0e2f5a2778c238c393144a1c 100644
--- a/matlab/utilities/doc/dynInfo.m
+++ b/matlab/utilities/doc/dynInfo.m
@@ -51,7 +51,7 @@ function dynInfo(fun)
 % Original author: stephane DOT adjemian AT univ DASH lemans DOT fr
 
 if isempty(strfind(fun,'@')) & (~isempty(strfind(fun,'/')) || ~isempty(strfind(fun,'\')) )
-    [pathstr1, name, ext] = fileparts(fun);
+    pathstr1 = fileparts(fun);
     addpath(pathstr1);
     rm_path = 1;
 else
@@ -80,4 +80,4 @@ end
 
 if rm_path
     rmpath(pathstr1)
-end
\ No newline at end of file
+end
diff --git a/matlab/utilities/general/clean_current_folder.m b/matlab/utilities/general/clean_current_folder.m
index 8243590128f3f9d7d6b4a8071e51488fe6b9d905..f458d6a641a1cc8459e31923ab2e80164f8ab242 100644
--- a/matlab/utilities/general/clean_current_folder.m
+++ b/matlab/utilities/general/clean_current_folder.m
@@ -21,7 +21,7 @@ a = dir('*.mod');
 
 
 for i = 1:length(a)
-    [~,basename,extension] = fileparts(a(i).name);
+    [~,basename] = fileparts(a(i).name);
     if exist([basename '.m'])
         delete([basename '.m']);
     end
@@ -38,4 +38,4 @@ for i = 1:length(a)
     if exist(['protect_' basename '_steadystate.m'])
         movefile(['protect_' basename '_steadystate.m'],[basename '_steadystate.m']);
     end
-end
\ No newline at end of file
+end
diff --git a/matlab/utilities/tests b/matlab/utilities/tests
index 3c4079a8e1e495be50be3cd81855eb37189fceab..a6b97581c379bde8e5022e230075c18421a8094a 160000
--- a/matlab/utilities/tests
+++ b/matlab/utilities/tests
@@ -1 +1 @@
-Subproject commit 3c4079a8e1e495be50be3cd81855eb37189fceab
+Subproject commit a6b97581c379bde8e5022e230075c18421a8094a
diff --git a/matlab/varlist_indices.m b/matlab/varlist_indices.m
index 3ea9a2fcd48b51790c2a0367cfa229ac60ec1217..7aee6877f04b8d5a76544e31ebfb59ef33f2ff83 100644
--- a/matlab/varlist_indices.m
+++ b/matlab/varlist_indices.m
@@ -56,7 +56,7 @@ if ~all(check)
 end
 
 nvar_present = length(i_var(check));
-[~, index_unique, ~] = unique(i_var, 'first');
+[~, index_unique] = unique(i_var, 'first');
 index_unique_present = index_unique(~ismember(index_unique,indices_not_present));
 index_unique_present = sort(index_unique_present);
 i_var_unique_present = i_var(index_unique_present);
diff --git a/tests/analytic_derivatives/BrockMirman_PertParamsDerivs.mod b/tests/analytic_derivatives/BrockMirman_PertParamsDerivs.mod
index 9c181990a94354a6f1bf84d2726af047040912d0..9a982fb000b2c6724bc1117e329deedffae874ec 100644
--- a/tests/analytic_derivatives/BrockMirman_PertParamsDerivs.mod
+++ b/tests/analytic_derivatives/BrockMirman_PertParamsDerivs.mod
@@ -131,7 +131,7 @@ identification(order=@{ORDER},nograph,no_identification_strength);
 indpmodel  = estim_params_.param_vals(:,1);
 indpstderr = estim_params_.var_exo(:,1);
 indpcorr   = estim_params_.corrx(:,1:2);
-[I,~]      = find(M_.lead_lag_incidence');
+[I, ~] = find(M_.lead_lag_incidence');
 
 %% Parameter derivatives of perturbation 
 @#if CREATE_SYMBOLIC == 1
diff --git a/tests/trend-component-and-var-models/vm4.mod b/tests/trend-component-and-var-models/vm4.mod
index 17c4ac3e64781a912bb71b1560ad14aab1c5c1aa..8bb16e73e517717dfd012e9d8dec2125af4d0988 100644
--- a/tests/trend-component-and-var-models/vm4.mod
+++ b/tests/trend-component-and-var-models/vm4.mod
@@ -43,10 +43,10 @@ shocks;
     var ez = 1.0;
 end;
 
-[~, b0, ~] = get_companion_matrix_legacy('toto');
+[~, b0] = get_companion_matrix_legacy('toto');
 
-[~, ~, b1, ~] = get_companion_matrix('toto');
+[~, ~, b1] = get_companion_matrix('toto');
 
 if any(abs(b0(:)-b1(:))>1e-9)
    error('get_companion_matrix and get_comapnion_matrix_legacy do not return the same AR matrices.')
-end
\ No newline at end of file
+end