In [9]:

    

import os
import sys
import numpy as np
import pandas as pd
from keras.utils import np_utils
from keras.optimizers import Adam
from keras.models import Sequential
from matplotlib import pyplot as plt
from sklearn.decomposition import PCA
from keras.callbacks import ModelCheckpoint
from sklearn.metrics import confusion_matrix
from sklearn.cross_validation import train_test_split
from keras.layers.pooling import GlobalAveragePooling1D
from keras.layers.normalization import BatchNormalization
from keras.layers import Dropout, Activation, Dense, Flatten
from keras.layers.convolutional import Convolution1D,AveragePooling1D,MaxPooling1D


    

In [2]:

    

Mn2_C = (pd.read_pickle('input_spectra/Mn2_Larger_Clean_Thin.pkl'))
Mn3_C = (pd.read_pickle('input_spectra/Mn3_Larger_Clean_Thin.pkl'))
Mn4_C = (pd.read_pickle('input_spectra/Mn4_Larger_Clean_Thin.pkl'))
Mn_All = (Mn2_C.append(Mn3_C, ignore_index=True)).append(Mn4_C, ignore_index=True)
Mn_All = np.array(Mn_All)

labels=[]
for i in range(0, len(Mn2_C)):
    labels.append(0)
for i in range(0, len(Mn3_C)):
    labels.append(1)
for i in range(0, len(Mn4_C)):
    labels.append(2)


    

In [14]:

    

#train-test split
X_train, X_test, y_train, y_test = train_test_split(Mn_All, labels, test_size=0.15, random_state=13)

#data augmentation using principal components
noise_aug = []
noise = np.copy(X_train)
mu = np.mean(noise, axis=0)
pca = PCA()
noise_model = pca.fit(noise)
nComp = 10
Xhat = np.dot(pca.transform(noise)[:,:nComp], pca.components_[:nComp,:])
noise_level = np.dot(pca.transform(noise)[:,nComp:], pca.components_[nComp:,:])
Xhat += mu
SNR = np.linspace(1,5,50)
for i  in range(len(SNR)):
    noise_aug.append(SNR[i]*noise_level + Xhat)
    j = 0
    for spectra in noise_aug[i]:
        noise_aug[i][j] = spectra/np.max(spectra)
        j += 1
X_train = np.array(noise_aug).reshape(50*2684,700)
y_train = [item for i in range(50) for item in y_train]

#cropping
X_train = X_train[:,100:600]
X_test = X_test[:,100:600]

#formatting for keras
X_train = np.array(X_train).astype('float32')
X_train = X_train.reshape(X_train.shape + (1,))
X_train -=  np.mean(X_train)
X_train /= np.max(X_train)
X_test = np.array(X_test).astype('float32')
X_test = X_test.reshape(X_test.shape + (1,))
X_test -= np.mean(X_test)   
X_test /= np.max(X_test)

y_train = np.array(y_train)
y_test = np.array(y_test)
y_train = np_utils.to_categorical(y_train)
y_test = np_utils.to_categorical(y_test)
num_classes = y_test.shape[1]

print("Total of "+str(num_classes)+" classes.")
print("Data mean-centered, normalized and hot-encoded.")
print("Total of "+str(len(X_train))+" training samples.")
plt.plot(X_test[0])
plt.plot(X_train[0])
plt.show()


    

    




Total of 3 classes.
Data mean-centered, normalized and hot-encoded.
Total of 134200 training samples.






    

In [5]:

    

model = Sequential()
activation = 'relu'
model.add(Convolution1D(2, 9, input_shape=(500,1), activation=activation))
model.add(AveragePooling1D())
model.add(BatchNormalization())

model.add(Convolution1D(2, 7, activation=activation))
model.add(AveragePooling1D())
model.add(BatchNormalization())

model.add(Convolution1D(4, 7, activation=activation))
model.add(AveragePooling1D())
model.add(BatchNormalization())

model.add(Convolution1D(4, 5, activation=activation))
model.add(AveragePooling1D())
model.add(BatchNormalization())

model.add(Convolution1D(8, 3, activation=activation))
model.add(AveragePooling1D())
model.add(BatchNormalization())

model.add(Dropout(0.10))
model.add(Convolution1D(3, 1))
model.add(GlobalAveragePooling1D())

