补2 · 补2 · KNN、机器学习词汇与 K-Means

# 下载数据与安装依赖
import os
import requests

os.makedirs("data", exist_ok=True)
data_files = {
    "data/penguins.csv":
        "https://raw.githubusercontent.com/mwaskom/seaborn-data/master/penguins.csv",
}
for path, url in data_files.items():
    if not os.path.exists(path):
        r = requests.get(url)
        r.raise_for_status()
        with open(path, "wb") as f:
            f.write(r.content)
        print(f"Downloaded {path}")
    else:
        print(f"Found {path}")

import numpy as np
import pandas as pd
import plotly.express as px

px.defaults.width = 800
pd.set_option("plotting.backend", "plotly")

os.makedirs("images", exist_ok=True)

penguins = pd.read_csv("data/penguins.csv")
penguins.head()

penguins.info()

# describe() 对数值列做汇总统计
penguins.describe()

print("shape:", penguins.shape)

penguins["species"].value_counts()

# iloc：按位置。第 0-4 行，前三列。
penguins.iloc[0:5, 0:3]

# loc：按标签。第 0-4 行（含 4！），指定列名。
penguins.loc[0:4, ["species", "island", "bill_length_mm"]]

# 单列是 Series；用列表取多列得到的是 DataFrame
print(type(penguins["bill_length_mm"]))
print(type(penguins[["bill_length_mm"]]))

penguins.isna().sum()

# 布尔过滤：哪些企鹅既重、鳍又长？
mask = (penguins["body_mass_g"] > 5000) & (penguins["flipper_length_mm"] > 220)
penguins[mask].head()

df = penguins.dropna()
print(len(penguins), "rows ->", len(df), "rows after dropping missing values")
df.head()

df.groupby("species")[["bill_length_mm", "flipper_length_mm", "body_mass_g"]].mean().round(1)

FEATURES = ["bill_length_mm", "flipper_length_mm"]
X = df[FEATURES].to_numpy()      # 特征矩阵，形状 (n_samples, n_features)
y = df["species"].to_numpy()     # 标签，形状 (n_samples,)

print("X shape:", X.shape, "| y shape:", y.shape)
print("X dtype:", X.dtype, "\n")
print("The first three rows of X:\n", X[:3])

# 向量化运算一次作用于整个数组，不需要写 Python 循环。
print("column means:", X.mean(axis=0).round(2))
print("column stds: ", X.std(axis=0).round(2))

fig = px.scatter(
    df, x="bill_length_mm", y="flipper_length_mm", color="species",
    title="Palmer Penguins: bill length vs flipper length",
    labels={
        "bill_length_mm": "Bill length (mm)",
        "flipper_length_mm": "Flipper length (mm)",
        "species": "Species",
    },
    height=520,
)
# Plotly 默认图例在右侧，在 notebook 里经常被裁掉，所以要手动挪位置。
fig.update_layout(
    showlegend=True,
    legend=dict(
        title_text="Species",
        x=0.01, y=0.99, xanchor="left", yanchor="top",
        bgcolor="white", bordercolor="lightgray", borderwidth=1,
    ),
)
fig.show()

from sklearn.neighbors import KNeighborsClassifier

knn1 = KNeighborsClassifier(n_neighbors=1)
knn1.fit(X, y)          # “训练”：把数据存下来
print("Model is trained.")

# 推断：预测一只新企鹅的物种
new_penguin = np.array([[45.0, 200.0]])   # 喙长 45mm，鳍长 200mm
print("Predicted species:", knn1.predict(new_penguin)[0])

flowchart LR
    A["新样本 x"] --> B["计算 x 到所有训练点的距离"]
    B --> C["取距离最近的 k 个点"]
    C --> D{"分类还是回归？"}
    D -- "分类" --> E["k 个邻居投票（多数票）"]
    D -- "回归" --> F["k 个邻居的目标值取平均"]

accuracy = knn1.score(X, y)
print(f"Accuracy on our data: {accuracy:.1%}")

# 为什么是 100%？看看离第一只企鹅最近的点是哪一个。
distances, indices = knn1.kneighbors(X[:1], n_neighbors=1)
print("Nearest neighbour of penguin 0 is penguin index:", indices[0][0])
print("at a distance of:", distances[0][0])

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(
    X, y, test_size=0.2, random_state=42, stratify=y
)
print(f"Training set: {X_train.shape[0]} penguins")
print(f"Test set:     {X_test.shape[0]} penguins")

