optimal way of defining a numerically stable sigmoid function for a list in python

You are right, you can do better by using np.where, the numpy equivalent of if:

def sigmoid(x):
    return np.where(x >= 0, 
                    1 / (1 + np.exp(-x)), 
                    np.exp(x) / (1 + np.exp(x)))

This function takes a numpy array x and returns a numpy array, too:

data = np.arange(-5,5)
sigmoid(data)
#array([0.00669285, 0.01798621, 0.04742587, 0.11920292, 0.26894142,
#       0.5       , 0.73105858, 0.88079708, 0.95257413, 0.98201379])

def sigmoid(x):
    """
    A numerically stable version of the logistic sigmoid function.
    """
    pos_mask = (x >= 0)
    neg_mask = (x < 0)
    z = np.zeros_like(x)
    z[pos_mask] = np.exp(-x[pos_mask])
    z[neg_mask] = np.exp(x[neg_mask])
    top = np.ones_like(x)
    top[neg_mask] = z[neg_mask]
    return top / (1 + z)

This piece of code comes from assignment3 of cs231n, I don't really understand why we should calculate it in this way, but I know this may be the code that you are looking for. Hope to be helpful.