# How much to log in a custom made pdf (for the linear ballistic accumulator)
**URL:** https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752
**Category:** Modeling
**Tags:** specification, cognitive-science
**Created:** [December 13, 2017, 3:12pm UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752 "2017-12-13T15:12:12Z")
**Posts on this page:** 7
**Page:** 1
### Author:  [@bnicenboim](https://discourse.mc-stan.org/u/bnicenboim)
#### Post date: [December 13, 2017, 3:12pm UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/1 "2017-12-13T15:12:12Z")
Hi,
I’m trying to implement a custom made for the linear ballistic accumulator with a drift from a normal+, to fit responses and reaction times. I’m taking the pdfs and cdf of each accumulator from here [single-lba.r](https://github.com/rtdists/rtdists/blob/master/R/single-lba.r) (dlba\_norm\_core and plba\_norm\_core), and then basically the PDF of the model is (for two accumulator) pdf\_r \* (1-cdf\_n), where r is the accumulator that finished first (and “gave” the answer.)
So this is (a simplification) of the code for the pdf of one accumulator in R :
eta\_i is real, b\>0, A \>b, s\_i\>0 , t\>0
```
pdf <- function(t=.3, A=1, b=1.4, eta_i=3.5, s_i=1){
Phi <- function(x) pnorm(x,0,1)
npdf <- function(x) dnorm(x,0,1)
denom <- Phi(eta_i/s_i)
zs <- t*s_i
zu <- t*eta_i
chiminuszu <- b-zu
chizu <- chiminuszu/zs
chizumax <- (chiminuszu-A)/zs
pdf <- (eta_i*(Phi(chizu)-Phi(chizumax)) + s_i*(npdf(chizumax)-npdf(chizu)))/(A*denom)
return(pdf)
}
```
I need to get the log(pdf), so my question is how “inside” should I have the logs. There are many things that can under/overflow very quickly, if eta\_i is very large positive or negative, if s\_i is very small or of b is very large. My first instinct was too transform everything into a sum of terms minus another sum of terms, and then
the lpdf becomes `log_diff_exp(log_sum_exp(logterms1),log_sum_exp(logterms2)) - log(A) - normal_lcdf(eta_i/s_i|0,1)`
But then the log of the Phi gets to 0 or to -Inf and the eta\_i can be negative during the sampling, and then I have a lot of if clauses to catch all these cases and avoid NaNs.
My lpdf of the accumulator ends up being horribly long, and even worst for the lccdf. **Am I on the right track?** My custom made pdfs and cdfs seem to give the same answers as the package single-lba.r, except for extreme values.
But I was wondering if what I did makes sense and if there is a way to optimize my code. (I want to make this hierarchical as a next step).
A second question is when (or if) I should use log(Phi\_approx(.)) instead of the normal\_lcdf (that are inside the lpdf and lccdf of the acummulator), (when should I use this function if at all?)
I’m attaching the code just in case. [tLBA.stan](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/a/a3f5b40089e74fabd8d4c111754b1423aa773004.stan) (6.8 KB)
And this is to test the model and to test the lpdf and lccdf: [test\_tnormLBA.R](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/f/f8cef5118792b1f1b35729ac364caf5afb3a22fa.R) (4.6 KB)
Bests,
Bruno
---
### Author:  [@bnicenboim](https://discourse.mc-stan.org/u/bnicenboim)
#### Post date: [December 14, 2017, 10:15am UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/2 "2017-12-14T10:15:29Z")
Less ugly (and working version) of the code here:[tLBA.stan](https://canada1.discourse-cdn.com/flex030/uploads/mc_stan/original/2X/f/f46e3df93e03f47f501d3d572813a9bf32561ad7.stan) (6.2 KB)
But all my questions are still relevant
---
### Author:  [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)
#### Post date: [January 9, 2018, 5:30am UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/3 "2018-01-09T05:30:14Z")
> [@bnicenboim](#):
>
> so my question is how “inside” should I have the logs.
You want to do as much as possible on the log scale, so I think your first instict is right.
We don’t have good examples in the code as we’re using standard CDF and CCDF algorithms and then just naively taking logs.
> [@bnicenboim](#):
>
> when (or if) I should use log(Phi\_approx(.)) instead of the normal\_lcdf
It uses a logistic approximation that’s much more efficient to calculate and pretty close. So I’d try it both ways and see if you get the same answer and then use `Phi_approx` where possible. That’d give you more opportunites to roll in some log terms.
---
### Author:  [@bnicenboim](https://discourse.mc-stan.org/u/bnicenboim)
#### Post date: [January 9, 2018, 8:16am UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/4 "2018-01-09T08:16:24Z")
Thanks for the answers.
> [@Bob\_Carpenter](#):
>
> It uses a logistic approximation that’s much more efficient to calculate and pretty close. So I’d try it both ways and see if you get the same answer and then use Phi\_approx where possible. That’d give you more opportunites to roll in some log terms.
Why would Phi\_approx allow me to roll in more log terms? normal\_lcdf is already logged, isn’t it?
And related to this, would Phi\_approx allow me to go more negative than -37.5 and more positive than 8.25 without under/overflowing? The manual is not very clear in this, and it seems not…
And since I have you here, is there a way to approximate log(Phi(b) - Phi(b-a)) to avoid over/underflow, this is the most problematic aspect of my pdf.
---
### Author:  [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)
#### Post date: [January 10, 2018, 8:36am UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/5 "2018-01-10T08:36:28Z")
> [@bnicenboim](#):
>
> Why would Phi\_approx allow me to roll in more log terms? normal\_lcdf is already logged, isn’t it?
The implementation of `normal_lcdf` just applies `log()` to `normal_cdf`, so it doesn’t operate on the log scale. The logistic form is simple enough that you could easily implement it stably on the log scale.
---
### Author:  [@bnicenboim](https://discourse.mc-stan.org/u/bnicenboim)
#### Post date: [January 10, 2018, 9:06am UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/6 "2018-01-10T09:06:51Z")
> [@Bob\_Carpenter](#):
>
> The implementation of normal\_lcdf just applies log() to normal\_cdf, so it doesn’t operate on the log scale. The logistic form is simple enough that you could easily implement it stably on the log scale.
Ok, I understand, and any insight about log(Phi(b) - Phi(b-a)) ? Or I can’t do more than log(Phi\_approx(b) - Phi\_approx(b-a))? (Or this is probably more for stackexchange than for here).
---
### Author:  [@Bob\_Carpenter](https://discourse.mc-stan.org/u/Bob_Carpenter)
#### Post date: [January 17, 2018, 6:47pm UTC](https://discourse.mc-stan.org/t/how-much-to-log-in-a-custom-made-pdf-for-the-linear-ballistic-accumulator/2752/7 "2018-01-17T18:47:02Z")
I’d try substituting the definition for `Phi_approx` into
`log(Phi_approx(b) - Phi_approx(b-a))`
and see if you can reduce it at all to keep everything on the log scale. We have a `log_diff_exp` function that can help, as this will be:
```
log_diff_exp(log(Phi_approx(b)), log(Phi_approx(b - a)))
```
so then it’s just a matter of figuring out if `log(Phi_approx(u))` can be reduced to something more stable.