Author(s): Greg Postalian-Yrausquin Originally published on Towards AI. Decision tree types of classification algorithms have the advantage that they produce results that are relatively easier to explain in terms of the impact of the predictors when compared to other supervised training algorithms, like neural networks, which is more like a closed box. To illustrate this, I use RandomForest and descriptive analytics see if it is possible to predict if a person earns more or less than 50k dollars import numpy as npimport pandas as pdimport matplotlib.pyplot as pltimport seaborn as snsimport klibimport dablfrom scipy.stats import chi2_contingencyfrom sklearn.preprocessing import OneHotEncoder, StandardScalerfrom sklearn.compose import ColumnTransformerfrom sklearn.ensemble import RandomForestClassifierfrom sklearn.model_selection import train_test_splitfrom sklearn.pipeline import make_pipelinefrom sklearn.metrics import confusion_matrixfrom sklearn import tree From the description of the dataset, this is what we know about the predictors: age: continuous. workclass: Private, Self-emp-not-inc, Self-emp-inc, Federal-gov, Local-gov, State-gov, Without-pay, Never-worked. fnlwgt: continuous. education: Bachelors, Some-college, 11th, HS-grad, Prof-school, Assoc-acdm, Assoc-voc, 9th, 7th-8th, 12th, Masters, 1st-4th, 10th, Doctorate, 5th-6th, Preschool. education-num: continuous. marital-status: Married-civ-spouse, Divorced, Never-married, Separated, Widowed, Married-spouse-absent, Married-AF-spouse. occupation: Tech-support, Craft-repair, Other-service, Sales, Exec-managerial, Prof-specialty, Handlers-cleaners, Machine-op-inspct, Adm-clerical, Farming-fishing, Transport-moving, Priv-house-serv, Protective-serv, Armed-Forces. relationship: Wife, Own-child, Husband, Not-in-family, Other-relative, Unmarried. race: White, Asian-Pac-Islander, Amer-Indian-Eskimo, Other, Black. sex: Female, Male. capital-gain: continuous. capital-loss: continuous. hours-per-week: continuous. native-country: United-States, Cambodia, England, Puerto-Rico, Canada, Germany, Outlying-US(Guam-USVI-etc), India, Japan, Greece, South, China, Cuba, Iran, Honduras, Philippines, Italy, Poland, Jamaica, Vietnam, Mexico, Portugal, Ireland, France, Dominican-Republic, Laos, Ecuador, Taiwan, Haiti, Columbia, Hungary, Guatemala, Nicaragua, Scotland, Thailand, Yugoslavia, El-Salvador, Trinadad&Tobago, Peru, Hong, Holand-Netherlands. maindataset = pd.read_csv("adultsalary.csv")print(maindataset)print(maindataset.columns) Right from the start, I can see that some columns are going to be correlated (marital status — relationship, education and education numerical) and some others are basically empty. The next step is a review of the numerical variables klib.dist_plot(maindataset['age'])klib.dist_plot(maindataset[' education-num'])klib.dist_plot(maindataset[' hours-per-week']) For the categorical variables: #I am going to group the countries of origin to reduce the number of values of the variableus = [' United-States']latin = [' Cuba',' South',' Puerto-Rico',' Honduras',' Columbia',' Ecuador',' Dominican-Republic',' El-Salvador',' Guatemala',' Peru',' Nicaragua']carib = [' Jamaica',' Haiti',' Outlying-US(Guam-USVI-etc)',' Trinadad&Tobago']mexico = [' Mexico']canada = [' Canada']westeur = [' England',' Germany',' France',' Scotland',' Ireland',' Holand-Netherlands']easteur = [' Poland',' Yugoslavia',' Hungary']southeur = [' Italy',' Portugal',' Greece']iran = [' Iran']india = [' India']seasia = [' Philippines',' Cambodia',' Thailand',' Laos',' Vietnam']china = [' China']easia = [' Taiwan',' Japan',' Hong']maindataset['region'] = np.nanmaindataset['region'] = np.where(maindataset[' native-country'].isin(us),'US',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(latin),'Other LA',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(carib),'Caribbean',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(mexico),'Mexico',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(canada),'Canada',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(westeur),'West Euro',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(easteur),'East Euro',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(southeur),'South Euro',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(iran),'Iran',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(india),'India',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(seasia),'SE Asia',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(china),'China',maindataset['region'])maindataset['region'] = np.where(maindataset[' native-country'].isin(easia),'E Asia',maindataset['region'])maindataset[' salary'] = np.where(maindataset[' salary']==' <=50K.',' <=50K',maindataset[' salary'])maindataset[' salary'] = np.where(maindataset[' salary']==' >50K.',' >50K',maindataset[' salary'])klib.cat_plot(maindataset[[' workclass', ' occupation', ' relationship', ' race', ' sex', 'region', ' salary']], top=5, bottom=5) Now I apply the filters and clean the data maindatasetF = maindataset[['age',' education-num',' hours-per-week',' workclass', ' occupation', ' relationship', ' race', ' sex', 'region', ' salary']]maindatasetF = maindatasetF[maindatasetF[' workclass']!=' ?']maindatasetF = maindatasetF[maindatasetF[' occupation']!=' ?']maindatasetF = maindatasetF.dropna()maindatasetF Exploring the correlations. For numerical data: Which are not strong, positive, correlations For the categorical data two tests are important: Chi-square: to test for independence between the categorical variables. ANOVA: to determine if between the values of the categorical variables there is a difference in the outcome (salary range). For the Chi-square test: maindatasetF = maindatasetF.rename(columns={' salary':'salary', ' workclass':'workclass', ' occupation':'occupation', ' relationship':'relationship', ' race':'race', ' sex':'sex'})catvar = ['workclass', 'occupation', 'relationship', 'race', 'sex', 'region', 'salary']xtabs = []pairs = []for i in catvar: for j in catvar: if i!=j: pair = pd.DataFrame([[i,j]], columns=['i','j']) data_xtab = pd.crosstab(maindatasetF[i],maindatasetF[j],margins=False) xtabs.append(data_xtab) pairs.append(pair)pairs = pd.concat(pairs, ignore_index=True, axis=0)ps = []for i in xtabs: stat, p, dof, expected = chi2_contingency(i) ps.append(p)pairs['p values'] = pspairs These very low values of p mean that there is dependency between the features. ANOVA test import statsmodels.api as smfrom statsmodels.formula.api import olsmodel = ols('num_salary ~ region + workclass + occupation + relationship + race + sex', data=maindatasetF).fit()aov_table = sm.stats.anova_lm(model, typ=2)aov_table These low values of P mean that we can reject the NULL hypothesis, specifying that there is a difference between the members of the categorical features, and so, there is predictive power in them. Next, it’s a random forest classification, and get the variable importance. num_columns = ['age',' education-num',' hours-per-week']cat_columns = ['workclass', 'occupation', 'relationship', 'race', 'sex', 'region']cat_preprocessor = OneHotEncoder(handle_unknown="ignore")num_preprocessor = StandardScaler()preprocessor = ColumnTransformer( [ ("one-hot-encoder", cat_preprocessor, cat_columns), ("standard_scaler", num_preprocessor, num_columns), ])train, test = train_test_split(maindatasetF, train_size=0.8)train = train.dropna()X_train = train[['age',' education-num',' hours-per-week','workclass', 'occupation', 'relationship', 'race', 'sex', 'region']]Y_train = train[['num_salary']]X = pd.DataFrame.sparse.from_spmatrix(preprocessor.fit_transform(Xtrain))catnames = preprocessor.transformers[0][1].get_feature_names_out(catcolumns).tolist()numnames = preprocessor.transformers[1][1].get_feature_names_out(num_columns).tolist()featnames = catnames + numnamesrf = RandomForestClassifier(n_estimators=100)rf.fit(X, Y_train) Plot of the importance of the variables imp = rf.feature_importances_imp = pd.Series(imp, index=featnames)std = pd.Series(np.std([tree.featureimportances for tree in rf.estimators_], axis=0), index=featnames)fig, ax = plt.subplots()imp.plot(kind='barh', yerr=std, ax=ax, figsize=(15,15))ax.set_title("Feature importances using MDI")ax.set_ylabel("Mean decrease in impurity")fig.tight_layout()plt.show() Age, education, hours worked and the fact that the subject is a husband dominate in the outcome of the salary, which is in line with what we know of the society we live in. Now, I am performing a test with a confusion matrix to see how accurately we can predict the salary from the Random Forest model. X_test = test[['age',' education-num',' hours-per-week','workclass', 'occupation', 'relationship', 'race', 'sex', 'region']]Y_test = test[['num_salary']]Xt = pd.DataFrame.sparse.from_spmatrix(preprocessor.fit_transform(X_test))Yt = rf.predict(Xt)cm = confusion_matrix(Y_test['num_salary'], Yt)sns.heatmap(cm, annot=True, cmap='PuBu', fmt='g') Which can be improved, but it is a good start. At least we know that the model used for variable importance, which is what we want to understand, makes sense. I would like to draw a decision tree as well, so see the variables act in a graph. tr = tree.DecisionTreeClassifier(max_depth=4)tr = tr.fit(X, Y_train)#I am limiting the depth to threefig, ax1 = plt.subplots(figsize=(25,15))tree.plot_tree(tr, ax=ax1, feature_names=featnames, proportion=True, filled=True, fontsize=7)plt.show() Blue indicates higher tendency towards high salaries. Note that the numeric variables are standardized (0 is the mean). “Yes’s” are to the left, and “No’s” to the right. So here is clear to see that 80% of: husbands, with higher education than the mean, that work more than the mean, and older than the mean, earn more than 50k. This group is only 9% […]