Capstone: Comparative Machine Learning Approaches for Predicting Scale Score Averages¶
In [ ]:
Training ML¶
Important: I added comments within the code and/or before/after the code for added explanation.¶
In [1]:
# Import packages
import pandas as pd
import numpy as np
import os
import joblib
import seaborn as sns
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split, GridSearchCV, cross_val_score, RandomizedSearchCV
from sklearn.linear_model import Ridge
from sklearn.linear_model import Lasso
from sklearn.metrics import r2_score, mean_squared_error, mean_absolute_error
from scipy.stats import pearsonr
from sklearn.tree import DecisionTreeRegressor
from sklearn.preprocessing import StandardScaler
from xgboost import XGBRegressor
# Suppress all warnings
import warnings
warnings.filterwarnings('ignore')
# producing visualizations in notebook
%matplotlib inline
# Ensure all columns are displayed
pd.set_option('display.max_columns', None) # Show all columns
pd.set_option('display.width', None) # Adjust display width to avoid wrapping
pd.set_option('display.max_colwidth', None) # Show full content of each column (if needed)
# Function to print bold headers
def print_bold(text):
print(f"\033[1m{text}\033[0m")
In [ ]:
Data¶
The datasets used in this study is taken from Delaware Open Data Council. Two datasets are used in this study, ‘student assessment performance’ dataset and ‘student discipline’ dataset.
Student Assessment Performance dataset (assessment.csv): https://data.delaware.gov/Education/Student-Assessment-Performance/ms6b-mt82/about_data Student Discipline (discipline.csv): https://data.delaware.gov/Education/Student-Discipline/yr4w-jdi4/about_data
In [ ]:
Research Question¶
How do traditional regression models compare to tree-based models in predicting school-level Scale Score Averages based on predictive accuracy and model performance metrics?
In [ ]:
Hypothesis¶
Null Hypothesis (H₀): There is no significant difference in predictive accuracy between tree-based models and traditional regression models in predicting school-level Scale Score Averages.
Alternative Hypothesis (H₁): Tree-based models demonstrate significantly higher predictive accuracy than traditional regression models in predicting school-level Scale Score Averages.
In [ ]:
Data storing: All datasets used here are saved in my local data store for a seamless process.¶
In [ ]:
Split Data into Training and Testing Sets¶
In [2]:
def split_data(filepath='4_ready_training.csv', test_size=0.2, random_state=42):
"""
Load the dataset, convert specific columns to categorical, dynamically select predictor columns,
and split the data into training and testing sets.
Parameters:
filepath (str): Path to the CSV file.
test_size (float): Proportion of the dataset to include in the test split.
random_state (int): Random seed for reproducibility.
Returns:
X_train (DataFrame): Training predictors.
X_test (DataFrame): Testing predictors.
y_train (Series): Training target variable.
y_test (Series): Testing target variable.
"""
# Load the dataset
df = pd.read_csv(filepath)
# Convert relevant columns to 'category' data type
df['School Year'] = df['School Year'].astype('category')
df['Grade_Numeric'] = df['Grade_Numeric'].astype('category')
# Find relevant columns
grade_numeric_cols = [col for col in df.columns if col.startswith('Grade_Numeric_')]
grade_group_cols = [col for col in df.columns if col.startswith('Grade_Group_')]
category_cols = [col for col in df.columns if col.lower().startswith('category_')]
content_area_cols = [col for col in df.columns if col.startswith('ContentArea_')]
special_demo_cols = [col for col in df.columns if col.startswith('SpecialDemo_')]
district_code_cols = [col for col in df.columns if col.startswith('DistrictCode_')]
school_year_cols = [col for col in df.columns if col.startswith('SchoolYear_')]
assessment_name_cols = [col for col in df.columns if col.startswith('AssessmentName_')]
# Define Dependent and Independent Variables
y = df['ScaleScoreAvg']
# Independent Variables (Predictors)
independent_vars = (
['Incidents', 'AvgDuration', 'PctEnrollment'] +
grade_numeric_cols + grade_group_cols + category_cols +
content_area_cols + special_demo_cols + district_code_cols +
school_year_cols + assessment_name_cols
)
X = df[independent_vars]
# Split Data into Training and Testing Sets
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=test_size, random_state=random_state
)
print("\n\nShape of X_train:", X_train.shape)
print("Shape of X_test:", X_test.shape)
print("Process successfully!\n")
return X_train, X_test, y_train, y_test
In [ ]:
Ridge Regression Model¶
In [4]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [ ]:
In [19]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Fit a Ridge Regression Model
ridge_model = Ridge()
ridge_model.fit(X_train, y_train)
# Evaluate the Model
y_pred_ridge_nottuned = ridge_model.predict(X_test)
r2 = r2_score(y_test, y_pred_ridge_nottuned)
mse = mean_squared_error(y_test, y_pred_ridge_nottuned)
mae = mean_absolute_error(y_test, y_pred_ridge_nottuned)
correlation, p_value = pearsonr(y_test, y_pred_ridge_nottuned)
print_bold("\n\nRidge Regression Analysis Results")
print(f"R² Score: {r2:.6f}")
print(f"Mean Squared Error (MSE): {mse:.6f}")
print(f"Mean Absolute Error (MAE): {mae:.6f}\n")
print(f"Pearson Correlation: {correlation:.3f}, P-value: {p_value:.3f}")
# Calculate Residual Metrics
residuals = y_test - y_pred_ridge_nottuned
mean_residual = np.mean(residuals)
std_residual = np.std(residuals)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
# Interpretation
print("\nInterpretation: R² Score")
print("- R² reflects the proportion of variance in the target variable explained by the predictors.")
print(" R² Score > 0.7: Strong model.")
print(" 0.4 < R² Score ≤ 0.7: Moderate model.")
print(" R² Score ≤ 0.4: Weak model.\n\n")
# Cross-Validation Scores
cv_scores = cross_val_score(ridge_model, X_train, y_train, cv=3, scoring='r2')
print("Cross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Visualize Residuals
# (Residuals already computed above)
plt.figure(figsize=(8, 4))
sns.histplot(residuals, kde=True, color='blue', bins=30)
plt.title("Residual Distribution", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Feature Importances (Coefficients)
feature_importances = pd.Series(ridge_model.coef_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
print("\nFeature Importances (Ridge Regression Coefficients):")
print(feature_importances.head(10))
# Visualize Feature Importances
plt.figure(figsize=(10, 6))
sns.barplot(x=feature_importances.values, y=feature_importances.index, palette='viridis')
plt.title('Ridge Regression Coefficients', fontsize=16)
plt.xlabel('Coefficient Value', fontsize=12)
plt.ylabel('Feature', fontsize=12)
plt.tight_layout()
plt.show()
# Save model
model_filename = f"{model_dir}/ridge_reg.joblib"
joblib.dump(ridge_model, model_filename)
print(f"Model saved: {model_filename}")
Ridge Regression Analysis Results
R² Score: 0.292510
Mean Squared Error (MSE): 0.279764
Mean Absolute Error (MAE): 0.401849
Pearson Correlation: 0.541, P-value: 0.000
Residual Analysis:
- Mean of residuals: 0.003585
- Standard deviation of residuals: 0.528915
- The mean residual is close to 0, indicating a good fit.
Interpretation: R² Score
- R² reflects the proportion of variance in the target variable explained by the predictors.
R² Score > 0.7: Strong model.
0.4 < R² Score ≤ 0.7: Moderate model.
R² Score ≤ 0.4: Weak model.
Cross-Validation R² Scores: [0.29512694 0.29289429 0.2958661 ]
Mean CV R² Score: 0.29462911180547957
Feature Importances (Ridge Regression Coefficients): PctEnrollment 0.283685 AssessmentName_Smarter Balanced Summative Assessment 0.247759 DistrictCode_9606 0.132770 DistrictCode_9605 0.119697 DistrictCode_88 0.087946 DistrictCode_82 0.087017 DistrictCode_87 0.073664 DistrictCode_9612 0.071901 DistrictCode_38 0.070775 DistrictCode_85 0.064318 dtype: float64
Model saved: saved_models/ridge_reg.joblib
In [ ]:
Ridge Regression - Hyperparameter Tuning with RandomizedSearchCV¶
In [14]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [15]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Hyperparameter grid for Ridge Regression
param_grid = {
'alpha': [0.001, 0.01, 0.1, 1, 10, 100],
}
# GridSearchCV setup
grid_search_ridge = GridSearchCV(
Ridge(),
param_grid=param_grid,
cv=3,
scoring='r2',
n_jobs=-1,
verbose=1
)
# Fit GridSearchCV
grid_search_ridge.fit(X_train, y_train)
# Predict using best model
best_ridge = grid_search_ridge.best_estimator_
y_pred_ridge = best_ridge.predict(X_test)
# Evaluation Metrics
r2_ridge = r2_score(y_test, y_pred_ridge)
mse_ridge = mean_squared_error(y_test, y_pred_ridge)
mae_ridge = mean_absolute_error(y_test, y_pred_ridge)
correlation_ridge, p_value_ridge = pearsonr(y_test, y_pred_ridge)
print("\n=== Best Hyperparameters ===")
print(grid_search_ridge.best_params_)
print("\n=== Ridge Regression Results ===")
print(f"R² Score: {r2_ridge:.3f}")
print(f"MSE: {mse_ridge:.6f}")
print(f"MAE: {mae_ridge:.6f}")
print(f"Pearson Correlation: {correlation_ridge:.3f}, P-value: {p_value_ridge:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(best_ridge, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Visualize Residuals
residuals_ridge = y_test - y_pred_ridge
plt.figure(figsize=(8, 4))
sns.histplot(residuals_ridge, kde=True, bins=30, color='orange')
plt.title("Residual Distribution (Tuned Ridge Regression)", fontsize=16)
plt.xlabel("Residuals")
plt.ylabel("Frequency")
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_ridge)
std_residual = np.std(residuals_ridge)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
# Feature Importances (Coefficients)
feature_importances = pd.Series(best_ridge.coef_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
print("\nFeature Importances (Ridge Regression Coefficients):")
print(feature_importances.head(10))
# Visualize Feature Importances
plt.figure(figsize=(10, 6))
sns.barplot(x=feature_importances.values, y=feature_importances.index, palette='viridis')
plt.title('Ridge Regression Coefficients', fontsize=16)
plt.xlabel('Coefficient Value', fontsize=12)
plt.ylabel('Feature', fontsize=12)
plt.tight_layout()
plt.show()
# Save Model
model_filename = f"{model_dir}/ridge_model_tuned.joblib"
joblib.dump(best_ridge, model_filename)
print(f"Model saved: {model_filename}")
Fitting 3 folds for each of 6 candidates, totalling 18 fits
=== Best Hyperparameters ===
{'alpha': 1}
=== Ridge Regression Results ===
R² Score: 0.293
MSE: 0.279764
MAE: 0.401849
Pearson Correlation: 0.541, P-value: 0.000
Cross-Validation R² Scores: [0.29512694 0.29289429 0.2958661 ]
Mean CV R² Score: 0.29462911180547957
Residual Analysis: - Mean of residuals: 0.003585 - Standard deviation of residuals: 0.528915 - The mean residual is close to 0, indicating a good fit. Feature Importances (Ridge Regression Coefficients): PctEnrollment 0.283685 AssessmentName_Smarter Balanced Summative Assessment 0.247759 DistrictCode_9606 0.132770 DistrictCode_9605 0.119697 DistrictCode_88 0.087946 DistrictCode_82 0.087017 DistrictCode_87 0.073664 DistrictCode_9612 0.071901 DistrictCode_38 0.070775 DistrictCode_85 0.064318 dtype: float64
Model saved: saved_models/ridge_model_tuned.joblib
In [ ]:
Lasso Regression Model¶
In [16]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [17]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Initialize and fit the Lasso Regression model
lasso_model = Lasso(alpha=1.0, random_state=42)
lasso_model.fit(X_train, y_train)
# Predict using Lasso Regression
y_pred_lasso = lasso_model.predict(X_test)
# Evaluate Lasso Model
r2_lasso = r2_score(y_test, y_pred_lasso)
mse_lasso = mean_squared_error(y_test, y_pred_lasso)
mae_lasso = mean_absolute_error(y_test, y_pred_lasso)
correlation_lasso, p_value_lasso = pearsonr(y_test, y_pred_lasso)
print("\n\n=== Lasso Regression Results ===")
print(f"R² Score: {r2_lasso:.6f}")
print(f"Mean Squared Error (MSE): {mse_lasso:.6f}")
print(f"Mean Absolute Error (MAE): {mae_lasso:.6f}")
print(f"Pearson Correlation: {correlation_lasso:.3f}, P-value: {p_value_lasso:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(lasso_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Residual Analysis for Lasso Regression
residuals_lasso = y_test - y_pred_lasso
plt.figure(figsize=(8, 4))
sns.histplot(residuals_lasso, kde=True, color='green', bins=30)
plt.title("Residual Distribution (Lasso Regression)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_lasso)
std_residual = np.std(residuals_lasso)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
# Feature Importance from Lasso Regression
# For Lasso, the model coefficients represent feature importance.
coefficients = pd.Series(lasso_model.coef_, index=X_train.columns)
# Sort features by the absolute value of coefficients in descending order
coefficients = coefficients.loc[coefficients.abs().sort_values(ascending=False).index]
print("\nFeature Importances (Lasso Regression Coefficients):")
print(coefficients.head(10))
plt.figure(figsize=(10, 6))
coefficients.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (Lasso Regression)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Coefficient Value', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the Lasso Regression Model
model_filename = f"{model_dir}/lasso_regression_model.joblib"
joblib.dump(lasso_model, model_filename)
print(f"\nModel saved: {model_filename}")
=== Lasso Regression Results === R² Score: -0.000032 Mean Squared Error (MSE): 0.395445 Mean Absolute Error (MAE): 0.506775 Pearson Correlation: nan, P-value: nan Cross-Validation R² Scores: [-2.23037546e-06 -7.59538655e-07 -3.89627643e-07] Mean CV R² Score: -1.1265139201466212e-06
Residual Analysis: - Mean of residuals: 0.003556 - Standard deviation of residuals: 0.628834 - The mean residual is close to 0, indicating a good fit. Feature Importances (Lasso Regression Coefficients): Incidents -0.0 DistrictCode_72 -0.0 DistrictCode_87 -0.0 DistrictCode_86 -0.0 DistrictCode_85 -0.0 DistrictCode_82 -0.0 DistrictCode_80 -0.0 DistrictCode_77 -0.0 DistrictCode_76 -0.0 DistrictCode_74 -0.0 dtype: float64
Model saved: saved_models/lasso_regression_model.joblib
In [ ]:
Ridge Regression Model: Hyperparameter Tuning with RandomizedSearchCV¶
In [5]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [19]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Hyperparameter Tuning with RandomizedSearchCV
print_bold("\nHyperparameter Tuning with RandomizedSearchCV (Lasso Regression)")
# Define the parameter grid for Lasso Regression
param_dist = {
'alpha': [0.001, 0.01, 0.1, 1.0, 10.0],
'max_iter': [1000, 5000, 10000]
}
# Initialize the Lasso Regression model
lasso_model = Lasso(random_state=42)
# Initialize RandomizedSearchCV
random_search = RandomizedSearchCV(
estimator=lasso_model,
param_distributions=param_dist,
n_iter=10,
scoring='r2',
cv=3,
random_state=42,
n_jobs=-1
)
# Fit the random search model
random_search.fit(X_train, y_train)
# Get the best Lasso model from the random search
best_lasso_model = random_search.best_estimator_
# Predict using the best model
y_pred_lasso_tuned = best_lasso_model.predict(X_test)
# Evaluate the best model
r2_lasso = r2_score(y_test, y_pred_lasso_tuned)
mse_lasso = mean_squared_error(y_test, y_pred_lasso_tuned)
mae_lasso = mean_absolute_error(y_test, y_pred_lasso_tuned)
correlation_lasso, p_value_lasso = pearsonr(y_test, y_pred_lasso_tuned)
print_bold("\n\nLasso Regression Results After Hyperparameter Tuning")
print(f"R² Score: {r2_lasso:.3f}")
print(f"Mean Squared Error (MSE): {mse_lasso:.6f}")
print(f"Mean Absolute Error (MAE): {mae_lasso:.6f}")
print(f"Pearson Correlation: {correlation_lasso:.3f}, P-value: {p_value_lasso:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(best_lasso_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Display the best hyperparameters found by RandomizedSearchCV
print("\nBest Hyperparameters found by RandomizedSearchCV:")
print(random_search.best_params_)
# Feature Importance from the Best Lasso Model
# For Lasso, the model coefficients serve as feature importances.
coefficients = pd.Series(best_lasso_model.coef_, index=X_train.columns)
# Sort by the absolute value of the coefficients in descending order
coefficients = coefficients.loc[coefficients.abs().sort_values(ascending=False).index]
# Residual Analysis for Lasso Regression
residuals_lasso_tuned = y_test - y_pred_lasso_tuned
plt.figure(figsize=(8, 4))
sns.histplot(residuals_lasso_tuned, kde=True, color='green', bins=30)
plt.title("Residual Distribution (Lasso Regression After Tuning)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_lasso_tuned)
std_residual = np.std(residuals_lasso_tuned)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
print("\nFeature Importances (Lasso Regression Coefficients):")
print(coefficients.head(10))
plt.figure(figsize=(10, 6))
coefficients.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (Lasso Regression)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Coefficient Value', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the Best Lasso Regression Model
model_filename = f"{model_dir}/lasso_regression_model_tuned.joblib"
joblib.dump(best_lasso_model, model_filename)
print(f"\nBest model saved: {model_filename}")
Hyperparameter Tuning with RandomizedSearchCV (Lasso Regression) Lasso Regression Results After Hyperparameter Tuning R² Score: 0.289 Mean Squared Error (MSE): 0.281102 Mean Absolute Error (MAE): 0.403611 Pearson Correlation: 0.538, P-value: 0.000 Cross-Validation R² Scores: [0.2915587 0.28941561 0.29208666] Mean CV R² Score: 0.29102032480315687 Best Hyperparameters found by RandomizedSearchCV: {'max_iter': 1000, 'alpha': 0.001}
Residual Analysis: - Mean of residuals: 0.003544 - Standard deviation of residuals: 0.530179 - The mean residual is close to 0, indicating a good fit. Feature Importances (Lasso Regression Coefficients): Incidents -1.511145 PctEnrollment 0.281220 AssessmentName_Smarter Balanced Summative Assessment 0.236320 AssessmentName_SAT School-Day (Spring) -0.163387 AssessmentName_Delaware System of Student Assessment (DeSSA) -0.123425 ContentArea_Essay -0.076836 ContentArea_Social Studies -0.047034 SchoolYear_2019 0.036298 category_Expulsion (Permanent Removal from School) 0.035242 SchoolYear_2022 -0.027942 dtype: float64
Best model saved: saved_models/lasso_regression_model_tuned.joblib
In [ ]:
Decision Tree Regressor¶
In [20]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [6]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Initialize and fit the Decision Tree Regressor
dt_model = DecisionTreeRegressor(max_depth=5, random_state=42)
dt_model.fit(X_train, y_train)
# Predict using Decision Tree Regressor
y_pred_dt = dt_model.predict(X_test)
# Evaluate Decision Tree Model
r2_dt = r2_score(y_test, y_pred_dt)
mse_dt = mean_squared_error(y_test, y_pred_dt)
mae_dt = mean_absolute_error(y_test, y_pred_dt)
correlation_dt, p_value_dt = pearsonr(y_test, y_pred_dt)
print("\n\n=== Decision Tree Regression Results ===")
print(f"R² Score: {r2_dt:.6f}")
print(f"Mean Squared Error (MSE): {mse_dt:.6f}")
print(f"Mean Absolute Error (MAE): {mae_dt:.6f}")
print(f"Pearson Correlation: {correlation_dt:.3f}, P-value: {p_value_dt:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(dt_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Residual Analysis for Decision Tree
residuals_dt = y_test - y_pred_dt
plt.figure(figsize=(8, 4))
sns.histplot(residuals_dt, kde=True, color='green', bins=30)
plt.title("Residual Distribution (Decision Tree)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_dt)
std_residual = np.std(residuals_dt)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
# Feature Importance from Decision Tree
feature_importances = pd.Series(dt_model.feature_importances_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
print("\nFeature Importances (Decision Tree):")
print(feature_importances.head(10))
plt.figure(figsize=(10, 6))
feature_importances.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (Decision Tree)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Importance', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the Decision Tree Model
model_filename = f"{model_dir}/decision_tree_model.joblib"
joblib.dump(dt_model, model_filename)
print(f"\nModel saved: {model_filename}")
=== Decision Tree Regression Results === R² Score: 0.741542 Mean Squared Error (MSE): 0.102202 Mean Absolute Error (MAE): 0.177782 Pearson Correlation: 0.861, P-value: 0.000 Cross-Validation R² Scores: [0.73905386 0.72939229 0.73929941] Mean CV R² Score: 0.7359151879355927
Residual Analysis: - Mean of residuals: 0.000940 - Standard deviation of residuals: 0.319690 - The mean residual is close to 0, indicating a good fit. Feature Importances (Decision Tree): Incidents 0.394420 AvgDuration 0.309712 PctEnrollment 0.264387 AssessmentName_Smarter Balanced Summative Assessment 0.030266 category_Suspended 0.001051 ContentArea_Computer Science 0.000164 Grade_Group_Elementary 0.000000 DistrictCode_76 0.000000 DistrictCode_88 0.000000 DistrictCode_87 0.000000 dtype: float64
Model saved: saved_models/decision_tree_model.joblib
In [ ]:
Decision Tree Regressor: Hyperparameter Tuning with RandomizedSearchCV¶
In [18]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [12]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Hyperparameter Tuning with RandomizedSearchCV
print_bold("\nHyperparameter Tuning with RandomizedSearchCV (Decision Tree Regression)")
# Define the parameter grid
param_dist = {
'max_depth': [None, 3, 5, 7, 10],
'min_samples_split': [2, 5, 10],
'min_samples_leaf': [1, 2, 4],
'max_features': [None, 'auto', 'sqrt', 'log2']
}
# Initialize the Decision Tree Regressor
dt_model = DecisionTreeRegressor(random_state=42)
# Initialize RandomizedSearchCV
random_search = RandomizedSearchCV(
estimator=dt_model,
param_distributions=param_dist,
n_iter=10,
scoring='r2',
cv=3,
random_state=42,
n_jobs=-1
)
# Fit the random search model
random_search.fit(X_train, y_train)
# Get the best Decision Tree model from the random search
best_dt_model = random_search.best_estimator_
# Predict using the best model
y_pred_dt_tuned = best_dt_model.predict(X_test)
# Evaluate the best model
r2_dt = r2_score(y_test, y_pred_dt_tuned)
mse_dt = mean_squared_error(y_test, y_pred_dt_tuned)
mae_dt = mean_absolute_error(y_test, y_pred_dt_tuned)
correlation_dt, p_value_dt = pearsonr(y_test, y_pred_dt_tuned)
print_bold("\n\nDecision Tree Regression Results After Hyperparameter Tuning")
print(f"R² Score: {r2_dt:.3f}")
print(f"Mean Squared Error (MSE): {mse_dt:.6f}")
print(f"Mean Absolute Error (MAE): {mae_dt:.6f}")
print(f"Pearson Correlation: {correlation_dt:.3f}, P-value: {p_value_dt:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(best_dt_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Display the best hyperparameters found by RandomizedSearchCV
print("\nBest Hyperparameters found by RandomizedSearchCV:")
print(random_search.best_params_)
# Feature Importance from the Best Decision Tree Model
feature_importances = pd.Series(best_dt_model.feature_importances_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
# Residual Analysis for Decision Tree Regression
residuals_dt_tuned = y_test - y_pred_dt_tuned
plt.figure(figsize=(8, 4))
sns.histplot(residuals_dt_tuned, kde=True, color='green', bins=30)
plt.title("Residual Distribution (Decision Tree After Tuning)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_dt_tuned)
std_residual = np.std(residuals_dt_tuned)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
print("\nFeature Importances (Best Decision Tree Model):")
print(feature_importances.head(10))
plt.figure(figsize=(10, 6))
feature_importances.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (Best Decision Tree Model)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Importance', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the Best Decision Tree Model
model_filename = f"{model_dir}/decision_tree_model_tuned.joblib"
joblib.dump(best_dt_model, model_filename)
print(f"\nBest model saved: {model_filename}")
Hyperparameter Tuning with RandomizedSearchCV (Decision Tree Regression) Decision Tree Regression Results After Hyperparameter Tuning R² Score: 0.852 Mean Squared Error (MSE): 0.058364 Mean Absolute Error (MAE): 0.117435 Pearson Correlation: 0.923, P-value: 0.000 Cross-Validation R² Scores: [0.8477795 0.84962428 0.84856838] Mean CV R² Score: 0.8486573891425078 Best Hyperparameters found by RandomizedSearchCV: {'min_samples_split': 5, 'min_samples_leaf': 1, 'max_features': None, 'max_depth': 10}
Residual Analysis: - Mean of residuals: -0.000297 - Standard deviation of residuals: 0.241586 - The mean residual is close to 0, indicating a good fit. Feature Importances (Best Decision Tree Model): Incidents 0.385484 AvgDuration 0.312013 PctEnrollment 0.254819 AssessmentName_Smarter Balanced Summative Assessment 0.032142 ContentArea_Essay 0.003676 SchoolYear_2022 0.001944 category_Out-of-School Suspension_Merged 0.001755 AssessmentName_Delaware System of Student Assessment (DeSSA) 0.001464 AssessmentName_SAT School-Day (Spring) 0.001355 category_Suspended 0.000910 dtype: float64
Best model saved: saved_models/decision_tree_model_tuned.joblib
In [ ]:
XGBoost Regressor¶
In [11]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [56]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Initialize and fit the XGBoost model
xgb_model = XGBRegressor(n_estimators=100, learning_rate=0.1, max_depth=5, random_state=42)
xgb_model.fit(X_train, y_train)
# Predict using XGBoost
y_pred_xgb = xgb_model.predict(X_test)
# Evaluate XGBoost Model
r2_xgb = r2_score(y_test, y_pred_xgb)
mse_xgb = mean_squared_error(y_test, y_pred_xgb)
mae_xgb = mean_absolute_error(y_test, y_pred_xgb)
correlation_dt, p_value_dt = pearsonr(y_test, y_pred_xgb)
print("\n\n=== XGBoost Regression Results ===")
print(f"R² Score: {r2_xgb:.6f}")
print(f"Mean Squared Error (MSE): {mse_xgb:.6f}")
print(f"Mean Absolute Error (MAE): {mae_xgb:.6f}")
print(f"Pearson Correlation: {correlation_dt:.3f}, P-value: {p_value_dt:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(xgb_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Residual Analysis for XGBoost
residuals_xgb = y_test - y_pred_xgb
plt.figure(figsize=(8, 4))
sns.histplot(residuals_xgb, kde=True, color='green', bins=30)
plt.title("Residual Distribution (XGBoost)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_xgb)
std_residual = np.std(residuals_xgb)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
# Feature Importance from XGBoost
feature_importances = pd.Series(xgb_model.feature_importances_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
print("\nFeature Importances (XGBoost):")
print(feature_importances.head(10))
plt.figure(figsize=(10, 6))
feature_importances.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (XGBoost)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Importance', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the XGBoost Model
model_filename = f"{model_dir}/XGBoost_Regressor.joblib"
joblib.dump(xgb_model, model_filename)
print(f"\nModel saved: {model_filename}")
=== XGBoost Regression Results === R² Score: 0.842192 Mean Squared Error (MSE): 0.062402 Mean Absolute Error (MAE): 0.135715 Pearson Correlation: 0.918, P-value: 0.000 Cross-Validation R² Scores: [0.83951415 0.83915408 0.8385428 ] Mean CV R² Score: 0.8390703440787001
Residual Analysis: - Mean of residuals: 0.000586 - Standard deviation of residuals: 0.249804 - The mean residual is close to 0, indicating a good fit. Feature Importances (XGBoost): PctEnrollment 0.212727 Incidents 0.184378 AvgDuration 0.164079 AssessmentName_Smarter Balanced Summative Assessment 0.152896 category_Suspended 0.043352 ContentArea_Essay 0.030131 AssessmentName_Delaware System of Student Assessment (DeSSA) 0.029917 AssessmentName_SAT School-Day (Spring) 0.027834 AssessmentName_DeSSA Alternate Assessment 0.022732 SchoolYear_2022 0.014596 dtype: float32
Model saved: saved_models/XGBoost_Regressor.joblib
In [ ]:
XGBoost Regressor: Hyperparameter Tuning with RandomizedSearchCV¶
In [19]:
# Call function: split_data
X_train, X_test, y_train, y_test = split_data()
Shape of X_train: (519869, 84) Shape of X_test: (129968, 84) Process successfully!
In [16]:
# Create a directory for saved models if it doesn't exist
model_dir = "saved_models"
os.makedirs(model_dir, exist_ok=True)
# Hyperparameter Tuning with RandomizedSearchCV
print_bold("\nHyperparameter Tuning with RandomizedSearchCV (XGBoost)")
# Define the parameter grid for RandomizedSearchCV
param_dist = {
'n_estimators': [50, 100, 200],
'learning_rate': [0.01, 0.1, 0.2],
'max_depth': [3, 5, 7],
'subsample': [0.6, 0.8, 1.0],
'colsample_bytree': [0.6, 0.8, 1.0],
'gamma': [0, 0.1, 0.2],
'min_child_weight': [1, 3, 5]
}
# Initialize the XGBoost Regressor
xgb_model = XGBRegressor(random_state=42)
# Initialize RandomizedSearchCV with XGBoost Regressor
random_search = RandomizedSearchCV(
estimator=xgb_model,
param_distributions=param_dist,
n_iter=10,
scoring='r2',
cv=3,
random_state=42,
n_jobs=-1
)
# Fit the random search model
random_search.fit(X_train, y_train)
# Get the best XGBoost model from the random search
best_xgb_model = random_search.best_estimator_
# Predict using the best model
y_pred_xgb_tuned = best_xgb_model.predict(X_test)
# Evaluate the best model
r2_xgb = r2_score(y_test, y_pred_xgb_tuned)
mse_xgb = mean_squared_error(y_test, y_pred_xgb_tuned)
mae_xgb = mean_absolute_error(y_test, y_pred_xgb_tuned)
correlation_xgb, p_value_xgb = pearsonr(y_test, y_pred_xgb_tuned)
print_bold("\n\nXGBoost Regression Results After Hyperparameter Tuning")
print(f"R² Score: {r2_xgb:.3f}")
print(f"Mean Squared Error (MSE): {mse_xgb:.6f}")
print(f"Mean Absolute Error (MAE): {mae_xgb:.6f}")
print(f"Pearson Correlation: {correlation_xgb:.3f}, P-value: {p_value_xgb:.3f}")
# Cross-Validation Scores
cv_scores = cross_val_score(best_xgb_model, X_train, y_train, cv=3, scoring='r2')
print("\nCross-Validation R² Scores:", cv_scores)
print("Mean CV R² Score:", np.mean(cv_scores))
# Display the best hyperparameters found by RandomizedSearchCV
print("\nBest Hyperparameters found by RandomizedSearchCV:")
print(random_search.best_params_)
# Feature Importance from the Best XGBoost Model
feature_importances = pd.Series(best_xgb_model.feature_importances_, index=X_train.columns)
feature_importances = feature_importances.sort_values(ascending=False)
# Residual Analysis for XGBoost
residuals_xgb_tuned = y_test - y_pred_xgb_tuned
plt.figure(figsize=(8, 4))
sns.histplot(residuals_xgb_tuned, kde=True, color='green', bins=30)
plt.title("Residual Distribution (XGBoost After Tuning)", fontsize=16)
plt.xlabel("Residuals", fontsize=12)
plt.ylabel("Frequency", fontsize=12)
plt.tight_layout()
plt.show()
# Calculate Residual Metrics
mean_residual = np.mean(residuals_xgb_tuned)
std_residual = np.std(residuals_xgb_tuned)
print("\nResidual Analysis:")
print(f"- Mean of residuals: {mean_residual:.6f}")
print(f"- Standard deviation of residuals: {std_residual:.6f}")
if np.isclose(mean_residual, 0, atol=1e-2):
print("- The mean residual is close to 0, indicating a good fit.")
else:
print("- The mean residual is not close to 0, suggesting potential model improvement.")
print("\nFeature Importances (Best XGBoost Model):")
print(feature_importances.head(10))
plt.figure(figsize=(10, 6))
feature_importances.head(10).plot(kind='bar', color='skyblue')
plt.title('Top 10 Feature Importances (Best XGBoost Model)', fontsize=16)
plt.xlabel('Features', fontsize=12)
plt.ylabel('Importance', fontsize=12)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()
# Save the Best XGBoost Model
model_filename = f"{model_dir}/XGBoost_Regressor_Tuned.joblib"
joblib.dump(best_xgb_model, model_filename)
print(f"\nBest model saved: {model_filename}")
Hyperparameter Tuning with RandomizedSearchCV (XGBoost) XGBoost Regression Results After Hyperparameter Tuning R² Score: 0.865 Mean Squared Error (MSE): 0.053433 Mean Absolute Error (MAE): 0.118371 Pearson Correlation: 0.930, P-value: 0.000 Cross-Validation R² Scores: [0.86180015 0.86157849 0.86200565] Mean CV R² Score: 0.8617947644426805 Best Hyperparameters found by RandomizedSearchCV: {'subsample': 0.8, 'n_estimators': 200, 'min_child_weight': 1, 'max_depth': 7, 'learning_rate': 0.2, 'gamma': 0.1, 'colsample_bytree': 0.6}
Residual Analysis: - Mean of residuals: 0.000088 - Standard deviation of residuals: 0.231156 - The mean residual is close to 0, indicating a good fit. Feature Importances (Best XGBoost Model): AvgDuration 0.243002 Incidents 0.191493 AssessmentName_Smarter Balanced Summative Assessment 0.107642 PctEnrollment 0.042667 AssessmentName_Delaware System of Student Assessment (DeSSA) 0.039766 AssessmentName_SAT School-Day (Spring) 0.033250 ContentArea_Mathematics 0.031046 category_Suspended 0.021201 AssessmentName_DeSSA Alternate Assessment 0.018355 ContentArea_English Language Arts 0.017583 dtype: float32
Best model saved: saved_models/XGBoost_Regressor_Tuned.joblib
In [ ]: