fitting latent class multinomial logit model in rstan

[malai...@gmail.com]

Hi All,

Thank you Stephen for the script. It is very useful. One question though;

Is there a line in you're script you specified the dependent variable?

Because in my data, i have one dependent variable which is injury severity (0-no injury, 1- severe injury) and also i have 40 predictor variables, some with categorical and some continuous. The sample size of my data is 846 crashes.

Due to the nature of my data i thought fitting the data with latent class logit model will help in explaining the unobserved heterogeneity and also yielding robust coefficient.

Please assist me in clarifying this situation.

Thanks.

data {
    int N;
    int J;
    int K;
    vector[N] x;
    matrix[N,J] predictor;

}
parameters{
    ordered[K] means;
    vector<lower=0>[K] sds;

    vector<lower=0>[K] pi_logit[N];
    vector<lower=0>[K] beta0;

    matrix[J,K] beta;
}
transformed parameters{
    simplex[K] pi[N];
//    vector[K] pi_logit[N];

    for(n in 1:N){
        pi[n] = softmax(pi_logit[n]);
    }
}
model{
    vector[K] lp_nk[N];
    vector[N] lp_n;
    
    means ~ normal(0,4);
    sds ~ cauchy(0,3);

    beta0 ~ logistic(0,3);
    for(k in 1:K){
        beta[1:J,k] ~ logistic(0,3);

    }
    for(n in 1:N){
        for(k in 1:K){

            pi_logit[n,k] ~ logistic(beta0[k] + predictor[n]*beta[1:J,k],1);

        }
    }
    
    for(n in 1:N){
        for(k in 1:K){
            lp_nk[n,k] = log(pi[n,k]) + normal_lpdf(x[n] | means[k], sds[k]);
        }
        lp_n[n] = log_sum_exp(lp_nk[n]);
    }
target += sum(lp_n);
}

[Stephen Martin]

In my model, it's a multinomial (softmax) regression predicting latent group membership to one of K clusters, the clusters which define separate gaussian distributions.
The outcome variable in my case is x[N], where x is assumed to come from one of K gaussian distributions. The softmax regression predicts one's log-odds of belonging to one gaussian over another.
You could modify this to predict one's log-odds of belonging to one normal-theory regression over another, where the K components differ in the regression coefficients and intercept.
If you have binary outcome data, then you need to modify the likelihood, such that you assume each bernoulli trial comes from one logistic model as opposed to K other logistic models.

That said, if you have that many predictors, I don't see this going well. I'm guessing it will be highly multimodal and you'll have quite a lot of label-switching. Even with, e.g., an ordered intercept constraint, that many predictors may cause the sampler to switch modes within a chain (and between chains, which is less of an issue than within-chain), because the constraint on an intercept may not be enough to identify the model; it'll still likely have multiple modes it can pass through.

[malai...@gmail.com]

Thank you Stephen for the clarification,

considering the large number of predictor variables that i have in my data, do you have any idea of a suitable model i can use to fit my model?

Thank you all once again.

[Stephen Martin]

It would really depend on what your model actually is; as in, what are your predictors?

Personally, I wouldn't really bother with such a massive list of predictors; if you can combine predictors into something more meaningful (whether that's via summary stats or dimension reduction [latent variables, PLS combination]), that would probably be better.
If these are 41 fairly unrelated variables, I'd just have to ask why you are using 41 predictors.

> --
> You received this message because you are subscribed to a topic in the
> Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> stan-users+...@googlegroups.com
> <mailto:stan-users+...@googlegroups.com>.
> To post to this group, send email to stan-...@googlegroups.com
> <mailto:stan-...@googlegroups.com>.
> For more options, visit https://groups.google.com/d/optout.

[malai...@gmail.com]

Hi all,

The predictor variable in my model are highly correlated with one another, and its important for me to use each of them in order to capture the relationship of each of the predictor variable with the dependent variable.

I need to capture the contribution of each of the independent variable that will be significant.

> <mailto:stan-users+unsub...@googlegroups.com>.

[Stephen Martin]

If many of your predictors are highly correlated, I would first ask *why* they are correlated, try reducing some of the 41 down to a smaller number of constructs. You can do this via PLS-like combined variables or CFA-like latent variables, but I would certainly hesitate to throw 41 predictors into a model, ESPECIALLY if they are highly correlated. Aside from the label-switching problem inherent in bayesian mixtures, you're going to hit some multicollinearity problems.

On 10/14/2016 11:01 AM, malai...@gmail.com wrote:
> Hi all,
>
> The predictor variable in my model are highly correlated with one
> another, and its important for me to use each of them in order to
> capture the relationship of each of the predictor variable with the
> dependent variable.
> I need to capture the contribution of each of the independent variable
> that will be significant.
>
>
> On Thursday, October 13, 2016 at 3:48:43 PM UTC-4, Stephen Martin wrote:
>
> It would really depend on what your model actually is; as in, what
> are your predictors?
>
> Personally, I wouldn't really bother with such a massive list of
> predictors; if you can combine predictors into something more
> meaningful (whether that's via summary stats or dimension reduction
> [latent variables, PLS combination]), that would probably be better.
> If these are 41 fairly unrelated variables, I'd just have to ask why
> you are using 41 predictors.

