10. 特征工程实战

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

df = pd.read_csv('data.csv')

# 基本信息
print(df.info())
print(df.describe())
print(df.isnull().sum() / len(df) * 100)  # 缺失率

# 分布分析
fig, axes = plt.subplots(2, 3, figsize=(15, 10))
for idx, col in enumerate(df.select_dtypes(include=[np.number]).columns[:6]):
    ax = axes[idx // 3, idx % 3]
    sns.histplot(df[col], kde=True, ax=ax)
    ax.set_title(f'Distribution of {col}')

# 相关性热力图
corr = df.corr(numeric_only=True)
plt.figure(figsize=(12, 10))
sns.heatmap(corr, annot=True, cmap='coolwarm', center=0)

# IQR 方法
def detect_outliers_iqr(df, column):
    Q1 = df[column].quantile(0.25)
    Q3 = df[column].quantile(0.75)
    IQR = Q3 - Q1
    lower = Q1 - 1.5 * IQR
    upper = Q3 + 1.5 * IQR
    return df[(df[column] < lower) | (df[column] > upper)]

# Z-Score 方法
from scipy import stats
z_scores = np.abs(stats.zscore(df.select_dtypes(include=[np.number])))
outliers = (z_scores > 3).any(axis=1)

# Isolation Forest
from sklearn.ensemble import IsolationForest
iso = IsolationForest(contamination=0.1, random_state=42)
outliers = iso.fit_predict(df.select_dtypes(include=[np.number]))

from sklearn.preprocessing import StandardScaler, MinMaxScaler, RobustScaler

# Z-Score 标准化（均值为 0，方差为 1）
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# Min-Max 归一化（缩放到 [0, 1]）
scaler = MinMaxScaler(feature_range=(0, 1))
X_normalized = scaler.fit_transform(X)

# RobustScaler（使用中位数和 IQR，对异常值鲁棒）
scaler = RobustScaler()
X_robust = scaler.fit_transform(X)

# 对数变换（处理右偏分布）
df['log_income'] = np.log1p(df['income'])

# Box-Cox 变换（自动找最优 λ）
from scipy.stats import boxcox
df['boxcox_income'], lambda_param = boxcox(df['income'] + 1)

# 平方根变换
df['sqrt_age'] = np.sqrt(df['age'])

# 分位数变换（映射到均匀/正态分布）
from sklearn.preprocessing import QuantileTransformer
qt = QuantileTransformer(output_distribution='normal')
X_transformed = qt.fit_transform(X)

# 等宽分箱
df['age_bin'] = pd.cut(df['age'], bins=5, labels=['very_young', 'young', 'middle', 'old', 'very_old'])

# 等频分箱（每个区间样本数相同）
df['income_quintile'] = pd.qcut(df['income'], q=5, labels=['Q1', 'Q2', 'Q3', 'Q4', 'Q5'])

# 基于业务规则分箱
def age_group(age):
    if age < 18: return '未成年'
    elif age < 35: return '青年'
    elif age < 50: return '中年'
    elif age < 65: return '中老年'
    else: return '老年'

df['age_group'] = df['age'].apply(age_group)

# One-Hot 编码（低基数类别）
df_encoded = pd.get_dummies(df, columns=['gender', 'city_level'], prefix=['gender', 'level'])

# Label Encoding（高基数或有序类别）
from sklearn.preprocessing import LabelEncoder
le = LabelEncoder()
df['brand_encoded'] = le.fit_transform(df['brand'])

# 目标编码 / 均值编码（高基数，用目标变量统计值编码）
def target_encode(df, col, target, smoothing=1.0):
    global_mean = df[target].mean()
    agg = df.groupby(col)[target].agg(['mean', 'count'])
    counts = agg['count']
    means = agg['mean']
    smooth = (counts * means + smoothing * global_mean) / (counts + smoothing)
    return df[col].map(smooth)

df['city_target_encoded'] = target_encode(df, 'city', 'target', smoothing=10)

# 频率编码（用出现次数编码）
freq_map = df['category'].value_counts().to_dict()
df['category_freq'] = df['category'].map(freq_map)

# 聚类编码（对高基数类别聚类）
from sklearn.cluster import KMeans

# 用其他特征聚类
category_features = df.groupby('category').agg({
    'price': ['mean', 'std'],
    'sales': ['mean', 'sum'],
    'rating': 'mean'
}).fillna(0)

kmeans = KMeans(n_clusters=20, random_state=42)
df['category_cluster'] = kmeans.fit_predict(category_features).reindex(df['category']).values

df['datetime'] = pd.to_datetime(df['timestamp'])

df['year'] = df['datetime'].dt.year
df['month'] = df['datetime'].dt.month
df['day'] = df['datetime'].dt.day
df['hour'] = df['datetime'].dt.hour
df['dayofweek'] = df['datetime'].dt.dayofweek
df['dayofyear'] = df['datetime'].dt.dayofyear
df['weekofyear'] = df['datetime'].dt.isocalendar().week
df['quarter'] = df['datetime'].dt.quarter

# 周期性编码（用 sin/cos 保留周期信息）
df['hour_sin'] = np.sin(2 * np.pi * df['hour'] / 24)
df['hour_cos'] = np.cos(2 * np.pi * df['hour'] / 24)
df['month_sin'] = np.sin(2 * np.pi * df['month'] / 12)
df['month_cos'] = np.cos(2 * np.pi * df['month'] / 12)

# 按用户分组，构造滞后特征
df = df.sort_values(['user_id', 'datetime'])

# 滞后特征
df['price_lag_1'] = df.groupby('user_id')['price'].shift(1)
df['price_lag_7'] = df.groupby('user_id')['price'].shift(7)

# 滑动窗口统计
for window in [7, 14, 30]:
    df[f'price_mean_{window}d'] = df.groupby('user_id')['price'].transform(
        lambda x: x.rolling(window=window, min_periods=1).mean()
    )
    df[f'price_std_{window}d'] = df.groupby('user_id')['price'].transform(
        lambda x: x.rolling(window=window, min_periods=1).std()
    )
    df[f'price_max_{window}d'] = df.groupby('user_id')['price'].transform(
        lambda x: x.rolling(window=window, min_periods=1).max()
    )

# 差分特征
df['price_diff_1'] = df.groupby('user_id')['price'].diff(1)
df['price_pct_change'] = df.groupby('user_id')['price'].pct_change()

# 历史累计特征
df['cumulative_spend'] = df.groupby('user_id')['amount'].cumsum()
df['order_count'] = df.groupby('user_id').cumcount() + 1

from scipy import signal

# 趋势分解
from statsmodels.tsa.seasonal import seasonal_decompose

decomposition = seasonal_decompose(df.set_index('datetime')['sales'], 
                                   model='additive', period=7)
df['trend'] = decomposition.trend.values
df['seasonal'] = decomposition.seasonal.values

# 统计特征：偏度、峰度
df['sales_skew_30d'] = df.groupby('user_id')['sales'].transform(
    lambda x: x.rolling(30, min_periods=5).skew()
)

# 多项式特征
from sklearn.preprocessing import PolynomialFeatures
poly = PolynomialFeatures(degree=2, interaction_only=False, include_bias=False)
X_poly = poly.fit_transform(X[['age', 'income']])

# 比率特征
df['income_per_age'] = df['income'] / df['age']
df['spend_ratio'] = df['spend'] / df['income']

# 加减交互
df['income_minus_avg'] = df['income'] - df.groupby('city')['income'].transform('mean')

# 分组统计特征
df['avg_price_by_category'] = df.groupby('category')['price'].transform('mean')
df['price_diff_from_category_avg'] = df['price'] - df['avg_price_by_category']

# 分组排名
df['price_rank_in_category'] = df.groupby('category')['price'].rank(method='dense')

df['text_length'] = df['text'].str.len()
df['word_count'] = df['text'].str.split().str.len()
df['avg_word_length'] = df['text'].apply(lambda x: np.mean([len(w) for w in x.split()]))
df['exclamation_count'] = df['text'].str.count('!')
df['question_count'] = df['text'].str.count('\?')
df['uppercase_ratio'] = df['text'].apply(lambda x: sum(1 for c in x if c.isupper()) / len(x))

from sklearn.feature_extraction.text import TfidfVectorizer, CountVectorizer

# TF-IDF
tfidf = TfidfVectorizer(
    max_features=1000,
    ngram_range=(1, 2),
    stop_words='english',
    min_df=2,
    max_df=0.8
)
tfidf_matrix = tfidf.fit_transform(df['text'])

# 主题特征 (LDA)
from sklearn.decomposition import LatentDirichletAllocation
lda = LatentDirichletAllocation(n_components=10, random_state=42)
topic_features = lda.fit_transform(tfidf_matrix)

# 经纬度距离
def haversine_distance(lat1, lon1, lat2, lon2):
    R = 6371  # 地球半径 km
    lat1, lon1, lat2, lon2 = map(np.radians, [lat1, lon1, lat2, lon2])
    dlat = lat2 - lat1
    dlon = lon2 - lon1
    a = np.sin(dlat/2)**2 + np.cos(lat1) * np.cos(lat2) * np.sin(dlon/2)**2
    return 2 * R * np.arcsin(np.sqrt(a))

df['distance_to_center'] = haversine_distance(
    df['lat'], df['lon'], 39.9042, 116.4074  # 北京中心
)

# 网格编码
df['lat_grid'] = (df['lat'] * 100).astype(int)
df['lon_grid'] = (df['lon'] * 100).astype(int)
df['grid_id'] = df['lat_grid'].astype(str) + '_' + df['lon_grid'].astype(str)

from sklearn.feature_selection import SelectKBest, f_classif, mutual_info_classif, chi2

# 方差分析 F 值
selector = SelectKBest(score_func=f_classif, k=20)
X_selected = selector.fit_transform(X, y)

# 互信息
selector = SelectKBest(score_func=mutual_info_classif, k=20)
X_selected = selector.fit_transform(X, y)

# 卡方检验（非负特征）
selector = SelectKBest(score_func=chi2, k=20)
X_selected = selector.fit_transform(X, y)

# 查看选中特征
selected_features = X.columns[selector.get_support()]
print(f"选中特征: {selected_features.tolist()}")

from sklearn.feature_selection import RFE, RFECV
from sklearn.ensemble import RandomForestClassifier

# 递归特征消除
estimator = RandomForestClassifier(n_estimators=100, random_state=42)
rfe = RFE(estimator, n_features_to_select=15, step=1)
X_rfe = rfe.fit_transform(X, y)

# 带交叉验证的 RFE
rfecv = RFECV(estimator, step=1, cv=5, scoring='f1_macro')
X_rfecv = rfecv.fit_transform(X, y)
print(f"最优特征数: {rfecv.n_features_}")

from sklearn.linear_model import Lasso, Ridge
from sklearn.ensemble import RandomForestClassifier

# Lasso 自动特征选择
lasso = Lasso(alpha=0.01)
lasso.fit(X, y)
selected = X.columns[np.abs(lasso.coef_) > 0]

# 随机森林特征重要性
rf = RandomForestClassifier(n_estimators=200, random_state=42)
rf.fit(X, y)
importances = pd.DataFrame({
    'feature': X.columns,
    'importance': rf.feature_importances_
}).sort_values('importance', ascending=False)

from sklearn.linear_model import RandomizedLasso

rlasso = RandomizedLasso(alpha=0.025)
rlasso.fit(X, y)
selected = X.columns[rlasso.scores_ > 0.5]

import featuretools as ft

# 定义实体集
es = ft.EntitySet(id='customer_data')

# 添加实体
es = es.add_dataframe(
    dataframe_name='transactions',
    dataframe=transactions_df,
    index='transaction_id',
    time_index='timestamp'
)

# 添加用户实体
es = es.add_dataframe(
    dataframe_name='users',
    dataframe=users_df,
    index='user_id'
)

# 定义关系
relationship = ft.Relationship(
    es["users"]["user_id"],
    es["transactions"]["user_id"]
)
es = es.add_relationship(relationship)

# 自动特征合成
feature_matrix, feature_defs = ft.dfs(
    entityset=es,
    target_dataframe_name='users',
    agg_primitives=['mean', 'sum', 'max', 'min', 'std', 'count', 'trend'],
    trans_primitives=['day', 'month', 'diff', 'absolute'],
    max_depth=2
)

print(f"生成 {len(feature_defs)} 个特征")

from autofeat import AutoFeatRegressor

model = AutoFeatRegressor(
    feateng_steps=2,
    featsel_runs=3,
    max_gb=16
)
X_new = model.fit_transform(X, y)

# Feast：特征存储
from feast import Entity, Feature, FeatureView, ValueType
from feast import FileSource

# 定义实体
user = Entity(name="user_id", value_type=ValueType.INT64)

# 定义特征视图
user_stats_view = FeatureView(
    name="user_stats",
    entities=["user_id"],
    ttl=timedelta(days=1),
    features=[
        Feature(name="avg_purchase_amount", dtype=ValueType.FLOAT),
        Feature(name="total_orders", dtype=ValueType.INT64),
        Feature(name="days_since_last_order", dtype=ValueType.INT64),
    ],
    online=True,
    source=FileSource(
        path="data/user_stats.parquet",
        event_timestamp_column="timestamp"
    ),
)

# 在线获取特征（低延迟）
from feast import FeatureStore
store = FeatureStore(repo_path=".")
features = store.get_online_features(
    features=["user_stats:avg_purchase_amount", "user_stats:total_orders"],
    entity_rows=[{"user_id": 123}, {"user_id": 456}]
).to_dict()