中高级算法实战教程

# ========== 决策树分类 ==========
from sklearn.tree import DecisionTreeClassifier, DecisionTreeRegressor, plot_tree
from sklearn.model_selection import train_test_split, GridSearchCV
from sklearn.metrics import accuracy_score, classification_report
import matplotlib.pyplot as plt

# 1. 数据准备
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# 2. 基础模型训练
dt = DecisionTreeClassifier(
    criterion='gini',      # 分裂标准：gini / entropy
    max_depth=5,            # 最大深度（核心防过拟合参数）
    min_samples_split=10,   # 内部节点最少样本数
    min_samples_leaf=5,     # 叶节点最少样本数
    max_features='sqrt',    # 每次分裂考虑的特征数
    random_state=42
)
dt.fit(X_train, y_train)

# 3. 评估
y_pred = dt.predict(X_test)
print(f'准确率: {accuracy_score(y_test, y_pred):.4f}')
print(classification_report(y_test, y_pred))

# 4. 特征重要性
importances = dt.feature_importances_
for name, imp in zip(feature_names, importances):
    print(f'{name}: {imp:.4f}')

# 5. 可视化决策树
plt.figure(figsize=(20, 10))
plot_tree(dt, feature_names=feature_names, class_names=class_names, 
          filled=True, rounded=True, fontsize=10)
plt.show()

# 6. 超参数调优（防止过拟合）
param_grid = {
    'max_depth': [3, 5, 7, 10, None],
    'min_samples_split': [2, 5, 10, 20],
    'min_samples_leaf': [1, 3, 5, 10],
    'criterion': ['gini', 'entropy']
}
grid = GridSearchCV(DecisionTreeClassifier(random_state=42), 
                    param_grid, cv=5, scoring='accuracy')
grid.fit(X_train, y_train)
print(f'最优参数: {grid.best_params_}')
print(f'最优CV分数: {grid.best_score_:.4f}')

# ========== 决策树回归 ==========
dt_reg = DecisionTreeRegressor(max_depth=5, min_samples_leaf=10, random_state=42)
dt_reg.fit(X_train, y_train)
print(f'回归 R²: {dt_reg.score(X_test, y_test):.4f}')

from sklearn.ensemble import RandomForestClassifier, RandomForestRegressor
from sklearn.model_selection import cross_val_score
import numpy as np

# ========== 随机森林分类 ==========
rf = RandomForestClassifier(
    n_estimators=100,        # 树的数量（越大越稳定，但边际收益递减）
    max_depth=10,            # 树的最大深度
    min_samples_split=5,     # 内部节点最小样本数
    min_samples_leaf=2,      # 叶节点最小样本数
    max_features='sqrt',     # 每次分裂考虑的特征数（分类默认sqrt，回归默认所有）
    bootstrap=True,          # 是否使用Bootstrap抽样
    oob_score=True,          # 是否计算OOB评分
    n_jobs=-1,               # 并行计算
    random_state=42
)
rf.fit(X_train, y_train)

# OOB评分（无需验证集）
print(f'OOB Score: {rf.oob_score_:.4f}')

# 特征重要性
importances = rf.feature_importances_
indices = np.argsort(importances)[::-1]
for i in indices[:10]:
    print(f'{feature_names[i]}: {importances[i]:.4f}')

# ========== 关键参数调优 ==========
from sklearn.model_selection import RandomizedSearchCV

param_dist = {
    'n_estimators': [50, 100, 200, 300],
    'max_depth': [5, 10, 15, 20, None],
    'min_samples_split': [2, 5, 10],
    'min_samples_leaf': [1, 2, 4],
    'max_features': ['sqrt', 'log2', None]
}
random_search = RandomizedSearchCV(
    RandomForestClassifier(random_state=42),
    param_dist, n_iter=30, cv=5, scoring='roc_auc', n_jobs=-1
)
random_search.fit(X_train, y_train)
print(f'最优参数: {random_search.best_params_}')

# ========== 学习曲线（判断是否需要更多树）==========
train_scores, test_scores = [], []
n_estimators_range = range(10, 200, 10)
for n in n_estimators_range:
    rf_temp = RandomForestClassifier(n_estimators=n, random_state=42)
    rf_temp.fit(X_train, y_train)
    train_scores.append(rf_temp.score(X_train, y_train))
    test_scores.append(rf_temp.score(X_test, y_test))
# 画出 n_estimators vs 准确率曲线，找到拐点

