机器学习模型评估与改进

import numpy as np
%matplotlib widget
import matplotlib.pyplot as plt
from sklearn.linear_model import LinearRegression, Ridge
from sklearn.preprocessing import StandardScaler, PolynomialFeatures
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error
import tensorflow as tf
from tensorflow.keras.models import Sequential
from tensorflow.keras.layers import Dense
from tensorflow.keras.activations import relu,linear
from tensorflow.keras.losses import SparseCategoricalCrossentropy
from tensorflow.keras.optimizers import Adam

import logging
logging.getLogger("tensorflow").setLevel(logging.ERROR)

tf.keras.backend.set_floatx('float64')
from assignment_utils import *

tf.autograph.set_verbosity(0)

# 生成一些数据
X, y, x_ideal, y_ideal = gen_data(18, 2, 0.7)
print("X.shape", X.shape, "y.shape", y.shape)

# 使用 sklearn 例程分割数据
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.33, random_state=1)
print("X_train.shape", X_train.shape, "y_train.shape", y_train.shape)
print("X_test.shape", X_test.shape, "y_test.shape", y_test.shape)

fig, ax = plt.subplots(1,1,figsize=(4,4))
ax.plot(x_ideal, y_ideal, "--", color = "orangered", label="y_ideal", lw=1)
ax.set_title("Training, Test",fontsize = 14)
ax.set_xlabel("x")
ax.set_ylabel("y")

ax.scatter(X_train, y_train, color = "red",           label="train")
ax.scatter(X_test, y_test,   color = dlc["dlblue"],   label="test")
ax.legend(loc='upper left')
plt.show()

def eval_mse(y, yhat):
    """
    计算数据集上的均方误差。
    参数:
        y (ndarray): 形状 (m,) 或 (m,1)，每个样本的目标值
        yhat (ndarray): 形状 (m,) 或 (m,1)，每个样本的预测值
    返回:
        err (scalar): 均方误差
    """
    m = len(y)
    err = 0.0
    for i in range(m):
        err_i = (yhat[i] - y[i]) ** 2
        err += err_i
    err = err / (2 * m)
    return err

# 在 sklearn 中创建一个模型，并在训练数据上训练
degree = 10
lmodel = lin_model(degree)
lmodel.fit(X_train, y_train)

# 预测训练数据，计算训练误差
yhat = lmodel.predict(X_train)
err_train = lmodel.mse(y_train, yhat)

# 预测测试数据，计算误差
yhat = lmodel.predict(X_test)
err_test = lmodel.mse(y_test, yhat)

print(f"training err {err_train:0.2f}, test err {err_test:0.2f}")

# 绘制预测曲线
x = np.linspace(0, int(X.max()), 100)
y_pred = lmodel.predict(x).reshape(-1, 1)

plt_train_test(X_train, y_train, X_test, y_test, x, y_pred, x_ideal, y_ideal, degree)

# Generate  data
X,y, x_ideal,y_ideal = gen_data(40, 5, 0.7)
print("X.shape", X.shape, "y.shape", y.shape)

#split the data using sklearn routine
X_train, X_, y_train, y_ = train_test_split(X,y,test_size=0.40, random_state=1)
X_cv, X_test, y_cv, y_test = train_test_split(X_,y_,test_size=0.50, random_state=1)
print("X_train.shape", X_train.shape, "y_train.shape", y_train.shape)
print("X_cv.shape", X_cv.shape, "y_cv.shape", y_cv.shape)
print("X_test.shape", X_test.shape, "y_test.shape", y_test.shape)

fig, ax = plt.subplots(1, 1, figsize=(4, 4))
ax.plot(x_ideal, y_ideal, "--", color="orangered", label="y_ideal", lw=1)
ax.set_title("Training, CV, Test", fontsize=14)
ax.set_xlabel("x")
ax.set_ylabel("y")

ax.scatter(X_train, y_train, color="red", label="train")
ax.scatter(X_cv, y_cv, color="orange", label="cv")
ax.scatter(X_test, y_test, color="blue", label="test")
ax.legend(loc='upper left')
plt.show()

