DOI: 10.5281/zenodo.18478615
title: Robust Design of Polypropylene Waste Modified Concrete Using Physics Guided Machine Learning

description: # =============================================================================
# TITLE: Thermo-Mechanical Performance of SWPP-Modified Concrete Using ML
# AUTHOR: 
# LICENSE: MIT / CC-BY 4.0
# DESCRIPTION: 
# This notebook performs data preprocessing, piecewise-aware augmentation,
# and machine learning modeling (XGBoost) as described in the paper.
# =============================================================================

import os
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
from sklearn.metrics import r2_score, mean_squared_error
from xgboost import XGBRegressor

# ---------------------------------------------------------
# CONFIGURATION & PATHS (Portable)
# ---------------------------------------------------------
# Paths are set relative to the 'src' folder to work on any machine
DATA_PATH = os.path.join("..", "data", "Data.xlsx")
OUT_DIR   = os.path.join("..", "outputs")

# Create output directory if it doesn't exist
os.makedirs(OUT_DIR, exist_ok=True)

print(f"Data Source: {DATA_PATH}")
print(f"Output Dir : {OUT_DIR}")

# =============================================================================
# STEP 1: DATA LOADING & PREPROCESSING (Wide -> Long)
# =============================================================================

# 1.1 Load Data
if not os.path.exists(DATA_PATH):
    raise FileNotFoundError(f"Data file not found at {DATA_PATH}. Please check folder structure.")

df_raw = pd.read_excel(DATA_PATH)
print("Data Loaded. Shape:", df_raw.shape)

# 1.2 Column Mapping
COL_SAMPLE = "Sample no"
COL_SWPP   = "SWPP (Shredded woven polypropylene)"

# Helper to clean numeric columns (handles comma decimals if any)
def numify(x):
    s = pd.Series(x) if not isinstance(x, pd.Series) else x
    s = s.astype(str).str.replace(",", ".", regex=False).str.strip()
    return pd.to_numeric(s, errors="coerce")

# Clean SWPP ratio column
df_raw["SWPP_ratio_pct"] = df_raw[COL_SWPP].astype(str).str.replace("%", "").str.replace(",", ".").astype(float)

# 1.3 Transform WIDE to LONG (Tidy) format
# Mapping: Temp -> (CS col, FS col, STS col)
temp_map = [
    (24,  "Compressive Strength-28 day, 24 C",       "Flexural Strength-28 day, 24 C", "Splitting Tensile Strength-28 day, 24 C"),
    (200, "Compressive Strength, after 200 C",       "Flexural Strength, after 200 C", "Splitting Tensile Strength, after 200 C"),
    (300, "Compressive Strength, after 300 C",       "Felxural Strength, after 300 C", "Splitting Tensile Strength, after 300 C"),
    (400, "Compressive Strength, after 400 C",       "Felxural Strength, after 400 C", "Splitting Tensile Strength, after 400 C"),
    (600, "Compressive Strength, after 600 C",       "Flexural Strength, after 600 C", "Splitting Tensile Strength, after 600 C"),
]

records = []
for _, r in df_raw.iterrows():
    for T, c_cs, c_fs, c_sts in temp_map:
        if c_cs in df_raw.columns:
            records.append({
                "Sample": str(r[COL_SAMPLE]).strip(),
                "SWPP_ratio_pct": r["SWPP_ratio_pct"],
                "Temperature_C": T,
                "CS_MPa":  float(numify(r[c_cs]).iloc[0]),
                "FS_MPa":  float(numify(r[c_fs]).iloc[0]),
                "STS_MPa": float(numify(r[c_sts]).iloc[0]),
            })

df = pd.DataFrame(records)

# 1.4 Define Exposure Regimes (Bins)
bins   = [-np.inf, 200, 300, 400, np.inf]
labels = ["RT-200", "200-300", "300-400", "400-600"]
df["Exposure_level"] = pd.cut(df["Temperature_C"], bins=bins, labels=labels, right=False)

# Save cleaned data
cleaned_path = os.path.join(OUT_DIR, "cleaned_data.xlsx")
df.to_excel(cleaned_path, index=False)
print(f"Data cleaned and saved to {cleaned_path}")