# ========== 回归版本 ==========
rf_reg = RandomForestRegressor(
    n_estimators=100, max_depth=10, 
    min_samples_leaf=5, n_jobs=-1, random_state=42
)
rf_reg.fit(X_train, y_train)
print(f'回归 R²: {rf_reg.score(X_test, y_test):.4f}')

import xgboost as xgb
from sklearn.model_selection import train_test_split, cross_val_score
from sklearn.metrics import accuracy_score, roc_auc_score, mean_squared_error
import numpy as np

# ========== XGBoost 分类 ==========
xgb_clf = xgb.XGBClassifier(
    n_estimators=100,        # 树的数量
    max_depth=6,             # 树的最大深度（典型值3-10）
    learning_rate=0.1,       # 学习率（典型值0.01-0.3）
    subsample=0.8,           # 每棵树用的样本比例（防止过拟合）
    colsample_bytree=0.8,    # 每棵树用的特征比例
    reg_alpha=0.1,           # L1正则化（稀疏特征时有用）
    reg_lambda=1.0,          # L2正则化（默认1）
    gamma=0,                 # 分裂所需的最小损失减少量
    min_child_weight=1,      # 叶节点最小样本权重和
    scale_pos_weight=1,      # 正负样本不平衡时的权重调整
    eval_metric='logloss',   # 评估指标
    early_stopping_rounds=20, # 早停轮数
    random_state=42,
    n_jobs=-1
)

# 训练（使用 eval_set 做早停）
xgb_clf.fit(
    X_train, y_train,
    eval_set=[(X_test, y_test)],
    verbose=False
)

# 评估
y_pred = xgb_clf.predict(X_test)
y_proba = xgb_clf.predict_proba(X_test)[:, 1]
print(f'准确率: {accuracy_score(y_test, y_pred):.4f}')
print(f'AUC: {roc_auc_score(y_test, y_proba):.4f}')

# 特征重要性（三种方式）
# weight: 特征被用于分裂的次数
# gain: 特征在所有分裂中的平均增益
# cover: 特征覆盖的样本数
importance = xgb_clf.get_booster().get_score(importance_type='gain')

# ========== XGBoost 回归 ==========
xgb_reg = xgb.XGBRegressor(
    n_estimators=200, max_depth=7, learning_rate=0.05,
    subsample=0.8, colsample_bytree=0.8,
    reg_alpha=0.1, reg_lambda=1.0,
    early_stopping_rounds=30, random_state=42
)
xgb_reg.fit(X_train, y_train, eval_set=[(X_test, y_test)], verbose=False)
y_pred_reg = xgb_reg.predict(X_test)
print(f'RMSE: {np.sqrt(mean_squared_error(y_test, y_pred_reg)):.4f}')

# ========== 超参数调优（关键参数）==========
# 调参顺序建议：n_estimators→max_depth+min_child_weight→gamma→subsample+colsample→reg_alpha+reg_lambda→learning_rate
from sklearn.model_selection import GridSearchCV

param_grid = {
    'max_depth': [3, 5, 7],
    'min_child_weight': [1, 3, 5],
    'gamma': [0, 0.1, 0.2],
    'subsample': [0.6, 0.8, 1.0],
    'colsample_bytree': [0.6, 0.8, 1.0]
}
grid = GridSearchCV(xgb.XGBClassifier(n_estimators=100, learning_rate=0.1, random_state=42),
                    param_grid, cv=5, scoring='roc_auc', n_jobs=-1)
grid.fit(X_train, y_train)
print(f'最优参数: {grid.best_params_}')

# ========== 自定义损失函数（示例：Focal Loss）==========
def focal_loss(y_pred, dtrain, alpha=0.25, gamma=2.0):
    y_true = dtrain.get_label()
    grad = ...  # 自定义一阶梯度
    hess = ...  # 自定义二阶梯度
    return grad, hess
# xgb_clf.fit(X_train, y_train, obj=focal_loss)

import lightgbm as lgb
from sklearn.metrics import accuracy_score, mean_squared_error
import numpy as np

# ========== LightGBM 分类 ==========
lgb_clf = lgb.LGBMClassifier(
    n_estimators=200,
    max_depth=7,
    num_leaves=31,           # 叶子数（核心参数！≈2^max_depth）
    learning_rate=0.05,
    min_child_samples=20,    # 叶节点最少样本数
    subsample=0.8,
    colsample_bytree=0.8,
    reg_alpha=0.1,
    reg_lambda=1.0,
    boosting_type='gbdt',    # gbdt / dart / goss
    objective='binary',      # 目标函数
    metric='auc',
    random_state=42,
    n_jobs=-1,
    verbose=-1
)