model.add(Activation('softmax', name='loss'))

model.compile(loss='categorical_crossentropy', optimizer='adam', metrics=['accuracy'])

print(model.summary())
print("CNN Model created.")


    

    




____________________________________________________________________________________________________
Layer (type)                     Output Shape          Param #     Connected to                     
====================================================================================================
convolution1d_1 (Convolution1D)  (None, 492, 2)        20          convolution1d_input_1[0][0]      
____________________________________________________________________________________________________
averagepooling1d_1 (AveragePooli (None, 246, 2)        0           convolution1d_1[0][0]            
____________________________________________________________________________________________________
batchnormalization_1 (BatchNorma (None, 246, 2)        8           averagepooling1d_1[0][0]         
____________________________________________________________________________________________________
convolution1d_2 (Convolution1D)  (None, 240, 2)        30          batchnormalization_1[0][0]       
____________________________________________________________________________________________________
averagepooling1d_2 (AveragePooli (None, 120, 2)        0           convolution1d_2[0][0]            
____________________________________________________________________________________________________
batchnormalization_2 (BatchNorma (None, 120, 2)        8           averagepooling1d_2[0][0]         
____________________________________________________________________________________________________
convolution1d_3 (Convolution1D)  (None, 114, 4)        60          batchnormalization_2[0][0]       
____________________________________________________________________________________________________
averagepooling1d_3 (AveragePooli (None, 57, 4)         0           convolution1d_3[0][0]            
____________________________________________________________________________________________________
batchnormalization_3 (BatchNorma (None, 57, 4)         16          averagepooling1d_3[0][0]         
____________________________________________________________________________________________________
convolution1d_4 (Convolution1D)  (None, 53, 4)         84          batchnormalization_3[0][0]       
____________________________________________________________________________________________________
averagepooling1d_4 (AveragePooli (None, 26, 4)         0           convolution1d_4[0][0]            
____________________________________________________________________________________________________
batchnormalization_4 (BatchNorma (None, 26, 4)         16          averagepooling1d_4[0][0]         
____________________________________________________________________________________________________
convolution1d_5 (Convolution1D)  (None, 24, 8)         104         batchnormalization_4[0][0]       
____________________________________________________________________________________________________
averagepooling1d_5 (AveragePooli (None, 12, 8)         0           convolution1d_5[0][0]            
____________________________________________________________________________________________________
batchnormalization_5 (BatchNorma (None, 12, 8)         32          averagepooling1d_5[0][0]         
____________________________________________________________________________________________________
dropout_1 (Dropout)              (None, 12, 8)         0           batchnormalization_5[0][0]       
____________________________________________________________________________________________________
convolution1d_6 (Convolution1D)  (None, 12, 3)         27          dropout_1[0][0]                  
____________________________________________________________________________________________________
globalaveragepooling1d_1 (Global (None, 3)             0           convolution1d_6[0][0]            
____________________________________________________________________________________________________
loss (Activation)                (None, 3)             0           globalaveragepooling1d_1[0][0]   
====================================================================================================
Total params: 405
Trainable params: 365
Non-trainable params: 40
____________________________________________________________________________________________________
None
CNN Model created.

In [7]:

    

#Params
epochs = 1000
batch_size = 512
seed = 7

# fit and run our model
np.random.seed(seed)
best_model_file = "weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch{epoch:02d}-val_acc{val_acc:.3f}_cnn.h5"
best_model = ModelCheckpoint(best_model_file, monitor='val_acc', verbose = 1, save_best_only = True)
hist = model.fit(X_train,
                 y_train,
                 validation_data=(X_test, y_test),
                 nb_epoch=epochs,
                 batch_size=batch_size,
                 callbacks = [best_model],
                 shuffle = True,
                 verbose=1)
print("done")


    

    




Train on 134200 samples, validate on 474 samples
Epoch 1/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.8008 - acc: 0.6255     - ETA: 1003s - loss: 1.1087 - acc: 0.3622 - ETA: 333s - loss: 1.0384 - acc: 0.4689 - ETA: 7s - loss: 0.8247 - acc: 0.6100Epoch 00000: val_acc improved from -inf to 0.33122, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch00-val_acc0.331_cnn.h5
134200/134200 [==============================] - 81s - loss: 0.8007 - acc: 0.6256 - val_loss: 1.2209 - val_acc: 0.3312
Epoch 2/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.3814 - acc: 0.8910 - ETA: 1s - loss: 0.4338 - acc: 0.8729 - ETA: 0s - loss: 0.4141 - acc: 0.8806Epoch 00001: val_acc improved from 0.33122 to 0.53586, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch01-val_acc0.536_cnn.h5
134200/134200 [==============================] - 3s - loss: 0.3813 - acc: 0.8910 - val_loss: 0.8770 - val_acc: 0.5359
Epoch 3/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.2453 - acc: 0.9259 - ETA: 2s - loss: 0.2750 - acc: 0.9190 - ETA: 1s - loss: 0.2677 - acc: 0.9204 - ETA: 1s - loss: 0.2613 - acc: 0.9218 - ETA: 0s - loss: 0.2489 - acc: 0.9249Epoch 00002: val_acc improved from 0.53586 to 0.89451, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch02-val_acc0.895_cnn.h5
134200/134200 [==============================] - 2s - loss: 0.2452 - acc: 0.9259 - val_loss: 0.2923 - val_acc: 0.8945
Epoch 4/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1908 - acc: 0.9405Epoch 00003: val_acc improved from 0.89451 to 0.94937, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch03-val_acc0.949_cnn.h5
134200/134200 [==============================] - 3s - loss: 0.1908 - acc: 0.9405 - val_loss: 0.1606 - val_acc: 0.9494
Epoch 5/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1633 - acc: 0.9472 - ETA: 0s - loss: 0.1660 - acc: 0.9464 - ETA: 0s - loss: 0.1660 - acc: 0.9463Epoch 00004: val_acc improved from 0.94937 to 0.97046, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch04-val_acc0.970_cnn.h5
134200/134200 [==============================] - 3s - loss: 0.1633 - acc: 0.9472 - val_loss: 0.1023 - val_acc: 0.9705
Epoch 6/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1487 - acc: 0.9500 - ETA: 0s - loss: 0.1508 - acc: 0.9493Epoch 00005: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1487 - acc: 0.9500 - val_loss: 0.1225 - val_acc: 0.9599
Epoch 7/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1399 - acc: 0.9523Epoch 00006: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1399 - acc: 0.9523 - val_loss: 0.0851 - val_acc: 0.9684
Epoch 8/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1322 - acc: 0.9549 - ETA: 1s - loss: 0.1332 - acc: 0.9546 - ETA: 1s - loss: 0.1333 - acc: 0.9545 - ETA: 1s - loss: 0.1327 - acc: 0.9549 - ETA: 1s - loss: 0.1323 - acc: 0.9551 - ETA: 0s - loss: 0.1321 - acc: 0.9551Epoch 00007: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1322 - acc: 0.9549 - val_loss: 0.0920 - val_acc: 0.9684
Epoch 9/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1266 - acc: 0.9564 - ETA: 2s - loss: 0.1280 - acc: 0.9562Epoch 00008: val_acc improved from 0.97046 to 0.98312, saving model to weights/CNN_Noise_DataAug/highest_val_acc_weights_epoch08-val_acc0.983_cnn.h5
134200/134200 [==============================] - 3s - loss: 0.1266 - acc: 0.9563 - val_loss: 0.0668 - val_acc: 0.9831
Epoch 10/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1231 - acc: 0.9575 - ETA: 1s - loss: 0.1208 - acc: 0.9592Epoch 00009: val_acc did not improve
134200/134200 [==============================] - 2s - loss: 0.1232 - acc: 0.9575 - val_loss: 0.0787 - val_acc: 0.9789
Epoch 11/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1188 - acc: 0.9594 - ETA: 1s - loss: 0.1188 - acc: 0.9596Epoch 00010: val_acc did not improve
134200/134200 [==============================] - 2s - loss: 0.1188 - acc: 0.9594 - val_loss: 0.0972 - val_acc: 0.9662
Epoch 12/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1133 - acc: 0.9611 - ETA: 0s - loss: 0.1134 - acc: 0.9613 - ETA: 0s - loss: 0.1136 - acc: 0.9614Epoch 00011: val_acc did not improve
134200/134200 [==============================] - 2s - loss: 0.1133 - acc: 0.9611 - val_loss: 0.0846 - val_acc: 0.9768
Epoch 13/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1100 - acc: 0.9624 - ETA: 0s - loss: 0.1100 - acc: 0.9624Epoch 00012: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1100 - acc: 0.9624 - val_loss: 0.0673 - val_acc: 0.9789
Epoch 14/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1070 - acc: 0.9633 - ETA: 2s - loss: 0.1053 - acc: 0.9648 - ETA: 2s - loss: 0.1064 - acc: 0.9643 - ETA: 1s - loss: 0.1071 - acc: 0.9630 - ETA: 0s - loss: 0.1071 - acc: 0.9631Epoch 00013: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1070 - acc: 0.9633 - val_loss: 0.0837 - val_acc: 0.9684
Epoch 15/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1053 - acc: 0.9633 - ETA: 3s - loss: 0.1055 - acc: 0.9628 - ETA: 2s - loss: 0.1032 - acc: 0.9637 - ETA: 1s - loss: 0.1048 - acc: 0.9635Epoch 00014: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1053 - acc: 0.9633 - val_loss: 0.0764 - val_acc: 0.9662
Epoch 16/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1026 - acc: 0.9638 - ETA: 2s - loss: 0.1064 - acc: 0.9628 - ETA: 1s - loss: 0.1037 - acc: 0.9641 - ETA: 0s - loss: 0.1029 - acc: 0.9640Epoch 00015: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1027 - acc: 0.9638 - val_loss: 0.0818 - val_acc: 0.9726
Epoch 17/1000
134144/134200 [============================>.] - ETA: 0s - loss: 0.1008 - acc: 0.9646 - ETA: 3s - loss: 0.1050 - acc: 0.9632 - ETA: 0s - loss: 0.1007 - acc: 0.9649 - ETA: 0s - loss: 0.1008 - acc: 0.9648 - ETA: 0s - loss: 0.1006 - acc: 0.9649 - ETA: 0s - loss: 0.1007 - acc: 0.9648Epoch 00016: val_acc did not improve
134200/134200 [==============================] - 3s - loss: 0.1008 - acc: 0.9646 - val_loss: 0.0795 - val_acc: 0.9705
Epoch 18/1000
 90112/134200 [===================>..........] - ETA: 0s - loss: 0.0977 - acc: 0.9651 - ETA: 1s - loss: 0.0979 - acc: 0.9649 - ETA: 1s - loss: 0.0976 - acc: 0.9651





    




---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
<ipython-input-7-63ad0b4f6090> in <module>()
     15                  callbacks = [best_model],
     16                  shuffle = True,
---> 17                  verbose=1)
     18 
     19 print("done")

/home/mike/anaconda3/lib/python3.5/site-packages/keras/models.py in fit(self, x, y, batch_size, nb_epoch, verbose, callbacks, validation_split, validation_data, shuffle, class_weight, sample_weight, initial_epoch, **kwargs)
    670                               class_weight=class_weight,
    671                               sample_weight=sample_weight,
--> 672                               initial_epoch=initial_epoch)
    673 
    674     def evaluate(self, x, y, batch_size=32, verbose=1,

/home/mike/anaconda3/lib/python3.5/site-packages/keras/engine/training.py in fit(self, x, y, batch_size, nb_epoch, verbose, callbacks, validation_split, validation_data, shuffle, class_weight, sample_weight, initial_epoch)
   1190                               val_f=val_f, val_ins=val_ins, shuffle=shuffle,
   1191                               callback_metrics=callback_metrics,
-> 1192                               initial_epoch=initial_epoch)
   1193 
   1194     def evaluate(self, x, y, batch_size=32, verbose=1, sample_weight=None):

/home/mike/anaconda3/lib/python3.5/site-packages/keras/engine/training.py in _fit_loop(self, f, ins, out_labels, batch_size, nb_epoch, verbose, callbacks, val_f, val_ins, shuffle, callback_metrics, initial_epoch)
    890                 batch_logs['size'] = len(batch_ids)
    891                 callbacks.on_batch_begin(batch_index, batch_logs)
--> 892                 outs = f(ins_batch)
    893                 if not isinstance(outs, list):
    894                     outs = [outs]

/home/mike/anaconda3/lib/python3.5/site-packages/keras/backend/tensorflow_backend.py in __call__(self, inputs)
   1898         session = get_session()
   1899         updated = session.run(self.outputs + [self.updates_op],
-> 1900                               feed_dict=feed_dict)
   1901         return updated[:len(self.outputs)]
   1902 

/home/mike/.local/lib/python3.5/site-packages/tensorflow/python/client/session.py in run(self, fetches, feed_dict, options, run_metadata)
    764     try:
    765       result = self._run(None, fetches, feed_dict, options_ptr,
--> 766                          run_metadata_ptr)
    767       if run_metadata:
    768         proto_data = tf_session.TF_GetBuffer(run_metadata_ptr)

/home/mike/.local/lib/python3.5/site-packages/tensorflow/python/client/session.py in _run(self, handle, fetches, feed_dict, options, run_metadata)
    962     if final_fetches or final_targets:
    963       results = self._do_run(handle, final_targets, final_fetches,
--> 964                              feed_dict_string, options, run_metadata)
    965     else:
    966       results = []

/home/mike/.local/lib/python3.5/site-packages/tensorflow/python/client/session.py in _do_run(self, handle, target_list, fetch_list, feed_dict, options, run_metadata)
   1012     if handle is None:
   1013       return self._do_call(_run_fn, self._session, feed_dict, fetch_list,
-> 1014                            target_list, options, run_metadata)
   1015     else:
   1016       return self._do_call(_prun_fn, self._session, handle, feed_dict,

/home/mike/.local/lib/python3.5/site-packages/tensorflow/python/client/session.py in _do_call(self, fn, *args)
   1019   def _do_call(self, fn, *args):
   1020     try:
-> 1021       return fn(*args)
   1022     except errors.OpError as e:
   1023       message = compat.as_text(e.message)

/home/mike/.local/lib/python3.5/site-packages/tensorflow/python/client/session.py in _run_fn(session, feed_dict, fetch_list, target_list, options, run_metadata)
   1001         return tf_session.TF_Run(session, options,
   1002                                  feed_dict, fetch_list, target_list,
-> 1003                                  status, run_metadata)
   1004 
   1005     def _prun_fn(session, handle, feed_dict, fetch_list):

KeyboardInterrupt: 

In [19]:

    

#summarize history for accuracy
plt.figure(figsize=(15, 5))
plt.rcParams.update({'font.size': 16})

plt.subplot(1, 2, 1)
plt.plot(hist.history['acc'], linewidth = 3)
plt.title('Model Training Accuracy')
plt.ylabel('Training Accuracy')
plt.xlabel('Epoch')

# summarize history for loss
plt.subplot(1, 2, 2)
plt.plot(hist.history['loss'], linewidth = 3)
plt.title('Model Training Loss')
plt.ylabel('Cross Entropy Loss')
plt.xlabel('Epoch')
plt.savefig('figures/training_accuracy.png')
plt.show()

plt.figure(figsize=(10, 8))

plt.plot(hist.history['val_acc'], linewidth = 3)
plt.plot(hist.history['acc'], linewidth = 3)
plt.title('Model Accuracy')
plt.ylabel('Accuracy')
plt.xlabel('Epoch')
plt.legend(['Test', 'Train'], loc='lower right')
plt.savefig('figures/test_accuracy.png')

plt.show()


    

In [10]:

    

y_test_pred, y_test_labels=[], []
for i in range(0, len(X_test)):
    y_test_pred.append(np.argmax(model.predict(X_test[i:i+1])))
    y_test_labels.append(np.argmax(y_test[i]))
print("Confusion Matrix of Test Set")
conf_matrix = pd.DataFrame(confusion_matrix(y_pred=y_test_pred, y_true=y_test_labels))
conf_matrix.columns = ["Mn2+", "Mn3+", "Mn4+" ]
conf_matrix = pd.DataFrame.transpose(conf_matrix)
conf_matrix.columns = ["Mn2+", "Mn3+", "Mn4+" ]
conf_matrix = pd.DataFrame.transpose(conf_matrix)
print(conf_matrix)
scores = model.evaluate(X_test, y_test, verbose=0)
print("Accuracy: %.2f%%" % (scores[1]*100))


    

    




Confusion Matrix of Test Set
      Mn2+  Mn3+  Mn4+
Mn2+   148     8     1
Mn3+     0   136     5
Mn4+     0     2   174
Accuracy: 96.62%