Parameters constraint

[Marcelo Silva]

I have two vectors parameters: g[J] and s[J], where g[j]~beta(1,1) and s[j]~beta(1,1), with j=1,2,..,J.

Could someone help me put the following restriction:

0 <= g[j] + s[j] < 1, for all j=1,2,...,J

[Bob Carpenter]

Your constraint is convenient because you're really
asking for the first two components of the cumulative
distribution over the simplex.

data {
int<lower=0> J;
}
parameters {
simplex[3] theta[J];
}
transformed parameters {
vector[J] g;
vector[J] s;

g <- to_vector(theta[, 1]);
s <- g + to_vector(theta[, 2]);
}

Because beta(1, 1) is uniform, those priors on s and j won't
do anything --- Stan only needs the log density up to a constant.
So you can just drop them.

The resulting prior is uniform over the (g,s) pairs meeting
the constraint (the points in (0, 1) x (0, 1) above the diagonal).

- Bob

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

[Marcelo Silva]

Thanks for responding, Bob Carpenter. 

I apologize for particularized the parameters for a Beta distribution, but in fact I would like to have the restriction 0 <= g[j] + s[j] < 1, for all j=1,2,...,J for any parameter  a_g, b_g, a_s and b_s, where

g[j]~beta(a_g,b_g) and s[j]~beta(a_s,b_s), with j=1,2,..,J.

[Bob Carpenter]

I read that as g[j] < s[j] < 1 rather than g[j] + s[j] < 1.
You can still use a simplex, but now you don't need to add
and can just have:

> g <- to_vector(theta[, 1]);

> s <- to_vector(theta[, 2]);

You can add whatever prior you like on top of the implicit
uniform over values meeting the constraint that 0 < g[j] + s[j] < 1.

But if you try to write

> g[j]~beta(a_g,b_g) and s[j]~beta(a_s,b_s)

and mean the joint prior to be

p(g[j], s[j] | a_g, b_g, a_s, b_s)

= beta(g[j] | a_g, b_g) * beta(s[j] | a_s, b_s)

you have a problem because it doesn't integrate to 1 over the
values of g and s meeting the constraint that 0 < g + s < 1.
If a and b are constants, it'll only be off by a constant, but
if they are parameters, you have to figure out how to account
for them, like taking the second distribution over s to be
truncated. That is, take the prior to be

p(g[j], s[j] | a_g, b_g, a_s, b_s)

= beta(g[j] | a_g, b_g) * beta(s[j] | a_s, b_s) T[0, g[j]]

where the second truncation means to divide by the cdf of beta
up to g[j] with parameters a_s and b_s.

- Bob

[Marcelo Silva]

Thank you for your attention, Bob Carpenter. I write the model code I am studying in Stan code, but the code did not work. Can you help me? I am working with a model called DINA and I'm having problems to implement it in Stan. 

Below, write an abstract of this model:

A questionnaire with J questions is answered by $n$ individuals in order to classify each individual based upon mastery and nonmastery of the $K$ skills. We have 2^K possible attribute profiles. Example: For K=4: we have 2^4=16 differents profiles: alpha_1 = (0,0,0,0), alpha_2 = (1,0,0,0), alpha_3 = (0,1,0,0), ... , alpha_16 = (1,1,1,1). In case, alpha_2 is a profile that have skill 1 and don't have skill 2,3 and 4.

The skills needed to succeed on item j are stored in vector q_j=(q_j1,q_j2,...,q_jK), where q_jk is 1 if item j requires skill k and zero otherwise.

The random variable Y_ij|\alpha_i,\Omega_j ~ Bernoulli(P_ij) is one if the individual i answered correctly the item j and zero otherwise.

Omega_j = (s_j, g_j), where s_j and g_j are item parameters and P_ij = P(Y_ij=1|alpha_i,Omega_j) = ((1-s_j)^eta_ij)*(g_j^(1-eta_ij)) and 

eta_ij = (alpha_i1^q_ij) * (alpha_i2^q_ij) * (alpha_i3^q_ij) * ... * (alpha_iK^q_ij)

I have the feeling that what I did in Stan's wrong. Follows the code attached.

For large volumes of data, the runtime is very large (48 hours to run in 1000 iterations).

[Marcelo Silva]

I wish estimate g, s, alpha and pi.

[Bob Carpenter]

You're in luck, because we have a tutorial
for the DINA model (with independent attributes) on the
way from Seung Yeon Lee at Berkeley. I don't know when he's going to be willing to share it,
but I don't want to go and implement it again when
there's another one coming. I'm bcc-ing him in case
he's willing to share a preliminary draft before it
goes up on our case studies.

You have to marginalize the discrete parameters
out (those profiles you indicate). And whatever you
do, you shouldn't be converting reals to integers unless
you're building spline models.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

> <code_dina.R>

[Marcelo Silva]

Thanks for your tip, Bob.

I will talk to my doctoral supervisor before contacting Seung Yeon Lee.

Also, I will study the marginalization of discrete parameters.

Marcelo Silva