lgb_clf.fit(
    X_train, y_train,
    eval_set=[(X_test, y_test)],
    eval_metric='auc',
    callbacks=[lgb.early_stopping(50), lgb.log_evaluation(0)]
)

y_pred = lgb_clf.predict(X_test)
y_proba = lgb_clf.predict_proba(X_test)[:, 1]
print(f'AUC: {roc_auc_score(y_test, y_proba):.4f}')

# 特征重要性
lgb.plot_importance(lgb_clf, max_num_features=15, figsize=(10, 6))

# ========== LightGBM 回归 ==========
lgb_reg = lgb.LGBMRegressor(
    n_estimators=300, num_leaves=63, max_depth=8,
    learning_rate=0.03, min_child_samples=30,
    subsample=0.8, colsample_bytree=0.8,
    random_state=42, verbose=-1
)
lgb_reg.fit(X_train, y_train, eval_set=[(X_test, y_test)],
            eval_metric='rmse',
            callbacks=[lgb.early_stopping(50)])
print(f'RMSE: {np.sqrt(mean_squared_error(y_test, lgb_reg.predict(X_test))):.4f}')

# ========== 使用GOSS（更快）==========
lgb_goss = lgb.LGBMClassifier(
    boosting_type='goss',    # 梯度单边采样，更快
    n_estimators=200, num_leaves=31,
    top_rate=0.2,            # 保留大梯度样本比例
    other_rate=0.1,          # 小梯度样本采样比例
    random_state=42
)

# ========== 使用DART（防过拟合）==========
lgb_dart = lgb.LGBMClassifier(
    boosting_type='dart',    # Dropouts meet Multiple Additive Regression Trees
    n_estimators=200, num_leaves=31,
    drop_rate=0.1,           # 每轮丢弃树的比例
    random_state=42
)

from catboost import CatBoostClassifier, CatBoostRegressor, Pool
from sklearn.metrics import accuracy_score, roc_auc_score

# ========== CatBoost 分类 ==========
# 标识类别特征列索引
cat_features = [0, 2, 5]  # 城市、品类、渠道等列的索引

cb_clf = CatBoostClassifier(
    iterations=500,          # 迭代次数
    depth=6,                 # 树深度（典型值4-10）
    learning_rate=0.03,
    l2_leaf_reg=3,           # L2正则化系数
    border_count=128,        # 数值特征分箱数
    random_seed=42,
    eval_metric='AUC',
    early_stopping_rounds=50,
    verbose=100,             # 每100轮打印一次
    task_type='CPU'          # CPU / GPU
)

cb_clf.fit(
    X_train, y_train,
    cat_features=cat_features,
    eval_set=(X_test, y_test),
    plot=True               # 自动画训练曲线
)

y_pred = cb_clf.predict(X_test)
y_proba = cb_clf.predict_proba(X_test)[:, 1]
print(f'AUC: {roc_auc_score(y_test, y_proba):.4f}')

# 特征重要性
cb_clf.get_feature_importance(prettified=True)

# ========== CatBoost 回归 ==========
cb_reg = CatBoostRegressor(
    iterations=800, depth=8, learning_rate=0.03,
    l2_leaf_reg=5, random_seed=42,
    early_stopping_rounds=50, verbose=100
)
cb_reg.fit(X_train, y_train, cat_features=cat_features,
           eval_set=(X_test, y_test))

# ========== 使用 Pool 对象（推荐方式）==========
train_pool = Pool(X_train, y_train, cat_features=cat_features)
test_pool = Pool(X_test, y_test, cat_features=cat_features)

cb_model = CatBoostClassifier(iterations=500, depth=6, random_seed=42)
cb_model.fit(train_pool, eval_set=test_pool, early_stopping_rounds=50)

# ========== 模型分析 ==========
# SHAP分析
cb_model.get_feature_importance(train_pool, type='ShapValues')

# 模型保存与加载
cb_model.save_model('catboost_model.cbm')
# cb_model.load_model('catboost_model.cbm')

from sklearn.naive_bayes import GaussianNB, MultinomialNB, BernoulliNB
from sklearn.model_selection import cross_val_score
from sklearn.metrics import accuracy_score, classification_report
import numpy as np

