14 P4 프로젝트: 와인 품질 예측

import sys
import sklearn
from packaging.version import Version

assert sys.version_info >= (3, 10)
assert Version(sklearn.__version__) >= Version("1.6.1")

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

RANDOM_STATE = 42
plt.rc("font", size=12)
plt.rc("axes", labelsize=13, titlesize=13)
plt.rc("legend", fontsize=11)
plt.rc("xtick", labelsize=10)
plt.rc("ytick", labelsize=10)

print("scikit-learn:", sklearn.__version__)

from pathlib import Path
from urllib.request import urlopen
from zipfile import ZipFile
from io import BytesIO

csv_path = Path("winequality-white.csv")

if csv_path.exists():
    df_raw = pd.read_csv(csv_path, sep=";")
else:
    uci_url = "https://archive.ics.uci.edu/static/public/186/wine+quality.zip"
    try:
        with urlopen(uci_url, timeout=30) as response:
            content = response.read()
        with ZipFile(BytesIO(content)) as archive:
            csv_member = next(name for name in archive.namelist()
                              if name.endswith("winequality-white.csv"))
            with archive.open(csv_member) as file:
                df_raw = pd.read_csv(file, sep=";")
    except Exception as exc:
        raise RuntimeError(
            "UCI 데이터 다운로드에 실패했습니다. Wine Quality 데이터셋의 "
            "winequality-white.csv를 내려받아 노트북과 같은 폴더에 놓아주세요."
        ) from exc

print("원본 데이터 크기:", df_raw.shape)
df_raw.head()

expected_features = [
    "fixed acidity", "volatile acidity", "citric acid", "residual sugar",
    "chlorides", "free sulfur dioxide", "total sulfur dioxide",
    "density", "pH", "sulphates", "alcohol"
]
assert set(expected_features + ["quality"]) == set(df_raw.columns)
print("열별 자료형:\n", df_raw.dtypes)
print("\n결측값 개수:\n", df_raw.isna().sum())
print("\n완전히 동일한 행의 중복 개수:", df_raw.duplicated().sum())

df = df_raw.drop_duplicates().reset_index(drop=True)
print("중복 제거 전:", len(df_raw))
print("중복 제거 후:", len(df))
print("제거한 행:", len(df_raw) - len(df))

assert df[expected_features + ["quality"]].notna().all().all()

from sklearn.model_selection import train_test_split

binary_label = (df["quality"] >= 6).astype(int)
dev_df, test_df = train_test_split(
    df, test_size=0.20, random_state=RANDOM_STATE,
    stratify=binary_label
)
train_df, val_df = train_test_split(
    dev_df, test_size=0.25, random_state=RANDOM_STATE,
    stratify=(dev_df["quality"] >= 6).astype(int)
)
print("훈련:", train_df.shape, "검증:", val_df.shape,
      "최종 테스트:", test_df.shape)
print("※ 최종 테스트셋의 점수 분포·성능은 선택을 마칠 때까지 확인하지 않는다.")

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
train_df["quality"].value_counts().sort_index().plot.bar(ax=axes[0])
axes[0].set_title("Quality scores: training data")
axes[0].set_xlabel("Quality score")
axes[0].set_ylabel("Count")
(train_df["quality"] >= 6).value_counts(normalize=True).sort_index().plot.bar(ax=axes[1])
axes[1].set_title("Binary class ratio: training data")
axes[1].set_xlabel("0: < 6   /   1: >= 6")
axes[1].set_ylabel("Proportion")
plt.tight_layout()
plt.show()

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
train_df.boxplot(column="alcohol", by="quality", ax=axes[0])
train_df.boxplot(column="volatile acidity", by="quality", ax=axes[1])
axes[0].set_title("Alcohol by quality")
axes[1].set_title("Volatile acidity by quality")
plt.suptitle("")
plt.tight_layout()
plt.show()

# TODO: 위에서 정한 EDA 질문을 확인하는 코드를 작성한다.
# 예: train_df.groupby(...)[...].mean()

feature_cols = expected_features

