Location parameter is nan, but must be finite!

Modeling
1

Hi all, I am new to stan while I got the error message like this when initializing.
Chain 3: Rejecting initial value:
Chain 3: Error evaluating the log probability at the initial value.
Chain 3: Exception: normal_lpdf: Location parameter is nan, but must be finite! (in ‘model339674e1a53_f0754ff79b50e64febd5efec7c18a911’ at line 87)

Any suggestion would be great appreciated!

New_Stan.stan (2.1 KB) Rstan_NewStan.R (1.8 KB)

line 87 is :
y_ijk[(n-1)*5+k] ~ lognormal(mu[(n-1)*5+k], sigma);

my Stan file is:
//

// This Stan program defines a simple model, with a
// vector of values ‘y’ modeled as normally distributed
// with mean ‘mu’ and standard deviation ‘sigma’.
//
// Learn more about model development with Stan at:
//
// http://mc-stan.org/users/interfaces/rstan.html
// https://github.com/stan-dev/rstan/wiki/RStan-Getting-Started
//

// The input data.
data {
// Define variables in data
// Number of level-1 observations (an integer)
int<lower=1> N; //number of units
// Number of level-2 clusters
int<lower=1> K; //number of conditions
// RH from the data
real RH[NK];
// Temp from the data
int<lower=1> T[NK];
// Continuous outcome
vector[NK] y_ijk;//number of the outputs
// Continuous predictor
real t_ijk[NK];
}

//Prepocessing of the data
transformed data{
real x1[NK];
real x2[NK];
x1[NK] = log(RH[NK]);
x2[NK] = 11605/(T[NK]+273.15);
}

// The parameters accepted by the model.
parameters {
// Define parameters to estimate
// Random effect
matrix[N,2] mat;
// Level-1
real<lower=1> sigma;
// Hyperparameters
vector[2] mu_mat;
real<lower=100> A[NK];
real<lower=1> B[NK];
real<lower=1> delta_H[N*K];
corr_matrix[2] sigma_mat;
corr_matrix[2] sigma_for_prior;
}

//Parameters processing before the postier is computed
transformed parameters{
// Random effect
real beta0[N];
real gamma[N];
row_vector[NK] beta1;
row_vector[NK] mu;
for(n in 1:N){
beta0[n] = mat[n,1];
gamma[n] = mat[n,2];
}
// Population slope
for(k in 1:K){
for(n in 1:N){
beta1[(n-1)5+k] = exp(log(A[nk]) + B[n*k]*x1[(n-1)5+k] + delta_H[nk]*x2[(n-1)*5+k]);
mu[(n-1)*5+k] = beta0[n] + beta1[(n-1)*5+k] * (t_ijk[(n-1)*5+k]^gamma[n]);
}
}
}

model {
sigma ~ gamma(1e-3,1e-3);
mu_mat[1] ~ normal(1e-3,1e-3);
mu_mat[2] ~ normal(1e-3,1e-3);
A ~ normal(0,1e-3);
B ~ normal(0, 1e-3);
delta_H ~ normal(0,1e-3);
mat[2] ~ multi_normal(mu_mat,sigma_mat);
sigma_mat ~ wishart(3,sigma_for_prior);
sigma_for_prior ~ lkj_corr(1);
for (k in 1:K){
for (n in 1:N){
y_ijk[(n-1)*5+k] ~ lognormal(mu[(n-1)*5+k], sigma);
}
}

}

in the R:
#load libraries
library(ggplot2)
library(StanHeaders)
library(rstan)
library(RColorBrewer)
#where the STAN model is saved
setwd("~/GitHub/Optical_media/")

Data generation

set.seed(3117)
#set up the model data
N <- 90
K <- 5

t_ijk <- ISO_data$Hours # Predictors

t_ijk <- rep(seq(100,2500,500),90)
y_ijk<- ISO_data$Deg. data
T <- ISO_data$Temp.
RH <- ISO_data$HR
logRH <- log(RH,base=exp(1))
#set the transformed data x1 and x2
x1_stan <- (logRH-log(40))/(log(85)-log(40))
head(x1_stan)
x2_stan <- (11605/(T+273.15)-11605/(15+273.15))/(11605/(85+273.15)-11605/(15+273.15))
head(x2_stan)

Set up the initial value of parameters

