machine learning

# import fundamental modules and prepare data
# for references of sklearn modules, see http://scikit-learn.org/stable/modules/classes.html
import numpy as np
import pandas as pd
from sklearn import datasets
from sklearn.cross_validation import train_test_split
# "cross_validation" module is outdated and deprecated
# train_test_split is now implemented in "model_selection"
# from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import accuracy_score

# load Iris and split into training & test sets
iris = datasets.load_iris()
X = iris.data[:, [2, 3]]
y = iris.target
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size = 0.3, random_state=0)

# Standardization
sc = StandardScaler()
sc.fit(X_train)
X_train_std = sc.transform(X_train)
X_test_std = sc.transform(X_test)

X_combined = np.vstack((X_train, X_test))
X_combined_std = np.vstack((X_train_std, X_test_std))
y_combined = np.hstack((y_train, y_test))

# load Wine data from UCI
# for details, see https://archive.ics.uci.edu/ml/datasets/wine
#df_wine = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/wine/wine.data', header=None)
df_wine = pd.read_csv('wine.csv', header=None)
df_wine.columns = ['Class', 'Alcohol', 'Malic acid', 'Ash', 'Alcalinity of ash', 'Magnesium', 'Total phenols',
                   'Flavanoids', 'Nonflavanoid phenols', 'Proanthocyanins', 'Color intensity', 'Hue',
                   'OD280/OD315 of diluted wines', 'Proline']
wX, wy = df_wine.iloc[:, 1:].values, df_wine.iloc[:, 0].values
wX_train, wX_test, wy_train, wy_test = train_test_split(wX, wy, test_size=0.3, random_state=0)

# Standardization
wsc = StandardScaler()
wsc.fit(wX_train)
wX_train_std = wsc.transform(wX_train)
wX_test_std = wsc.transform(wX_test)

# load Housing data from UCI
# for details, see https://archive.ics.uci.edu/ml/datasets/Housing
#df_housing = pd.read_csv('https://archive.ics.uci.edu/ml/machine-learning-databases/housing/housing.data',
#                         header=None, sep='\s+')
df_housing = pd.read_csv('housing.txt', header=None, sep='\s+')
df_housing.columns = ['CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM', 'AGE', 'DIS', 'RAD', 'TAX', 'PTRATIO',
                      'B', 'LSTAT', 'MEDV']
hX = df_housing.iloc[:, :-1].values
hy = df_housing['MEDV'].values
hX_train, hX_test, hy_train, hy_test = train_test_split(hX, hy, test_size=0.3, random_state=0)
# plot for classification region
from matplotlib.colors import ListedColormap
import matplotlib.pyplot as plt
%matplotlib inline

def plot_decision_regions(X, y, classifier, test_idx=None, resolution=0.02):
    markers = ('s', 'x', 'o', '^', 'v')
    colors = ('red', 'blue', 'lightgreen', 'gray', 'cyan')
    cmap = ListedColormap(colors[:len(np.unique(y))])
    
    x1_min, x1_max = X[:, 0].min() - 1, X[:, 0].max() + 1
    x2_min, x2_max = X[:, 1].min() - 1, X[:, 1].max() + 1
    xx1, xx2 = np.meshgrid(np.arange(x1_min, x1_max, resolution),
                           np.arange(x2_min, x2_max, resolution))

    Z = classifier.predict(np.array([xx1.ravel(), xx2.ravel()]).T)
    Z = Z.reshape(xx1.shape)
    plt.contourf(xx1, xx2, Z, alpha=0.4, cmap=cmap)
    plt.xlim(xx1.min(), xx1.max())
    plt.ylim(xx2.min(), xx2.max())
    
    for idx, cl in enumerate(np.unique(y)):
        plt.scatter(x=X[y == cl, 0], y=X[y == cl, 1], alpha=0.8, c=cmap(idx), marker=markers[idx], label=cl)
    
    if test_idx:
        X_test, y_test = X[test_idx, :], y[test_idx]
        plt.scatter(X_test[:, 0], X_test[:, 1], c='', alpha=1.0, linewidths=1, marker='o', s=55, label='test set')
# Linear SVM using Iris
from sklearn.svm import SVC