# ========== 高斯朴素贝叶斯（数值特征）==========
gnb = GaussianNB(var_smoothing=1e-9)  # 方差平滑，防止零方差
gnb.fit(X_train, y_train)
y_pred = gnb.predict(X_test)
y_proba = gnb.predict_proba(X_test)  # 输出概率
print(f'准确率: {accuracy_score(y_test, y_pred):.4f}')

# ========== 多项式朴素贝叶斯（文本分类）==========
from sklearn.feature_extraction.text import CountVectorizer
# 文本→词频矩阵
vectorizer = CountVectorizer(max_features=5000)
X_counts = vectorizer.fit_transform(texts)
mnb = MultinomialNB(alpha=1.0)  # 拉普拉斯平滑
mnb.fit(X_counts, labels)

# 查看每个类别的Top关键词
feature_names = vectorizer.get_feature_names_out()
for i, class_label in enumerate(mnb.classes_):
    top_indices = np.argsort(mnb.feature_log_prob_[i])[-10:][::-1]
    top_words = [feature_names[idx] for idx in top_indices]
    print(f'类别 {class_label}: {top_words}')

# ========== 混合类型特征处理 ==========
# 对于既有数值又有类别的数据，可用Pipeline
from sklearn.pipeline import Pipeline
from sklearn.compose import ColumnTransformer
from sklearn.preprocessing import StandardScaler

# 先对数值特征做标准化，再喂给GaussianNB
preprocessor = ColumnTransformer([
    ('num', StandardScaler(), numerical_cols)
])
pipeline = Pipeline([
    ('preprocessor', preprocessor),
    ('classifier', GaussianNB())
])
pipeline.fit(X_train, y_train)

# ========== 半监督学习（部分数据无标签）==========
# 用已标注数据训练 → 预测未标注数据 → 将高置信度样本加入训练集 → 迭代
gnb = GaussianNB()
for iteration in range(5):
    gnb.fit(X_labeled, y_labeled)
    probas = gnb.predict_proba(X_unlabeled)
    confidence = np.max(probas, axis=1)
    high_conf_idx = confidence > 0.95
    # 将高置信度样本加入训练集
    X_labeled = np.vstack([X_labeled, X_unlabeled[high_conf_idx]])
    y_labeled = np.hstack([y_labeled, gnb.predict(X_unlabeled[high_conf_idx])])
    X_unlabeled = X_unlabeled[~high_conf_idx]

from sklearn.svm import SVC, SVR
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import GridSearchCV
from sklearn.pipeline import Pipeline

# ========== SVM分类（关键：特征必须标准化！）==========
pipeline = Pipeline([
    ('scaler', StandardScaler()),  # SVM对尺度敏感！
    ('svm', SVC(kernel='rbf', C=1.0, gamma='scale', 
                probability=True, random_state=42))
])
pipeline.fit(X_train, y_train)
y_pred = pipeline.predict(X_test)
y_proba = pipeline.predict_proba(X_test)

# ========== 不同核函数对比 ==========
kernels = ['linear', 'poly', 'rbf', 'sigmoid']
for kernel in kernels:
    svm = SVC(kernel=kernel, C=1.0, gamma='scale', random_state=42)
    svm.fit(X_train_scaled, y_train)
    score = svm.score(X_test_scaled, y_test)
    print(f'{kernel}: {score:.4f}')

# ========== 超参数调优（RBF核）==========
param_grid = {
    'C': [0.01, 0.1, 1, 10, 100],        # 惩罚系数
    'gamma': [0.001, 0.01, 0.1, 1, 'scale', 'auto'],  # RBF核宽度
    'kernel': ['rbf', 'linear']
}
grid = GridSearchCV(SVC(random_state=42), param_grid, 
                    cv=5, scoring='accuracy', n_jobs=-1)
grid.fit(X_train_scaled, y_train)
print(f'最优参数: {grid.best_params_}')
print(f'最优CV分数: {grid.best_score_:.4f}')

# ========== SVM回归 ==========
svr = SVR(kernel='rbf', C=1.0, epsilon=0.1, gamma='scale')
svr.fit(X_train_scaled, y_train)
print(f'R²: {svr.score(X_test_scaled, y_test):.4f}')

# ========== 大数据集下用LinearSVC（更快）==========
from sklearn.svm import LinearSVC
linear_svc = LinearSVC(C=1.0, max_iter=10000, dual=False, random_state=42)
linear_svc.fit(X_train_scaled, y_train)

