# Exception: categorical\_logit\_lpmf: log odds parameter\[2\] is nan, but must be finite!

**URL:** <https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443>\
**Category:** General\
**Created:** [August 31, 2024, 4:21pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443 "2024-08-31T16:21:23Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![dmp](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/dmp/32/16409_2.png) [@dmp](https://discourse.mc-stan.org/u/dmp)\
**Post date:** [August 31, 2024, 4:21pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/1 "2024-08-31T16:21:23Z")

</div>

My model returning the message:  
‘Exception: categorical\_logit\_lpmf: log odds parameter[2] is nan, but must be finite! (in ‘test0825 copy 3.stan’, line 93, column 16 to column 66)’  
here is the block that involves this error message.

```no-highlight
            if (conte_it==0){
                vector[4] util = rep_vector(0,4);
                real rho = inv_logit(mit);
                real ur;
                ur =mean(Initiated[(ob-min(itd,14)+1):ob]);
                vector[3] nuse;
                for(j in 1:3){
                    nuse[j] = ur*(1-rho^condur[j])/(1-rho);
                    util[j+1] = -exp(-r*nuse[j]*mit+(r*nuse[j]^2/2*vit)) + intcept[j] + (Constatus[ob]==j+1)*restpar[5];
                }

                target += categorical_logit_lpmf(Choice[ob]|util);
            }

```

The model was able to sample, but the some of the parameter estimates do not converge. Here is the traceplot. How to deal with the ‘NA’ issue? Is the nonconvergence due to this problem?

 ![image](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/3X/5/3/53fde253767b29324623a405eb49fae1c9d496a7.png)

---

<div class="post-metadata">

**Author:** ![Corey.Plate](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/corey.plate/32/20600_2.png) [@Corey.Plate](https://discourse.mc-stan.org/u/Corey.Plate)\
**Post date:** [September 3, 2024, 4:45pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/2 "2024-09-03T16:45:42Z")

</div>

Can you share the entire code?

---

<div class="post-metadata">

**Author:** ![dmp](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/dmp/32/16409_2.png) [@dmp](https://discourse.mc-stan.org/u/dmp)\
**Post date:** [September 3, 2024, 7:09pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/3 "2024-09-03T19:09:13Z")

</div>

Yes, sure. Here is the code.

```no-highlight
data{
    int ncust;
    int nob;
    int nupdt;

    array[nob] int Initiate;
    array[nob] int Complete;
    array[nob] int Choice;
    array[nob] int Itd;
    array[nob] int Conidx;
    array[nob] int Dur;
    array[nob] int Contd;
    array[nob] int Conte;
    array[nob] int Constatus;
    array[nob] int Updt;
    array[nob] int updt_idx;
    array[nob] int updt_idx_ini;

    array[nupdt] real <lower=1e-500> Vit;
    array[nupdt] int Itd_updt;
    array[nupdt] int K;
}

parameters {
   real M;
   real<lower=0.01> V;
   real<lower=0.001> r;
   array[3] real restpar;
   array[nupdt] real Ms;
   array[3] real intcept;
}

transformed parameters {
   array[nupdt] real Ms_bar;
   array[nupdt] real<lower= 0> Vs_bar;
   array[nupdt] real <lower= 0> Vs;
   array[nupdt] real <lower= 0> sqrt_Vbar;
   
   for (up in 1:nupdt){
    real m0 = (Itd_updt[up]==1)? M:Ms[up-1];
    real v0 = (Itd_updt[up]==1)? 10:Vs[up-1];
    
    Vs[up] = 1/(1/v0 + K[up]/V);
    Vs_bar[up] = K[up]*Vs[up]^2/V;
    sqrt_Vbar[up] = sqrt(Vs_bar[up]);
    Ms_bar[up] = Vs[up]/v0*m0 + K[up]*Vs[up]/V*M;
   }
}

model{
    M ~ normal(2,1);
    V ~ lognormal(log(12),0.5);
    r ~ lognormal(log(2),0.5);
    //real r = 2;
    Ms ~ normal(Ms_bar,sqrt_Vbar);
    //row_vector[3] intcept=[-1,-1,-2];
    intcept ~ normal(-1,1);
    restpar ~ normal(1,1);

    row_vector[3] condur = [1,3,12];

    for (ob in 1:nob){
        int itd = Itd[ob];
        int conte_it = Conte[ob];
        int contd_it = Contd[ob];
        int conidx_it = Conidx[ob];
        int updt_idx_it = updt_idx[ob];
        int dur_it = Dur[ob];

        if (itd>1){
            target += bernoulli_logit_lpmf(Initiate[ob]|Ms[updt_idx_ini[ob]]*(conidx_it==0)+(contd_it*1.0/dur_it<=0.3)*restpar[1]+(conte_it*1.0/dur_it<=0.3)*restpar[2]);

            if (Initiate[ob]==1){
                target += bernoulli_logit_lpmf(Complete[ob]|M);
            }

            if(conte_it==0){
                vector[4] util = rep_vector(0,4);
                real mit = Ms[updt_idx_it];
                real rho = inv_logit(mit);
                real ur;
                ur = mean(Complete[(ob-min(itd,14)+1):ob]);
                vector[3] nuse;
                for(j in 1:3){
                    nuse[j] = ur*(1-rho^condur[j])/(1-rho);
                    util[j+1] = nuse[j]*mit - r*nuse[j]^2*Vs[updt_idx_it] + intcept[j] + (Constatus[ob]==j+1)*restpar[3];
                }
                real maxutil = max(util);
                target += categorical_logit_lpmf(Choice[ob]|util-maxutil);
            }

        }

    }

}

```

actually, I fixed some bug in the code. and it seems to converge now although the error message still pops up.  
_Exception: categorical\_logit\_lpmf: log odds parameter[2] is nan, but must be finite! (in ‘test0902.stan’, line 105, column 16 to column 74)_  
_Exception: normal\_lpdf: Scale parameter[1] is 0, but must be positive! (in ‘test0902.stan’, line 70, column 4 to column 34)_  
_Exception: normal\_lpdf: Scale parameter[1] is 0, but must be positive! (in ‘test0902.stan’, line 70, column 4 to column 34)_

Here is the traceplot. For the parameter ‘V’ (on the bottom), however, it still seems to have some autocorrelation.

 ![image](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/3X/0/f/0f23714ec9afb0f6ec1ef19661af98690aa813e6.png)  
the summary of the estimates is as follows:  
 ![image](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/3X/f/9/f9210fc6ecc1ac5ed19b1f8c4bd7ac7ad342c6b7.png)

I am new to both Bayesian Estimation and Stan. Do you think the results have good mixing? I am running with 500 warmups, and 900 iterations just for trial run. the estimate of M is close enough but the true V is 12 which still has some distance from the estimates 13. I am increasing the number of warmups and iterations and running the model again. Hopefully, it could get closer to the true value of V.

---

<div class="post-metadata">

**Author:** ![Corey.Plate](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/corey.plate/32/20600_2.png) [@Corey.Plate](https://discourse.mc-stan.org/u/Corey.Plate)\
**Post date:** [September 3, 2024, 7:33pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/4 "2024-09-03T19:33:18Z")

</div>

Your MSCE for V is pretty high, but given the StdDev, I would not expect a V of 12 to be recovered with more iterations. But that could be because of the data set itself, not the model. There is some clear autocorrelation in V that isn’t present in M. If you want to really check for convergence sensitively, you should use rank plots. A code like this in MATLAB:

```no-highlight
for c = 1:length(C)
 chains((c-1)*size(C{c}(:,p),1)+1:c*size(C{c}(:,p),1)) = C{c}(:,p);
end

chains = tiedrank(chains);
chains = reshape(chains,length(C),[]);

set(figure(1),'OuterPosition',[400 100 1000 1000])

for c = 1:length(C)
 subplot(2,2,c)
 histogram(chains(c,:),'BinEdges',linspace(0,4*length(chains(c,:)),11))
 hold on
 yline(length(chains(c,:))/10,'LineStyle','--','LineWidth',1.5)
end

```

You can read more about rank plots here:  
[https://bayesiancomputationbook.com/markdown/chp\_02.html#exploratory-analysis-of-bayesian-models](https://bayesiancomputationbook.com/markdown/chp_02.html#exploratory-analysis-of-bayesian-models)

> **[Troubleshooting Bayesian cognitive models - PubMed](https://pubmed.ncbi.nlm.nih.gov/36972080/)**
>
> Using Bayesian methods to apply computational models of cognitive processes, or Bayesian cognitive modeling, is an important new trend in psychological research. The rise of Bayesian cognitive modeling has been accelerated by the introduction of...

Rank plots are a sensitive test of both convergence and mixing.

Additionally, you could generate ten data sets from the same generative parameters and see if the recovered parameter estimates fall around the true parameter set.

I would do a print statement on util, maxutil, and util-maxutil to see where this nan is coming from.

---

<div class="post-metadata">

**Author:** ![dmp](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/dmp/32/16409_2.png) [@dmp](https://discourse.mc-stan.org/u/dmp)\
**Post date:** [September 4, 2024, 5:20am UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/5 "2024-09-04T05:20:12Z")

</div>

thank you very much! let me look into it.

---

<div class="post-metadata">

**Author:** ![Bob\_Carpenter](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bob_carpenter/32/9230_2.png) [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)\
**Post date:** [September 8, 2024, 5:42pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/6 "2024-09-08T17:42:00Z")

</div>

> [@dmp](#):
>
> I am new to both Bayesian Estimation and Stan. Do you think the results have good mixing?

Was this for the default 1000 warmup/1000 sampling iterations in 4 chains (i.e., 4000 posterior draws). If so, it’s fine—it’s an integrated autocorrelation time of roughly 7, which is generally considered good mixing.

The way to evaluate if it’s mixing through the right region of the posterior is to simulate data then fit it and check whether the posterior has the right coverage (e.g., 50% intervals contain about 50% of the true values).

> [@dmp](#):
>
> How to deal with the ‘NA’ issue? Is the nonconvergence due to this problem?

This is a `NaN` issue, which is not-a-number, a fixed floating point value for errors that is part of the IEEE standard. `NA` is what R uses for undefined values. NaN values come up when arithmetic is incoherent, like `0 / 0`, or `inf - inf`. I’m guessing that what’s happening is that you are rounding `rho` to `1`, so that `nuse[j]` becomes NaN.

You can instrument your code with print statements to see where things are going wrong. For instance, print out the values of `rho` and `ur` and `nuse[j]` and `util[j + 1]` and you can see where `util` gets a NaN introduced.

Are you sure you want `-exp(...)` before applying `categorical_logit`? The categorical\_logit distribution applies softmax (exp(x) / sum(exp(x))) to the argument, so it’s getting doubly exponentiated. This is a quick route to overflow or underflow.

The thing to do to mitigate this problems is to work on the log scale for as long as you can. This may, unfortunately, require unfolding categorical\_logit and coding it up manually, because that can be part of the problem.

For example, we know that your `target +=` statement will amount to this.

```no-highlight
target += log(softmax(util))[Choice[ob]];

```

Now you can start unfolding the algebra back through the definition of `util`.

But given that things seem to be mixing OK, I wouldn’t bother with this. I would try to do the validation with simulated data, though.

> [@dmp](#):
>
> I am increasing the number of warmups and iterations and running the model again. Hopefully, it could get closer to the true value of V.

this is very unlikely given that the existing fit has an MCMC standard error (MCSE) of only 0.02. You wouldn’t expect to find results further away than 4 times the MCSE.

You have to be careful about what you mean by the “true” value of `V`. We’re trying to find the posterior mean of `V`, not recover the simulated value. Here’s an example. If I simulate a bunch of data i.i.d. via `V ~ normal(0, 1)`, it might have a mean of 0.3 and a standard deviation of 0.9, so that you don’t expect to recover `mu = 0` and `sigma = 1` by fitting `V ~ normal(mu, sigma)`. In the limit of infinite data, you will you get convergence to the true values with a well-behaved model.

---

<div class="post-metadata">

**Author:** ![dmp](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/dmp/32/16409_2.png) [@dmp](https://discourse.mc-stan.org/u/dmp)\
**Post date:** [September 10, 2024, 6:18pm UTC](https://discourse.mc-stan.org/t/exception-categorical-logit-lpmf-log-odds-parameter-2-is-nan-but-must-be-finite/36443/7 "2024-09-10T18:18:37Z")

</div>

Hi Bob,  
thank you very much for your such detailed explainations. To avoid NaN that might be produced due to rho → 1, I set rho = (rho,1.0-1e-6), and I unfolded the loglikelihood as you suggested using log\_sum\_exp. Now the error terms disapeared.

The full code is as follows.

```no-highlight

data{
    int ncust;
    int nob;
    int nupdt;

    array[nob] int Initiate;
    array[nob] int Complete;
    array[nob] int Choice;
    array[nob] int Itd;
    array[nob] int Conidx;
    array[nob] int Dur;
    array[nob] int Contd;
    array[nob] int Conte;
    array[nob] int Constatus;
    array[nob] int Updt;
    array[nob] int updt_idx;
    array[nob] int updt_idx_ini;

    array[nupdt] real <lower=1e-500> Vit;
    array[nupdt] int Itd_updt;
    array[nupdt] int K;
}

parameters {
   real M;
   real<lower=0.01> V;
   array[3] real restpar;
   array[nupdt] real Ms;
   array[3] real intcept;
}

transformed parameters {
   array[nupdt] real Ms_bar;
   array[nupdt] real<lower= 1e-10> Vs_bar;
   array[nupdt] real <lower= 1e-10> Vs;
   array[nupdt] real <lower= 1e-10> sqrt_Vbar;
   
   for (up in 1:nupdt){
    real m0 = (Itd_updt[up]==1)? M:Ms[up-1];
    real v0 = (Itd_updt[up]==1)? 10:Vs[up-1];
    
    Vs[up] = 1/(1/v0 + K[up]/V);
    Vs_bar[up] = K[up]*Vs[up]^2/V;
    sqrt_Vbar[up] = sqrt(Vs_bar[up]);
    Ms_bar[up] = Vs[up]/v0*m0 + K[up]*Vs[up]/V*M;
   }
}

model{
    M ~ normal(1,1);
    V ~ lognormal(log(4),0.5);
    //r ~ lognormal(log(1),0.5);
    real r = 1.5;
    Ms ~ normal(Ms_bar,sqrt_Vbar);
    //row_vector[3] intcept=[-1,-1,-2];
    intcept ~ normal(2,1);
    restpar ~ normal(1,1);

    row_vector[3] condur = [1,3,12];

    for (ob in 1:nob){
        int itd = Itd[ob];
        int conte_it = Conte[ob];
        int contd_it = Contd[ob];
        int conidx_it = Conidx[ob];
        int updt_idx_it = updt_idx[ob];
        int dur_it = Dur[ob];

        if (itd>1){
            target += bernoulli_logit_lpmf(Initiate[ob]|Ms[updt_idx_ini[ob]]*(conidx_it==0)+(contd_it*1.0/dur_it<=0.3)*restpar[1]+(conte_it*1.0/dur_it<=0.3)*restpar[2]);

            if (Initiate[ob]==1){
                target += bernoulli_logit_lpmf(Complete[ob]|M);
            }

            if(conte_it==0){
                vector[4] util = rep_vector(0,4);
                real mit = Ms[updt_idx_it];
                real rho = inv_logit(mit);
                rho = min(rho,1.0-1e-6);
                real ur;
                ur = mean(Complete[(ob-min(itd,7)+1):ob]);
                vector[3] nuse;
                for(j in 1:3){
                    nuse[j] = ur*(1-rho^condur[j])/(1-rho);
                    //util[j+1] = nuse[j]*mit - r*nuse[j]^2*Vs[updt_idx_it] + intcept[j] + (Constatus[ob]==j+1)*restpar[3];
                    util[j+1] = -exp(-r*nuse[j]*mit +(r*nuse[j])^2/2*Vs[updt_idx_it]) + intcept[j] + (Constatus[ob]==j+1)*restpar[3];
                }
                //real maxutil = max(util);
                //target += categorical_logit_lpmf(Choice[ob]|util-maxutil);
                real maxutil = max(util);
                vector[4] log_softmax_util;
                for (k in 1:4){
                    log_softmax_util[k] = util[k] - maxutil-log_sum_exp(util-maxutil);
                }

                target += log_softmax_util[Choice[ob]];
            }

        }

    }

}

```

However, the model still produces the following error terms:  
**error 1.**  
Exception: normal\_lpdf: Scale parameter[1] is 0, but must be positive! (in ‘/tmp/tmpr\_\_a\_zjj/tmpvebvcty9.stan’, line 74, column 4 to column 34)  
**error 2.**  
18:46:35 - cmdstanpy - WARNING - Some chains may have failed to converge.  
Chain 3 had 2 iterations at max treedepth (0.2%)  
Chain 4 had 15 iterations at max treedepth (1.9%)  
Use the “diagnose()” method on the CmdStanMCMC object to see further information.

the first error refers to the line  
Ms ~ normal(Ms\_bar,sqrt\_Vbar);  
I have already constrained sqrt\_Vbar with the lower bound 1e-10, why the scale parameter would still be 0?

for the error 2, would this message be a concern?

The traceplot is as follows. Do you think the traceplot for V looks acceptable? I am still trying to get a sense of what looks acceptable and what does not. Thank you very much in advance for your help.

 ![image](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/3X/e/c/ec59db4ce31b2119df9f2d5eb95699abae243913.jpeg)