svm_iris = SVC(kernel='linear', C=1.0, random_state=0)
svm_iris.fit(X_train_std, y_train)
svm_iris_pre = svm_iris.predict(X_test_std)
print('Accuracy: %.2f' % accuracy_score(y_test, svm_iris_pre) )
plot_decision_regions(X_combined_std, y_combined, classifier=svm_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.98
<matplotlib.legend.Legend at 0x21a8839ffd0>

# random XOR data

np.random.seed(0)
X_xor = np.random.randn(200, 2)
y_xor = np.logical_xor(X_xor[:, 0] > 0, X_xor[:, 1] > 0)
y_xor = np.where(y_xor, 1, -1)
plt.scatter(X_xor[y_xor==1, 0], X_xor[y_xor==1, 1], c='b', marker='x', label='1')
plt.scatter(X_xor[y_xor==-1, 0], X_xor[y_xor==-1, 1], c='r', marker='s', label='-1')
plt.xlim([-3, 3])
plt.ylim([-3, 3])
plt.legend(loc='best')
<matplotlib.legend.Legend at 0x21a8984a6a0>

# Linear SVM using random XOR data
from sklearn.svm import SVC

svm_xor = SVC(kernel='linear', C=1.0, random_state=0)
svm_xor.fit(X_xor, y_xor)
plot_decision_regions(X_xor, y_xor, classifier=svm_xor)
plt.legend(loc='upper left')
<matplotlib.legend.Legend at 0x21a895b7d30>

# RBF-kernel SVM using random XOR data
from sklearn.svm import SVC

svm_rbf_xor = SVC(kernel='rbf', random_state=0, gamma=0.10, C=10.0)
svm_rbf_xor.fit(X_xor, y_xor)
plot_decision_regions(X_xor, y_xor, classifier=svm_rbf_xor)
plt.legend(loc='upper left')
<matplotlib.legend.Legend at 0x21a89b9a4a8>

# RBF-kernel SVM using Iris
from sklearn.svm import SVC

# try various gamma parameter (such as 200.0)
svm_rbf_iris = SVC(kernel='rbf', random_state=0, gamma=0.2, C=1.0)
svm_rbf_iris.fit(X_train_std, y_train)
svm_rbf_iris_pre = svm_rbf_iris.predict(X_test_std)
print('Accuracy: %.2f' % accuracy_score(y_test, svm_rbf_iris_pre) )
plot_decision_regions(X_combined_std, y_combined, classifier=svm_rbf_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.98
<matplotlib.legend.Legend at 0x21a88467e80>

# Decision tree
from sklearn.tree import DecisionTreeClassifier

tree_iris = DecisionTreeClassifier(criterion='entropy', max_depth=3, random_state=0)
tree_iris.fit(X_train, y_train)
tree_iris_pre = tree_iris.predict(X_test)
print('Accuracy: %.2f' % accuracy_score(y_test, tree_iris_pre) )
plot_decision_regions(X_combined, y_combined, classifier=tree_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.98
<matplotlib.legend.Legend at 0x21a89cd98d0>

# Random forest
from sklearn.ensemble import RandomForestClassifier

forest_iris = RandomForestClassifier(criterion='entropy', n_estimators=10, random_state=1)
forest_iris.fit(X_train, y_train)
forest_iris_pre = forest_iris.predict(X_test)
print('Accuracy: %.2f' % accuracy_score(y_test, forest_iris_pre) )
plot_decision_regions(X_combined, y_combined, classifier=forest_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.96
<matplotlib.legend.Legend at 0x21a8a568b70>

# Perceptron
from sklearn.linear_model import Perceptron

ppn_iris = Perceptron(n_iter=40, eta0=0.1, random_state=0, shuffle=True)
ppn_iris.fit(X_train_std, y_train)
ppn_iris_pre = ppn_iris.predict(X_test_std)
print('Accuracy: %.2f' % accuracy_score(y_test, ppn_iris_pre) )
plot_decision_regions(X_combined_std, y_combined, classifier=ppn_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.91
<matplotlib.legend.Legend at 0x21a8a55a9b0>

# Logistic regression
from sklearn.linear_model import LogisticRegression

lr_iris = LogisticRegression(C=1000.0, random_state=0)
lr_iris.fit(X_train_std, y_train)
lr_iris_pre = lr_iris.predict(X_test_std)
print('Accuracy: %.2f' % accuracy_score(y_test, lr_iris_pre) )
plot_decision_regions(X_combined_std, y_combined, classifier=lr_iris, test_idx=range(105,150))
plt.xlabel('petal length [standardized]')
plt.ylabel('petal width [standardized]')
plt.legend(loc='upper left')
Accuracy: 0.98
<matplotlib.legend.Legend at 0x21a89b9e780>

# PCA on wine data
from sklearn.decomposition import PCA
from sklearn.linear_model import LogisticRegression

pca_wine = PCA(n_components=2)
wX_train_pca = pca_wine.fit_transform(wX_train_std)
wX_test_pca = pca_wine.transform(wX_test_std)

# try Logistic regression
# first, test without PCA
lr_wine = LogisticRegression()
lr_wine.fit(wX_train_std, wy_train)
lr_wine_pre = lr_wine.predict(wX_test_std)
print('Accuracy: %.2f' % accuracy_score(wy_test, lr_wine_pre) )

# then, test with PCA
lr_wine_pca = LogisticRegression()
lr_wine_pca.fit(wX_train_pca, wy_train)
lr_wine_pca_pre = lr_wine_pca.predict(wX_test_pca)
print('Accuracy: %.2f' % accuracy_score(wy_test, lr_wine_pca_pre) )
plot_decision_regions(wX_train_pca, wy_train, classifier=lr_wine_pca)
plt.xlabel('PC1')
plt.ylabel('PC2')
plt.legend(loc='lower left')
Accuracy: 1.00
Accuracy: 0.98
<matplotlib.legend.Legend at 0x21a884476a0>

# LDA on wine data
from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA
from sklearn.linear_model import LogisticRegression

lda_wine = LDA(n_components=2)
wX_train_lda = lda_wine.fit_transform(wX_train_std, wy_train)
wX_test_lda = lda_wine.transform(wX_test_std)

# try Logistic regression
lr_wine_lda = LogisticRegression()
lr_wine_lda.fit(wX_train_lda, wy_train)
lr_wine_lda_pre = lr_wine_pca.predict(wX_test_lda)
print('Accuracy: %.2f' % accuracy_score(wy_test, lr_wine_lda_pre) )
plot_decision_regions(wX_test_lda, wy_test, classifier=lr_wine_lda)
plt.xlabel('PC1')
plt.ylabel('PC2')
plt.legend(loc='lower left')
Accuracy: 1.00
<matplotlib.legend.Legend at 0x21a8a55c0f0>

# kernel PCA on wine data
from sklearn.decomposition import KernelPCA
from sklearn.linear_model import LogisticRegression

# using Gaussian kernel
kpca_wine = KernelPCA(n_components=2, kernel='rbf')
wX_train_kpca = kpca_wine.fit_transform(wX_train_std)
wX_test_kpca = kpca_wine.transform(wX_test_std)

# try Logistic regression
lr_wine_kpca = LogisticRegression()
lr_wine_kpca.fit(wX_train_kpca, wy_train)
lr_wine_kpca_pre = lr_wine_kpca.predict(wX_test_kpca)
print('Accuracy: %.2f' % accuracy_score(wy_test, lr_wine_kpca_pre) )
plot_decision_regions(wX_train_kpca, wy_train, classifier=lr_wine_kpca)
plt.xlabel('PC1')
plt.ylabel('PC2')
plt.legend(loc='lower left')
Accuracy: 1.00
<matplotlib.legend.Legend at 0x21a8b6e9f98>

# Feature importance by RandomForest
from sklearn.ensemble import RandomForestClassifier
from sklearn.feature_selection import SelectFromModel
from sklearn.linear_model import LogisticRegression
import matplotlib.pyplot as plt

feat_labels = df_wine.columns[1:]
feat_forest = RandomForestClassifier(n_estimators=10000, random_state=0)
feat_forest.fit(wX_train, wy_train)
feat_importances = feat_forest.feature_importances_
feat_indices = np.argsort(feat_importances)[::-1]

# show importance in descending order
for f in range(wX_train.shape[1]):
    print("%2d) %-*s %f" %
          (f + 1, 30, feat_labels[feat_indices[f]], feat_importances[feat_indices[f]]))

plt.bar(range(wX_train.shape[1]), feat_importances[feat_indices], color='lightblue', align='center')
plt.xticks(range(wX_train.shape[1]), feat_labels[feat_indices], rotation=90)
plt.xlim([-1, wX_train.shape[1]])
plt.title('Feature importances')

# select features and try logistic regression
sfm = SelectFromModel(feat_forest, prefit=True, threshold=0.15)
wX_train_rf = sfm.transform(wX_train)
wX_test_rf = sfm.transform(wX_test)
lr_wine_rf = LogisticRegression()
lr_wine_rf.fit(wX_train_rf, wy_train)
lr_wine_rf_pre = lr_wine_rf.predict(wX_test_rf)
print('Accuracy: %.2f' % accuracy_score(wy_test, lr_wine_rf_pre) )
 1) Color intensity                0.182483
 2) Proline                        0.158610
 3) Flavanoids                     0.150948
 4) OD280/OD315 of diluted wines   0.131987
 5) Alcohol                        0.106589
 6) Hue                            0.078243
 7) Total phenols                  0.060718
 8) Alcalinity of ash              0.032033
 9) Malic acid                     0.025400