from prophet import Prophet
import pandas as pd
import matplotlib.pyplot as plt

# ========== 数据准备 ==========
# Prophet要求两列：ds(日期)和y(数值)
df = pd.DataFrame({
    'ds': pd.date_range('2024-01-01', periods=365, freq='D'),
    'y': daily_gmv_values  # 每日GMV
})

# ========== 添加中国节假日 ==========
from prophet.make_holidays import make_holidays_df
# 获取中国节假日（春节、国庆等）
chinese_holidays = make_holidays_df(
    year_list=[2024, 2025, 2026],
    country='CN'
)

# 添加自定义大促日期
custom_events = pd.DataFrame({
    'holiday': '双11',
    'ds': pd.to_datetime(['2024-11-11', '2025-11-11']),
    'lower_window': -7,  # 提前7天预热
    'upper_window': 3,   # 延后3天
})

# ========== 模型训练 ==========
model = Prophet(
    growth='linear',              # 趋势类型：linear / logistic
    seasonality_mode='additive',  # 周期性模式：additive / multiplicative
    yearly_seasonality=True,      # 年周期（傅里叶级数阶数默认10）
    weekly_seasonality=True,      # 周周期（默认3阶）
    daily_seasonality=False,      # 日周期（小时级数据才开）
    holidays=pd.concat([chinese_holidays, custom_events]),
    changepoint_prior_scale=0.05, # 趋势灵活性（越大越灵活，越容易过拟合）
    seasonality_prior_scale=10.0, # 周期性强度
    holidays_prior_scale=10.0,    # 节假日效应强度
    interval_width=0.95,          # 置信区间宽度
)

model.fit(df)

# ========== 预测 ==========
future = model.make_future_dataframe(periods=90)  # 预测未来90天
forecast = model.predict(future)

# ========== 可视化 ==========
# 1. 整体预测图
model.plot(forecast)
plt.title('GMV 预测（含置信区间）')
plt.show()

# 2. 组件分解图
model.plot_components(forecast)
plt.show()

# 3. 趋势变点
from prophet.plot import add_changepoints_to_plot
fig = model.plot(forecast)
add_changepoints_to_plot(fig.gca(), model, forecast)

# ========== 评估 ==========
from prophet.diagnostics import cross_validation, performance_metrics

# 时间序列交叉验证
df_cv = cross_validation(
    model, initial='180 days',  # 初始训练窗口
    period='30 days',           # 每次向前滚动
    horizon='30 days'           # 预测视野
)

# 评估指标
df_metrics = performance_metrics(df_cv)
print(df_metrics[['horizon', 'rmse', 'mape', 'mae']].head())

# ========== 预测结果解读 ==========
print(forecast[['ds', 'yhat', 'yhat_lower', 'yhat_upper']].tail(7))

import tensorflow as tf
from tensorflow import keras
from tensorflow.keras import layers, models
from tensorflow.keras.preprocessing.image import ImageDataGenerator

# ========== 基础CNN模型 ==========
model = models.Sequential([
    # 卷积块1
    layers.Conv2D(32, (3,3), activation='relu', padding='same', input_shape=(224,224,3)),
    layers.BatchNormalization(),
    layers.Conv2D(32, (3,3), activation='relu', padding='same'),
    layers.BatchNormalization(),
    layers.MaxPooling2D((2,2)),
    layers.Dropout(0.25),
    
    # 卷积块2
    layers.Conv2D(64, (3,3), activation='relu', padding='same'),
    layers.BatchNormalization(),
    layers.Conv2D(64, (3,3), activation='relu', padding='same'),
    layers.BatchNormalization(),
    layers.MaxPooling2D((2,2)),
    layers.Dropout(0.25),
    
    # 卷积块3
    layers.Conv2D(128, (3,3), activation='relu', padding='same'),
    layers.BatchNormalization(),
    layers.MaxPooling2D((2,2)),
    layers.Dropout(0.25),
    
    # 全连接层
    layers.GlobalAveragePooling2D(),
    layers.Dense(256, activation='relu'),
    layers.BatchNormalization(),
    layers.Dropout(0.5),
    layers.Dense(num_classes, activation='softmax')
])

model.compile(
    optimizer=keras.optimizers.Adam(learning_rate=0.001),
    loss='categorical_crossentropy',
    metrics=['accuracy']
)
model.summary()