mu <- rnorm(450,1,1)
sigma <- rnorm(450,7,2)
mat <- matrix(runif(450,1,3),nrow=N, ncol=2)
mu_mat <- matrix(runif(450,1,3),nrow=N,ncol=2)
sigma_mat <- matrix(runif(450,0,2),nrow=2, ncol=2)
A <- rnorm(90,100,1)
B <- rnorm(90,10,3)
delta_H <- rnorm(90,10,1.5)
sigma_for_prior <- matrix(c(1,1,1,1),2,2)
for(n in 1:N){
for(k in 1:K){
beta1[(n-1)*5+k] <- exp(log(A[n]) + B[n]*x1_stan[(n-1)*5+k] + delta_H[n]*x2_stan[(n-1)*5+k])
}
}

beta0 <- mat[,1]
gamma <- mat[,2]

Set up the initial value of outcome

for(n in 1:N){
for(k in 1:K){
mu[(n-1)*5+k] <- beta0[n] + beta1[(n-1)*5+k] * (t_ijk[(n-1)*5+k]^gamma[n])
y_ijk[(n-1)*5+k] <- mu[(n-1)*5+k] + sigma[(n-1)*5+k]
}
}

Set model data

stan_data <- list(N=N,
K=K,
t_ijk=t_ijk,
y_ijk=y_ijk,
RH=RH,
T=T,
x1=x1_stan,
x2=x2_stan)

Load Stan file

fileName <- “New_Stan.stan”
stan_code <- readChar(fileName,file.info(fileName)$size)

cat(stan_code)

Run Stan

runStan <- stan(model_code=stan_code,data=stan_data,
chains = 3, iter = 3000, warmup = 500, thin = 10, init_r = .1)

print(runStan, pars=c(“mat”,“A”,“B”,“delta_H”,“sigma_mat”,“sigma”))

2

Hey @Amanda welcome to the forum! I think the problem is either with the mu or mu\_mat parameter. I’d say more likely mu without running the code. When the sampler starts up, unless you tell it otherwise, it’ll sample starting values for all parameters you haven’t given it a prior for to uniform(-2,2). The trouble is that for x<0, log(x)=NaN, which is the error you’re getting.

So to fix it, you could either give those parameters priors (try mu first) or change the range of the initial values (init= parameter). I’d suggest the former.

3

Actually, priors don’t do anything for initialization; they can only guide the sampler after it starts moving.
That init_r = .1 means all (unconstrained) parameters start are drawn from uniform(-0.1,0.1).
While mu indeed seems to be the problem I don’t see where the NaN values could come from.
beta0 and gamma should be initialized in the small range near zero and t_ijk is in the 100-2500 range. None of those seem extreme enough to cause numerical errors. beta1 is more complicated but just eyeballing the formula I don’t think it’s too bad.

There is something else that seems wrong to me here. In the R code A, B, and delta_H have size 90 but in the Stan code their size is N*K, not N.

And minor weirdness: the exception says normal_lpdf but line 87 has lognormal. Also that’s line 88 in the attached file.

1 Like
4

Hi @emiruz Thank you for the suggestion! I re-write the code immediately after seeing your reply. I tried give a prior to mu, but it still stuck there. Then I set all initial values as positive, still can’t see where the NaN values come from.
I was trying to use is_nan() function to test mu, but it didn’t run. besides, I don’t know how to print mu in RStan.

5

Hi Niko,
Thank you for the contributing. sorry I did post wrong error message…it’s at line 88, lognormal(mu,sigma).
I also changed the length of A, B and delta_H in R script, which should be 450, you’re right, thank you.
Here my model is logy_ijk = D_ijk + error_ijk.
where D_ijk = beta0_i + beta1_j * (t_ijk^gamma_i) + error_ijk
therefore, I wrote y_ijk ~ lognormal( mu, sigma)
where mu =beta0_i + beta1_j * (t_ijk^gamma_i)
Here beta0_i and gamma_i are random effect. beta0_i = mu0 + error0_i, and gamma_i = gamma0 + error1_i.
So I let mat = [beta_0, gamma] ~ multinormal( mu_mat, sigma_mat)

I uploaded the model as below, any thoughts about the mu? Thank you very much!

30%20PM
Screen Shot 2019-10-17 at 4.36.30 PM.png926×1396 162 KB

6

The indexing A[n*k] is also strange. Shouldn’t it be A[(n-1)*5+k], like beta1?

Does it fail when init_r = 0?

7

Yes, it didn’t run even with init_r =0 and I replaced all N*K to be (N-1)*5+K.
I’m stucked.

8

Print statements should work in RStan. Add something like print("mu: ", mu) and print("beta1: ", beta1) and maybe the NaN will turn up somewhere in the output.

4 Likes
9

Thank you Niko, my code runs after using the print statement to debug. the error was in x1 and x2, I forgot to add them into for loops, which cause the nan values.

2 Likes