10) Proanthocyanins                0.022351
11) Magnesium                      0.022078
12) Nonflavanoid phenols           0.014645
13) Ash                            0.013916
Accuracy: 0.94

# Regression analysis using Housing data
import matplotlib.pyplot as plt
import seaborn as sns

# plot scatter matrix for some attributes
sns.set(style='whitegrid', context='notebook')
cols = ['LSTAT', 'INDUS', 'NOX', 'RM', 'MEDV']
sns.pairplot(df_housing[cols])
<seaborn.axisgrid.PairGrid at 0x21a884c3358>

# plot correlation matrix for some attributes
housing_cm = np.corrcoef(df_housing[cols].values.T)
sns.set(font_scale=1.5)
housing_hm = sns.heatmap(housing_cm, cbar=True, annot=True, square=True, fmt='.2f', annot_kws={'size':15},
                         yticklabels=cols, xticklabels=cols)

# multivariate linear regression
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error
from sklearn.metrics import r2_score
import matplotlib.pyplot as plt

mlr_housing = LinearRegression()
mlr_housing.fit(hX_train, hy_train)
mlr_housing_train_pre = mlr_housing.predict(hX_train)
mlr_housing_test_pre = mlr_housing.predict(hX_test)

# residual plot
plt.scatter(mlr_housing_train_pre, mlr_housing_train_pre - hy_train, c='blue', marker='o', label='Training data')
plt.scatter(mlr_housing_test_pre, mlr_housing_test_pre - hy_test, c='lightgreen', marker='s', label='Test data')
plt.xlabel('Predicted values')
plt.ylabel('Residuals')
plt.legend(loc='upper left')
plt.hlines(y=0, xmin=-10, xmax=50, lw=2, color='red')
plt.xlim([-10, 50])
plt.ylim([-30, 30])

