1. How to deal with 0's when taking np.log? 2. scipy,optimize.fmin_tnc gives shape error even after taking transpose

Viewed 296

I have two questions concerning this python 3 implementation of logistic regression (Andrew Ng course):

  1. When I take alpha=0.01, I get two errors:

    a. zero value encountered while taking log

    b. matrix multiplication error

I knew the sigmoid function would only return values between (0,1) but on printing the hypothesis while running gradient descent, I realized some values were rounded off to 1 (making 1- hyp=0 and hence the error). So, I thought increasing the precision of theta to np.float128 would help but it didn't!

However, taking alpha as 0.001 doesn't give any errors but I have to increase the iterations to 1000000 to decrease the cost from 0.693 to 0.224.

  1. I also tried to use scipy's optimizer to get the optimum value of theta. However, it gives the error that I've attached with the code. Even on passing theta.T, i get the same error.
def sigmoid(z):
    return 1/(1+np.exp(-z)) 
data_set.insert(0,'Ones',1)
X= data_set.iloc[:,0:3]
Y=data_set.iloc[:,3]

#convert X and Y to numpy matrices
X= np.matrix(X.values)
Y= np.matrix(Y.values)
#intilize theta
theta= np.float128(np.zeros([1,3]))
theta= np.matrix(theta)
Y= Y.T
#now let's define our cost functio
def costfunction(theta,X,Y):
    m=len(Y)
    hypothesis= sigmoid(np.dot(X,theta.T))
    error= (np.multiply(-Y,np.log(hypothesis)) - np.multiply((1-Y),np.log(1-hypothesis)))
    return 1/m * np.sum(error)
#let's define our gradient descent function now
def gradientdescent(X,Y,theta,alpha,iters):
    parameters=3
    temp= np.matrix(np.zeros(theta.shape))
    cost= np.zeros(iters)
    m= len(Y)
    
    for i in range(iters):
        error= sigmoid(X*theta.T) - Y
        for j in range(parameters):
            term= np.multiply(error,X[:,j])
            temp[0,j]= theta[0,j] - ((alpha/m) * np.sum(term))
    
        theta=temp
        cost[i]= costfunction(theta,X,Y)
    
    return theta, cost
alpha=0.001
iters=1000000
param,cost= gradientdescent(X,Y,theta,alpha,iters)
#We can also the optimum values for theta using scipy's optimize function
#so, let's define a gradient function now
def gradient(theta,X,Y):
    parameters=3
    grad= np.zeros(parameters)
    m=len(Y)
    
    for i in range(parameters):
        error= sigmoid((X*theta.T)) -Y
        term= np.multiply(error,X[:,i])
        grad[i]= np.sum(term)/m
        
    
    return grad
#now let's use scipy
import scipy.optimize as opt
result= opt.fmin_tnc(func=costfunction,x0=theta, fprime= gradient, args=(X,Y))
costfunction(result[0],X,Y)
---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-226-ac3f2f801635> in <module>
      1 #now let's use scipy
      2 import scipy.optimize as opt
----> 3 result= opt.fmin_tnc(func=costfunction,x0=theta, fprime= gradient, args=(X,Y))
      4 costfunction(result[0],X,Y)

~/anaconda3/lib/python3.7/site-packages/scipy/optimize/tnc.py in fmin_tnc(func, x0, fprime, args, approx_grad, bounds, epsilon, scale, offset, messages, maxCGit, maxfun, eta, stepmx, accuracy, fmin, ftol, xtol, pgtol, rescale, disp, callback)
    273             'disp': False}
    274 
--> 275     res = _minimize_tnc(fun, x0, args, jac, bounds, callback=callback, **opts)
    276 
    277     return res['x'], res['nfev'], res['status']

~/anaconda3/lib/python3.7/site-packages/scipy/optimize/tnc.py in _minimize_tnc(fun, x0, args, jac, bounds, eps, scale, offset, mesg_num, maxCGit, maxiter, eta, stepmx, accuracy, minfev, ftol, xtol, gtol, rescale, disp, callback, **unknown_options)
    407                                         offset, messages, maxCGit, maxfun,
    408                                         eta, stepmx, accuracy, fmin, ftol,
--> 409                                         xtol, pgtol, rescale, callback)
    410 
    411     funv, jacv = func_and_grad(x)

~/anaconda3/lib/python3.7/site-packages/scipy/optimize/tnc.py in func_and_grad(x)
    370         def func_and_grad(x):
    371             f = fun(x, *args)
--> 372             g = jac(x, *args)
    373             return f, g
    374 

<ipython-input-225-ad5800c8116a> in gradient(theta, X, Y)
      7 
      8     for i in range(parameters):
----> 9         error= sigmoid((X*theta.T)) -Y
     10         term= np.multiply(error,X[:,i])
     11         grad[i]= np.sum(term)/m

~/anaconda3/lib/python3.7/site-packages/numpy/matrixlib/defmatrix.py in __mul__(self, other)
    218         if isinstance(other, (N.ndarray, list, tuple)) :
    219             # This promotes 1-D vectors to row vectors
--> 220             return N.dot(self, asmatrix(other))
    221         if isscalar(other) or not hasattr(other, '__rmul__') :
    222             return N.dot(self, other)

<__array_function__ internals> in dot(*args, **kwargs)

ValueError: shapes (100,3) and (1,3) not aligned: 3 (dim 1) != 1 (dim 0)
1 Answers

I'm no expert in scipy, however, if you want to make you sigmoid function never return exactly 0 or 1, you can use the numpy minimum and maximum:

def sigmoid(z):
    sig = 1 / (1 + np.exp(-z))     # Define sigmoid function
    sig = np.minimum(sig, 0.9999)  # Set upper bound
    sig = np.maximum(sig, 0.0001)  # Set lower bound
    return sig

However, your real problem isn't the rounding in the cost calculation (even Octave/MATLAB return nan for some of the thetas your code produces). Your real problem is that your gradient descent implementation is diverging unless the learning rate is very small. Using gradient descent instead of more advanced optimization algorithms (like 'fminunc' in Octave/MATLAB) forces you to choose small learning rates and do many iterations. It might help to do some sort of feature normalization/standardization if you have not already.

Related