knn1 = KNeighborsClassifier(n_neighbors=1).fit(X_train, y_train)

print(f"Training accuracy: {knn1.score(X_train, y_train):.1%}")
print(f"Test accuracy:     {knn1.score(X_test, y_test):.1%}")

flowchart TD
    D["全部数据（344 只企鹅）"] --> A["训练集（80%）"]
    D --> B["测试集（20%）"]
    A --> A1["训练集（训练用）"]
    A --> A2["验证集（调超参数用）"]
    B --> B1["只在最后用一次<br/>衡量泛化"]

import plotly.graph_objects as go
from sklearn.preprocessing import LabelEncoder

def decision_boundary_figure(k, X_tr, y_tr, resolution=250):
    le = LabelEncoder().fit(y_tr)
    model = KNeighborsClassifier(n_neighbors=k).fit(X_tr, le.transform(y_tr))

    # 在特征空间上打一张细密的网格，对每个格点做预测
    pad = 1.0
    xs = np.linspace(X_tr[:, 0].min() - pad, X_tr[:, 0].max() + pad, resolution)
    ys = np.linspace(X_tr[:, 1].min() - pad, X_tr[:, 1].max() + pad, resolution)
    xx, yy = np.meshgrid(xs, ys)
    zz = model.predict(np.c_[xx.ravel(), yy.ravel()]).reshape(xx.shape)

    fig = go.Figure()
    # 用热力图铺出类别区域
    fig.add_trace(go.Heatmap(x=xs, y=ys, z=zz, showscale=False, opacity=0.30,
                             colorscale=[[0, "#002675"], [0.5, "#FDB515"], [1, "#028842"]]))
    # 再把训练点叠上去
    for i, cls in enumerate(le.classes_):
        m = y_tr == cls
        fig.add_trace(go.Scatter(x=X_tr[m, 0], y=X_tr[m, 1], mode="markers", name=cls))
    fig.update_layout(title=f"KNN decision boundary, k = {k}",
                      xaxis_title="Bill length (mm)", yaxis_title="Flipper length (mm)",
                      width=800, height=520)
    return fig

decision_boundary_figure(1, X_train, y_train).show()

decision_boundary_figure(15, X_train, y_train).show()

decision_boundary_figure(100, X_train, y_train).show()

X_tr, X_val, y_tr, y_val = train_test_split(
    X_train, y_train, test_size=0.25, random_state=42, stratify=y_train
)
print(f"train {len(X_tr)} | validation {len(X_val)} | test {len(X_test)}")

ks = list(range(1, 101))
train_acc, val_acc = [], []
for k in ks:
    m = KNeighborsClassifier(n_neighbors=k).fit(X_tr, y_tr)
    train_acc.append(m.score(X_tr, y_tr))
    val_acc.append(m.score(X_val, y_val))

best_k = ks[int(np.argmax(val_acc))]
print(f"Best k on the validation set: {best_k}  (validation accuracy {max(val_acc):.1%})")

fig = go.Figure()
fig.add_trace(go.Scatter(x=ks, y=train_acc, name="Training accuracy"))
fig.add_trace(go.Scatter(x=ks, y=val_acc, name="Validation accuracy"))
fig.add_vline(x=best_k, line_dash="dash",
              annotation_text=f"best k = {best_k}")
fig.update_layout(title="Training and validation accuracy as k increases",
                  xaxis_title="k (number of neighbours)", yaxis_title="Accuracy",
                  width=800, height=520)
fig.show()

# 只有现在，在选好 k 之后，我们才碰测试集。而且只碰一次。
final = KNeighborsClassifier(n_neighbors=best_k).fit(X_train, y_train)
print(f"Final test accuracy with k={best_k}: {final.score(X_test, y_test):.1%}")

from sklearn.pipeline import make_pipeline
from sklearn.preprocessing import StandardScaler

scaled = make_pipeline(StandardScaler(), KNeighborsClassifier(n_neighbors=best_k))
scaled.fit(X_train, y_train)

print(f"Unscaled test accuracy: {final.score(X_test, y_test):.1%}")
print(f"Scaled test accuracy:   {scaled.score(X_test, y_test):.1%}")

from sklearn.neighbors import KNeighborsRegressor