max_degree = 9
err_train = np.zeros(max_degree)
err_cv = np.zeros(max_degree)
x = np.linspace(0, int(X.max()), 100)
y_pred = np.zeros((100, max_degree))

for degree in range(max_degree):
    lmodel = lin_model(degree + 1)
    lmodel.fit(X_train, y_train)
    yhat = lmodel.predict(X_train)
    err_train[degree] = lmodel.mse(y_train, yhat)
    yhat = lmodel.predict(X_cv)
    err_cv[degree] = lmodel.mse(y_cv, yhat)
    y_pred[:, degree] = lmodel.predict(x)

optimal_degree = np.argmin(err_cv) + 1

plt.close("all")
plt_optimal_degree(X_train, y_train, X_cv, y_cv, x, y_pred, x_ideal, y_ideal,
                   err_train, err_cv, optimal_degree, max_degree)

lambda_range = np.array([0.0, 1e-6, 1e-5, 1e-4, 1e-3, 1e-2, 1e-1, 1, 10, 100])
num_steps = len(lambda_range)
degree = 10
err_train = np.zeros(num_steps)
err_cv = np.zeros(num_steps)
x = np.linspace(0, int(X.max()), 100)
y_pred = np.zeros((100, num_steps))

for i in range(num_steps):
    lambda_ = lambda_range[i]
    lmodel = lin_model(degree, regularization=True, lambda_=lambda_)
    lmodel.fit(X_train, y_train)
    yhat = lmodel.predict(X_train)
    err_train[i] = lmodel.mse(y_train, yhat)
    yhat = lmodel.predict(X_cv)
    err_cv[i] = lmodel.mse(y_cv, yhat)
    y_pred[:, i] = lmodel.predict(x)

optimal_reg_idx = np.argmin(err_cv)

plt.close("all")
plt_tune_regularization(X_train, y_train, X_cv, y_cv, x, y_pred, err_train, err_cv, optimal_reg_idx, lambda_range)

X_train, y_train, X_cv, y_cv, x, y_pred, err_train, err_cv, m_range, degree = tune_m()
plt_tune_m(X_train, y_train, X_cv, y_cv, x, y_pred, err_train, err_cv, m_range, degree)

# 生成和分割数据集
X, y, centers, classes, std = gen_blobs()

# 分割数据。较大的 CV 群体以供演示
X_train, X_, y_train, y_ = train_test_split(X, y, test_size=0.50, random_state=1)
X_cv, X_test, y_cv, y_test = train_test_split(X_, y_, test_size=0.20, random_state=1)
print("X_train.shape:", X_train.shape, "X_cv.shape:", X_cv.shape, "X_test.shape:", X_test.shape)

plt_train_eq_dist(X_train, y_train, classes, X_cv, y_cv, centers, std)

def eval_cat_err(y, yhat):
    """
    计算分类误差。
    参数:
        y (ndarray): 形状 (m,) 或 (m,1)，每个样本的目标值
        yhat (ndarray): 形状 (m,) 或 (m,1)，每个样本的预测值
    返回:
        cerr (scalar): 分类误差
    """
    m = len(y)
    incorrect = 0
    for i in range(m):
        if yhat[i] != y[i]:
            incorrect += 1
    cerr = incorrect / m
    return cerr

y_hat = np.array([1, 2, 0])
y_tmp = np.array([1, 2, 3])
print(f"categorization error {np.squeeze(eval_cat_err(y_hat, y_tmp)):0.3f}, expected:0.333" )
y_hat = np.array([[1], [2], [0], [3]])
y_tmp = np.array([[1], [2], [1], [3]])
print(f"categorization error {np.squeeze(eval_cat_err(y_hat, y_tmp)):0.3f}, expected:0.250" )

test_eval_cat_err(eval_cat_err)

tf.random.set_seed(1234)
model = Sequential(
    [
        Dense(120, activation='relu', name="L1"),
        Dense(40, activation='relu', name="L2"),
        Dense(classes, activation='linear', name="L3")
    ], name="Complex"
)
model.compile(
    loss=tf.keras.losses.SparseCategoricalCrossentropy(from_logits=True),
    optimizer=tf.keras.optimizers.Adam(0.01),
)

