Linear Decision Boundary in Logistic Regression

Which line(s) of code are responsible for the additional linear decision boundaries, 15%, 30%, etc., in the following code?

import matplotlib.pyplot as plt
import numpy as np
from sklearn import datasets
from sklearn.linear_model import LogisticRegression

iris = datasets.load_iris()
print(list(iris.keys()))

X = iris['data'][:,3:] #petal width
y = (iris['target']==2).astype(np.int32) #1 if virginica, else 0

X = iris["data"][:, (2, 3)]  # petal length, petal width
y = (iris["target"] == 2).astype(np.int32)

log_reg = LogisticRegression(solver="lbfgs", C=10**10, random_state=42)
log_reg.fit(X, y)

x0, x1 = np.meshgrid(
        np.linspace(2.9, 7, 500).reshape(-1, 1),
        np.linspace(0.8, 2.7, 200).reshape(-1, 1),
    )

X_new = np.c_[x0.ravel(), x1.ravel()]

y_proba = log_reg.predict_proba(X_new)

plt.figure(figsize=(10, 4))
plt.plot(X[y==0, 0], X[y==0, 1], "bs")
plt.plot(X[y==1, 0], X[y==1, 1], "g^")

zz = y_proba[:, 1].reshape(x0.shape)
contour = plt.contour(x0, x1, zz, cmap=plt.cm.brg)

left_right = np.array([2.9, 7])
boundary = -(log_reg.coef_[0][0] * left_right + log_reg.intercept_[0]) / log_reg.coef_[0][1]

plt.clabel(contour, inline=1, fontsize=12)
plt.plot(left_right, boundary, "k--", linewidth=3)
plt.text(3.5, 1.5, "Not Iris virginica", fontsize=14, color="b", ha="center")
plt.text(6.5, 2.3, "Iris virginica", fontsize=14, color="g", ha="center")
plt.xlabel("Petal length", fontsize=14)
plt.ylabel("Petal width", fontsize=14)
plt.axis([2.9, 7, 0.8, 2.7])

plt.show()

I'm working though a ML textbook with undocumented code and this is the one feature I can't figure out the origin of.

It is from this line

contour = plt.contour(x0, x1, zz, cmap=plt.cm.brg)

You can double-check on this by commenting out that line from your code.

Then, the labels (30%, 60%,...) are drawn using contour object above:

plt.clabel(contour, inline=1, fontsize=12)

Confirmed by looking into contour object:

> {'axes': <Axes: xlabel='Petal length', ylabel='Petal width'>,
'levels': array([0. , 0.15, 0.3 , 0.45, 0.6 , 0.75, 0.9 , 1.05]),...

Basically, the idea is as follows.
First, create x0 and x1 using Numpy meshgrid, then plot a contour with the probability (i.e., zz).

For further reading on meshgrid, you can refer to this link: https://www.sharpsightlabs.com/blog/numpy-meshgrid/