Nonlinear growth model - initialization error

Interfaces RStan
1

I am running a rather complicated nonlinear growth model (code below) - and am having trouble getting it to initialize. I provide initial values for all parameters in the parameter block, but it still gives the lovely error suggesting that it is generating inits between (-2,2):

"SAMPLING FOR MODEL ‘STAN_VB_TempDep_Weight&Length_FixedD_EstA&B’ NOW (CHAIN 1).
Chain 1: Rejecting initial value:
Chain 1: Error evaluating the log probability at the initial value.
Chain 1: Exception: normal_lpdf: Location parameter is nan, but must be finite! (in ‘model2a8c62983a3a_STAN_VB_TempDep_Weight_Length_FixedD_EstA_B’ at line 131)

Chain 1:
Chain 1: Initialization between (-2, 2) failed after 1 attempts.
Chain 1: Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
[1] “Error in sampler$call_sampler(args_list[[i]]) : Initialization failed.”
[1] “error occurred during calling the sampler; sampling not done”"

Here is the R code and an uploaded STAN file. Thanks for the help.

	
# Bundle the data together for JAGS 
standata <- list(lengths = lengths, lengthages = length.ages, lengthN = length.N, lengthyears = length.years, 
	nage = nage, nyears = nyears, temp = as.matrix(temp), monthtemp = month.temp, alpha = as.matrix(alpha), beta = as.matrix(betas),
	weights = weights.obs, weightsN = weights.N, weightsages = weights.ages, weightsyears = weights.years, 
	n = 1, Qk = 2, StressTemp = 28, d = 0.66)

inits <- function(chain_id = 1 ) { 
	list(
	H = runif(1,1.6,1.8), 
	k = runif(1,0.15,0.2), 
	g = runif(1,0.05,0.09),
	cvlength = runif(1,0.05,0.1), 
	cvweight = runif(1, 0.05, 0.1),
	Qc = runif(1, 2.1, 2.4)
	)}

nc <- 1 

init_list <- lapply(1:nc, function(id) inits(chain_id = nc))	
	
params <- c('H','k', 'Qc', 'cvlength','cvweight','g')#,
	#'lengthage102','weightage102','weightloglik','lengthloglik')

ni <- 1000
na <- 1000 
nb <- ni/2

library(rstan)

rstan_options(auto_write=T)
Sys.setenv(LOCAL_CPPFLAGS = '-march=corei7 -mtune=corei7')


stanmod <- stan_model(file = 'STAN_VB_TempDep_Weight&Length_FixedD_EstA&B.stan')

starttime <- Sys.time() 

options(mc.cores = nc) 

fit <- rstan::sampling(object = stanmod,
	data = standata, pars = params, chains = nc, iter = ni, warmup = nb, init = init_list,  cores = 1, 
	sample_file = '//STANdiags//'))

STAN_VB_TempDep_Weight&Length_FixedD_EstA&B.stan (3.4 KB)
I am very confident that this is not operating system or R version specific - as I have tried multiple. Perhaps I am missing some initial value somewhere?

2