# ========== 数据增强（防过拟合+提升泛化）==========
train_datagen = ImageDataGenerator(
    rescale=1./255,
    rotation_range=20,      # 随机旋转
    width_shift_range=0.2,  # 水平平移
    height_shift_range=0.2, # 垂直平移
    shear_range=0.2,        # 剪切变换
    zoom_range=0.2,         # 随机缩放
    horizontal_flip=True,   # 水平翻转
    fill_mode='nearest'
)

# ========== 训练 ==========
history = model.fit(
    train_datagen.flow(X_train, y_train, batch_size=32),
    validation_data=(X_test, y_test),
    epochs=50,
    callbacks=[
        keras.callbacks.EarlyStopping(patience=10, restore_best_weights=True),
        keras.callbacks.ReduceLROnPlateau(factor=0.5, patience=5, min_lr=1e-6)
    ]
)

# ========== Grad-CAM 可视化（看模型关注哪里）==========
import cv2, numpy as np
def grad_cam(model, img, layer_name):
    grad_model = keras.Model(
        [model.inputs], 
        [model.get_layer(layer_name).output, model.output]
    )
    with tf.GradientTape() as tape:
        conv_outputs, predictions = grad_model(img)
        class_idx = tf.argmax(predictions[0])
        loss = predictions[:, class_idx]
    grads = tape.gradient(loss, conv_outputs)
    # ...生成热力图叠加到原图...

import tensorflow as tf
from tensorflow import keras

# ========== 使用预训练模型（ResNet50）==========
base_model = keras.applications.ResNet50(
    weights='imagenet',          # 加载ImageNet权重
    include_top=False,           # 去掉原始分类头
    input_shape=(224, 224, 3)
)

# 冻结底层（策略：先冻结全部训练分类头，再解冻部分微调）
base_model.trainable = False

# 构建新分类头
model = keras.Sequential([
    base_model,
    keras.layers.GlobalAveragePooling2D(),
    keras.layers.Dense(512, activation='relu'),
    keras.layers.BatchNormalization(),
    keras.layers.Dropout(0.5),
    keras.layers.Dense(num_classes, activation='softmax')
])

model.compile(
    optimizer=keras.optimizers.Adam(1e-3),
    loss='categorical_crossentropy',
    metrics=['accuracy']
)

# 第一阶段：只训练分类头
history1 = model.fit(train_data, epochs=10, validation_data=val_data)

# 第二阶段：解冻顶层进行微调
base_model.trainable = True
# 冻结前100层，只微调后面
for layer in base_model.layers[:100]:
    layer.trainable = False

model.compile(
    optimizer=keras.optimizers.Adam(1e-5),  # 极小学习率
    loss='categorical_crossentropy',
    metrics=['accuracy']
)

history2 = model.fit(train_data, epochs=20, validation_data=val_data,
                     callbacks=[keras.callbacks.EarlyStopping(patience=5)])

# ========== 常用预训练模型选择 ==========
# ResNet50/101/152：通用性强，最经典
# EfficientNetB0-B7：效率高，准确率好
# MobileNetV2/V3：轻量级，适合移动端
# ConvNeXt：2020s新架构，SOTA性能
# ViT (Vision Transformer)：适合大数据集

efficientnet = keras.applications.EfficientNetB3(
    weights='imagenet', include_top=False, input_shape=(300,300,3)
)

# ========== 数据量不同时的策略 ==========
# <1000张/类：只训练分类头（base_model.trainable=False）
# 1000-10000张/类：解冻最后几层微调
# >10000张/类：可以全量微调

import numpy as np
from tensorflow import keras
from tensorflow.keras import layers
from sklearn.preprocessing import MinMaxScaler

# ========== 时序数据准备（滑动窗口）==========
def create_sequences(data, seq_length=30):
    """将时序数据转为监督学习格式"""
    X, y = [], []
    for i in range(len(data) - seq_length):
        X.append(data[i:i+seq_length])
        y.append(data[i+seq_length])
    return np.array(X), np.array(y)

# 数据标准化
scaler = MinMaxScaler()
data_scaled = scaler.fit_transform(data.reshape(-1, 1))

# 创建序列
seq_length = 30  # 用过去30天预测下一天
X, y = create_sequences(data_scaled, seq_length)
X = X.reshape(X.shape[0], X.shape[1], 1)  # (样本数, 时间步, 特征数)