> <https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe>.

> > To unsubscribe from this group and all its topics, send an email to

> > stan-users+...@googlegroups.com <javascript:>
> > <mailto:stan-users+...@googlegroups.com <javascript:>>.

> > To post to this group, send email to stan-...@googlegroups.com

> <javascript:>
> > <mailto:stan-...@googlegroups.com <javascript:>>.

> > For more options, visit https://groups.google.com/d/optout

> <https://groups.google.com/d/optout>.

>
> --
> You received this message because you are subscribed to a topic in the
> Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> stan-users+...@googlegroups.com

> <mailto:stan-users+...@googlegroups.com>.

[stijn masschelein]

While Stephen's suggestion is more sensible than what I am going to suggest, there are other things you could do. I never haven been able to use it myself but Andrew published a paper which suggests a Bayesian approach to your problem (section 5.2.3). If you have multiple variables x1, x2, x3 that are measuring a similar construct e.g. economic performance of a country as measured by (GDP, employment, ...), you can model their contribution with the parameters a, b1, b2, b3 and a strong prior on b1, b2, b3.

y = a (b1 x1 + b2 x2 + b3 x3) 

b1, b2, b3 ~ normal(1/3, sigma)

and a suitable prior for a. a is than the contribution of economic performance to your outcome and b1, b2, b3 are the contribution of your variables to economic performance. The joint prior for b1, b2, and b3 should help with potential multicollinearity issues.

Stijn

>     > <mailto:stan-users+unsub...@googlegroups.com <javascript:>>.

>     > To post to this group, send email to stan-...@googlegroups.com
>     <javascript:>
>     > <mailto:stan-...@googlegroups.com <javascript:>>.
>     > For more options, visit https://groups.google.com/d/optout
>     <https://groups.google.com/d/optout>.
>
> --
> You received this message because you are subscribed to a topic in the
> Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> stan-users+...@googlegroups.com

> <mailto:stan-users+unsub...@googlegroups.com>.

[Stephen Martin]

That looks very similar to the pls-like approach actually.

>     > <mailto:stan-users+unsubscribe@googlegroups.com <javascript:>>.

>     > To post to this group, send email to stan-...@googlegroups.com
>     <javascript:>
>     > <mailto:stan-...@googlegroups.com <javascript:>>.
>     > For more options, visit https://groups.google.com/d/optout
>     <https://groups.google.com/d/optout>.
>
> --
> You received this message because you are subscribed to a topic in the
> Google Groups "Stan users mailing list" group.
> To unsubscribe from this topic, visit
> https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe.
> To unsubscribe from this group and all its topics, send an email to
> stan-users+...@googlegroups.com

> <mailto:stan-users+unsubscribe@googlegroups.com>.

> To post to this group, send email to stan-...@googlegroups.com
> <mailto:stan-...@googlegroups.com>.
> For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/cDnjAVCN90Y/unsubscribe.

To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.

[malai...@gmail.com]

Hi All,

Thank you for the useful suggestion. I will work on the proposed solutions.

[malai...@gmail.com]

Hi Stephen,

The reason predictor variables are correlated in my data is because they highly depend on each other. Example variables weather condition, and road surface condition i.e. mostly the when it rains the road surface will be wet.

Another variable is traffic volume and speed, i.e. normally roads with high traffic volume have high speed limits too.

On Thursday, October 13, 2016 at 9:05:41 AM UTC-4, malai...@gmail.com wrote:

[krishna mohan]

Hi

Another suggestion that I would make to help your modelling a little easier is to potentially use Shapley Value Regression on some of your predictors. Here are the steps I would recommend.

1. Download relaimpo R package that allows you to use Shapley Value regression with the metric labelled as LMG in the package (Prof Gromping's website has a link 

http://prof.beuth-hochschule.de/groemping/software/relaimpo-relative-importance-of-regressors/download-package-non-us-version/

2. Either via PLSR or CFA models or even exploratory components through a PCA, you can reduce the 41 predictors to a more manageable set. I prefer using PLSR as the variance components are extracted to maximize the varianbce in your dependent variable. 

3. Then take each of the top loading predictor and predict that IV as a function of other IV predictors using Shapley Value Regression. This will allows you to get a fairly accurate assessment of how much the different predictors are contributing to the main predictor vairables.

4. Then you can potentially build a hierarchy of predictors driving your DV.

5. Please note that running Shapley Value via relaimpo can drain your machine's memory as it does this analyses by testing all possible combinations and orders of variables entered into the analyses. However, it in my experience is the best way to do it.

Krishna