soil_incubation.stan

[will246]

Hello there,

When I try to run soil_incubation.stan with the example data from the SoilR package (https://mc-stan.org/users/documentation/case-studies/soil-knit.html), I get the following error:

`Error in stanc(file = file, model_code = model_code, model_name = model_name, :
0

Semantic error in 'string', line 93, column 13 to line 94, column 57:

Ill-typed arguments supplied to function 'integrate_ode'. Available signatures:
((real, real[], real[], data real[], data int[]) => real[], real[], real, real[], real[], data real[], data int[]) => real[,]
Instead supplied arguments of incompatible type: (real, real[], real[], real[], int[]) => real[], real[], real, real[], real[], real[], int[].`

Do you know of any possible reasons for this? Thank you.

[dlebauer]

Hi @will246, thanks for reporting this! I was able to generate the same error, but unfortunately don't know the solution. I do suspect that it has to do with the fact that a) the code is almost a decade old, and the underlying software has likely changed! e.g. the eCO2 dataset from SoilR is different and you can load the original like this:

download.file('https://github.com/cran/SoilR/raw/1.1-23/data/eCO2.rda', 
             'eCO2.rda')
load('eCO2.rda')

In addition, stan has changed and the function integrate_ode has been deprecated as described in https://mc-stan.org/docs/functions-reference/functions-old-ode-solver.html.

Hopefully the author and Stan developer @bob-carpenter has some insight.

[nsiccha]

@dlebauer is likely correct. I had also recently had a look at the case study, and also found some conceptual errors. I'd recommend to rework it a bit more thoroughly and then re-upload it.

[bob-carpenter]

Hi, @dlebauer. There are a bunch of things that have changed in our interfaces, so probably best to rewrite. Hi, @nsiccha --- what did you want to fix conceptually? I'm happy to review.

[ktoddbrown]

Hi folks. I have some students who might be interested in revising this if you don't mind a longer time horizon.

@nsiccha I'm curious what you meant by conceptual errors here.

[nsiccha]

I hope everyone had nice holidays.

@bob-carpenter and @ktoddbrown: The following are things I would change.

• I believe the transfer rates alpha_{...} have to be between 0 and 1. Mass conservation implies for values of alpha above 1 that you are sequestering carbon from the enviroment, but in a way that I do not believe you want to allow.
• I would discuss the need for/utility of/alternatives to the Beta(10,1) prior for gamma. AFAICT, it's there to break the symmetry intrinsic to the model, but as the the two compartments are "made up" anyways, maybe some discussion would be helpful? The model seems to be very flexible and able to reproduce the observations just as well with a flat prior on gamma, it's just that you "lose" identifiability with the flat prior.
• A bunch of output in the rendered case study is weird:
  • The trace plots and the pair plots look like the chains are not mixing well at all. I didn't follow up for my quick reimplementation, but I did observe as well that it was surprisingly difficult to sample from similar posteriors in that chains sometimes got stuck or just did not mix well.
  • For some spooky reason the fit_me objects has both 100 or 1000 draws per chain and uses both 2 or 4 chains? It gets initialized via fit_me <- stan("soil_incubation_measurement_err.stan", [...], chains=2, iter=100, seed=1234);, the output after that block twice shows Iteration: 100 / 100 [100%] (Sampling), but the summary after that prints ## 4 chains, each with iter=1000; warmup=500; thin=1; ? The trace plot after the summary than again shows 100 iterations per chain

I'd be happy to help with rewriting, improving and/or extending the case study, though I don't think I'll have time before February. What kind of time frame did you have in mind, @ktoddbrown?

[will246]

Thank you for your input everyone @dlebauer @ktoddbrown @nsiccha. I was eventually able to get it to run. I used soil_incubation.stan from https://github.com/soil-metamodel/stan/tree/master/soil-incubation. It seems that the syntax in that file is different to that shown on the example page (https://mc-stan.org/users/documentation/case-studies/soil-knit.html). @nsiccha, I noticed the issues with mixing as well. I can look into this further too.

[bob-carpenter]

• (0, 1) is a reasonable constraint on transfer rate and I'd be in favor of enforcing it. I have this discussion with @andrewgelman from time to time, who likes to leave things unconstrained to see what the data want to do with the fit. The same issue comes up in something like an AR(N) time series model where you might want to constrain to a stationary distribution.
• A better way to identify otherwise unidentified mixture models is to impose an ordering on the mixture components' parameters. @betanalpha has a nice case study: https://betanalpha.github.io/assets/case_studies/identifying_mixture_models.html
• The model should mix well enough. One problem we have with these models is stiffness under initialization. In those cases, using the stiff ODE solver is slower, but more robust.
• I'd just rewrite everything in CmdStanPy if it were me. Or CmdStanR if we want to keep an R version.

[will246]

Thank you for your input everyone @dlebauer @ktoddbrown @nsiccha. I was eventually able to get it to run. I used soil_incubation.stan from https://github.com/soil-metamodel/stan/tree/master/soil-incubation. It seems that the syntax in that file is different to that shown on the example page (https://mc-stan.org/users/documentation/case-studies/soil-knit.html). @nsiccha, I noticed the issues with mixing as well. I can look into this further too.

Actually, I used soil_incubation.stan from here: https://github.com/stan-dev/example-models/blob/master/knitr/soil-carbon/soil_incubation.stan.