Xr = df[["flipper_length_mm"]].to_numpy()
yr = df["body_mass_g"].to_numpy()

Xr_train, Xr_test, yr_train, yr_test = train_test_split(Xr, yr, test_size=0.2, random_state=42)
print("Predicting body mass (g) from flipper length (mm)")

grid = np.linspace(Xr.min(), Xr.max(), 400).reshape(-1, 1)

fig = go.Figure()
fig.add_trace(go.Scatter(x=Xr_train.ravel(), y=yr_train, mode="markers", name="Training data"))
for k in [1, 25, 200]:
    reg = KNeighborsRegressor(n_neighbors=k).fit(Xr_train, yr_train)
    fig.add_trace(go.Scatter(x=grid.ravel(), y=reg.predict(grid), mode="lines", name=f"k = {k}"))
fig.update_layout(title="KNN regression: body mass vs flipper length",
                  xaxis_title="Flipper length (mm)", yaxis_title="Body mass (g)",
                  width=800, height=520)
fig.show()

for k in [1, 25, 200]:
    reg = KNeighborsRegressor(n_neighbors=k).fit(Xr_train, yr_train)
    print(f"k={k:>3}  train R^2 = {reg.score(Xr_train, yr_train):.3f}   "
          f"test R^2 = {reg.score(Xr_test, yr_test):.3f}")

# 把标签扔掉。
fig = px.scatter(df, x="bill_length_mm", y="flipper_length_mm",
                 title="The same penguins, with no labels",
                 labels={"bill_length_mm": "Bill length (mm)",
                         "flipper_length_mm": "Flipper length (mm)"},
                 height=520)
fig.update_traces(marker=dict(color="gray", size=7))
fig.show()

from sklearn.cluster import KMeans

Xs = StandardScaler().fit_transform(X)   # 又是基于距离的方法，所以要先标准化

kmeans = KMeans(n_clusters=3, n_init=10, random_state=42)
cluster = kmeans.fit_predict(Xs)

df_c = df.copy()
df_c["cluster"] = cluster.astype(str)
px.scatter(df_c, x="bill_length_mm", y="flipper_length_mm", color="cluster",
           title="K-means with K = 3",
           labels={"bill_length_mm": "Bill length (mm)",
                   "flipper_length_mm": "Flipper length (mm)"},
           height=520).show()

def lloyd_steps(Xs, K=3, n_iter=6, seed=1):
    rng = np.random.default_rng(seed)
    centres = Xs[rng.choice(len(Xs), K, replace=False)].copy()
    history = []
    for _ in range(n_iter):
        # 用广播一次算出所有点到所有质心的平方距离，形状 (n_points, K)
        d = ((Xs[:, None, :] - centres[None, :, :]) ** 2).sum(axis=2)
        assign = d.argmin(axis=1)                       # 分配步
        history.append((centres.copy(), assign.copy()))
        for j in range(K):                              # 更新步
            if (assign == j).any():
                centres[j] = Xs[assign == j].mean(axis=0)
    return history

history = lloyd_steps(Xs)
print(f"{len(history)} iterations recorded")

# 把每一轮迭代做成动画的一帧
fig = go.Figure()
frames = []
for step, (centres, assign) in enumerate(history):
    data = []
    for j in range(3):
        m = assign == j
        data.append(go.Scatter(x=Xs[m, 0], y=Xs[m, 1], mode="markers",
                               name=f"cluster {j}"))
    # 质心用黑色叉号标出
    data.append(go.Scatter(x=centres[:, 0], y=centres[:, 1], mode="markers",
                           name="centroids"))
    frames.append(go.Frame(data=data, name=str(step)))

fig.add_traces(frames[0].data)
fig.frames = frames
fig.update_layout(title="Lloyd's algorithm: assign, then update",
                  xaxis_title="Bill length (standardized)",
                  yaxis_title="Flipper length (standardized)",
                  width=800, height=560)
fig.show()

flowchart TD
    A["初始化：随机挑 K 个点当质心"] --> B["分配步<br/>每个点归到离它最近的质心"]
    B --> C["更新步<br/>每个质心移到本簇所有点的均值"]
    C --> D{"质心还在移动吗？"}
    D -- "是（未收敛）" --> B
    D -- "否（收敛或达到迭代上限）" --> E["输出每个点的簇标签"]