model.fit(X_train, y_train, epochs=1000)

model.summary()
model_test(model, classes, X_train.shape[1])

Model: "Complex"
_________________________________________________________________
 Layer (type)                Output Shape              Param #
=================================================================
 L1 (Dense)                  (None, 120)               360

 L2 (Dense)                  (None, 40)                4840

 L3 (Dense)                  (None, 6)                 246

=================================================================
Total params: 5,446
Trainable params: 5,446
Non-trainable params: 0
_________________________________________________________________

training_cerr_complex = eval_cat_err(y_train, model_predict(X_train))
cv_cerr_complex = eval_cat_err(y_cv, model_predict(X_cv))
print(f"分类误差，训练，复杂模型: {training_cerr_complex:0.3f}")
print(f"分类误差，CV，复杂模型: {cv_cerr_complex:0.3f}")

tf.random.set_seed(1234)
model_s = Sequential(
    [
        Dense(6, activation='relu', name="L1"),
        Dense(classes, activation='linear', name="L2")
    ], name="Simple"
)
model_s.compile(
    loss=tf.keras.losses.SparseCategoricalCrossentropy(from_logits=True),
    optimizer=tf.keras.optimizers.Adam(0.01),
)

model_s.fit(X_train, y_train, epochs=1000)

model_s.summary()
model_s_test(model_s, classes, X_train.shape[1])

#make a model for plotting routines to call
model_predict_s = lambda Xl: np.argmax(tf.nn.softmax(model_s.predict(Xl)).numpy(),axis=1)
plt_nn(model_predict_s,X_train,y_train, classes, X_cv, y_cv, suptitle="Simple Model")

Model: "Simple"
_________________________________________________________________
 Layer (type)                Output Shape              Param #
=================================================================
 L1 (Dense)                  (None, 6)                 18

 L2 (Dense)                  (None, 6)                 42

=================================================================
Total params: 60
Trainable params: 60
Non-trainable params: 0
_________________________________________________________________

training_cerr_simple = eval_cat_err(y_train, model_predict_s(X_train))
cv_cerr_simple = eval_cat_err(y_cv, model_predict_s(X_cv))
print(f"分类误差，训练，简单模型: {training_cerr_simple:0.3f}, 复杂模型: {training_cerr_complex:0.3f}")
print(f"分类误差，CV，简单模型: {cv_cerr_simple:0.3f}, 复杂模型: {cv_cerr_complex:0.3f}")

tf.random.set_seed(1234)
model_r = Sequential(
    [
        Dense(120, activation='relu', kernel_regularizer=tf.keras.regularizers.l2(0.1), name="L1"),
        Dense(40, activation='relu', kernel_regularizer=tf.keras.regularizers.l2(0.1), name="L2"),
        Dense(classes, activation='linear', name="L3")
    ], name="ComplexRegularized"
)
model_r.compile(
    loss=tf.keras.losses.SparseCategoricalCrossentropy(from_logits=True),
    optimizer=tf.keras.optimizers.Adam(0.01),
)

model_r.fit(X_train, y_train, epochs=1000)

model_r.summary()
model_r_test(model_r, classes, X_train.shape[1])

Model: "sequential"
_________________________________________________________________
 Layer (type)                Output Shape              Param #
=================================================================
 L1 (Dense)                  (None, 120)               360

 L2 (Dense)                  (None, 40)                4840

 L3 (Dense)                  (None, 6)                 246

=================================================================
Total params: 5,446
Trainable params: 5,446
Non-trainable params: 0
_________________________________________________________________

#make a model for plotting routines to call
model_predict_r = lambda Xl: np.argmax(tf.nn.softmax(model_r.predict(Xl)).numpy(),axis=1)
plt_nn(model_predict_r, X_train,y_train, classes, X_cv, y_cv, suptitle="Regularized")