# mean squared error (MSE)
print('MSE train: %.3f, test: %.3f' % (mean_squared_error(hy_train, mlr_housing_train_pre),
                                       mean_squared_error(hy_test, mlr_housing_test_pre)))

# R^2
print('R^2 train: %.3f, test: %.3f' % (r2_score(hy_train, mlr_housing_train_pre),
                                       r2_score(hy_test, mlr_housing_test_pre)))
MSE train: 19.958, test: 27.196
R^2 train: 0.765, test: 0.673

# RandomForest regression
from sklearn.ensemble import RandomForestRegressor

rfr_housing = RandomForestRegressor(n_estimators=1000, criterion='mse', random_state=1)
rfr_housing.fit(hX_train, hy_train)
rfr_housing_train_pre = rfr_housing.predict(hX_train)
rfr_housing_test_pre = rfr_housing.predict(hX_test)

print('MSE train: %.3f, test: %.3f' % (mean_squared_error(hy_train, rfr_housing_train_pre),
                                       mean_squared_error(hy_test, rfr_housing_test_pre)))
print('R^2 train: %.3f, test: %.3f' % (r2_score(hy_train, rfr_housing_train_pre),
                                       r2_score(hy_test, rfr_housing_test_pre)))

# residual plot
plt.scatter(rfr_housing_train_pre, rfr_housing_train_pre - hy_train, c='black', marker='o', label='Training data')
plt.scatter(rfr_housing_test_pre, rfr_housing_test_pre - hy_test, c='lightgreen', marker='s', label='Test data')
plt.xlabel('Predicted values')
plt.ylabel('Residuals')
plt.legend(loc='upper left')
plt.hlines(y=0, xmin=-10, xmax=50, lw=2, color='red')
plt.xlim([-10, 50])
plt.ylim([-30, 30])
MSE train: 1.409, test: 14.354
R^2 train: 0.983, test: 0.828
(-30, 30)