Hello,
I need to simulate data from additive Bayesian networks (ABN) with parameters that will be given to Stan with an adjacency matrix. All variables in ABNs are modeled with an appropriate regression model.
Considering the network A → B → C where A and C are normally distributed variables and B is binary, the model would look as follows:
A ~ Normal(μ_A, σ^2_A), with μ_A = β_{A,0}
B ~ Bernoulli(π), with logit(π) = β_{B,0} + β_{B,1} * A
C ~ Normal(μ_C, σ^2_C), with μ_C = β_{C,0} + β_{C,1} * B
I got this to work with the following Stan code
data{
int<lower = 1> N;
matrix[N, N] X;
vector[N] sigma;
}
model{
}
generated quantities {
real A;
int B;
real C;
real piB;
A = normal_rng(X[1, 1], sigma[1]);
piB = 1/(1+exp(-(X[2,2]+X[1,2]*A)));
B = bernoulli_rng(piB);
C = normal_rng(X[3,3] + X[2,3]*B, sigma[3]);
}
where N = 3, sigma=(1, 0, 1.5) and X=
(0, 0.8, 0
0, 2, 1.5
0, 0, 0)
are given as data to Stan.
However in the general case, this doesn’t seem to work where we need some form of matrix-vector multiplication. In a purely normal case, it would be something along the lines of this:
functions {
vector delete_element(vector x, int m) {
int N;
vector[num_elements(x)-1] y;
N=num_elements(x);
for (i in 1:(m-1)) {
y[i]=x[i];
}
for (i in (m+1):N) {
y[i-1]=x[i];
}
return y;
}
}
data{
int<lower = 1> N;
matrix[N, N] X;
vector[N] sigma;
}
model {}
generated quantities {
vector[N] nodes;
for (i in 1:N) {
nodes[i] = normal_rng(X[i, i] + dot_product(delete_element(to_vector(X[i, ]), i),
delete_element(nodes, i)), sigma[i]);
}
This doesn’t work, since it includes undefined variables. Even when they are multiplied by zero, Stan gives:
normal_rng: Location parameter is nan, but must be finite!
I have also tried things with specified parameter and model block but could not find a solution so far. Any help would be much appreciated.