1
$\begingroup$

I am fitting a Weibull distribution to some data in Stan. I am trying to reproduce some published values of parameters from a paper. However I am running into problems because I believe the normalizing constant does not match. The paper gives the pdf equation as follows:

$p(x|\mu,\nu) = \frac{1}{K} (\frac{x}{\nu})^{\mu - 1} exp(-(\frac{x}{\nu})^{\mu})$

However the Weibull pdf in Stan is:

$p(x|\mu,\nu) = \frac{\mu}{\nu} (\frac{x}{\nu})^{\mu - 1} exp(-(\frac{x}{\nu})^{\mu})$

When I fit the distribution to the same data as the paper, I get different fitted values for the shape and scale parameters ($\mu$ and $\nu$) than the ones in the paper, but the paper gives no indication of how to find the normalizing constant K. Is there a way to determine the correct value for the constant so that I can get the correct values of shape and scale parameters?

Here is the (very simple) Stan model I fit:

data {
    int<lower=0> N;
    vector<lower=0>[N] x;
}

parameters {
    // Weibull density
    real<lower=0> mu;
    real<lower=0> nu;
}

model {
    // Priors: Weibull density
    mu ~ lognormal(1, 1);
    nu ~ lognormal(1, 1);

    // Likelihood: Weibull density
    x ~ weibull(mu, nu);
}
|cite|improve this question
$\endgroup$
  • 1
    $\begingroup$ I'm also interested, in a more general sense, in whether improperly specifying the normalizing constant would cause different parameter values to be fit. $\endgroup$ – qdread Feb 21 '18 at 12:55
  • 1
    $\begingroup$ Do you have a link to the paper? $\endgroup$ – Jarle Tufto Feb 23 '18 at 12:52
  • 1
    $\begingroup$ @JarleTufto Unfortunately it's behind a paywall but here is the link (the relevant thing is equation 4b): onlinelibrary.wiley.com/doi/10.1111/j.1461-0248.2006.00915.x/… $\endgroup$ – qdread Feb 23 '18 at 13:03
  • 1
    $\begingroup$ As far as I understand, what they are fitting are truncated weibull distributions with a lower truncation point $D_0>0$ being the size (diameter) at recruitment. So this leads to a different normalizing constant (some function of $\nu$, $\mu$ and $D_0$ that they derive in Appendix S1), also different from what is used by Stan which only has the non-truncated distribution built in I guess. $\endgroup$ – Jarle Tufto Feb 23 '18 at 13:43
  • $\begingroup$ @ Jarle Tufto is it possible to specify the log pdf of truncated Weibull manually in Stan? $\endgroup$ – qdread Feb 23 '18 at 14:11
3
+25
$\begingroup$

Perhaps the answer you are looking for is to change your likelihood statement to

for (i in 1:N) x[i] ~ weibull(mu, nu) T[L, ];

where L is the lower truncation point. But you should not actually do it that way. Stan arrives at that answer by subtracting the logarithm of the complimentary CDF evaluated at L. Thus, it is equivalent (and faster) to write

target += weibull_lpdf(x | mu, nu) - N * weibull_lccdf(L | mu, nu);

presuming that L is a scalar truncation point that applies to all N observations. Chapter 12 of the Stan User Manual has more examples.

If you actually want to see the math, you can use the Mathematica syntax

PowerExpand[Log[FullSimplify[ PDF[WeibullDistribution[\[Mu],\[Nu]], x] / (1 - CDF[WeibullDistribution[\[Mu],\[Nu]], L]), Assumptions->{x>0, L>0}]]]

which evaluates to

(L^\[Mu]-x^\[Mu]) \[Nu]^-\[Mu]-Log[x]+Log[\[Mu]]+\[Mu] (Log[x]-Log[\[Nu]])
|cite|improve this answer
$\endgroup$
  • $\begingroup$ can you specify a lower and upper truncation like this: target += weibull_lpdf(x | mu, nu) - N * weibull_lccdf(L | mu, nu) - N * weibull_lccdf(U | mu, nu); @Ben Goodrich $\endgroup$ – qdread Feb 26 '18 at 12:11
  • $\begingroup$ For two-sided truncation, you need to divide the PDF by the probability of being on the (L, U) interval, which is the CDF evaluated at U minus the CDF evaluated at L. And then take the logarithm because Stan uses the log-kernel. So, it would be target += weibull_lpdf(x | mu, nu) - N * log(weibull_cdf(U | mu, nu) - weibull_cdf(L | mu, nu)). Again, that is assuming U and L are scalars. If not, then you have to do the correction one observation at a time or analytically work out the log of the difference in the CDFs to hopefully obtain an expression that can be vectorized. $\endgroup$ – Ben Goodrich Feb 27 '18 at 15:57
  • $\begingroup$ This is the way that ended up working for me. After fitting the truncated Weibull in Stan with both a lower limit ll and upper limit ul and outputting predicted values in R, I divided the predicted values by sum(1 - pweibull(q = c(ll,ul), shape, scale)) which rescales the area under the truncated Weibull to 1. @Ben Goodrich apologies for not officially awarding the bounty, I was out of town and missed the deadline. $\endgroup$ – qdread Mar 22 '18 at 14:33
1
$\begingroup$

If you can get easy access to it the book The Weibull Distribution: A Handbook by Horst Rinne covers the general three parameter Weibull down to the one. There is also a section on the truncated Weibull distributions. You may want to consult it for how different parameterizations are derived.

|cite|improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.