X_train = train_df[feature_cols]
X_val = val_df[feature_cols]
y_train_reg = train_df["quality"]
y_val_reg = val_df["quality"]
y_train_clf = (y_train_reg >= 6).astype(int)
y_val_clf = (y_val_reg >= 6).astype(int)

assert "quality" not in X_train.columns
assert X_train.columns.tolist() == X_val.columns.tolist()
print("입력 특성 수:", X_train.shape[1])
print("훈련/검증 샘플 수:", len(X_train), len(X_val))

from sklearn.dummy import DummyRegressor
from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import (
    BaggingRegressor, RandomForestRegressor, GradientBoostingRegressor
)
from sklearn.metrics import mean_absolute_error, root_mean_squared_error, r2_score
from sklearn.base import clone

try:
    from xgboost import XGBRegressor, XGBClassifier
    HAS_XGBOOST = True
except ImportError:
    HAS_XGBOOST = False
    print("XGBoost 미설치: 해당 모델을 건너뜁니다. 필요하면 %pip install xgboost 실행 후 다시 시작하세요.")

reg_models = {
    "Dummy": DummyRegressor(strategy="mean"),
    "Tree (depth=2)": DecisionTreeRegressor(max_depth=2, random_state=RANDOM_STATE),
    "Tree (depth=3)": DecisionTreeRegressor(max_depth=3, random_state=RANDOM_STATE),
    "Tree (unlimited)": DecisionTreeRegressor(random_state=RANDOM_STATE),
    "Bagging": BaggingRegressor(
        estimator=DecisionTreeRegressor(random_state=RANDOM_STATE),
        n_estimators=100, max_samples=0.8, bootstrap=True,
        oob_score=True, random_state=RANDOM_STATE, n_jobs=-1
    ),
    "Random Forest": RandomForestRegressor(
        n_estimators=100, min_samples_leaf=2, oob_score=True,
        random_state=RANDOM_STATE, n_jobs=-1
    ),
    "Gradient Boosting": GradientBoostingRegressor(
        n_estimators=100, learning_rate=0.1, max_depth=2,
        random_state=RANDOM_STATE
    ),
    "Gradient Boosting (slow)": GradientBoostingRegressor(
        n_estimators=150, learning_rate=0.05, max_depth=2,
        random_state=RANDOM_STATE
    )
}
if HAS_XGBOOST:
    reg_models["XGBoost"] = XGBRegressor(
        n_estimators=100, learning_rate=0.1, max_depth=3,
        objective="reg:squarederror", tree_method="hist",
        random_state=RANDOM_STATE, n_jobs=2
    )
print("비교할 회귀 모델:", list(reg_models))

reg_rows = []
for name, model in reg_models.items():
    model.fit(X_train, y_train_reg)
    pred_train = model.predict(X_train)
    pred_val = model.predict(X_val)
    reg_rows.append({
        "모델": name,
        "훈련 RMSE": root_mean_squared_error(y_train_reg, pred_train),
        "검증 RMSE": root_mean_squared_error(y_val_reg, pred_val),
        "검증 MAE": mean_absolute_error(y_val_reg, pred_val),
        "검증 R2": r2_score(y_val_reg, pred_val)
    })
reg_results = pd.DataFrame(reg_rows).set_index("모델").sort_values("검증 RMSE")
reg_results.round(3)

ax = reg_results[["훈련 RMSE", "검증 RMSE"]].sort_values("검증 RMSE").plot.bar(figsize=(11, 4))
ax.set_ylabel("RMSE (quality score)")
ax.set_title("Regression: train vs validation RMSE")
ax.set_xlabel("Model")
plt.xticks(rotation=35, ha="right")
plt.tight_layout()
plt.show()

print("Bagging OOB R2:", round(reg_models["Bagging"].oob_score_, 3))
print("Random Forest OOB R2:", round(reg_models["Random Forest"].oob_score_, 3))
print("검증 R2와 OOB R2는 다른 평가 방식의 결과입니다.")