Looks like the backticks: ``` ate the uploaded Stan file/model. Could you reupload the Stan code / edit your original post?

3

The error is misleading. It tried your provided initial values and gave up. The fact it failed after 1 attempt shows that all initial values were fully specified; if there were any room for randomness it would make 100 attempts.

2 Likes
4

Thanks. I have updated the post.

5

OK, so it is actually trying my initial values and choking, rather than attempting the (-2,2). This is very misleading if that is the case. I assumed that I might have somehow left out an initial value, as a related model gave this error whenever I did not provide inits for the hierarchical effects (since STAN requires inits at all levels and not just for the hyperpriors). A related question - is it possible to get STAN to try 100 attempts from specified inits?

6

If all inits are specified then the result is the same for each attempt; doing the calculation again would just fail again. If there are some missing inits Stan will automatically try 100 times to fill in some random values.

Anyway, line 131 is this

lengths[l] ~ normal(lengthmu[lengthages[l],lengthyears[l]], sdlength[lengthages[l], lengthyears[l]]);

so the error is saying that lengthmu is somehow NaN. I recommend you try a smaller dataset, ideally with lengthN=1 and see if there’s a single datapoint that causes the failure. You can also insert print() statements in the model to see what intermediate values are.
lengthmu is sum of dLs so some dL so most likely this expression

dL[a,y] =
 ha[monthtemp[a],y]*lengthmu[a,y]^da[monthtemp[a],y]*FcT[monthtemp[a],y]
  - ka[monthtemp[a],y]*lengthmu[a,y]^na[monthtemp[a],y]*FkT[monthtemp[a],y];

is somehow ill-defined. For example, if da and na are large this may overflow and give infinity minus infinity, which is NaN.

7

Thank you for the tips! I will try some things. These inits work in JAGS - with the same formula, so I can’t quite understand why it wouldn’t work here. I’ll take a good look at the code to make sure there are no errors, and then check with a smaller dataset. And yes, the entire nonlinear model is “poorly” defined, but there is no way to fix that without setting constraints. Or re-parameterizing, but making great changes to this published model is above my ability.

8

OK - I have now fixed the initial values to constants and I have verified in R that these initial values do not produce nonsensical predictions (for example, mean lengths of NA). However, STAN is still throwing the error. Shouldn’t the initial values at least allow the program to start if they give reasonable results in R? Or do I misunderstand something that is happening in the background here?

9

Yes, it should initialize.
Add this on line 129, or right before the loop that fails

print("na: ", na);
print("da: ", da);
print("dL: ", dL);
print("lengthmu: ", lengthmu);

and compare the printed numbers to the values calculated in R.

1 Like
10

Thanks very much for the code to help troubleshoot in STAN - I didn’t know these print commands would work. I have done this, and the outputs do not match (STAN produces NAs for dL and lengthmu as expected).
As I was beginning to suspect, there must be something wrong in my STAN script. I am a JAGS user - and the models run fine there. As I understand, STAN is similar to JAGS and the order of operations does not matter. Is that the case?

For a more specific example, does the order of these lines below matter?

for(y in 1:nyears) { 	

	for(a in 1:(nage)) { 

		dL[a,y] = ha[monthtemp[a],y]*lengthmu[a,y]^da[monthtemp[a],y]*FcT[monthtemp[a],y] - ka[monthtemp[a],y]*lengthmu[a,y]^na[monthtemp[a],y]*FkT[monthtemp[a],y];
		
		} 
	}


for(y in 1:nyears) { 
	lengthmu[1,y] = 2; 
}

for(y in 1:nyears) { 

	for(a in 2:nage) { 

		lengthmu[a,y] = lengthmu[a-1,y] + dL[a-1,y]; 
		
		}
	}
11

OK, I have answered my own previous question. STAN now runs when I changed these loops to the following, which is more similar to what R would expect. So I guess within a loop format STAN is expecting variables to be available?

for(y in 1:nyears) {

}

for(y in 1:nyears) { 	

	lengthmu[1,y] = 2; 
	weightmu[1,y] = 0.1; 

	for(a in 1:(nage-1)) { 

		dL[a,y] = ha[monthtemp[a],y]*lengthmu[a,y]^da[monthtemp[a],y]*FcT[monthtemp[a],y] - ka[monthtemp[a],y]*lengthmu[a,y]^na[monthtemp[a],y]*FkT[monthtemp[a],y];

		dW[a,y] = H*weightmu[a,y]^d*FcT[monthtemp[a],y] - k*weightmu[a,y]^n*FkT[monthtemp[a],y]; 
		
		lengthmu[a+1,y] = lengthmu[a,y] + dL[a,y]; 
		
		weightmu[a+1,y] = weightmu[a,y] + dW[a,y];
		
		sdlength[a+1,y] = cvlength*lengthmu[a,y]; 
	
		sdweight[a+1,y] = cvweight*weightmu[a,y];
		} 
	}
1 Like