censored multivariate normal

554 views
Skip to first unread message

Keith Lau

unread,
Dec 5, 2013, 7:16:07 AM12/5/13
to stan-...@googlegroups.com
Dear Stan developers,

I looked up Stan manual for censored data. However, only one example is given for univariate normal data. Is there any function for multivariate case?
For example, I have Y1, Y2, and Y3 distributed from multivariate normal distribution with a parameterized covariance matrix, given y1=1 if Y1>0, y2=0 if Y2<=0, and y3=1 if Y3>0.

Bob Carpenter

unread,
Dec 5, 2013, 2:47:00 PM12/5/13
to stan-...@googlegroups.com
I'm afraid this has come up before, and we don't
have any kind of cumulative distribution functions for
multivariate normals. It's on our to-do list, but not
prioritized for the immediate future. (We do accept
contributed code!)

Is there some way you can rotate the multinormal and
turn it into a product of univariate normals? If so,
you might be able to use a product of univariate cdfs.

Ben probably knows the answer to this one.

- Bob

Daniel Lee

unread,
Jan 4, 2014, 1:16:05 PM1/4/14
to stan-...@googlegroups.com
Hi Keith,

Sorry about not responding earlier. We don't have any truncated multivariate distribution functions as of v2.1.0, but I don't think that's what you're asking about.

If you haven't gotten an answer to this, can you reword it or clarify it?




Daniel




--
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.
For more options, visit https://groups.google.com/groups/opt_out.

Keith Lau

unread,
Jan 7, 2014, 1:42:09 AM1/7/14
to

Thanks for your positive attention to my question. I am encountering the estimation issue of IRT models similar to multivariate probit model, which can be expressed like

y1* = link(theta-mu1) + e1 for item one, where y1* is latent continuous variable, link() is a linking function (linear or nonlinear), theta is a random variable following N(0,1), mu1 is an intercept, and e1 is error term following N(0,sig). In the same manner, y2* = link(theta-mu2) + e2 is for second item, and so on.

Now, for example, the observable responses are y1 = 1, y2=0, y3=0,…, yn=1, which are determined if y1*>=0, y2*<0, y3*<0,…, yn*>0. Additionally, the model assumes the error terms are correlated with parameterized covariance matrix due to the experimental design. Therefore, I have to estimate theta, mu, and a covariance matrix of error term. However, to calculate the multivariate cumulative normal probability is a horrible drag and impractical. So, as far as I know, the best answer may regard the problem in the perspective of censoring mechanism. That is, the y* samples can be drawn from truncated multivariate normal distribution.

By the way, I have tried the JAGS to run this model by using one of its function called dinterval. So, the y* may be sampled from truncated multivariate normal dist (I’m not sure if it indeed works or not). However, I did not get satisfactory parameter recovery in a simulation study.

Hope above description is clear enough for you.


Daniel Lee於 2014年1月5日星期日UTC+8上午2時16分05秒寫道:

Daniel Lee

unread,
Jan 7, 2014, 11:21:59 AM1/7/14
to stan-...@googlegroups.com
Hi Keith,

Just based on your response, it seems like you're not specifying the model correctly.

In Stan, the sampling statements don't actually draw random variables from a distribution. You're specifying a joint model over all the parameters. If that doesn't make sense, try rereading chapter 1 of the manual.

Have you taken a look at the IRT models in the Stan repository? Maybe you can start there.




Daniel





On Tue, Jan 7, 2014 at 7:39 AM, Keith Lau <genw...@gmail.com> wrote:

Thanks for your positive attention to my question. I am encountering the estimation issue of IRT models similar to multivariate probit model, which can be expressed like

y1* = link(theta-mu1) + e1 for item one, where y1* is latent continuous variable, link() is a linking function (linear or nonlinear), theta is a random variable following N(0,1), mu1 is an intercept, and e1 is error term following N(0,sig). Like same manner, y2* = link(theta-mu2) + e2 is for second item, and so on.

