1% STK_PARAM_ESTIM_OPTIMIZE [STK internal] 2% 3% INTERNAL FUNCTION WARNING: 4% 5% This function is currently considered as internal. API-breaking 6% changes are very likely to happen in future releases. 7 8% Copyright Notice 9% 10% Copyright (C) 2015-2019 CentraleSupelec 11% Copyright (C) 2014 Ashwin Ravisankar 12% Copyright (C) 2011-2014 SUPELEC 13% 14% Authors: Julien Bect <julien.bect@centralesupelec.fr> 15% Emmanuel Vazquez <emmanuel.vazquez@centralesupelec.fr> 16% Ashwin Ravisankar <ashwinr1993@gmail.com> 17 18% Copying Permission Statement 19% 20% This file is part of 21% 22% STK: a Small (Matlab/Octave) Toolbox for Kriging 23% (http://sourceforge.net/projects/kriging) 24% 25% STK is free software: you can redistribute it and/or modify it under 26% the terms of the GNU General Public License as published by the Free 27% Software Foundation, either version 3 of the License, or (at your 28% option) any later version. 29% 30% STK is distributed in the hope that it will be useful, but WITHOUT 31% ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 32% or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public 33% License for more details. 34% 35% You should have received a copy of the GNU General Public License 36% along with STK. If not, see <http://www.gnu.org/licenses/>. 37 38function [model_opt, info] = stk_param_estim_optim ... 39 (model0, xi, zi, criterion, covparam_select, noiseparam_select) 40 41select = [covparam_select; noiseparam_select]; 42 43% Starting point 44v0 = stk_get_optimizable_parameters (model0); 45w0 = v0(select); 46 47% Bounds 48% FIXME: this could (should) be implemented directly for models 49[covparam_lb, covparam_ub] = stk_param_getdefaultbounds (model0.covariance_type, model0.param, xi, zi); 50[covparam_lb, covparam_ub] = select_bounds (covparam_lb, covparam_ub, covparam_select); 51[noiseparam_lb, noiseparam_ub] = stk_param_getdefaultbounds_lnv (model0, model0.lognoisevariance, xi, zi); 52[noiseparam_lb, noiseparam_ub] = select_bounds (noiseparam_lb, noiseparam_ub, noiseparam_select); 53lb = [covparam_lb; noiseparam_lb]; 54ub = [covparam_ub; noiseparam_ub]; 55 56% Define objective function 57f = @(v)(crit_wrapper (model0, v, xi, zi, criterion, covparam_select, noiseparam_select)); 58 59bounds_available = (~ isempty (lb)) && (~ isempty (ub)); 60 61% Sanity check (part 1) 62crit0 = f (w0); 63if ~ (isscalar (crit0) && isfinite (crit0)) 64 errmsg = '*** PANIC: crit0 is not a finite scalar value. ***'; 65 stk_error (errmsg, 'OptimizationFailure'); 66end 67 68if bounds_available 69 A = stk_options_get ('stk_param_estim', 'minimize_box'); 70 [w_opt, crit_opt] = stk_minimize_boxconstrained (A, f, w0, lb, ub); 71else 72 A = stk_options_get ('stk_param_estim', 'minimize_unc'); 73 [w_opt, crit_opt] = stk_minimize_unconstrained (A, f, w0); 74end 75 76% Sanity check (part 2) 77if crit0 < crit_opt 78 s1 = '*** PANIC: Something went SERIOUSLY WRONG during the optimization ***'; 79 s2 = sprintf ('crit0 = %f, crit_opt = %f: crit0 < crit_opt', crit0, crit_opt); 80 stk_error (sprintf ('%s\n%s\n', s1, s2), 'OptimizationFailure'); 81end 82 83% Create outputs 84v_opt = v0; 85v_opt(select) = w_opt; 86model_opt = stk_set_optimizable_parameters (model0, v_opt); 87 88% Create 'info' structure, if requested 89if nargout > 1 90 info.criterion = criterion; 91 info.crit_opt = crit_opt; 92 info.starting_point = w0; 93 info.final_point = w_opt; 94 info.lower_bounds = lb; 95 info.upper_bounds = ub; 96 info.param_select = covparam_select; 97 info.noiseparam_select = noiseparam_select; 98end 99 100end % function 101 102%#ok<*CTCH,*LERR,*SPWRN,*WNTAG> 103 104 105function [C, dC] = crit_wrapper ... 106 (model, w, xi, zi, criterion, covparam_select, noise_select) 107 108v = stk_get_optimizable_parameters (model); 109v([covparam_select; noise_select]) = w; 110model = stk_set_optimizable_parameters (model, v); 111 112if nargout == 1 113 114 % Compute only the value of the criterion 115 C = criterion (model, xi, zi); 116 117elseif any (noise_select) 118 119 % Compute the value of the criterion and the gradients 120 % FIXME: We might be computing a lot of derivatives that we don't really need... 121 [C, dC_param, dC_lnv] = criterion (model, xi, zi); 122 123 dC = [dC_param(covparam_select); dC_lnv(noise_select)]; 124 125else 126 127 % Compute the value of the criterion and the gradients 128 % FIXME: We might be computing a lot of derivatives that we don't really need... 129 [C, dC_param] = criterion (model, xi, zi); 130 131 dC = dC_param(covparam_select); 132 133end 134 135end % function 136 137 138function [lb, ub] = select_bounds (lb, ub, select) 139 140if ~ isempty (lb) 141 lb = lb(select); 142end 143 144if ~ isempty (ub) 145 ub = ub(select); 146end 147 148end % function 149