Require real return type for probability functions error

Modeling
, ,
1

I am trying to use stan with user-defined function, here is my code,

functions{
  real copulagaussian_log(real[,,] y, real[,,] Z, real[,,] X,int p,real[,] v1,real[,] v2,real[] nbeta1,real[] nbeta2,real[,] mu1,real[,] mu2,int N,int n,real ramda,real sigma, real rtilde,int q,real w1,real w2){
 real lprob;
 lprob = 0;
 for(i in 1:N){
    for(k in 1:n){
       lprob=lprob-0.5*log(1-rtilde^2)-0.5/(1-rtilde^2)*(rtilde^2*v1[i,k]*v1[i,k]+rtilde^2*v2[i,k]*v2[i,k]-2*rtilde*v1[i,k]*v2[i,k])-0.5*log(2*pi())-log(sigma)-0.5*(y[1,i,k]-mu1[k,i])*(y[1,i,k]-mu1[k,i])/sigma^2+log(ramda)-ramda*y[2,i,k];
       }
    }
    for(i in 1:p)
    {
      lprob=lprob+2*log(1/(sqrt(2*pi())*10))-0.5*(nbeta1[i]-1)^2/100-0.5*(nbeta2[i]-2)^2/100;
    }
    lprob=lprob+log(0.5)-(2*q+2*q+1)*0.5*log(w1)-0.5*w2-2*q^2*log(2)-(q^2-0.5*q)*log(pi());
    for(i in 1:(2*q))
    {
      lprob=lprob-log(tgamma(q+0.5*(1-i)));
    }
// return the log likelihood value
 return lprob;
 }
}

data{
 int N; //number of observations
 int n; //number of branches
 int p; //dim of X
 int q; //dim of Z
 matrix[(2*q),(2*q)] I;
 vector[(2*q)] bmu;
 real y[2,N,n]; //outcomes array
 real X[N,n,p]; //branch-varying variables array
 real Z[N,n,q]; //option-varying variables array
}

parameters{
 real<lower=-1, upper=1> rtilde; //scale-parameters
 vector[p] nbeta1; //option-varying parameters
 vector[p] nbeta2; //option-varying parameters
 matrix[(2*q),N] B;
 real<lower=0> ramda;
 real<lower=0> sigma;
 cov_matrix[(2*q)] sigmab;
}

transformed parameters{
  real mu1[n,N];
  real<lower=0> mu2[n,N];
  real v1[N,n];
  real v2[N,n];
  real<lower=0> w1;
  real w2;
  w1=determinant(sigmab);
  w2=trace(I*inverse(sigmab));
  for(i in 1:N)
  {
    for(k in 1:n)
    {
      mu1[k,i]=dot_product(to_row_vector(X[i,k,]),nbeta1)+dot_product(to_row_vector(Z[i,k,]),B[1:3,i]);
      mu2[k,i]=exp(dot_product(to_row_vector(X[i,k,]),nbeta2)+dot_product(to_row_vector(Z[i,k,]),B[4:6,i]));
      v1[i,k]=y[1,i,k];
      v2[i,k]=inv_Phi(1-exp(-ramda*y[2,i,k]));
    }
  }
}
model{
//specify priors
 for(i in 1:p) 
 {
   nbeta1[i]~normal(1,10);
   nbeta2[i]~normal(2,10);
 }
 ramda~exponential(0.5);
 sigma~normal(2,10);
 sigmab~inv_wishart((2*q), I);
 rtilde~uniform(-1,1);
 for(i in 1:N){
    B[,i]~multi_normal(bmu,sigmab);
 }

           y~copulagaussian(Z,X,p,v1,v2,nbeta1,nbeta2,mu1,mu2,N,n,ramda,sigma,rtilde,q,w1,w2);
  
}
2

What error do you get?

1 Like
3

In the last line, it indicates error of requiring real return type for probability functions. It has tilde under w2, so I guess w2 is the problem.

4

If you juxtapose the signature call and the available signature, you get

 real[ , , ] ~ copulagaussian(real[ , , ], real[ , , ], int, real[ , ], real[ , ], vector, vector, real[ , ], real[ , ], int, int, real, real, real, int, real, real)
 real[ , , ] ~ copulagaussian(real[ , , ], real[ , , ], int, real[ , ], real[ , ], real[ ], real[ ], real[ , ], real[ , ], int, int, real, real, real, int, real, real)

So there are a couple things Stan expects to be real arrays but are actually vectors.

Changing the function call to

real copulagaussian_log(real[,,] y, real[,,] Z, real[,,] X,int p,real[,] v1,real[,] v2,vector nbeta1, vector nbeta2,real[,] mu1,real[,] mu2,int N,int n,real ramda,real sigma, real rtilde,int q,real w1,real w2){
...
}

Seems to solve the issue.

1 Like
5

Yes, it works now. Thank you so much.