Now, for example, the observable responses are y1 = 1, y2=0, y3=0,…, yn=1, which are determined if y1*>=0, y2*<0, y3*<0,…, yn*>0. Additionally, the model assumes the error terms are correlated with parameterized covariance matrix due to the experimental design. Therefore, I have to estimate theta, mu, and a covariance matrix of error term. However, to calculate the multivariate cumulative normal probability is a horrible drag and impractical. So, as far as I know, the best answer may to regard the problem in the perspective of censoring mechanism. That is, the y* samples can be drawn from truncated multivariate normal distribution.

By the way, I have tried the JAGS to run this model by using one of its function called dinterval. So, the y* may be sampled from truncated multivariate normal dist (I’m not sure if it indeed works or not). However, I did not get satisfactory parameter recovery in a simulation study.

Hope above description is clear enough for you.


Daniel Lee於 2014年1月5日星期日UTC+8上午2時16分05秒寫道:
Hi Keith,

Keith Lau

unread,
Jan 8, 2014, 4:05:12 AM1/8/14
to stan-...@googlegroups.com
Sorry for my poor mcmc language. I  mean that the the location of y* is only known in a restricted region; that is, observable data y = 1 if y* >=0, else y=0 if y* <0.

I have read the IRT section, but traditional IRT model like Rasch assumes uncorrelated error between items. That is not what I am doing with.
I note that my model is something similar to multivariate probit model except that I included the person's ability parameter within model.
I also googled this forum about multivariate probit model and found this https://groups.google.com/forum/#!topic/stan-users/1QzlcGem-2k
It seems there was still no clear solution to the model at that time. Is Stan 2.1.0 capable of it?


Daniel Lee於 2014年1月8日星期三UTC+8上午12時21分59秒寫道:

Bob Carpenter

unread,
Jan 9, 2014, 6:06:17 AM1/9/14
to stan-...@googlegroups.com
You can define constraints on your parameters,
but I don't understand what you're trying to say about
y and y*. Is y* a parameter?

We have truncation in the language in the same way
as BUGS/JAGS, but you need to make sure the constraint
on your variables match the truncation in order for
sampling to be at all efficient.