# ========== LSTM模型 ==========
model = keras.Sequential([
    layers.LSTM(64, return_sequences=True, input_shape=(seq_length, 1)),
    layers.Dropout(0.2),
    layers.LSTM(32, return_sequences=False),
    layers.Dropout(0.2),
    layers.Dense(16, activation='relu'),
    layers.Dense(1)
])

model.compile(
    optimizer=keras.optimizers.Adam(learning_rate=0.001),
    loss='mse',
    metrics=['mae']
)

history = model.fit(
    X_train, y_train,
    validation_data=(X_test, y_test),
    epochs=100,
    batch_size=32,
    callbacks=[
        keras.callbacks.EarlyStopping(patience=15, restore_best_weights=True),
        keras.callbacks.ReduceLROnPlateau(factor=0.5, patience=5)
    ]
)

# ========== 多步预测（预测未来N天）==========
def multi_step_forecast(model, last_sequence, steps=7):
    """用最后已知序列，逐步预测未来steps天"""
    predictions = []
    current_seq = last_sequence.copy()
    for _ in range(steps):
        pred = model.predict(current_seq.reshape(1, seq_length, 1), verbose=0)
        predictions.append(pred[0,0])
        # 将预测值加入序列，滑动窗口
        current_seq = np.append(current_seq[1:], pred)
    return scaler.inverse_transform(np.array(predictions).reshape(-1,1))

# ========== BiLSTM（双向LSTM，更好捕捉上下文）==========
from tensorflow.keras.layers import Bidirectional

model_bilstm = keras.Sequential([
    Bidirectional(layers.LSTM(64, return_sequences=True), input_shape=(seq_length, 1)),
    layers.Dropout(0.2),
    Bidirectional(layers.LSTM(32)),
    layers.Dropout(0.2),
    layers.Dense(16, activation='relu'),
    layers.Dense(1)
])

import pandas as pd
import numpy as np
from statsmodels.tsa.arima.model import ARIMA
from statsmodels.tsa.statespace.sarimax import SARIMAX
from statsmodels.graphics.tsaplots import plot_acf, plot_pacf
from statsmodels.tsa.stattools import adfuller
import pmdarima as pm

# ========== 平稳性检验 ==========
def check_stationarity(series):
    result = adfuller(series.dropna())
    print(f'ADF统计量: {result[0]:.4f}')
    print(f'p值: {result[1]:.4f}')
    return result[1] < 0.05

# ========== ACF/PACF 确定 p,q ==========
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
plot_acf(series, lags=40, ax=axes[0])
plot_pacf(series, lags=40, ax=axes[1])
plt.show()

# ========== 自动搜索最优参数 ==========
auto_model = pm.auto_arima(
    train_series,
    seasonal=True, m=7,          # 周周期=7天
    start_p=0, max_p=5,
    start_q=0, max_q=5,
    start_P=0, max_P=3,
    start_Q=0, max_Q=3,
    trace=True,                  # 打印搜索过程
    error_action='ignore',
    suppress_warnings=True,
    stepwise=True
)
print(auto_model.summary())

# ========== SARIMAX（含外生变量）==========
model = SARIMAX(
    train_series,
    exog=exog_train,            # 外生变量：节假日、天气等
    order=(2, 1, 2),            # ARIMA(p,d,q)
    seasonal_order=(1, 1, 1, 7), # (P,D,Q,s) s=7表示周周期
    enforce_stationarity=False,
    enforce_invertibility=False
)
results = model.fit(disp=False)
print(results.summary())

# 预测
forecast = results.get_forecast(steps=30, exog=exog_forecast)
forecast_mean = forecast.predicted_mean
forecast_ci = forecast.conf_int(alpha=0.05)

# ========== 残差诊断 ==========
residuals = results.resid
# Ljung-Box检验：残差是否为白噪声
from statsmodels.stats.diagnostic import acorr_ljungbox
lb_test = acorr_ljungbox(residuals, lags=[10, 20, 30])
print(lb_test)  # p>0.05说明残差是白噪声，模型充分提取了信息

from sklearn.ensemble import IsolationForest
import numpy as np

# ========== 孤立森林异常检测 ==========
iso_forest = IsolationForest(
    n_estimators=100,         # 树的数量
    contamination=0.05,       # 预期的异常比例（核心参数）
    max_samples='auto',       # 每棵树用的样本数
    max_features=1.0,         # 每棵树用的特征比例
    bootstrap=False,          # 是否用Bootstrap
    random_state=42,
    n_jobs=-1
)

