Initialization Error for Censored Model

Modeling
1

Hi,

I’m trying to fit a model with censored data. The example on in the manual works great for that case, but I’d like to make it bit more complicated and eventually do something like a regression where the outcome is censored. Therefore, I tried to add a transformed parameter block in order to simulate preparing a set of parameters for the model block, however I keep running into an initialization error…

Compiling /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.pyx because it changed.
[1/1] Cythonizing /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.pyx
building 'stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421' extension
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T
creating /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2
gcc -Wno-unused-result -Wsign-compare -Wunreachable-code -DNDEBUG -g -fwrapv -O3 -Wall -Wstrict-prototypes -I/Users/thauck/miniconda3/envs/p3/include -arch x86_64 -I/Users/thauck/miniconda3/envs/p3/include -arch x86_64 -DBOOST_RESULT_OF_USE_TR1 -DBOOST_NO_DECLTYPE -DBOOST_DISABLE_ASSERTS -I/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2 -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/stan/src -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/stan/lib/stan_math_2.15.0 -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/stan/lib/stan_math_2.15.0/lib/eigen_3.2.9 -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/stan/lib/stan_math_2.15.0/lib/boost_1.62.0 -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/stan/lib/stan_math_2.15.0/lib/cvodes_2.9.0/include -I/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/numpy/core/include -I/Users/thauck/miniconda3/envs/p3/include/python3.6m -c /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.cpp -o /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.o -O2 -ftemplate-depth-256 -Wno-unused-function -Wno-uninitialized
g++ -bundle -undefined dynamic_lookup -L/Users/thauck/miniconda3/envs/p3/lib -L/Users/thauck/miniconda3/envs/p3/lib -arch x86_64 /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.o -L/Users/thauck/miniconda3/envs/p3/lib -o /var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.cpython-36m-darwin.so
---------------------------------------------------------------------------
RemoteTraceback                           Traceback (most recent call last)
RemoteTraceback: 
"""
Traceback (most recent call last):
  File "/Users/thauck/miniconda3/envs/p3/lib/python3.6/multiprocessing/pool.py", line 119, in worker
    result = (True, func(*args, **kwds))
  File "/Users/thauck/miniconda3/envs/p3/lib/python3.6/multiprocessing/pool.py", line 44, in mapstar
    return list(map(*args))
  File "stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.pyx", line 368, in stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421._call_sampler_star (/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.cpp:8321)
  File "stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.pyx", line 401, in stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421._call_sampler (/var/folders/xj/9f5d9t5d23lc008bxndjw05r0000gn/T/tmpwtnm23l2/stanfit4anon_model_4fd699559b463f37644fa7a9130b4121_6572719757130037421.cpp:8973)
RuntimeError: Initialization failed.
"""

The above exception was the direct cause of the following exception:

RuntimeError                              Traceback (most recent call last)
<ipython-input-32-00cd3dec3631> in <module>()
     50 }
     51 
---> 52 fit = pystan.stan(model_code=code, data=data, verbose=True)

/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/api.py in stan(file, model_name, model_code, fit, data, pars, chains, iter, warmup, thin, init, seed, algorithm, control, sample_file, diagnostic_file, verbose, boost_lib, eigen_lib, n_jobs, **kwargs)
    393                      sample_file=sample_file, diagnostic_file=diagnostic_file,
    394                      verbose=verbose, algorithm=algorithm, control=control,
--> 395                      n_jobs=n_jobs, **kwargs)
    396     return fit

/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/model.py in sampling(self, data, pars, chains, iter, warmup, thin, seed, init, sample_file, diagnostic_file, verbose, algorithm, control, n_jobs, **kwargs)
    719         call_sampler_args = izip(itertools.repeat(data), args_list, itertools.repeat(pars))
    720         call_sampler_star = self.module._call_sampler_star
--> 721         ret_and_samples = _map_parallel(call_sampler_star, call_sampler_args, n_jobs)
    722         samples = [smpl for _, smpl in ret_and_samples]
    723 

