0
$\begingroup$

I've been trying to estimate the two-parameter (a,b) Weibull distribution (loc. param. = 0).

$$f(t;a,b)=\frac{b}{a}\left(\frac{t}{a}\right)^{b-1}\exp(- \left(\frac{t}{a}\right)^b) $$

To find the estimates I am using maximum likelihood and NR method.

t = np.array(list(range(1, 1000))) #data input
def gradient(a,b): 
    grad_a = -10*b/a + b/a*np.sum((t/a)**b)
    grad_b = 10/b - 10*(np.log(a)) + np.sum(np.log(t)) - np.sum(((t/a)**b)*np.log(t/a))   
    grad_matrix = np.array([grad_a, grad_b])
    return grad_matrix
    
def hessian(a,b): 
    hess_a = 10*b/a**2 + (b*(b+1)/a**2)*np.sum((t/a)**b)
    hess_b = 10/(b**2) + np.sum((t/a)**b* (np.log(t/a))**2)
    hessians = np.array([hess_a, hess_b]) 
    return hessians  

iters = 0     
a0, b0 = 2,6 #initial condition

while iters < 4:  
    if hessian(a0,b0).any() == 0.0:
        print('Divide by zero error!') 
    else:
        a = a0 - gradient(a0,b0)[0]/hessian(a0,b0)[0]
        b = b0 - gradient(a0,b0)[1]/hessian(a0,b0)[1]    
        print('Iteration-%d, a = %0.6f, b= %0.6f, e1 = %0.6f, e2 = %0.6f' % (iters, a,b,a-a0,b-b0))
    if gradient(a,b)[0]==0 and  gradient(a,b)[1]==0: 
        break
    else:
        a0,b0 = a,b
        iters = iters +1
print(a,b)
print(iters)

Output:
Iteration-0, a = 1.714391, b= 6.656262, e1 = -0.285609, e2 = 0.656262
Iteration-1, a = 1.490482, b= 7.249362, e1 = -0.223909, e2 = 0.593099
Iteration-2, a = 1.309805, b= 7.794966, e1 = -0.180677, e2 = 0.545604
Iteration-3, a = 1.160878, b= 8.303167, e1 = -0.148927, e2 = 0.508201

a goes to zero, b to infinity, neither of which is a correct estimation. It feels like I am doing something wrong in defining gradients and hessians, though these seem to be correct first and second derivatives of the log-likelihood function.

$\endgroup$
2
  • $\begingroup$ I'm curious what the factors of 10 in your gradient and Hessian are for. $\endgroup$ Mar 24 at 23:56
  • $\begingroup$ @Aruralreader, this is just the size of t (number of data points), because I am taking gradient and hessian from the log-likelihood function $\endgroup$
    – kittycat
    Mar 25 at 0:07

Your Answer

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

Browse other questions tagged or ask your own question.