Hi, thanks for your feedback. You have probably solved the problem you were facing already as I am way too late in my reply! For what it's worth, I have corrected an error in the aux-equation for estimating alpha in my article. If you use the corrected version of the equation, you might get a better value for alpha.

Alternately, you could try using the builtin support for the NB model in statsmodels. Statsmodels will estimate the alpha for you as part of the model fitting process. See sample code below:

import statsmodels.discrete.discrete_model as dm

#build the model
nb_model = dm.NegativeBinomial(y_train, X_train, loglike_method='nb2')

#print out the estimated alpha. The second param printed in the output is the alpha estimated by statsmodels:
print(nb_model._get_start_params_null())

#fit the model
nb2_training_results = nb_model.fit(maxiter=100)

'best