diff --git a/matlab/check_posterior_sampler_options.m b/matlab/check_posterior_sampler_options.m
index d518abedbeb5c56e94bc7b87974a11f7c025408f..13e9307a2794965e3af23a7b69820524baf3c8af 100644
--- a/matlab/check_posterior_sampler_options.m
+++ b/matlab/check_posterior_sampler_options.m
@@ -301,6 +301,12 @@ if init
                   case 'save_tmp_file'
                     posterior_sampler_options.save_tmp_file = options_list{i,2};
 
+                  case 'maximize'
+                    posterior_sampler_options.maximize = options_list{i,2};
+                    
+                  case 'mode_compute'
+                    posterior_sampler_options.mode_compute = options_list{i,2};
+                    
                   otherwise
                     warning(['slice_sampler: Unknown option (' options_list{i,1} ')!'])
                 end
diff --git a/matlab/default_option_values.m b/matlab/default_option_values.m
index 0ebbe0e54a52406c02e91366656c9ad2fd5dddf8..bc55e9d222a62f3e8aa929cf0afad6f585fe3144 100644
--- a/matlab/default_option_values.m
+++ b/matlab/default_option_values.m
@@ -470,6 +470,8 @@ options_.posterior_sampler_options.slice.use_mh_covariance_matrix=0;
 options_.posterior_sampler_options.slice.WR=[];
 options_.posterior_sampler_options.slice.mode_files=[];
 options_.posterior_sampler_options.slice.mode=[];
+options_.posterior_sampler_options.slice.maximize=0;
+options_.posterior_sampler_options.slice.mode_compute=5;
 options_.posterior_sampler_options.slice.initial_step_size=0.8;
 options_.posterior_sampler_options.slice.save_tmp_file=1;
 % Independent Metropolis-Hastings
diff --git a/matlab/posterior_sampler_iteration.m b/matlab/posterior_sampler_iteration.m
index 185e990ddf696bff8c51a2d24350ba17a363858f..ae1f36990905481de25a5849ae05db49cadfc72b 100644
--- a/matlab/posterior_sampler_iteration.m
+++ b/matlab/posterior_sampler_iteration.m
@@ -50,8 +50,51 @@ mh_bounds = sampler_options.bounds;
 switch posterior_sampling_method
   case 'slice'
 
+      if  ~sampler_options.maximize
+          
     [par, logpost, neval] = slice_sampler(TargetFun,last_draw, [mh_bounds.lb mh_bounds.ub], sampler_options,varargin{:});
     accepted = 1;
+          
+      else
+          options_=varargin{3};
+          bayestopt_=varargin{6};
+          options_.mode_compute = sampler_options.mode_compute;
+          if options_.mode_compute==5
+              %get whether outer product Hessian is requested
+              newratflag=[];
+              if ~isempty(options_.optim_opt)
+                  options_list = read_key_value_string(options_.optim_opt);
+                  for i=1:rows(options_list)
+                      if strcmp(options_list{i,1},'Hessian')
+                          newratflag=options_list{i,2};
+                      end
+                  end
+              end
+              if options_.analytic_derivation
+                  options_analytic_derivation_old = options_.analytic_derivation;
+                  options_.analytic_derivation = -1;
+                  if ~isempty(newratflag) && newratflag~=0 %numerical hessian explicitly specified
+                      error('newrat: analytic_derivation is incompatible with numerical Hessian.')
+                  else %use default
+                      newratflag=0; %exclude DYNARE numerical hessian
+                  end
+              elseif ~options_.analytic_derivation
+                  if isempty(newratflag)
+                      newratflag=options_.newrat.hess; %use default numerical dynare hessian
+                  end
+              end
+          end
+          options_.newrat.Save_files = sampler_options.curr_block;
+          [par, fval, exitflag, hess_mat_optimizer, options_, Scale, new_rat_hess_info] = ...
+              dynare_minimize_objective(TargetFun,last_draw(:),options_.mode_compute,options_,[mh_bounds.lb mh_bounds.ub],bayestopt_.name,bayestopt_,[], ...
+                  varargin{:}); %inputs for objective
+
+          logpost = -fval;
+          accepted = 1;
+          neval = 1;
+          
+      end
+  
   case 'random_walk_metropolis_hastings'
     neval = 1;
     ProposalFun = sampler_options.proposal_distribution;