import numpy as np
import plotnine as pn
import pandas as pd
import statsmodels.api as sm
from sklearn.metrics import confusion_matrixLogistic Regression
# carnivora
df = pd.read_csv("https://raw.githubusercontent.com/roualdes/data/refs/heads/master/carnivora.csv")
df = df.dropna(subset = ["SuperFamily", "LY"])df["dog"] = (df["SuperFamily"] == "Caniformia").astype(np.int64)fit = sm.GLM.from_formula("dog ~ LY",
data = df,
family=sm.families.Binomial() ).fit()fit.summary()| Dep. Variable: | dog | No. Observations: | 49 |
| Model: | GLM | Df Residuals: | 47 |
| Model Family: | Binomial | Df Model: | 1 |
| Link Function: | Logit | Scale: | 1.0000 |
| Method: | IRLS | Log-Likelihood: | -33.400 |
| Date: | Tue, 28 Apr 2026 | Deviance: | 66.799 |
| Time: | 13:11:54 | Pearson chi2: | 48.8 |
| No. Iterations: | 4 | Pseudo R-squ. (CS): | 0.02237 |
| Covariance Type: | nonrobust |
| coef | std err | z | P>|z| | [0.025 | 0.975] | |
| Intercept | -0.8570 | 0.842 | -1.017 | 0.309 | -2.508 | 0.794 |
| LY | 0.0045 | 0.004 | 1.025 | 0.306 | -0.004 | 0.013 |
x = np.linspace(np.min(df["LY"]), np.max(df["LY"]), 101)
ndf = pd.DataFrame({"LY": x})
ndf["phat"] = fit.predict(ndf)pn.ggplot() + \
pn.geom_point(df, pn.aes("LY", "dog")) + \
pn.geom_line(ndf, pn.aes("LY", "phat"))df["p_dog"] = fit.predict() > 0.5#np.mean((df["dog"] == 1) & (df["p_dog"] == 1)) # True Positive
#np.mean((df["dog"] == 1) & (df["p_dog"] == 0)) # False Negative
np.mean((df["dog"] == 0) & (df["p_dog"] == 1)) # False Positive
#np.mean((df["dog"] == 0) & (df["p_dog"] == 0)) # True Negative0.14285714285714285
confusion_matrix(df["dog"], df["p_dog"], normalize = "all")array([[0.36734694, 0.14285714],
[0.34693878, 0.14285714]])
| Predicted 0 | Predicted 1 | |
|---|---|---|
| Actual 0 | TN 😀 | FP 🙁 |
| Actual 1 | FN 🙁 | TP 😀 |
df = pd.read_csv("http://roualdes.sfo3.digitaloceanspaces.com/data/maize.csv")
df = df.dropna(subset = ["yield", "plantheight"])df["lots"] = (df["yield"] > np.median(df["yield"])).astype(np.int64)pn.ggplot() + \
pn.geom_point(df, pn.aes("plantheight", "lots"))fit = sm.GLM.from_formula("lots ~ plantheight",
data = df,
family=sm.families.Binomial() ).fit()df["p_lots"] = fit.predict()pn.ggplot() + \
pn.geom_point(df, pn.aes("plantheight", "lots")) + \
pn.geom_line(df, pn.aes("plantheight", "p_lots"))def bootstrap(arr, T, R = 1_000):
N = np.shape(arr)[0]
Ts = np.zeros(R)
rng = np.random.default_rng()
for r in range(R):
idx = rng.integers(N, size = N)
if type(arr) is np.ndarray:
Ts[r] = T(arr[idx])
else:
Ts[r] = T(arr.iloc[idx])
return Tsm = np.mean(df["plantheight"])
def logistic_ci(_data):
fit = sm.GLM.from_formula("lots ~ plantheight",
data = _data,
family=sm.families.Binomial() ).fit()
ndf = pd.DataFrame({"plantheight": [m, m + 1]})
return np.diff(fit.predict(ndf))[0]b = bootstrap(df, logistic_ci)np.quantile(b, [0.025, 0.975])array([0.00948544, 0.01030603])