# 训练（注意：无监督，y不需要！）
iso_forest.fit(X)

# 预测：1=正常，-1=异常
anomaly_labels = iso_forest.predict(X)
anomaly_scores = iso_forest.decision_function(X)  # 负值越大=越异常

# 获取异常样本
anomalies = X[anomaly_labels == -1]
print(f'异常样本数: {len(anomalies)} / {len(X)} ({len(anomalies)/len(X)*100:.1f}%)')

# ========== 确定最优contamination ==========
# 用业务规则交叉验证：已知某天是正常/异常的
# 调整contamination使误报率和漏报率达到平衡
from sklearn.metrics import precision_recall_curve

# ========== 多维度异常监控 ==========
# 对每个维度独立检测，再综合判断
dimensions = ['GMV', '订单量', '核销率', '客单价']
all_anomalies = pd.DataFrame()
for dim in dimensions:
    iso = IsolationForest(contamination=0.05, random_state=42)
    scores = iso.fit_predict(X[[dim]])
    all_anomalies[dim+'_anomaly'] = (scores == -1).astype(int)

# 综合判断：>=2个维度异常才告警
all_anomalies['alert'] = all_anomalies.sum(axis=1) >= 2

# ========== 特征重要性（SHAP解释）==========
import shap
explainer = shap.TreeExplainer(iso_forest)
shap_values = explainer.shap_values(X)
shap.summary_plot(shap_values, X)

from sklearn.cluster import DBSCAN
from sklearn.preprocessing import StandardScaler
from sklearn.neighbors import NearestNeighbors
import numpy as np

# ========== 确定最优 eps（K距离图）==========
def find_optimal_eps(X, min_pts=5):
    neighbors = NearestNeighbors(n_neighbors=min_pts)
    neighbors_fit = neighbors.fit(X)
    distances, _ = neighbors_fit.kneighbors(X)
    # 取每个点到第min_pts近邻的距离，排序
    k_distances = np.sort(distances[:, -1])
    # 画K距离图，找拐点
    plt.plot(k_distances)
    plt.xlabel('Points sorted by distance')
    plt.ylabel(f'{min_pts}-th nearest distance')
    plt.title('K-Distance Graph (找肘部)')
    plt.show()
    return k_distances

# ========== DBSCAN 聚类 ==========
# 数据标准化（DBSCAN对尺度敏感！）
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

dbscan = DBSCAN(
    eps=0.5,           # 邻域半径（根据K距离图确定）
    min_samples=5,     # 核心点最少邻居数
    metric='euclidean',
    n_jobs=-1
)
labels = dbscan.fit_predict(X_scaled)

# labels=-1 表示噪声点
n_clusters = len(set(labels)) - (1 if -1 in labels else 0)
n_noise = list(labels).count(-1)
print(f'聚类数: {n_clusters}, 噪声点: {n_noise} ({n_noise/len(labels)*100:.1f}%)')

# ========== 聚类可视化（PCA降维）==========
from sklearn.decomposition import PCA
pca = PCA(n_components=2)
X_pca = pca.fit_transform(X_scaled)

plt.figure(figsize=(10, 8))
colors = plt.cm.tab20(np.linspace(0, 1, n_clusters))
for i in range(n_clusters):
    plt.scatter(X_pca[labels==i, 0], X_pca[labels==i, 1], 
                c=[colors[i]], label=f'Cluster {i}', s=20, alpha=0.7)
plt.scatter(X_pca[labels==-1, 0], X_pca[labels==-1, 1], 
            c='red', marker='x', label='Noise', s=40)
plt.legend()
plt.title(f'DBSCAN: {n_clusters} clusters, {n_noise} noise points')
plt.show()

# ========== 噪声分析（噪声点有什么特征？）==========
noise_df = pd.DataFrame(X, columns=feature_names)
noise_df['label'] = labels
noise_only = noise_df[noise_df['label'] == -1]
print(noise_only.describe())  # 噪声点的统计特征

# ========== HDBSCAN（DBSCAN的进阶版，自动选eps）==========
# pip install hdbscan
import hdbscan
hdb = hdbscan.HDBSCAN(
    min_cluster_size=10,     # 最小簇大小
    min_samples=5,           # 核心点参数
    metric='euclidean',
    cluster_selection_method='eom'  # eom/leaf
)
hdb_labels = hdb.fit_predict(X_scaled)
print(f'HDBSCAN 聚类数: {len(set(hdb_labels))-1}')