0
$\begingroup$

I want to fit some samples to the right tail of a Weibull distribution. To fix the notation:

  • the samples are $\{X_i\}_{i=1,\ldots,n}$,

  • all samples are greater than a fixed threshold $L$,

  • the $X_i$'s are distributed according to a Weibull distribution left truncated at $L$: $$ X_i \sim F(x)=1-\mathrm e ^{-\frac{x^\alpha-L^\alpha}{\beta^\alpha}}. $$

For simplicity in the following I assume $L=1$ and replace $X_i$ with $Y_i=X_i/L$. I want to estimate $\alpha$ and $\beta$ with MLE. The log-likelihood function is: $$ \ell(\alpha,\beta;Y_i) = \sum_i\left(-\frac{Y_i^\alpha-1}{\beta^\alpha}+\log\alpha-\alpha\log\beta+(\alpha-1)\log Y_i\right). $$

By setting $\partial\ell/\partial\beta = 0$ I obtain $$ \beta^\alpha = \frac{\sum_i (Y_i^\alpha-1)}{n}. $$

Now to find $\alpha$ I can proceed in two different ways. I can set $\partial\ell/\partial\alpha = 0$ and then substitute $\beta$ in the resulting equation. This way I get: $$ 0= \frac{1}{\alpha} + \frac{\sum_i\log Y_i}n - \frac{\sum_i Y_i^\alpha\log Y_i}{\sum_i (Y_i^\alpha-1)}. $$

I can also define $\tilde\ell(\alpha)$ by substituting the expression for $\beta$ in $\ell$: $$ \tilde\ell(\alpha)=-n+n\log\alpha-n\log\left(\frac{\sum_i (Y_i^\alpha-1)}{n}\right) +(\alpha-1)\frac{\sum_i\log Y_i}n, $$ and setting $\partial\tilde\ell(\alpha)/\partial\alpha=0$, reaching the same result as above.

The equation for $\alpha$ can then be solved numerically.

EDIT I found a couple of mistakes in the derivations. Now everything works as expected.

Leaving the previous text below:

Now I have two problems:

Problem 1. The two equations for $\alpha$ are different and not equivalent. I suspect this is one or more calculation mistake on my side, but I can't find it.

Problem 2. When I apply the formulas to some simulated data, the estimates for $\alpha$ are wildly wrong.

$\endgroup$

1 Answer 1

1
$\begingroup$

Without reading through your derivation, I imagine there was an issue with your implementation of the optimization.

If you wanted to implement it numerically, here's a quick example (in R).

set.seed(10101)

n = 1000
a = 3
b = 0.1

x = rweibull(n, a, b)
L = quantile(x, 0.2)

x = x[x>=L]

lik = function(pars, x, L) {
  
  n = length(x)
  
  -(n*(log(pars[1]) - log(pars[2])) + (pars[1]-1)*sum(log(x/pars[2])) - sum((x/pars[2])^pars[1]) + n*(L/pars[2])^pars[1])
  
}

o = optim(par = c(1, 2),
          fn = lik,
          method = "BFGS",
          L = L,
          x = x)

o$par
$\endgroup$
2
  • $\begingroup$ Thanks for the R code. Unfortunately I have to eventually implement this in VBA :( and I need the intermediate steps. $\endgroup$
    – AndreA
    Feb 10, 2023 at 12:42
  • $\begingroup$ I have managed to find my mistake. Thanks anyways! $\endgroup$
    – AndreA
    Feb 10, 2023 at 15:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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