Local individual and population parameters

Modeling
, ,
1

Hi All,
I am having issues with local variables implemented within a model.
Some parameters are for individual specimens, some are for the entire population.

Problem Outline #modeling

  • I want to implement individual local parameters which are not reported. These include noise1, Epsilon.
  • I have the model running where noise1 and Epsilon are declared as parameters, but I don’t want this as it slows down the process with my non-test data (and I have others to ultimately include)
  • Error obtained: Chain 1: Initialization between (-2, 2) failed after 100 attempts.
  • I understand I can not put bounds on local variables, I have queried moving things to Transformed parameters block with no luck (and am not even sure if that is the way to go). I have also tried pulling in initial conditions to push past initialization…

I am grateful for any pearls of wisdom.
Kind regards,
Steve

STAN CODE

data {
  int<lower=1> N;
  int<lower=50> s2[N];
  int<lower=50> s1[N];
  int<lower=1> time[N];
}

parameters {
  real<lower = 0.05, upper=0.4> GrowthConstant;
  real<lower = 160, upper = 250> SizeMax;
  real<lower=0.01, upper=0.3> Epsilonsd;
}

model {
  //Declare local variables
  real Epsilon[N];
  real noise1[N];
  
  //Priors
  GrowthConstant ~ normal(0.2, 0.05);
  SizeMax ~ normal(190, 10);
  Epsilonsd ~ gamma(2, 10);
  
  for (n in 1:N){
    noise1[n] ~ normal(0, 0.98);
    Epsilon[n] ~ normal(1, Epsilonsd);
  }

  //Likelihood
  for (n in 1:N) {
    s2[n] ~ normal(s1[n] - noise1[n] + Epsilon[n]*(SizeMax - s1[n] + noise1[n])*(1-exp(-GrowthConstant * time[n])), 0.98);
  }

}

SAMPLE DATA CREATION

n <- 2000
df <- data.frame(time = numeric(n), size1= numeric(n), size2= numeric(n))

df$size1 <- rnorm(2000, 140, 20)
df$size1[df$size1>190] <- 190

df$time <- rgamma(2000, 1.5, 1) %>% ceiling()
maxSize <- 185
GrowthConstant <- 0.22

df$size2 <- df$size1 + rnorm(n=dim(df)[1], mean=1, sd=0.05)(maxSize - df$size1)(1-exp(-GrowthConstant*(df$time)))

df$s1 <- (df$size1 + rnorm(dim(df)[1],0, 0.98)) %>% round(0)
df$s2 <- (df$size2 + rnorm(dim(df)[1],0, 0.98)) %>% round(0)

df <- df[, c(“time”, “s1”, “s2”)]

data <- list(N=dim(df)[1]
, s2=as.numeric(df$s2)
, s1=as.numeric(df$s1)
, time=as.numeric(df$time)
)

2

You have to specify parameters in the parameters block.
You can declare them as vectors in parameters block:

vector[N] Epsilon;
vector[N] noise1;

and omit the loops

noise1 ~ normal(0, 0.98);
Epsilon ~ normal(1, Epsilonsd);

Further you can exclude them in the reporting:

sfit <- sampling(smod
        , data = sdata
        , pars = c("noise1", "Epsilon")
        , include = FALSE
        , iter = 1000
        , chains = 4
)

or exclusively include by include = TRUE.

3

Thank you Andre. I am very grateful!