training_cerr_reg = eval_cat_err(y_train, model_predict_r(X_train))
cv_cerr_reg = eval_cat_err(y_cv, model_predict_r(X_cv))
print(f"分类误差，训练，正则化模型: {training_cerr_reg:0.3f}, 简单模型: {training_cerr_simple:0.3f}, 复杂模型: {training_cerr_complex:0.3f}")
print(f"分类误差，CV，正则化模型: {cv_cerr_reg:0.3f}, 简单模型: {cv_cerr_simple:0.3f}, 复杂模型: {cv_cerr_complex:0.3f}")

tf.random.set_seed(1234)
lambdas = [0.0, 0.001, 0.01, 0.05, 0.1, 0.2, 0.3]
models = [None] * len(lambdas)

for i in range(len(lambdas)):
    lambda_ = lambdas[i]
    models[i] = Sequential(
        [
            Dense(120, activation='relu', kernel_regularizer=tf.keras.regularizers.l2(lambda_)),
            Dense(40, activation='relu', kernel_regularizer=tf.keras.regularizers.l2(lambda_)),
            Dense(classes, activation='linear')
        ]
    )
    models[i].compile(
        loss=tf.keras.losses.SparseCategoricalCrossentropy(from_logits=True),
        optimizer=tf.keras.optimizers.Adam(0.01),
    )
    models[i].fit(X_train, y_train, epochs=1000)
    print(f"完成 lambda = {lambda_}")

plot_iterate(lambdas, models, X_train, y_train, X_cv, y_cv)

plt_compare(X_test, y_test, classes, model_predict_s, model_predict_r, centers)

def bias_variance_decomposition(model, X_train, y_train, X_test, y_test, n_iterations=100):
    """
    通过Bootstrap采样估计偏差和方差
    """
    predictions = []

    for i in range(n_iterations):
        # Bootstrap采样
        indices = np.random.choice(len(X_train), len(X_train), replace=True)
        X_boot = X_train[indices]
        y_boot = y_train[indices]

        # 训练模型
        model_copy = clone(model)
        model_copy.fit(X_boot, y_boot)
        predictions.append(model_copy.predict(X_test))

    predictions = np.array(predictions)

    # 计算偏差和方差
    mean_prediction = predictions.mean(axis=0)
    bias_squared = np.mean((mean_prediction - y_test) ** 2)
    variance = np.mean(predictions.var(axis=0))

    return bias_squared, variance

# 使用示例
from sklearn.base import clone
bias_sq, var = bias_variance_decomposition(model, X_train, y_train, X_cv, y_cv)
print(f"偏差²: {bias_sq:.4f}, 方差: {var:.4f}")

# 目标函数
J(w) = MSE(w) + λ||w||²

# 等价于约束优化
minimize MSE(w)
subject to ||w||² ≤ t

# 目标函数
J(w) = MSE(w) + λ||w||₁

# 等价于约束优化
minimize MSE(w)
subject to ||w||₁ ≤ t

from sklearn.linear_model import Ridge, Lasso

# Ridge回归
ridge = Ridge(alpha=0.1)
ridge.fit(X_train, y_train)
print("Ridge权重:", ridge.coef_)

# Lasso回归
lasso = Lasso(alpha=0.1)
lasso.fit(X_train, y_train)
print("Lasso权重:", lasso.coef_)
print("非零权重数量:", np.sum(lasso.coef_ != 0))

from sklearn.model_selection import KFold, cross_val_score

# 5折交叉验证
kfold = KFold(n_splits=5, shuffle=True, random_state=42)
scores = cross_val_score(model, X, y, cv=kfold, scoring='neg_mean_squared_error')

print(f"CV MSE: {-scores.mean():.4f} (+/- {scores.std():.4f})")

from sklearn.model_selection import StratifiedKFold

# 保持类别比例的5折交叉验证
skfold = StratifiedKFold(n_splits=5, shuffle=True, random_state=42)
scores = cross_val_score(model, X, y, cv=skfold, scoring='accuracy')

from sklearn.model_selection import TimeSeriesSplit

# 时间序列专用交叉验证（避免数据泄露）
tscv = TimeSeriesSplit(n_splits=5)
for train_idx, test_idx in tscv.split(X):
    X_train_fold, X_test_fold = X[train_idx], X[test_idx]
    y_train_fold, y_test_fold = y[train_idx], y[test_idx]
    # 训练和评估