# =============================================================================
# STEP 2: LEAKAGE-SAFE SPLITTING & AUGMENTATION
# =============================================================================

# 2.1 Split by Sample to prevent data leakage
train_samples, test_samples = train_test_split(df["Sample"].unique(), test_size=0.30, random_state=42)

train_df = df[df["Sample"].isin(train_samples)].copy()
test_df  = df[df["Sample"].isin(test_samples)].copy()

print(f"Train samples: {len(train_samples)} | Test samples: {len(test_samples)}")

# 2.2 Piecewise-Aware Augmentation Function
rng = np.random.default_rng(42)

def piecewise_interpolate_regime(df_regime, n_synth, jitter=True):
    if len(df_regime) < 2: return pd.DataFrame()
    
    synth_rows = []
    X_swpp = df_regime["SWPP_ratio_pct"].values
    X_temp = df_regime["Temperature_C"].values
    Y_mat  = df_regime[["CS_MPa", "FS_MPa", "STS_MPa"]].values
    
    for _ in range(n_synth):
        idx = rng.choice(len(df_regime), 2, replace=False)
        lam = rng.uniform(0.05, 0.95)
        
        new_swpp = lam * X_swpp[idx[0]] + (1 - lam) * X_swpp[idx[1]]
        new_temp = lam * X_temp[idx[0]] + (1 - lam) * X_temp[idx[1]]
        new_y    = lam * Y_mat[idx[0]]  + (1 - lam) * Y_mat[idx[1]]
        
        if jitter:
            new_swpp += rng.normal(0, 0.5) 
            new_y    *= rng.normal(1.0, 0.02, size=3) 
            
        row = df_regime.iloc[idx[0]].copy()
        row["Sample"] = f"SYN_{int(rng.random()*10000)}"
        row["SWPP_ratio_pct"] = np.clip(new_swpp, 0, 100)
        row["Temperature_C"] = new_temp
        row[["CS_MPa", "FS_MPa", "STS_MPa"]] = np.maximum(new_y, 0)
        row["is_synthetic"] = True
        synth_rows.append(row)
        
    return pd.DataFrame(synth_rows)

# 2.3 Apply Augmentation
aug_list = [train_df.assign(is_synthetic=False)]
MULTIPLIER = 2

for regime in labels:
    regime_data = train_df[train_df["Exposure_level"] == regime]
    n_synth = int(len(regime_data) * MULTIPLIER)
    if n_synth > 0:
        synth_df = piecewise_interpolate_regime(regime_data, n_synth)
        aug_list.append(synth_df)

train_aug = pd.concat(aug_list, ignore_index=True)
aug_path = os.path.join(OUT_DIR, "train_augmented.xlsx")
train_aug.to_excel(aug_path, index=False)
print(f"Augmentation complete. Saved to {aug_path}")

# =============================================================================
# STEP 3: ML MODELING (XGBoost)
# =============================================================================

features = ["SWPP_ratio_pct", "Temperature_C"]
targets  = ["CS_MPa", "FS_MPa", "STS_MPa"]

model = XGBRegressor(n_estimators=300, max_depth=3, learning_rate=0.05, random_state=42)

print("\n--- Model Performance on Test Set ---")
for target in targets:
    model.fit(train_aug[features], train_aug[target])
    y_pred = model.predict(test_df[features])
    y_true = test_df[target]
    
    r2 = r2_score(y_true, y_pred)
    rmse = mean_squared_error(y_true, y_pred, squared=False)
    
    print(f"{target} -> R2: {r2:.4f} | RMSE: {rmse:.4f}")
    
    # Save Plot
    plt.figure(figsize=(5,5))
    plt.scatter(y_true, y_pred, alpha=0.7)
    plt.plot([y_true.min(), y_true.max()], [y_true.min(), y_true.max()], 'r--')
    plt.xlabel("Experimental")
    plt.ylabel("Predicted")
    plt.title(f"{target} Prediction")
    plt.tight_layout()
    plt.savefig(os.path.join(OUT_DIR, f"Parity_{target}.png"), dpi=300)
    plt.close()

print(f"\nAll outputs saved to: {OUT_DIR}")