How do I manually `predict_proba` from logistic regression model in scikit-learn?-CodePudding

I am trying to manually predict a logistic regression model using the coefficient and intercept outputs from a scikit-learn model. However, I can't match up my probability predictions with the predict_proba method from the classifier.

I have tried:

from sklearn.datasets import load_iris
from scipy.special import expit
import numpy as np

X, y = load_iris(return_X_y=True)
clf = LogisticRegression(random_state=0).fit(X, y)
# use sklearn's predict_proba function
sk_probas = clf.predict_proba(X[:1, :])
# and attempting manually (using scipy's inverse logit)
manual_probas = expit(np.dot(X[:1], clf.coef_.T) clf.intercept_)
# with a completely manual inverse logit
full_manual_probas = 1/(1 np.exp(-(np.dot(iris_test, iris_coef.T) clf.intercept_)))

outputs:

>>> sk_probas
array([[9.81815067e-01, 1.81849190e-02, 1.44120963e-08]])
>>> manual_probas
array([[9.99352591e-01, 9.66205386e-01, 2.26583306e-05]])
>>> full_manual_probas
array([[9.99352591e-01, 9.66205386e-01, 2.26583306e-05]])

I do seem to get the classes to match (using np.argmax), but the probabilities are different. What am I missing?

I've looked at this and this but haven't managed to figure it out yet.

CodePudding user response：

The documentation states that

For a multi_class problem, if multi_class is set to be “multinomial” the softmax function is used to find the predicted probability of each class

That is, in order to get the same values as sklearn you have to normalize using softmax, like this:

from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression
import numpy as np

X, y = load_iris(return_X_y=True)
clf = LogisticRegression(random_state=0, max_iter=1000).fit(X, y)
decision = np.dot(X[:1], clf.coef_.T) clf.intercept_
print(clf.predict_proba(X[:1]))
print(np.exp(decision) / np.exp(decision).sum())

To use sigmoids instead you can do it like this:

from sklearn.datasets import load_iris
from sklearn.linear_model import LogisticRegression
import numpy as np

X, y = load_iris(return_X_y=True)
clf = LogisticRegression(random_state=0, max_iter=1000, multi_class='ovr').fit(X, y) # Notice the extra argument
full_manual_probas = 1/(1 np.exp(-(np.dot(X[:1], clf.coef_.T) clf.intercept_)))
print(clf.predict_proba(X[:1]))
print(full_manual_probas / full_manual_probas.sum())