How to declare the type of a variable inside the stan::math:hessian functional

General
1

Hello. I’m trying to get the Hessian with respect to some parameters of a posterior density, and I’m using stan::math::hessian function to do so. My point is that I don’t know exactly how to declare the type of some variables that are functions of the parameters I want to differentiate.

Just to exemplify, for the vector of parameters “f” (that will be differentiated) and the fixed matrix “K” let’s suppose that we want to get the Hessian for g(f,K)=f^{T}Kf.

The code to perform it can be

// [[Rcpp::export]]
Eigen::MatrixXd Hf(Eigen::VectorXd f, Eigen::VectorXd K){
  double fx;
  Eigen::VectorXd grad_fx;
  Eigen::MatrixXd H;
  stan::math::hessian([&K](auto f){return(f.transpose()*K*f);},
    f, fx, grad_fx, H);
  return H;
}

and it will work fine. I can rewrite it like this

// [[Rcpp::export]]
Eigen::MatrixXd Hf(Eigen::VectorXd f, Eigen::VectorXd K){
  double fx;
  Eigen::VectorXd grad_fx;
  Eigen::MatrixXd H;
  stan::math::hessian([&K](auto f){
  auto aux = f.transpose()*K*f; //****
  return(aux);},
    f, fx, grad_fx, H);
  return H;
}

and it will work fine too, but I wouldn’t know how to declare the type of the variable “aux” in line marked with **** in specific way (not “auto”). How could I declare it?

For the functional stan::math::gradient I declared as stan::math::var and it worked fine, but for the Hessian it didn’t.

Thank’s in advance.

1 Like
2

Hi, sorry for not getting to you eariler.

For Hessian, I think this is nested reverse mode in forward mode autodiff, so you would need either Eigen::Matrix<fvar<var>, Eigen::Dynamic, 1> or maybe just fvar<var> (this is how it is called from stan::math::hessian: math/hessian.hpp at 92075708b1d1796eb82e3b284cd11e544433518e · stan-dev/math · GitHub)

Hope that helps at least a bit!

1 Like
3

Thank you to reply Martin! As you suggested, I declared the variable aux as stan::math::fvar<stan::math::var> but it was given an error message that I attatched in the file “error-example.txt”.

I’ll also show the full “.cpp” file just for you to confirm if I’m missing to include any package, etc.

// [[Rcpp::depends(BH)]]
// [[Rcpp::depends(RcppEigen)]]
// [[Rcpp::depends(RcppParallel)]]
// [[Rcpp::depends(StanHeaders)]]

#include <stan/math/rev.hpp>
#include <stan/math/fwd.hpp>
#include <stan/math/mix.hpp> // then stuff from mix/ must come next
#include <stan/math.hpp>                         // finally pull in everything from rev/ & prim/
#include <Rcpp.h>
#include <RcppEigen.h>

// [[Rcpp::plugins(cpp14)]]

// [[Rcpp::export]]
auto Hf(Eigen::VectorXd f, Eigen::MatrixXd K) { // Hessian by AD using Stan
  double fx;
  Eigen::VectorXd grad_fx;
  Eigen::MatrixXd H;
  stan::math::hessian([&K](auto f) {
    stan::math::fvar<stan::math::var> aux = stan::math::quad_form(K,f);
    return aux; },f, fx, grad_fx, H);
  return H;
}

Thanks Martin!
error-example.txt (1.2 MB)