/Users/thauck/miniconda3/envs/p3/lib/python3.6/site-packages/pystan/model.py in _map_parallel(function, args, n_jobs)
     78             n_jobs = None
     79         with multiprocessing.Pool(processes=n_jobs) as pool:
---> 80             map_result = pool.map(function, args)
     81     else:
     82         map_result = list(map(function, args))

/Users/thauck/miniconda3/envs/p3/lib/python3.6/multiprocessing/pool.py in map(self, func, iterable, chunksize)
    258         in a list that is returned.
    259         '''
--> 260         return self._map_async(func, iterable, mapstar, chunksize).get()
    261 
    262     def starmap(self, func, iterable, chunksize=None):

/Users/thauck/miniconda3/envs/p3/lib/python3.6/multiprocessing/pool.py in get(self, timeout)
    606             return self._value
    607         else:
--> 608             raise self._value
    609 
    610     def _set(self, i, obj):

RuntimeError: Initialization failed.

The python code (w/ stan model as a string) is…

weights = pd.Series(np.random.normal(99, 1, 100))
print(weights.mean())

weights[weights > 100] = 100

code = """
  data {
    int<lower=0> N_obs;
    int<lower=0> N_cens;
    real U;
    real<upper=U> y[N_obs];
  }

  parameters {
    real<lower=U> y_cens[N_cens];
    real mu;
    real<lower=0> sigma;
  }
  
  transformed parameters {
      real<lower=U> y_hat_cens[N_obs];
      real<upper=U> y_hat_obs[N_cens];
      
      for (i in 1:N_obs) {
          y_hat_obs[i] = mu;
      };
      
      for (i in 1:N_cens) {
          y_hat_cens[i] = mu;
      };
  }

  model {
 
      
      for (i in 1:N_obs) {
          y[i] ~ normal(y_hat_obs[i], sigma);
      };
      
      for (i in 1:N_cens) {
          y_cens[i] ~ normal(y_hat_cens[i], sigma);
      };
  }
"""

data = {
    "N_obs": (weights != 100).sum(),
    "y": weights[~(weights == 100)].tolist(),
    "U": 100,
    "N_cens": (weights == 100).sum()
}

fit = pystan.stan(model_code=code, data=data, verbose=True)

I suspect it’s due to the conflicting upper/lower bounds in the for the parameters, and the transformed parameters block setting those values to mu.

Thanks very much, happy to provide more info if I can.

2

Yeah I think the error makes sense.

There’s an upper bound, U, on y_hat_obs and a lower bound U on y_hat_cens. It looks like they’re both being set equal to mu, which means only mu == U is ever going to be valid (so the sampler is not gonna be able to initialize itself).

I think you’re looking at section 11.3 in the manual? What direction are you trying to take this model? I’m not super familiar with this stuff but maybe someone else could give you a pointer on where to head.

Ben

3

Hi Ben,

Thanks for your response.

I am looking at section 11.3. Following the scale example from the manual (a scale that only weights to 300), my goal is have a model where weight ~ gender + height, so I started playing with the transformed parameters in order to do the initial transformation of the parameters before passing them to the model block.

Trent

4

Just to follow up on this, it was a simple as removing the constraints on the transformed parameters.

5

Yeah, if you’re talking about doing a regression with gender + height as predictors and you know the scales don’t go over 300, I think it makes sense to just throw the regression in the 11.3 model.

Hope it worked out!

1 Like
6

The efficient way to program your model would drop the transformed model block and replace the body of the model block with:

y ~ normal(mu, sigma);
y_cens ~ normal(mu, sigma);

Maybe you were getting it ready to add some predictors or something? If so, you often run into problems with bounds like these.

Also, if the censoring is this simple, you can often just use the CCDF, removing the y_cens parameters and replacing their sampling statement in the model body with

target += N_cens * normal_lccdf(U | mu, sigma);
1 Like