Ks = range(1, 11)
inertia = [KMeans(n_clusters=k, n_init=10, random_state=42).fit(Xs).inertia_ for k in Ks]

fig = px.line(x=list(Ks), y=inertia, markers=True,
              title="Elbow method: inertia vs K",
              labels={"x": "K (number of clusters)", "y": "Within-cluster sum of squares"})
fig.update_layout(width=800, height=480)
fig.show()

comparison = pd.crosstab(df_c["cluster"], df_c["species"])
comparison

fig = px.scatter(df_c, x="bill_length_mm", y="flipper_length_mm",
                 color="species", symbol="cluster",
                 title="True species (colour) vs discovered clusters (symbol)",
                 labels={"bill_length_mm": "Bill length (mm)",
                         "flipper_length_mm": "Flipper length (mm)"},
                 height=560)
fig.show()

s = pd.Series([3750, 3800, 3250], index=["p0", "p1", "p2"], name="body_mass_g")
print(s)
print("\nindex:", s.index.tolist())
print("values:", s.values)

frame = pd.DataFrame({
    "species": ["Adelie", "Gentoo", "Chinstrap"],
    "bill_length_mm": [39.1, 46.1, 46.5],
    "island": ["Torgersen", "Biscoe", "Dream"],
})
frame

print(penguins.head(3))          # 前几行
print(penguins.tail(3))          # 后几行
print(penguins.sample(3))        # 随机几行
print(penguins.columns.tolist()) # 列名
print(penguins.dtypes)           # 各列类型
print(penguins["island"].unique())

print(penguins.iloc[0])              # 第 0 行，作为一个 Series
print(penguins.iloc[0:3])            # 第 0、1、2 行
print(penguins.iloc[:, 2])           # 位置为 2 的那一列
print(penguins.iloc[[0, 5, 10], [0, 2]])   # 任意行与任意列

print(penguins.loc[0:2])                          # 标签为 0、1、2 的行（含 2）
print(penguins.loc[:, "species":"bill_length_mm"]) # 按列名切片
print(penguins.loc[penguins["species"] == "Gentoo", "body_mass_g"].mean())

tmp = df.copy()
tmp["mass_kg"] = tmp["body_mass_g"] / 1000              # 新增一列
tmp["bill_ratio"] = tmp["bill_length_mm"] / tmp["bill_depth_mm"]
tmp = tmp.drop(columns=["bill_ratio"])                  # 删除一列
tmp.sort_values("body_mass_g", ascending=False).head()

print(df.groupby("species")["body_mass_g"].agg(["mean", "std", "count"]).round(1))
print()
print(df.groupby(["species", "island"])["flipper_length_mm"].mean().round(1))
print()
print(df.pivot_table(index="species", columns="island",
                     values="body_mass_g", aggfunc="mean").round(0))

islands = pd.DataFrame({
    "island": ["Torgersen", "Biscoe", "Dream"],
    "latitude": [-64.77, -65.43, -64.73],
})
merged = df.merge(islands, on="island", how="left")   # 也可以试 how="inner" 和 how="outer"
merged[["species", "island", "latitude"]].head()

a = np.arange(12).reshape(3, 4)
print(a)
print("shape:", a.shape, "| ndim:", a.ndim)
print("row sums:", a.sum(axis=1))
print("col means:", a.mean(axis=0))
print("boolean mask:", a[a > 6])
print("broadcasting:", (a - a.mean(axis=0)).round(2))

# 手算欧氏距离，这正是 KNN 计算的东西
p, q = X[0], X[1]
print("manual :", np.sqrt(((p - q) ** 2).sum()))
print("numpy  :", np.linalg.norm(p - q))

import matplotlib.pyplot as plt
import seaborn as sns

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
axes[0].hist(df["body_mass_g"], bins=25, color="#002675")
axes[0].set_title("Body mass"); axes[0].set_xlabel("g")
sns.scatterplot(data=df, x="bill_length_mm", y="flipper_length_mm",
                hue="species", ax=axes[1])
axes[1].set_title("By species")
plt.tight_layout(); plt.show()

# Plotly：直方图、箱线图，以及按岛屿分面的散点图
px.histogram(df, x="body_mass_g", color="species", nbins=30, height=420).show()
px.box(df, x="species", y="flipper_length_mm", color="species", height=420).show()
px.scatter(df, x="bill_length_mm", y="flipper_length_mm",
           color="species", facet_col="island", height=420).show()