rf_reg = reg_models["Random Forest"]
rf_importance = pd.Series(rf_reg.feature_importances_, index=feature_cols).sort_values(ascending=False)
rf_importance.sort_values().plot.barh(figsize=(9, 5))
plt.xlabel("Impurity-based feature importance")
plt.title("Random Forest Regressor")
plt.tight_layout()
plt.show()
rf_importance.head(5).round(3)

boost_candidates = [
    {"n_estimators": 50, "learning_rate": 0.1},
    {"n_estimators": 100, "learning_rate": 0.1},
    {"n_estimators": 150, "learning_rate": 0.05},
]
boost_rows = []
for params in boost_candidates:
    model = GradientBoostingRegressor(max_depth=2, random_state=RANDOM_STATE, **params)
    model.fit(X_train, y_train_reg)
    boost_rows.append({
        **params,
        "훈련 RMSE": root_mean_squared_error(y_train_reg, model.predict(X_train)),
        "검증 RMSE": root_mean_squared_error(y_val_reg, model.predict(X_val))
    })
boost_results = pd.DataFrame(boost_rows).sort_values("검증 RMSE")
boost_results.round(3)

from sklearn.dummy import DummyClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.ensemble import BaggingClassifier, RandomForestClassifier, GradientBoostingClassifier
from sklearn.metrics import (
    accuracy_score, precision_score, recall_score, f1_score,
    roc_auc_score, ConfusionMatrixDisplay
)

clf_models = {
    "Dummy": DummyClassifier(strategy="most_frequent"),
    "Tree (depth=2)": DecisionTreeClassifier(max_depth=2, random_state=RANDOM_STATE),
    "Tree (depth=3)": DecisionTreeClassifier(max_depth=3, random_state=RANDOM_STATE),
    "Tree (unlimited)": DecisionTreeClassifier(random_state=RANDOM_STATE),
    "Bagging": BaggingClassifier(
        estimator=DecisionTreeClassifier(random_state=RANDOM_STATE),
        n_estimators=100, max_samples=0.8, max_features=1.0,
        bootstrap=True, oob_score=True,
        random_state=RANDOM_STATE, n_jobs=-1
    ),
    "Random Forest": RandomForestClassifier(
        n_estimators=100, min_samples_leaf=2,
        random_state=RANDOM_STATE, n_jobs=-1
    ),
    "Gradient Boosting": GradientBoostingClassifier(
        n_estimators=100, learning_rate=0.1, max_depth=2,
        random_state=RANDOM_STATE
    )
}
if HAS_XGBOOST:
    clf_models["XGBoost"] = XGBClassifier(
        n_estimators=100, learning_rate=0.1, max_depth=3,
        objective="binary:logistic", eval_metric="logloss",
        tree_method="hist", random_state=RANDOM_STATE, n_jobs=2
    )
print("비교할 분류 모델:", list(clf_models))

clf_rows = []
for name, model in clf_models.items():
    model.fit(X_train, y_train_clf)
    train_pred = model.predict(X_train)
    val_pred = model.predict(X_val)
    val_proba = model.predict_proba(X_val)[:, 1]
    clf_rows.append({
        "모델": name,
        "훈련 Accuracy": accuracy_score(y_train_clf, train_pred),
        "검증 Accuracy": accuracy_score(y_val_clf, val_pred),
        "검증 Precision": precision_score(y_val_clf, val_pred, zero_division=0),
        "검증 Recall": recall_score(y_val_clf, val_pred, zero_division=0),
        "검증 F1": f1_score(y_val_clf, val_pred, zero_division=0),
        "검증 ROC-AUC": roc_auc_score(y_val_clf, val_proba)
    })
clf_results = pd.DataFrame(clf_rows).set_index("모델").sort_values("검증 ROC-AUC", ascending=False)
clf_results.round(3)

ax = clf_results[["훈련 Accuracy", "검증 Accuracy"]].sort_values("검증 Accuracy").plot.bar(figsize=(11, 4))
ax.set_ylim(0, 1.05)
ax.set_ylabel("Accuracy")
ax.set_title("Classification: train vs validation accuracy")
ax.set_xlabel("Model")
plt.xticks(rotation=35, ha="right")
plt.tight_layout()
plt.show()

