JohannesBuchner
commented
10 years ago The problem is that the output column values are actually written last index fastest, while the header claims first index fastest indexing!
Closed JohannesBuchner closed 10 years ago
JohannesBuchner
commented
10 years ago The problem is that the output column values are actually written last index fastest, while the header claims first index fastest indexing!
bob-carpenter
commented
10 years ago I'm not quite following the problem here. I tried to create a slightly simpler example.
parameters {
real y[2,2];
}
model {
print("y=",y);
y[1,1] ~ normal(1,1);
y[2,1] ~ normal(100,1);
y[1,2] ~ normal(10000,1);
y[2,2] ~ normal(1000000,1);
}I then create the init file in R:
> y <- matrix(NA,2,2)
> y[1,1] <- 1
> y[2,1] <- 100
> y[1,2] <- 10000
> y[2,2] <- 1000000
> y
[,1] [,2]
[1,] 1 1e+04
[2,] 100 1e+06
> stan_rdump(c('y'),"init.init.R")which produces a file init.init.R with content
y <-
structure(c(1, 100, 10000, 1e+06),
.Dim = c(2, 2))And then I compile and run from CmdStan, turning warmup off:
~/temp2/init$ ./init sample num_warmup=0 num_samples=5 init=init.init.Rand what I get is what I'd expect with no warmup:
Iteration: 1 / 10 [ 10%] (Sampling)
y=[[1,10000],[100,1e+06]]
y=[[1.93513,10001.1],[100.572,1e+06]]
y=[[0.0648746,9998.93],[99.4284,1e+06]]
y=[[0.0648746,9998.93],[99.4284,1e+06]]
y=[[1.93513,10001.1],[100.572,1e+06]]
y=[[2.09643,9999.64],[99.1611,1e+06]]and so is the output
~/temp2/init$ cat output.csv
...
lp__,accept_stat__,stepsize__,treedepth__,n_leapfrog__,n_divergent__,y.1.1,y.2.1,y.1.2,y.2.2
...
0,0.725398,1,2,3,0,1,100,10000,1e+06
-0.669897,0.845798,1,2,3,0,0.16092,99.3665,10000.5,1e+06
-1.52923,0.807111,1,2,3,0,0.495168,101.343,9999.68,999999
-2.36202,0.75204,1,2,3,0,0.724809,101.193,9998.79,1e+06
-2.36202,0.458734,1,2,3,0,1.27519,98.8069,10001.2,999999That is, I get y.1.1 equal to 1, y.2.1 equal to 100, y.1.2 equal to 10000, and y.2.2. equal to 1000000.
mbrubake
commented
10 years ago I can confirm this with CmdStan on current develop. It also is only in issue with arrays. If mymean is a matrix, the ordering is correct and consistent. The bug is in how the write_array/write_csv/etc orders its loops for arrays which is inconsistent with write_csv_header.
@bob-carpenter @syclik I think we should just to everything column major since that's also how things are read in. Thoughts?
bob-carpenter
commented
10 years ago I guess I'm still being too dense to see what the problem is. @Marcus --- can you create a simpler example for me to follow?
Or is it somehow just a problem with optimization?
mbrubake
commented
10 years ago @bob-carpenter Look in the generated cpp file for the original model. You'll find this:
void write_csv_header(std::ostream& o__) const {
stan::io::csv_writer writer__(o__);
for (int k_1__ = 1; k_1__ <= 3; ++k_1__) {
for (int k_0__ = 1; k_0__ <= 2; ++k_0__) {
writer__.comma();
o__ << "mymean" << '.' << k_0__ << '.' << k_1__;
}
}
writer__.newline();
}
template <typename RNG>
void write_csv(RNG& base_rng__,
std::vector<double>& params_r__,
std::vector<int>& params_i__,
std::ostream& o__,
std::ostream* pstream__ = 0) const {
...
size_t dim_mymean_0__ = 2;
mymean.resize(dim_mymean_0__);
for (size_t k_0__ = 0; k_0__ < dim_mymean_0__; ++k_0__) {
size_t dim_mymean_1__ = 3;
for (size_t k_1__ = 0; k_1__ < dim_mymean_1__; ++k_1__) {
mymean[k_0__].push_back(in__.scalar_constrain());
writer__.write(mymean[k_0__][k_1__]);
}
}
...Note that the loop ordering is different between write_csv_header() and write_csv().
JohannesBuchner
commented
10 years ago The problem seems to only appear when running "optimize", not "sample". I think.
bob-carpenter
commented
10 years ago Thanks for the clarification --- I see the problem. I'd like to get rid of the dependencies on write_csv and write_csv_header so that the generated model code doesn't need to do anything with formatting.
@mbrubake Would you mind changing optimization so that it uses constrained_param_names() and write_array(), both of which are last-index fastest?
@maverickg Does RStan use write_csv or write_csv_header?
@ariddell Does PyStan use write_csv or write_csv_header?
mbrubake
commented
10 years ago Ok, changing now. When did those functions become deprecated? It's also somewhat of a shame as we now have duplicated code between command.hpp and in mcmc_writer.hpp because we're not using functions like write_csv or write_csv_header. Perhaps we should simply rewrite the write_csv* functions in terms of write_array?
JohannesBuchner
commented
10 years ago Could you please help me: How can I print out the current model parameters? I would like to insert some debugging code into command.hpp after the initialisation is loaded. In another test code, I get -inf values in cont_params after loading. I don't know the actual model parameter values. Can I dump them somehow?
betanalpha
commented
10 years ago write_csv is only used in the optimization — it’s in the process of being deprecated, a process which was put on hold due to the upcoming refactor of command.
On Feb 4, 2014, at 5:52 PM, Marcus Brubaker notifications@github.com wrote:
Ok, changing now. When did those functions become deprecated? It's also somewhat of a shame as we now have duplicated code between command.hpp and in mcmc_writer.hpp because we're not using functions like write_csv or write_csv_header. Perhaps we should simply rewrite the write_csv* functions in terms of write_array?
— Reply to this email directly or view it on GitHub.
bob-carpenter
commented
10 years ago @JohannesBuchner Current in what sense? And on the constrained or unconstrained scale? You can print the constrained versions using a print statement in the model. Those will print during each leapfrog step during sampling. I'm not sure what happens during optimization.
We're not going to add general printing to the Stan C++ API.
bob-carpenter
commented
10 years ago I'd rather not have the compiled model have anything built in to do printing.
I'd be OK with a write_csv-like utility function on the outside to use in CmdStan. PyStan and RStan both need to get the values as C++ data structures.
On 2/4/14, 6:52 PM, Marcus Brubaker wrote:
Ok, changing now. When did those functions become deprecated? It's also somewhat of a shame as we now have duplicated code between command.hpp and in mcmc_writer.hpp because we're not using functions like write_csv or write_csv_header. Perhaps we should simply rewrite the write_csv* functions in terms of write_array?
— Reply to this email directly or view it on GitHub https://github.com/stan-dev/stan/issues/532#issuecomment-34087102.
JohannesBuchner
commented
10 years ago Oh this worked:
std::cout << "initial log-probability: " << init_log_prob << std::endl;
std::cout << "values: " << std::endl;
for (unsigned int i = 0; i < model.num_params_r(); i++) {
std::cout << " " << cont_params[i] << std::endl;
}
model.write_csv(base_rng, cont_params, disc_params, std::cout, &std::cout);You're right, using print in the model is much simpler.
There is a problem with the way Stan reads the initialisation values for array parameters.
I prepared a minimal example code:
I generated test data from this model:
I optimize and the file written correctly contains:
I create the following initialisation file:
If I use this in R, it gives:
This is exactly what the output file describes.
Now I use this as a initialisation file for sampling.
The problem is that init_log_prob is then extremely low (-1.64101e+17). (For debugging, I inserted a line into command.hpp to print out init_log_prob). The correct log_probability (-303.011) is only given if I re-order the array in the initialisation:
Such that
Either the order the array is read or written incorrectly. I use Stan commit 65ceb15f5f1b.