We do not have anything built in for truncated multivariates
(which you'd need to define censoring).

There's a chapter in the manual on censored and truncated
data models.

- Bob

On 1/8/14, 10:05 AM, Keith Lau wrote:
> Sorry for my poor mcmc language. I mean that the the location of y* is only known in a restricted region; that is,
> observable data y = 1 if y* >=0, else y=0 if y* <0.
>
> I have read the IRT section, but traditional IRT model like Rasch assumes uncorrelated error between items. That is not
> what I am doing with.
> I note that my model is something similar to multivariate probit model except that I included the person's ability
> parameter within model.
> I also googled this forum about multivariate probit model and found this
> https://groups.google.com/forum/#!topic/stan-users/1QzlcGem-2k
> It seems there was still no clear solution to the model at that time. Is Stan 2.1.0 capable of it?
>
>
> Daniel Lee嚙踝蕭 2014嚙羯1嚙踝蕭8嚙踝蕭P嚙踝蕭嚙確UTC+8嚙磕嚙踝蕭12嚙踝蕭21嚙踝蕭59嚙踝蕭g嚙瘩嚙瘦
>
> Hi Keith,
>
> Just based on your response, it seems like you're not specifying the model correctly.
>
> In Stan, the sampling statements don't actually draw random variables from a distribution. You're specifying a joint
> model over all the parameters. If that doesn't make sense, try rereading chapter 1 of the manual.
>
> Have you taken a look at the IRT models in the Stan repository? Maybe you can start there.
>
>
>
>
> Daniel
>
>
>
>
>
> On Tue, Jan 7, 2014 at 7:39 AM, Keith Lau <genw...@gmail.com <javascript:>> wrote:
>
> Thanks for your positive attention to my question. I am encountering the estimation issue of IRT models similar
> to multivariate probit model, which can be expressed like
>
> y1* = link(theta-mu1) + e1 for item one, where y1* is latent continuous variable, link() is a linking function
> (linear or nonlinear), theta is a random variable following N(0,1), mu1 is an intercept, and e1 is error term
> following N(0,sig). Like same manner, y2* = link(theta-mu2) + e2 is for second item, and so on.
>
> Now, for example, the observable responses are y1 = 1, y2=0, y3=0,嚙皺, yn=1, which are determined if y1*>=0,
> y2*<0, y3*<0,嚙皺, yn*>0. Additionally, the model assumes the error terms are correlated with parameterized
> covariance matrix due to the experimental design. Therefore, I have to estimate theta, mu, and a covariance
> matrix of error term. However, to calculate the multivariate cumulative normal probability is a horrible drag
> and impractical. So, as far as I know, the best answer may to regard the problem in the perspective of censoring
> mechanism. That is, the y* samples can be drawn from truncated multivariate normal distribution.
>
> By the way, I have tried the JAGS to run this model by using one of its function called dinterval. So, the y*
> may be sampled from truncated multivariate normal dist (I嚙踝蕭m not sure if it indeed works or not). However, I did
> not get satisfactory parameter recovery in a simulation study.
>
> Hope above description is clear enough for you.
>
>

Keith Lau

unread,
Jan 9, 2014, 11:14:05 PM1/9/14
to
In binary data, we assume there a continuous latent variable y* that decide whether observable categorical response y will be 1 or 0. Thus, y* is parameter, while y is known data.
Take y = 1 for example, we only know that the y* is above 0.
In Stan manual, there are two methods to deal with this problem. The first one is to constrain the region of y*, the other is to calculate the univariate cumulative probability.
However, in multivariate case, I found it is very inefficient in MCMC by using former method.
I found there are a bit papers tackling this multivariate probit model like this one: a multivariate probit latent variable model for analyzing dichotomous ....
The model (3) is what I am especially interested in his paper, but they adopted MCEM algorithm to estimate parameters.
Is it possible to estimate the model in Stan without EM algorithm?

Thanks for both your patience for my problem.


Bob Carpenter於 2014年1月9日星期四UTC+8下午7時06分17秒寫道:

Bob Carpenter

unread,
Jan 10, 2014, 5:42:12 AM1/10/14
to stan-...@googlegroups.com


On 1/10/14, 5:13 AM, Keith Lau wrote:
> In binary data, we assume there a continuous latent variable y* that decide whether observable categorical response y
> will be 1 or 0. Thus, y* is parameter, while y is known data.
> Take y = 1 for example, we only know that the y* is below 0.
> In Stan manual, there are two methods to deal with this problem. The first one is to constrain the region of y*, the
> other is to calculate the univariate cumulative probability.

You always need to declare a variable with constraints such that
the support of the function is unconstrained (i.e., on all of R^n).

We haven't used the Albert and Chib-type coding for regressions in
Stan because there's nothing to be gained from conjugacy in HMC (other
than perhaps slightly faster iterations). Might it be possible to
just define the regression directly as we do in our logistic regression
examples?

> However, in multivariate case, I found it is very inefficient in MCMC by using former method.
> I found there are a bit papers tackling this multivariate probit model like this one: a multivariate probit latent
> variable model for analyzing dichotomous ...
> <http://www.google.com.hk/url?sa=t&rct=j&q=&esrc=s&source=web&cd=39&ved=0CIcBEBYwCDge&url=http%3A%2F%2Fwww3.stat.sinica.edu.tw%2Fstatistica%2Foldpdf%2FA15n34.pdf&ei=unDPUt-OLIeUiQeV2YHQCw&usg=AFQjCNH3UN9SqqUwopLfndK7JjzQMNLLwQ&sig2=UwhwVCde3HbG-4NdWj-ggw>.
> The model (3) is what I am especially interested in his paper, but they adopted MCEM algorithm to estimate parameters.
> Is it possible to estimate the model in Stan without EM algorithm?

Thanks for the reference. We really haven't thought about
multivariate cumulative distributions, nor do we have a very
sophisticated constraint language.

I'm afraid I, at least, don't have time to follow up with
all the paper references we're getting about models.

- Bob
Reply all
Reply to author
Forward
0 new messages