suggested_reg_name = reg_results["검증 RMSE"].idxmin()
suggested_clf_name = clf_results["검증 ROC-AUC"].idxmax()

chosen_reg_name = suggested_reg_name  # TODO: 다른 모델을 선택하면 직접 수정
chosen_clf_name = suggested_clf_name  # TODO: 다른 모델을 선택하면 직접 수정

assert chosen_reg_name in reg_models
assert chosen_clf_name in clf_models
print("최종 회귀 모델:", chosen_reg_name)
print("최종 분류 모델:", chosen_clf_name)

X_dev = dev_df[feature_cols]
y_dev_reg = dev_df["quality"]
y_dev_clf = (y_dev_reg >= 6).astype(int)

X_test = test_df[feature_cols]
y_test_reg = test_df["quality"]
y_test_clf = (y_test_reg >= 6).astype(int)

final_reg_model = clone(reg_models[chosen_reg_name])
final_clf_model = clone(clf_models[chosen_clf_name])
final_reg_model.fit(X_dev, y_dev_reg)
final_clf_model.fit(X_dev, y_dev_clf)

reg_test_pred = final_reg_model.predict(X_test)
clf_test_pred = final_clf_model.predict(X_test)
clf_test_proba = final_clf_model.predict_proba(X_test)[:, 1]

reg_test_scores = pd.Series({
    "RMSE": root_mean_squared_error(y_test_reg, reg_test_pred),
    "MAE": mean_absolute_error(y_test_reg, reg_test_pred),
    "R2": r2_score(y_test_reg, reg_test_pred)
}, name="회귀 테스트 성능")
clf_test_scores = pd.Series({
    "Accuracy": accuracy_score(y_test_clf, clf_test_pred),
    "Precision": precision_score(y_test_clf, clf_test_pred, zero_division=0),
    "Recall": recall_score(y_test_clf, clf_test_pred, zero_division=0),
    "F1": f1_score(y_test_clf, clf_test_pred, zero_division=0),
    "ROC-AUC": roc_auc_score(y_test_clf, clf_test_proba)
}, name="분류 테스트 성능")

print("[회귀]", chosen_reg_name)
print(reg_test_scores.round(3))
print("\n[분류]", chosen_clf_name)
print(clf_test_scores.round(3))

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
axes[0].scatter(y_test_reg, reg_test_pred, alpha=0.25)
lo = min(float(y_test_reg.min()), float(reg_test_pred.min()))
hi = max(float(y_test_reg.max()), float(reg_test_pred.max()))
axes[0].plot([lo, hi], [lo, hi], "k--")
axes[0].set_xlabel("Actual quality")
axes[0].set_ylabel("Predicted quality")
axes[0].set_title("Regression: actual vs predicted")

reg_errors = pd.DataFrame({
    "quality": y_test_reg.to_numpy(),
    "absolute_error": np.abs(y_test_reg.to_numpy() - reg_test_pred)
})
reg_errors.groupby("quality")["absolute_error"].mean().plot.bar(ax=axes[1])
axes[1].set_title("Mean absolute error by quality")
axes[1].set_ylabel("MAE")
plt.tight_layout()
plt.show()

ConfusionMatrixDisplay.from_predictions(
    y_test_clf, clf_test_pred, display_labels=["< 6", ">= 6"],
    cmap="Blues", values_format="d"
)
plt.title("Binary quality classification: test set")
plt.show()

test_examples = X_test.copy()
test_examples["quality"] = y_test_reg.to_numpy()
test_examples["reg_pred"] = reg_test_pred
test_examples["reg_abs_error"] = np.abs(y_test_reg.to_numpy() - reg_test_pred)
test_examples["clf_proba_ge6"] = clf_test_proba
test_examples["clf_pred_ge6"] = clf_test_pred

test_examples.sort_values("reg_abs_error", ascending=False).head(10)

# TODO: 선택 활동용 코드
# 예: 학습용 데이터에서 min_samples_leaf를 변경하고 검증 RMSE를 비교한다.