-
Notifications
You must be signed in to change notification settings - Fork 0
Expand file tree
/
Copy pathomds_functions.py
More file actions
544 lines (461 loc) · 21.1 KB
/
omds_functions.py
File metadata and controls
544 lines (461 loc) · 21.1 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import seaborn as sns
from scipy import stats
import importlib
def show_missing_columns(df, lower_bound, upper_bound):
missing_percent = (df.isnull().sum() / len(df)) * 100
filtered_missing = missing_percent[(missing_percent > lower_bound) & (missing_percent <= upper_bound)]
count = len(filtered_missing)
# Generate markdown table
table_df = filtered_missing.reset_index()
table_df.columns = ['Column', 'Missing %']
table_df['Missing %'] = table_df['Missing %'].round(2)
print(table_df.to_markdown(index=False))
print(f"There are \033[1m{count}\033[0m columns with missing values between {lower_bound}% and {upper_bound}% in this dataset.")
return filtered_missing, count
def find_missing(df):
missing_summary = pd.DataFrame({
'Column': df.columns,
'Missing_Count': df.isnull().sum().values,
'Missing_Percent': (df.isnull().sum() / len(df) * 100).values
})
missing_summary = missing_summary.sort_values('Missing_Percent', ascending=False)
print(missing_summary)
return missing_summary
def find_outliers(dataframe):
df = dataframe.select_dtypes(include=[np.number])
for column in df.columns:
q1 = df[column].quantile(0.25)
q3 = df[column].quantile(0.75)
iqr = q3 - q1
lower_bound = q1 - 1.5 * iqr
upper_bound = q3 + 1.5 * iqr
outliers = df[(df[column] < lower_bound) | (df[column] > upper_bound)]
print(f"Outliers in column '{column}':")
print(outliers[[column]])
def calculate_r2_for_datasets(datasets, target_map, test_size=0.2, random_state=42):
"""Calculate test-set R2 for each dataset in a dict.
Args:
datasets: dict[str, pd.DataFrame]
target_map: dict[str, str] mapping dataset name to target column
test_size: fraction of rows for the test split
random_state: split seed for reproducibility
Returns:
pd.DataFrame with columns: dataset, r2, note
"""
ColumnTransformer = importlib.import_module("sklearn.compose").ColumnTransformer
SimpleImputer = importlib.import_module("sklearn.impute").SimpleImputer
LinearRegression = importlib.import_module("sklearn.linear_model").LinearRegression
r2_score = importlib.import_module("sklearn.metrics").r2_score
train_test_split = importlib.import_module("sklearn.model_selection").train_test_split
Pipeline = importlib.import_module("sklearn.pipeline").Pipeline
OneHotEncoder = importlib.import_module("sklearn.preprocessing").OneHotEncoder
results = []
for name, df in datasets.items():
target_col = target_map.get(name)
if target_col is None:
results.append({"dataset": name, "r2": None, "note": "No target in target_map"})
continue
if target_col not in df.columns:
results.append({"dataset": name, "r2": None, "note": f"Target '{target_col}' not found"})
continue
data = df.copy().dropna(subset=[target_col])
X = data.drop(columns=[target_col])
y = data[target_col]
if len(data) < 3:
results.append({"dataset": name, "r2": None, "note": "Not enough rows"})
continue
numeric_cols = X.select_dtypes(include=["number", "bool"]).columns.tolist()
categorical_cols = X.select_dtypes(exclude=["number", "bool"]).columns.tolist()
preprocessor = ColumnTransformer(
transformers=[
(
"num",
Pipeline(steps=[("imputer", SimpleImputer(strategy="median"))]),
numeric_cols,
),
(
"cat",
Pipeline(
steps=[
("imputer", SimpleImputer(strategy="most_frequent")),
("encoder", OneHotEncoder(handle_unknown="ignore")),
]
),
categorical_cols,
),
],
remainder="drop",
)
model = Pipeline(
steps=[
("preprocessor", preprocessor),
("regressor", LinearRegression()),
]
)
X_train, X_test, y_train, y_test = train_test_split(
X, y, test_size=test_size, random_state=random_state
)
model.fit(X_train, y_train)
y_pred = model.predict(X_test)
r2 = r2_score(y_test, y_pred)
results.append({"dataset": name, "r2": float(r2), "note": "ok"})
return pd.DataFrame(results).sort_values("r2", ascending=False, na_position="last")
def regplotter(df, feature1, feature1_title, feature2, feature2_title, feature3, feature3_title):
# Accept either a string column name or a one-item list like ['col_name'].
f1 = feature1[0] if isinstance(feature1, (list, tuple)) else feature1
f2 = feature2[0] if isinstance(feature2, (list, tuple)) else feature2
f3 = feature3[0] if isinstance(feature3, (list, tuple)) else feature3
featurelist = [f1, f2, f3]
df_clean = df.dropna(subset=featurelist)
# Set style
sns.set(style="whitegrid")
plt.rcParams['figure.figsize'] = (10, 8)
plt.figure(figsize=(10, 8))
# Create scatter plot
scatter = sns.scatterplot(
data=df_clean,
x=f1,
y=f2,
hue=feature3,
palette='viridis',
alpha=0.7,
s=60,
edgecolor='k',
legend=False
)
# Add regression line (using all data points, not colored by state)
reg_line = sns.regplot(
data=df_clean,
x=f1,
y=f2,
scatter=False, # Don't show the scatter points again
color='red',
line_kws={'linewidth': 2.5, 'label': 'Regression Line'},
ci=95, # Show 95% confidence interval
)
# Calculate and display regression statistics.
x_values = df_clean[f1].to_numpy(dtype=float)
y_values = df_clean[f2].to_numpy(dtype=float)
slope, intercept = np.polyfit(x_values, y_values, 1)
r_value = np.corrcoef(x_values, y_values)[0, 1]
r_squared = float(r_value ** 2)
p_value = float("nan")
# Add text annotation with regression statistics
text_str = f'Regression Statistics:\nSlope: {slope:.2f}\nR²: {r_squared:.3f}\nP-value: {p_value:.4f}'
plt.text(0.80, 0.15, text_str, transform=plt.gca().transAxes,
fontsize=11, verticalalignment='top',
bbox=dict(boxstyle='round', facecolor='wheat', alpha=0.5))
plt.title(f'Relationship Between {feature1_title} and {feature2_title} (with Regression Analysis)', fontsize=16)
plt.xlabel(f'{feature1_title} ({f1})', fontsize=12)
plt.ylabel(f'{feature2_title} ({f2})', fontsize=12)
plt.axhline(0, color='darkgray', linestyle='--', linewidth=1.5, label='Break-even Point')
plt.tight_layout()
plt.show()
# Optional: Print detailed regression output
print("=" * 60)
print("REGRESSION ANALYSIS SUMMARY")
print("=" * 60)
print(f"Dependent Variable: {feature2_title} ({f2})")
print(f"Independent Variable: {feature1_title} ({f1})")
print(f"\nRegression Equation: y = {intercept:.2f} + ({slope:.2f})x")
print(f"R-squared: {r_squared:.3f}")
print(f"P-value: {p_value:.4f}")
print(f"\nInterpretation:")
print(f"- For every 1-unit increase in {f1}, {f2} changes by {slope:.2f}")
print(f"- R² of {r_squared:.3f} indicates {'strong' if r_squared > 0.5 else 'moderate' if r_squared > 0.2 else 'weak'} correlation")
print(f"- P-value {'< 0.05 (statistically significant)' if p_value < 0.05 else '> 0.05 (not statistically significant)'}")
print("=" * 60)
return slope, intercept, r_squared, p_value
def tree_compare(df, target_col, feature_cols, test_size=0.2, random_state=42):
import numpy as np
from sklearn.ensemble import RandomForestClassifier, RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.metrics import root_mean_squared_error
X_cmp = df[feature_cols].dropna()
y_cmp = (df.loc[X_cmp.index, target_col] > 0).astype(int)
X_train_cmp, X_test_cmp, y_train_cmp, y_test_cmp = train_test_split(
X_cmp, y_cmp, test_size=0.2, random_state=42, stratify=y_cmp
)
rf_clf_cmp = RandomForestClassifier(n_estimators=200, max_depth=10, random_state=42)
rf_clf_cmp.fit(X_train_cmp, y_train_cmp)
y_pred_clf_cmp = rf_clf_cmp.predict(X_test_cmp)
rmse_clf = root_mean_squared_error(y_test_cmp, y_pred_clf_cmp)
rf_reg_cmp = RandomForestRegressor(n_estimators=200, max_depth=10, random_state=42)
rf_reg_cmp.fit(X_train_cmp, y_train_cmp)
y_pred_reg_cmp = rf_reg_cmp.predict(X_test_cmp)
rmse_reg = root_mean_squared_error(y_test_cmp, y_pred_reg_cmp)
print(f"RandomForestClassifier RMSE: {rmse_clf:.4f}")
print(f"RandomForestRegressor RMSE: {rmse_reg:.4f}")
print("Lower RMSE:", "Classifier" if rmse_clf < rmse_reg else "Regressor")
def regplottter(df, feature1, feature1_title, feature2, feature2_title, feature3, feature3_title):
"""Backward-compatible wrapper for the common misspelling of regplotter."""
return regplotter(df, feature1, feature1_title, feature2, feature2_title, feature3, feature3_title)
def compare_trees_cal_housing_data(metric_choice="rmse", single_tree_params=None, bagging_params=None, rf_params=None):
# Compare single DecisionTree, Bagged trees, and RandomForest on California housing
# - Dataset: California housing (sklearn)
# - Models:
# * Single DecisionTreeRegressor
# * BaggingRegressor with DecisionTreeRegressor as estimator (bagging only)
# * RandomForestRegressor (bagging + per-split feature subsampling)
# - For each model we print train/test metrics and feature importances.
#
# User-selectable places:
# - metric_choice: choose "mse", "rmse", "mae", "r2", or "explained_variance"
# - single_tree_params: dict for DecisionTreeRegressor (e.g., {"max_depth": 8})
# - bagging_params: dict for BaggingRegressor (e.g., {"n_estimators": 50})
# - rf_params: dict for RandomForestRegressor (e.g., {"n_estimators": 100, "max_features": "sqrt"})
#
# Note: sklearn APIs vary by version. Modern sklearn uses `estimator=` for BaggingRegressor;
# older versions used `base_estimator=`. This script uses `estimator=`.
import numpy as np
import pandas as pd
from sklearn.datasets import fetch_california_housing
from sklearn.tree import DecisionTreeRegressor
from sklearn.ensemble import BaggingRegressor, RandomForestRegressor
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error, mean_absolute_error, r2_score, explained_variance_score
# ------------------ User-selectable choices ------------------
# Choose performance metric: "mse", "rmse", "mae", "r2", "explained_variance"
metric_choice = "rmse" # change to "mse", "mae", "r2", or "explained_variance"
# Single decision tree hyperparameters
single_tree_params = {
"criterion": "squared_error", # "squared_error" (MSE) in sklearn >=1.0
"max_depth": 8, # set to None to grow fully
"min_samples_leaf": 1,
"random_state": 0
}
# Bagging parameters (estimator + bagging settings)
bagging_params = {
"estimator": DecisionTreeRegressor(criterion="squared_error", max_depth=None, min_samples_leaf=1, random_state=0),
"n_estimators": 50, # number of trees in the bag
"max_samples": 1.0, # fraction or int (bootstrap sample size)
"max_features": 1.0, # fraction or int (features per base estimator)
"bootstrap": True,
"bootstrap_features": False,
"n_jobs": -1,
"random_state": 0
}
# Random forest parameters
rf_params = {
"n_estimators": 100,
"criterion": "squared_error",
"max_depth": None,
"min_samples_leaf": 1,
"max_features": "sqrt", # per-split feature subsampling; change as desired
"random_state": 0,
"n_jobs": -1
}
# --------------------------------------------------------------
# Metric wrapper
def compute_metric(y_true, y_pred, choice="rmse"):
if choice == "mse":
return mean_squared_error(y_true, y_pred)
elif choice == "rmse":
return np.sqrt(mean_squared_error(y_true, y_pred))
elif choice == "mae":
return mean_absolute_error(y_true, y_pred)
elif choice == "r2":
return r2_score(y_true, y_pred)
elif choice == "explained_variance":
return explained_variance_score(y_true, y_pred)
else:
raise ValueError("Unknown metric_choice: " + str(choice))
# Load data
cal = fetch_california_housing()
X = pd.DataFrame(cal.data, columns=cal.feature_names)
y = pd.Series(cal.target, name="MedHouseVal")
feature_names = list(X.columns)
# Train/test split (user can adjust test_size and random_state here if desired)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=0)
# Helper to format and print feature importances
def print_feature_importances(importances, names, title):
s = pd.Series(importances, index=names).sort_values(ascending=False)
print(title)
print(s.to_string())
print()
# 1) Single Decision Tree
dt = DecisionTreeRegressor(**single_tree_params)
dt.fit(X_train, y_train)
pred_train_dt = dt.predict(X_train)
pred_test_dt = dt.predict(X_test)
metric_train_dt = compute_metric(y_train, pred_train_dt, metric_choice)
metric_test_dt = compute_metric(y_test, pred_test_dt, metric_choice)
print("=== Single Decision Tree ===")
print("Params:", single_tree_params)
print(f"Train {metric_choice}: {metric_train_dt:.6f}")
print(f"Test {metric_choice}: {metric_test_dt:.6f}")
print_feature_importances(dt.feature_importances_, feature_names, "Feature importances (single tree):")
# 2) Bagged trees (BaggingRegressor)
bag = BaggingRegressor(**bagging_params)
bag.fit(X_train, y_train)
pred_train_bag = bag.predict(X_train)
pred_test_bag = bag.predict(X_test)
metric_train_bag = compute_metric(y_train, pred_train_bag, metric_choice)
metric_test_bag = compute_metric(y_test, pred_test_bag, metric_choice)
# Compute averaged feature importances across base estimators if they expose feature_importances_
base_importances = []
for est in bag.estimators_:
# in BaggingRegressor, estimators_ are fitted clones of the provided estimator
if hasattr(est, "feature_importances_"):
base_importances.append(est.feature_importances_)
if len(base_importances) > 0:
avg_importances = np.mean(base_importances, axis=0)
else:
avg_importances = np.zeros(len(feature_names))
print("=== Bagged Trees (BaggingRegressor) ===")
print("Params:", {k: v for k, v in bagging_params.items() if k != "estimator"})
print(f"n_estimators: {bagging_params['n_estimators']}")
print(f"Train {metric_choice}: {metric_train_bag:.6f}")
print(f"Test {metric_choice}: {metric_test_bag:.6f}")
print_feature_importances(avg_importances, feature_names, "Averaged feature importances (bagged trees):")
# 3) Random Forest
rf = RandomForestRegressor(**rf_params)
rf.fit(X_train, y_train)
pred_train_rf = rf.predict(X_train)
pred_test_rf = rf.predict(X_test)
metric_train_rf = compute_metric(y_train, pred_train_rf, metric_choice)
metric_test_rf = compute_metric(y_test, pred_test_rf, metric_choice)
print("=== Random Forest ===")
print("Params:", rf_params)
print(f"Train {metric_choice}: {metric_train_rf:.6f}")
print(f"Test {metric_choice}: {metric_test_rf:.6f}")
print_feature_importances(rf.feature_importances_, feature_names, "Feature importances (random forest):")
summary = pd.DataFrame({
"model": ["DecisionTree", "BaggedTrees", "RandomForest"],
"train_" + metric_choice: [metric_train_dt, metric_train_bag, metric_train_rf],
"test_" + metric_choice: [metric_test_dt, metric_test_bag, metric_test_rf]
})
print("=== Summary ===")
print(summary.to_string(index=False))
def compare_random_forest(
df,
target_col,
feature_cols=None,
test_size=0.2,
random_state=42,
threshold=0.0,
n_estimators=200,
max_depth=10,
):
"""Compare RandomForestClassifier vs RandomForestRegressor for DCD 2025 data.
Parameters
----------
df : pandas.DataFrame
Input dataset.
target_col : str
Continuous target column (for example, a delta temperature metric).
feature_cols : list[str] | None
Feature columns to use. If None, all columns except target_col are used.
test_size : float
Fraction of rows reserved for testing.
random_state : int
Seed for reproducibility.
threshold : float
Threshold used to create binary labels for the classifier: y > threshold.
n_estimators : int
Number of trees for both RF models.
max_depth : int | None
Maximum tree depth for both RF models.
Returns
-------
dict
RMSE values and the better model name.
"""
import pandas as pd
from sklearn.compose import ColumnTransformer
from sklearn.ensemble import RandomForestClassifier, RandomForestRegressor
from sklearn.impute import SimpleImputer
from sklearn.metrics import root_mean_squared_error
from sklearn.model_selection import train_test_split
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import OneHotEncoder
if target_col not in df.columns:
raise ValueError(f"target_col '{target_col}' is not present in the dataframe.")
if feature_cols is None:
feature_cols = [c for c in df.columns if c != target_col]
missing_features = [c for c in feature_cols if c not in df.columns]
if missing_features:
raise ValueError(f"feature_cols contain missing columns: {missing_features}")
work_df = df[feature_cols + [target_col]].copy()
work_df = work_df.dropna(subset=[target_col])
X = work_df[feature_cols]
y_reg = pd.to_numeric(work_df[target_col], errors="coerce")
valid_rows = y_reg.notna()
X = X.loc[valid_rows]
y_reg = y_reg.loc[valid_rows]
if len(X) < 10:
raise ValueError("Not enough valid rows after cleaning to train/test split.")
y_clf = (y_reg > threshold).astype(int)
stratify_vec = y_clf if y_clf.nunique() > 1 else None
X_train, X_test, y_train_reg, y_test_reg, y_train_clf, y_test_clf = train_test_split(
X, y_reg, y_clf,
test_size=test_size,
random_state=random_state,
stratify=stratify_vec,
)
numeric_cols = X.select_dtypes(include=["number", "bool"]).columns.tolist()
categorical_cols = [c for c in X.columns if c not in numeric_cols]
preprocessor = ColumnTransformer(
transformers=[
(
"num",
Pipeline(steps=[("imputer", SimpleImputer(strategy="median"))]),
numeric_cols,
),
(
"cat",
Pipeline(
steps=[
("imputer", SimpleImputer(strategy="most_frequent")),
("onehot", OneHotEncoder(handle_unknown="ignore")),
]
),
categorical_cols,
),
],
remainder="drop",
)
clf_model = Pipeline(
steps=[
("prep", preprocessor),
(
"rf",
RandomForestClassifier(
n_estimators=n_estimators,
max_depth=max_depth,
random_state=random_state,
),
),
]
)
reg_model = Pipeline(
steps=[
("prep", preprocessor),
(
"rf",
RandomForestRegressor(
n_estimators=n_estimators,
max_depth=max_depth,
random_state=random_state,
),
),
]
)
clf_model.fit(X_train, y_train_clf)
y_pred_clf = clf_model.predict(X_test)
rmse_clf = root_mean_squared_error(y_test_clf, y_pred_clf)
reg_model.fit(X_train, y_train_reg)
y_pred_reg = reg_model.predict(X_test)
rmse_reg = root_mean_squared_error(y_test_reg, y_pred_reg)
better_model = "Classifier" if rmse_clf < rmse_reg else "Regressor"
print(f"RandomForestClassifier RMSE: {rmse_clf:.4f}")
print(f"RandomForestRegressor RMSE: {rmse_reg:.4f}")
print(f"Lower RMSE: {better_model}")
return {
"rmse_classifier": rmse_clf,
"rmse_regressor": rmse_reg,
"better_model": better_model,
}