<jupyter_start><jupyter_text>GLM VS ML - Decision Tree Deconstruction and MC Simulations ![](https://datascience.foundation/img/pdf_images/understanding_decision_trees_with_python_decision_tree.png) **Imports**<jupyter_code>import importlib
import numpy as np
import pandas as pd
import scipy.stats as ss
import seaborn as sns
import matplotlib.pyplot as plt
from sklearn.model_selection import train_test_split
import statsmodels.api as sm
from sklearn import ensemble
from sklearn import tree
from sklearn import metrics
plt.style.use('seaborn')<jupyter_output><empty_output><jupyter_text>Auto Claims Data**[Dataset](https://www.kaggle.com/xiaomengsun/car-insurance-claim-data)**|VARIABLE NAME|DEFINITION                              |THEORETICAL EFFECT                                                                               ||-------------|----------------------------------------|-------------------------------------------------------------------------------------------------||INDEX        |Identification Variable (do not use)    |None                                                                                             ||TARGET FLAG  |Was Car in a crash? 1=YES 0=NO          |None                                                                                             ||TARGET AMT   |If car was in a crash, what was the cost|None                                                                                             ||AGE          |Age of Driver                           |Very young people tend to be risky. Maybe very old people also.                                  ||BLUEBOOK     |Value of Vehicle                        |Unknown effect on probability of collision, but probably effect the payout if there is a crash   ||CAR AGE      |Vehicle Age                             |Unknown effect on probability of collision, but probably effect the payout if there is a crash   ||CAR TYPE     |Type of Car                             |Unknown effect on probability of collision, but probably effect the payout if there is a crash   ||CAR USE      |Vehicle Use                             |Commercial vehicles are driven more, so might increase probability of collision                  ||CLM FREQ     | Claims (Past 5 Years)                 |The more claims you filed in the past, the more you are likely to file in the future             ||EDUCATION    |Max Education Level                     |Unknown effect, but in theory more educated people tend to drive more safely                     ||HOMEKIDS     | Children at Home                      |Unknown effect                                                                                   ||HOME VAL     |Home Value                              |In theory, home owners tend to drive more responsibly                                            ||INCOME       |Income                                  |In theory, rich people tend to get into fewer crashes                                            ||JOB          |Job Category                            |In theory, white collar jobs tend to be safer                                                    ||KIDSDRIV     | Driving Children                      |When teenagers drive your car, you are more likely to get into crashes                           ||MSTATUS      |Marital Status                          |In theory, married people drive more safely                                                      ||MVR PTS      |Motor Vehicle Record Points             |If you get lots of traffic tickets, you tend to get into more crashes                            ||OLDCLAIM     |Total Claims (Past 5 Years)             |If your total payout over the past five years was high, this suggests future payouts will be high||PARENT1      |Single Parent                           |Unknown effect                                                                                   ||RED CAR      |A Red Car                               |Urban legend says that red cars (especially red sports cars) are more risky. Is that true?       ||REVOKED      |License Revoked (Past 7 Years)          |If your license was revoked in the past 7 years, you probably are a more risky driver.           ||SEX          |Gender                                  |Urban legend says that women have less crashes then men. Is that true?                           ||TIF          |Time in Force                           |People who have been customers for a long time are usually more safe.                            ||TRAVTIME     |Distance to Work                        |Long drives to work usually suggest greater risk                                                 ||URBANICITY   |Home/Work Area                          |Unknown                                                                                          ||YOJ          |Years on Job                            |People who stay at a job for a long time are usually more safe                                   | **Load Data and Clean Data**<jupyter_code># Load Data
df = pd.read_csv('~/Downloads/car_insurance_claim.csv')

# make columns lowercase
df.columns = df.columns.str.lower()

# drop useless columns
df = df.drop(['kidsdriv','parent1','revoked','mvr_pts','travtime','id','birth'],axis=1)

# clean money amounts
df[['home_val','bluebook','oldclaim','clm_amt','income']] = df[['home_val','bluebook','oldclaim','clm_amt','income']].apply(lambda x: x.str.replace('$','',regex=False).str.replace(',','',regex=False)).astype(float)

# clean values from columns
to_clean = ['education','occupation','mstatus','gender','car_type']
for col in to_clean:
    df[col] = df[col].str.replace('z_','',regex=False).str.replace('<','',regex=False)

df['urbanicity'] = df['urbanicity'].str.split('/ ',expand=True)[1]

to_clean = ['mstatus','red_car']
for col in to_clean:
    df[col] = df[col].str.lower().replace({ 'yes': True, 'no': False}).astype(int)
    
df = df.drop(['car_age','occupation','home_val','income','yoj'],axis=1).dropna()
  
df[:3]<jupyter_output><empty_output><jupyter_text>**Select Features and Targets**<jupyter_code>features = ['age','gender','car_type','red_car','tif','car_use','bluebook','oldclaim']
binary_target = 'claim_flag'
severity_target = 'clm_amt'
processed_data = df.copy()
processed_data = pd.get_dummies(processed_data)
processed_data.columns = processed_data.columns.str.replace(' ','_',regex=False)
processed_data_features = list(set(processed_data.columns.tolist()) - set([binary_target, severity_target]))
processed_data
train, test = train_test_split(processed_data, test_size=0.45, random_state=42, stratify=processed_data[binary_target])

train_severity = train.loc[train[severity_target] > 0].copy()<jupyter_output><empty_output><jupyter_text>GLMs **GLM Binomial Model**<jupyter_code>binomial_train_X = sm.add_constant(train[processed_data_features])
binomial = sm.GLM(train[binary_target],binomial_train_X, family=sm.families.Binomial())
binomial_fitted = binomial.fit()
binomial_fitted.summary()
gamma_train_X = sm.add_constant(train_severity[processed_data_features])
gamma = sm.GLM(train_severity[severity_target], gamma_train_X, family=sm.families.Gamma(sm.families.links.log()))
gamma_fitted = gamma.fit()
gamma_fitted.summary()
from modules import dt_deconstruction
importlib.reload(dt_deconstruction)
from modules.dt_deconstruction import DTRuler,DTRulerAvg<jupyter_output><empty_output><jupyter_text>Tree Based Models<jupyter_code>size = 20
seed = 42
np.random.seed(seed)
dummy_data = pd.DataFrame({
    'car_type': np.random.choice(['Audi','BMW','Porsche'], size=size),
    'weight': np.random.normal(2500,size=size),
    'hp': np.random.randint(75,200,size=size),
    'speed': np.random.gamma(200,scale=np.random.choice([0.3, 0.5, 0.7, 1.2],size=size),size=size)
})
dummy_data = pd.get_dummies(dummy_data)
dummy_features = dummy_data.drop(['speed'],axis=1).columns.tolist()
dummy_data[:10]
rf_reg = ensemble.RandomForestRegressor(n_estimators=3, random_state=seed)
rf_reg.fit(dummy_data[dummy_features], dummy_data['speed'])
rf_reg.estimators_
fig,axs = plt.subplots(3,1,figsize=(25,25))
for i, es in enumerate(rf_reg.estimators_):
    
    tree.plot_tree(es,feature_names=dummy_features, ax=axs[i],fontsize=12);
    axs[i].set_title(f'Decision Trees in the Estimator[{i}]',fontsize=25)
    
plt.tight_layout()<jupyter_output><empty_output><jupyter_text>**Decision Tree Classifier Model**<jupyter_code>clf = ensemble.RandomForestClassifier(n_estimators=10,random_state=seed)
clf.fit(train[processed_data_features],train[binary_target])<jupyter_output><empty_output><jupyter_text>**Decision Tree Regressor Model**<jupyter_code>reg = ensemble.RandomForestRegressor(n_estimators=10, random_state=seed)
reg.fit(train_severity[processed_data_features], train_severity[severity_target])<jupyter_output><empty_output><jupyter_text>Tree Deconstruction DTRulerAvg will deconstruct the decision trees and extract the datafiltering rules<jupyter_code>%%time
clf_dtruler = dt_deconstruction.DTRulerAvg(model=clf,training_data=train,features=processed_data_features)
%%time
reg_dtruler = dt_deconstruction.DTRulerAvg(model=reg,training_data=train_severity,features=processed_data_features)
def custom_func(d):
    return {
        'mean': d.mean(),
        'std': d.std(),
        'var': d.var()
    }<jupyter_output><empty_output><jupyter_text>After we have our instances for DT rulers we can call the predict method which will concatenate all the samples from training data that matches the decision path for the given sample, in the same time we pass our custom aggregation function which will calculate in the currentc case mean, std and variance for the target variable in each data group/cluster<jupyter_code>%%time
clf_data_by_avg_trees = clf_dtruler.predict(test.reset_index(), aggfunc=lambda x: custom_func(x[binary_target]))
%%time
reg_data_by_avg_trees = reg_dtruler.predict(test.reset_index(), aggfunc=lambda x: custom_func(x[severity_target]))<jupyter_output>CPU times: user 7min 1s, sys: 542 ms, total: 7min 1s
Wall time: 7min 1s<jupyter_text>Then we want to validate our distributions from each group/cluster against the original distribution and to do that we'll use the next statistica test.  Kolmogorov-Smirnov testPerforms the (one-sample or two-sample) Kolmogorov-Smirnov test for goodness of fit.The one-sample test compares the underlying distribution F(x) of a sample against a given distribution G(x). The two-sample test compares the underlying distributions of two independent samples. Both tests are valid only for continuous distributions.<jupyter_code>reg_data_by_avg_trees_copy = reg_data_by_avg_trees.copy()
clf_data_by_avg_trees_copy = clf_data_by_avg_trees.copy()<jupyter_output><empty_output><jupyter_text>Compute the statistical test for each group<jupyter_code>%%time
reg_data_by_avg_trees_copy_kstest = reg_data_by_avg_trees_copy.all_trees_data.apply(lambda x: ss.kstest(x[severity_target], test.loc[test[severity_target]>0,severity_target]))
clf_data_by_avg_trees_copy_kstest = clf_data_by_avg_trees_copy.all_trees_data.apply(lambda x: ss.kstest(x[binary_target], test[binary_target]))<jupyter_output>CPU times: user 1min 21s, sys: 158 ms, total: 1min 22s
Wall time: 1min 21s<jupyter_text>Extract the statisict score and the P-values<jupyter_code>reg_data_by_avg_trees_copy[['kstest_statistic','kstest_pvalue']] = pd.DataFrame.from_dict(reg_data_by_avg_trees_copy_kstest.apply(lambda x: dict(kstest_statistic=x.statistic, kstest_pvalue=x.pvalue)).values.tolist()).values
reg_data_by_avg_trees_copy['statistical_significant'] = reg_data_by_avg_trees_copy['kstest_pvalue'].between(-0.05,0.05)
clf_data_by_avg_trees_copy[['kstest_statistic','kstest_pvalue']] = pd.DataFrame.from_dict(clf_data_by_avg_trees_copy_kstest.apply(lambda x: dict(kstest_statistic=x.statistic, kstest_pvalue=x.pvalue)).values.tolist()).values
clf_data_by_avg_trees_copy['statistical_significant'] = clf_data_by_avg_trees_copy['kstest_pvalue'].between(-0.05,0.05)
fig,axs = plt.subplots(1,2,figsize=(15,5),sharey=True)
sns.scatterplot(data=clf_data_by_avg_trees_copy, x='kstest_pvalue',y='kstest_statistic',hue='statistical_significant',ax=axs[0]);
sns.scatterplot(data=reg_data_by_avg_trees_copy, x='kstest_pvalue',y='kstest_statistic',hue='statistical_significant',ax=axs[1]);
axs[0].set_title('[Kolmogorov-Smirnov-test](Binomial) Leaf Data vs Original Distribution')
axs[1].set_title('[Kolmogorov-Smirnov-test](Gamma) Leaf Data vs Original Distribution')
plt.tight_layout()<jupyter_output><empty_output><jupyter_text>Monte Carlo Simulations Now that we validated our distributions, we want to apply Monte Carlo Simiulations using the extracted parameters from DT ruler and compare the results aginast the GLM models Add constant to the test dataframe<jupyter_code>X_test = sm.add_constant(test[processed_data_features])<jupyter_output>/opt/homebrew/Caskroom/miniforge/base/envs/ml/lib/python3.9/site-packages/statsmodels/tsa/tsatools.py:142: FutureWarning: In a future version of pandas all arguments of concat except for the argument 'objs' will be keyword-only
  x = pd.concat(x[::order], 1)<jupyter_text>Estimate the gama scale for each group/cluster<jupyter_code>estimated_scale = 1 / (reg_data_by_avg_trees_copy['mean'] / reg_data_by_avg_trees_copy['std']) ** 2
reg_data_by_avg_trees_copy['estimated_scale'] = estimated_scale<jupyter_output><empty_output><jupyter_text>Simulate frequency and severity for policies<jupyter_code>years = 1000

test_copy = test[[binary_target,severity_target]].copy().reset_index().drop(['index'],axis=1)

test_copy['glm_binomial_prediction'] = binomial_fitted.predict(X_test).values
test_copy['tree_binomial_prediction'] = clf.predict_proba(test[processed_data_features])[:,1]

test_copy['glm_gamma_prediction'] = gamma_fitted.predict(X_test).values
test_copy['tree_gamma_prediction'] = reg.predict(test[processed_data_features])


test_copy['glm_binomial_simulation'] = pd.Series(binomial_fitted.get_prediction(X_test).se_mean).apply(lambda x: np.random.binomial(1, p=x,size=years).mean()).values
test_copy['tree_binomial_simulation'] = clf_data_by_avg_trees_copy['mean'].apply(lambda x: np.random.binomial(1, p=x,size=years).mean()).values

test_copy['glm_gamma_simulation'] = pd.Series(gamma_fitted.get_prediction(X_test).se_mean).apply(lambda x: np.random.gamma(x, scale=gamma_fitted.scale,size=years).mean()).values
test_copy['tree_gamma_simulation'] = reg_data_by_avg_trees_copy.apply(lambda x: np.random.gamma(x['mean'], scale=x['estimated_scale'],size=years).mean() ,axis=1).values

test_copy<jupyter_output><empty_output><jupyter_text>**Actual Claim Frequency distribution against the Predicted and Simulated**<jupyter_code>test_copy[['claim_flag','glm_binomial_prediction','tree_binomial_prediction','glm_binomial_simulation','tree_binomial_simulation',]].applymap(
lambda x: 1 if x >=0.5 else 0
).hist(layout=(1,5),figsize=(25,3),sharey=True);
tmp = test_copy[['claim_flag','glm_binomial_prediction','tree_binomial_prediction','glm_binomial_simulation','tree_binomial_simulation',]].applymap(
lambda x: 1 if x >=0.5 else 0)

for x in ['glm_binomial_prediction','tree_binomial_prediction','glm_binomial_simulation','tree_binomial_simulation']:
    print(f'\nactual claims frequency vs {x}\n[accuracy_score={metrics.accuracy_score(tmp["claim_flag"],tmp[x])}, recall_score={metrics.recall_score(tmp["claim_flag"],tmp[x])}, f1_score={metrics.f1_score(tmp["claim_flag"],tmp[x])}]')<jupyter_output>actual claims frequency vs glm_binomial_prediction
[accuracy_score=0.7750917332182171, recall_score=0.37550689375506896, f1_score=0.47052845528455284]

actual claims frequency vs tree_binomial_prediction
[accuracy_score=0.7299805741420246, recall_score=0.4468775344687753, f1_score=0.4683382915427114]

actual claims frequency vs glm_binomial_simulation
[accuracy_score=0.7338657457371034, recall_score=0.0, f1_score=0.0]

actual claims frequency vs tree_binomial_simulation
[accuracy_score=0.7617094755018347, recall_score=0.32116788321167883, f1_score=0.4177215189873418]<jupyter_text>**Actual Claim Severity distribution against the Predicted and Simulated**<jupyter_code>test_copy.loc[test_copy.clm_amt>0][['clm_amt','glm_gamma_prediction','tree_gamma_prediction','glm_gamma_simulation','tree_gamma_simulation',]].hist(layout=(1,5),figsize=(25,3),sharey=True);
tmp = test_copy.loc[test_copy.clm_amt>0][['clm_amt','glm_gamma_prediction','tree_gamma_prediction','glm_gamma_simulation','tree_gamma_simulation',]]

for x in ['glm_gamma_prediction','tree_gamma_prediction','glm_gamma_simulation','tree_gamma_simulation',]:
    print(f'\nactual claims severity vs {x}\n[RMSE={metrics.mean_squared_error(tmp["clm_amt"],tmp[x], squared=False)}]')<jupyter_output>actual claims severity vs glm_gamma_prediction
[RMSE=7258.709605379099]

actual claims severity vs tree_gamma_prediction
[RMSE=7801.673419812962]

actual claims severity vs glm_gamma_simulation
[RMSE=8381.243766533747]

actual claims severity vs tree_gamma_simulation
[RMSE=13148.207709796603]
{'fraction_non_alphanumeric': 0.11483984719365266, 'fraction_numerical': 0.022803408756979136, 'mean_word_length': 3.828603859250851, 'pattern_counts': {'":': 0, '<': 65, '<?xml version=': 0, '>': 70, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 46, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '17', 'dir': 'jupyter-structured-clean-dedup', 'id': '560707', 'license': 'RSA-MD', 'n_tokens_mistral': 5595, 'n_tokens_neox': 5166, 'n_words': 1091, 'path': 'notebooks/2021-07-07-DTrees-Deconstruction.ipynb', 'repo_name': 'cristianexer/trainings-and-tutorials', 'size': '17015'}