add a new distribution similar the basic distributions in Stan

[Bistoon]

Dear All,
Hi.

I use the Kumaraswamy distribution instead of the Beta distribution with below pdf:
f(x; a,b) = a * b * x^(a - 1)* (1 - x^a)^(b - 1)  for a,b > 0 and 0 < x < 1

In below code I defined kumaraswamy_lpdf and use as a new distribution for "x" observations.
The below model work correctly:

// beginning model
functions
    {
    real kumaraswamy_lpdf(real x, real a, real b)
        {
        return log(a) + log(b) + (a-1)*log(x) + (b-1)*log(1-pow(x,a)) ;       
        }
    }
data
    {
    int<lower=1> n;
    vector<lower=0,upper=1>[n] x;
    }
parameters
    {
    real<lower=0> a;
    real<lower=0> b;
    }
model
    {
    a ~ gamma(.001,.001);
    b ~ gamma(.001,.001);
    for(i in 1:n) x[i] ~ kumaraswamy(a,b);
    }
// end model

I expect to be able to use "x ~ kumaraswamy(a,b)" instead of "for(i in 1:n) x[i] ~ kumaraswamy(a,b)" but it does not work. How can I define kumaraswamy_lpdf to use "x ~ kumaraswamy(a,b)" simillar "x ~ beta(a,b)"? This is so useful in the new distribution research.
Thanks in advance.

Best,

Bistoon

[Daniel Lee]

Hi Bistoon,

The first argument if your function is real x. Change that to vector x and work through the math to return the sum of the lpdf evaluated at each x. 

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

[Ben Goodrich]

On Tuesday, September 6, 2016 at 4:17:37 PM UTC-4, Bistoon wrote:

model
    {
    a ~ gamma(.001,.001);
    b ~ gamma(.001,.001);
    for(i in 1:n) x[i] ~ kumaraswamy(a,b);
    }

You don't want gamma priors with tiny shapes and rates. Also, you don't want that parameterization of the Kumaraswamy distribution. Try this one

https://scholar.google.com/scholar?cluster=18346740915892174762&hl=en&as_sdt=0,33&sciodt=1,33

Ben

[Bistoon]

Dear Lee and Goodrich

Hi,

If I change _lpdf to vector, I will have below error:

Require real return type for probability functions ending in _log, _lpdf, _lpmf, _lcdf, or _lccdf.

Also, the Kumaraswamy distribution and Gamma priors was example to explain main question in here. I need to add a new distribution in the core of Stan similar WBDev in BUGS for have a speed sampling.

Best Regards,

Bistoon

[Ben Goodrich]

On Tuesday, September 6, 2016 at 11:40:29 PM UTC-4, Bistoon wrote:

If I change _lpdf to vector, I will have below error:

Require real return type for probability functions ending in _log, _lpdf, _lpmf, _lcdf, or _lccdf.

The return type needs to be real rather than a vector or anything else. The arguments can be vectors. In your case, it would be

real kumaraswamy_lpdf(vector x, real a, real b) {
  int N;
  real lpdf;
  N = rows(x);
  lpdf = 0;
  for (n in 1:N) lpdf = lpdf + log1m(x[n] ^ a);
  lpdf = lpdf * (b - 1) + (a - 1) * sum(log(x)) + N * (log(a) + log(b));
  return lpdf;
}

It would be better to pass in sum(log(x)) as an argument. But you still don't want to use that parameterization.

Ben

 

[Bistoon]

Thank you very muck Ben. This is better and faster code. I do not know about parameterization and I want to use similar code for my new distribution in the future.

Bistoon

[Bob Carpenter]

There are various ways to parameterize distributions; for instance,
the normal can use scale (sigma), variance (sigma^2), or precision (1 / sigma^2).
negative_binomial has many parameterizations as does the gamma; check
out the Wikipedia articles.

Using _lpdf functions with sampling notation won't work until Stan 2.12,
which is due out any day now (it's already been tagged and is just
winding through the interfaces). Until then, you'd have to use the old
_log syntax (see the deprecated functions section of the manual) if
you want to use ~.

- Bob

[Bistoon]

Thank you very much Bob. I am eager to see Stan 2.12. Also Thanks about explanation about parameterization.

Let me here to thank Stan Group for good and fast response. The Bayesian analysis is comfortable with you and your software.

With Best Wishes,
Bistoon