from sklearn.model_selection import learning_curve

train_sizes, train_scores, val_scores = learning_curve(
    model, X, y, cv=5,
    train_sizes=np.linspace(0.1, 1.0, 10),
    scoring='neg_mean_squared_error'
)

plt.plot(train_sizes, -train_scores.mean(axis=1), label='训练误差')
plt.plot(train_sizes, -val_scores.mean(axis=1), label='验证误差')
plt.xlabel('训练集大小')
plt.ylabel('MSE')
plt.legend()
plt.show()

from sklearn.model_selection import GridSearchCV

param_grid = {'alpha': np.logspace(-4, 4, 20)}
grid_search = GridSearchCV(Ridge(), param_grid, cv=5,
                           scoring='neg_mean_squared_error')
grid_search.fit(X_train, y_train)

print(f"最优λ: {grid_search.best_params_['alpha']:.4f}")
print(f"最优CV MSE: {-grid_search.best_score_:.4f}")

from skopt import BayesSearchCV

bayes_search = BayesSearchCV(
    Ridge(),
    {'alpha': (1e-4, 1e4, 'log-uniform')},
    n_iter=32,
    cv=5,
    scoring='neg_mean_squared_error'
)
bayes_search.fit(X_train, y_train)

from tensorflow.keras.layers import Dropout, BatchNormalization
from tensorflow.keras.regularizers import l2
from tensorflow.keras.callbacks import EarlyStopping, ReduceLROnPlateau

model = Sequential([
    Dense(128, activation='relu', kernel_regularizer=l2(0.01)),
    BatchNormalization(),  # 批归一化稳定训练
    Dropout(0.3),          # 30%神经元随机失活

    Dense(64, activation='relu', kernel_regularizer=l2(0.01)),
    BatchNormalization(),
    Dropout(0.2),

    Dense(10, activation='softmax')
])

# 早停法
early_stop = EarlyStopping(
    monitor='val_loss',
    patience=10,
    restore_best_weights=True
)

# 学习率衰减
lr_schedule = ReduceLROnPlateau(
    monitor='val_loss',
    factor=0.5,
    patience=5,
    min_lr=1e-6
)

model.compile(optimizer=Adam(learning_rate=0.001),
              loss='sparse_categorical_crossentropy',
              metrics=['accuracy'])

history = model.fit(X_train, y_train,
                    validation_data=(X_cv, y_cv),
                    epochs=100,
                    callbacks=[early_stop, lr_schedule])

from sklearn.feature_selection import SelectKBest, f_regression

# 生成多项式特征
poly = PolynomialFeatures(degree=5)
X_poly = poly.fit_transform(X)

# 选择最重要的K个特征
selector = SelectKBest(f_regression, k=20)
X_selected = selector.fit_transform(X_poly, y)

print(f"原始特征数: {X_poly.shape[1]}")
print(f"选择后特征数: {X_selected.shape[1]}")

from sklearn.decomposition import PCA

pca = PCA(n_components=0.95)  # 保留95%方差
X_pca = pca.fit_transform(X_poly)

print(f"降维后特征数: {X_pca.shape[1]}")
print(f"解释方差比: {pca.explained_variance_ratio_.sum():.2%}")

# ✅ 数据分割检查
assert len(X_train) + len(X_cv) + len(X_test) == len(X), "数据分割错误"
assert 0.6 <= len(X_train)/len(X) <= 0.8, "训练集比例异常"

# ✅ 数据泄露检查
train_indices = set(range(len(X_train)))
cv_indices = set(range(len(X_train), len(X_train)+len(X_cv)))
assert train_indices.isdisjoint(cv_indices), "训练集和验证集有重叠！"

# ✅ 特征缩放检查
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_cv_scaled = scaler.transform(X_cv)  # 注意：只用transform，不用fit
assert np.abs(X_train_scaled.mean()) < 0.1, "特征未正确标准化"

# ✅ 过拟合检查
train_error = mean_squared_error(y_train, model.predict(X_train))
cv_error = mean_squared_error(y_cv, model.predict(X_cv))
if cv_error > 2 * train_error:
    print("⚠️ 警告：模型可能过拟合！")