<jupyter_start><jupyter_text>Features Extraction & Classification [Small Samples]In this notebook we will 1. **Features extraction**2. **CNN on Augmented data**3. **TSNE Representation**The afored mentionded steps are applied on the 800 ms records samples.<jupyter_code>import helpers 
import os
import pandas as pd
import IPython as IP
import struct
import matplotlib.pyplot as plt
import numpy as np
import math
import pickle
import sys
import IPython as IP
import numpy as np
import librosa
import librosa.display
import random
import warnings
import tensorflow as tf
import soundfile as sf
import tqdm
from   matplotlib                    import cm
from   tqdm                          import tqdm
from   pydub                         import AudioSegment
from   pydub.utils                   import make_chunks
from   pydub.silence                 import split_on_silencefrom   
from   numba                         import cuda
from   datetime                      import datetime
from   include                       import helpers
from   tensorflow                    import keras
from   sklearn.manifold              import TSNE
from   sklearn.preprocessing         import LabelEncoder
from   sklearn.model_selection       import train_test_split
from   tensorflow.keras.utils        import to_categorical
from   tensorflow.keras              import backend as keras_backend
from   tensorflow.keras.models       import Sequential, load_model
from   tensorflow.keras.layers       import Dense, SpatialDropout2D, Activation, Conv2D, MaxPooling2D, BatchNormalization, GlobalAveragePooling2D, LeakyReLU
from   tensorflow.keras.optimizers   import Adam
from   tensorflow.keras.utils        import to_categorical
from   tensorflow.keras.callbacks    import ModelCheckpoint 
from   tensorflow.keras.regularizers import l2
from   sklearn.preprocessing         import LabelEncoder
from   sklearn.metrics               import classification_report, confusion_matrix
from   helpers                       import *
# Tricks to solved :Failed to get convolution algorithm. 
# This is probably because cuDNN failed to initialize
# See https://stackoverflow.com/questions/43147983/could-not-create-cudnn-handle-cudnn-status-internal-error
gpus= tf.config.experimental.list_physical_devices('GPU')
tf.config.experimental.set_memory_growth(gpus[0], True)
warnings.filterwarnings('ignore')
%matplotlib inline
%reload_ext autoreload
%autoreload 2
path_to_output_sounds = r"sounds_output"<jupyter_output><empty_output><jupyter_text>Part1. Features extraction<jupyter_code>def get_mfcc(file_path, mfcc_max_padding=0, n_mfcc=40, robots_noise = None, noise_amp = 1):
    try:
        # Load audio file
        y, sr = librosa.load(file_path,sr=None)
        
        if robots_noise != None : 
            y_n, _ = librosa.load(robots_noise)
            y = (y + noise_amp * y_n)/(noise_amp + 1)
            
        # Normalize audio data between -1 and 1
        normalized_y = librosa.util.normalize(y)

        # Compute MFCC coefficients
        mfcc = librosa.feature.mfcc(y=normalized_y, sr=sr, n_mfcc=n_mfcc)

        # Normalize MFCC between -1 and 1
        normalized_mfcc = librosa.util.normalize(mfcc)

        # Should we require padding
        shape = normalized_mfcc.shape[1]

        if (shape < mfcc_max_padding):
            pad_width = mfcc_max_padding - shape
            normalized_mfcc = np.pad(normalized_mfcc, 
                        pad_width=((0, 0), (0, pad_width)), 
                        mode     ='constant', 
                        constant_values=(0,))

    except Exception as e:
        print("Error parsing wavefile: ", e)
        return None 
    return normalized_mfcc
path_augmented     = "sounds_output/SMALLER_SEQUENCE/augmented" 
list_manip         = ["Vertical", "Horizontal", "180-left","Random"]
path_save_features = "sounds_output/SMALLER_SEQUENCE/Features"
# Water have been removed
objects      = ['Rice', 'Corn_Flakes', 'Caramel', 'Vitamines','No-Sounds']
features      = []
labels        = []
frames_max    = 0
n_mfcc        = 40

for object_ in objects :
    # Loop over manip
    print("Processing object_ : ",object_)
    class_label = object_
    if object_ == "No-Sounds" : 
        object_manip_path = os.path.join(path_augmented,object_)
        for idx1, file_sample in enumerate(tqdm(os.listdir(object_manip_path))) :
            input_path  = os.path.join(object_manip_path,file_sample)     
            # Extract MFCCs (do not add padding)
            mfccs = get_mfcc(input_path, 0, n_mfcc)

            # Save current frame count
            num_frames = mfccs.shape[1]

            # Add row (feature / label)
            features.append(mfccs)
            labels.append(class_label)

            # Update frames maximum
            if (num_frames > frames_max):
                frames_max = num_frames
    else : 
        for manip in list_manip : 
            print("Processing manip : ",manip)
            object_manip_path = os.path.join(path_augmented,object_,manip)

            for idx1, file_sample in enumerate(tqdm(os.listdir(object_manip_path))) :
                input_path  = os.path.join(object_manip_path,file_sample)     
                # Extract MFCCs (do not add padding)
                mfccs = get_mfcc(input_path, 0, n_mfcc)

                # Save current frame count
                num_frames = mfccs.shape[1]

                # Add row (feature / label)
                features.append(mfccs)
                labels.append(class_label)

                # Update frames maximum
                if (num_frames > frames_max):
                    frames_max = num_frames
# Verify that all the size are egal
padded = []
mfcc_max_padding = frames_max
count = 0 
for i in range(len(features)):
    size = len(features[i][0])
    if (size < mfcc_max_padding):
        print(size)
        count += 1
print(count)    
# Convert features (X) and labels (y) to Numpy arrays
X = np.array(features)
y = np.array(labels)
os.path.join(path_save_features,"Features.npy")
os.path.join(path_save_features,"labels.npy")
# Optionally save the features to disk
np.save(os.path.join(path_save_features,"Features.npy"), X)
np.save(os.path.join(path_save_features,"labels.npy"), y)<jupyter_output><empty_output><jupyter_text>PART 2 CNN on Augmented data<jupyter_code>indexes = []
total = X.shape[0]
indexes = list(range(0, total))

# Randomize indexes
random.shuffle(indexes)

# Divide the indexes into Train and Test
test_split_pct = 20
split_offset   = math.floor(test_split_pct * total / 100)

# Split the metadata
test_split_idx  = indexes[0:split_offset]
train_split_idx = indexes[split_offset:total]


# Split the features with the same indexes
X_test  = np.take(X, test_split_idx, axis=0)
y_test  = np.take(y, test_split_idx, axis=0)
X_train = np.take(X, train_split_idx, axis=0)
y_train = np.take(y, train_split_idx, axis=0)

# Print status
print("X test shape: {} \t X train shape: {}".format(X_test.shape, X_train.shape))
print("y test shape: {} \t\t y train shape: {}".format(y_test.shape, y_train.shape))
# One hot encode labels
le = LabelEncoder()
y_test_encoded = to_categorical(le.fit_transform(y_test))
y_train_encoded = to_categorical(le.fit_transform(y_train))
# Reshape data
num_rows     = 40 # MFCC
num_columns  = 69 # frames_max
num_channels = 1

# Reshape to fit the network input (channel last!)
X_train = X_train.reshape(X_train.shape[0], num_rows, num_columns, num_channels)
X_test  = X_test.reshape(X_test.shape[0], num_rows, num_columns, num_channels)

# Total number of labels to predict (equal to the network output nodes)
num_labels = y_train_encoded.shape[1]
print("X test shape: {} \t X train shape: {}".format(X_test.shape, X_train.shape))
print("y test shape: {} \t\t y train shape: {}".format(y_test.shape, y_train.shape))
def create_model(spatial_dropout_rate_1=0, spatial_dropout_rate_2=0, l2_rate=0):
    
    # Create a secquential object
    model = Sequential()

    # Conv 1
    model.add(Conv2D(filters=32, 
                     kernel_size=(3, 3), 
                     kernel_regularizer=l2(l2_rate), 
                     input_shape=(num_rows, num_columns, num_channels)))
    model.add(LeakyReLU(alpha=0.1))
    model.add(BatchNormalization())

    model.add(SpatialDropout2D(spatial_dropout_rate_1))
    model.add(Conv2D(filters=32, 
                     kernel_size=(3, 3), 
                     kernel_regularizer=l2(l2_rate)))
    model.add(LeakyReLU(alpha=0.1))
    model.add(BatchNormalization())


    # Max Pooling #1
    model.add(MaxPooling2D(pool_size=(2, 2)))

    model.add(SpatialDropout2D(spatial_dropout_rate_1))
    model.add(Conv2D(filters=64, 
                     kernel_size=(3, 3), 
                     kernel_regularizer=l2(l2_rate)))
    model.add(LeakyReLU(alpha=0.1))
    model.add(BatchNormalization())

    model.add(SpatialDropout2D(spatial_dropout_rate_2))
    model.add(Conv2D(filters=64, 
                     kernel_size=(3,3), 
                     kernel_regularizer=l2(l2_rate)))
    model.add(LeakyReLU(alpha=0.1))
    model.add(BatchNormalization())
    
   
    # Reduces each h×w feature map to a single number by taking the average of all h,w values.
    model.add(GlobalAveragePooling2D())


    # Softmax output
    model.add(Dense(num_labels, activation='softmax'))
    
    return model

# Regularization rates
spatial_dropout_rate_1 = 0.07
spatial_dropout_rate_2 = 0.14
l2_rate = 0.0005

model = create_model(spatial_dropout_rate_1, spatial_dropout_rate_2, l2_rate)
adam = Adam(lr=1e-4, beta_1=0.99, beta_2=0.999)
model.compile(
    loss='categorical_crossentropy', 
    metrics=['accuracy'], 
    optimizer=adam)

model.summary()
models_path = "sounds_output/SMALLER_SEQUENCE/Model"
model_file = 'one_second_train-augmented.hdf5'
model_path = os.path.join(models_path, model_file)
num_epochs = 250
num_batch_size = 128

# Save checkpoints
checkpointer = ModelCheckpoint(filepath=model_path, 
                               verbose=1, 
                               save_best_only=True)
start = datetime.now()
history = model.fit(X_train, 
                    y_train_encoded, 
                    batch_size=num_batch_size, 
                    epochs=num_epochs, 
                    validation_split=1/12.,
                    callbacks=[checkpointer], 
                    verbose=1)

duration = datetime.now() - start
print("Training completed in time: ", duration)
def evaluate_model(model, X_train, y_train, X_test, y_test):
    train_score = model.evaluate(X_train, y_train, verbose=0)
    test_score = model.evaluate(X_test, y_test, verbose=0)
    return train_score, test_score


def model_evaluation_report(model, X_train, y_train, X_test, y_test, calc_normal=True):
    dash = '-' * 38

    # Compute scores
    train_score, test_score = evaluate_model(model, X_train, y_train, X_test, y_test)

    # Pint Train vs Test report
    print('{:<10s}{:>14s}{:>14s}'.format("", "LOSS", "ACCURACY"))
    print(dash)
    print('{:<10s}{:>14.4f}{:>14.4f}'.format( "Training:", train_score[0], 100 * train_score[1]))
    print('{:<10s}{:>14.4f}{:>14.4f}'.format( "Test:", test_score[0], 100 * test_score[1]))


    # Calculate and report normalized error difference?
    if (calc_normal):
        max_err = max(train_score[0], test_score[0])
        error_diff = max_err - min(train_score[0], test_score[0])
        normal_diff = error_diff * 100 / max_err
        print('{:<10s}{:>13.2f}{:>1s}'.format("Normal diff ", normal_diff, ""))
        
def plot_train_history(history, x_ticks_vertical=False):
    history = history.history

    # min loss / max accs
    min_loss = min(history['loss'])
    min_val_loss = min(history['val_loss'])
    max_accuracy = max(history['accuracy'])
    max_val_accuracy = max(history['val_accuracy'])

    # x pos for loss / acc min/max
    min_loss_x = history['loss'].index(min_loss)
    min_val_loss_x = history['val_loss'].index(min_val_loss)
    max_accuracy_x = history['accuracy'].index(max_accuracy)
    max_val_accuracy_x = history['val_accuracy'].index(max_val_accuracy)

    # summarize history for loss, display min
    plt.figure(figsize=(16,8))
    plt.plot(history['loss'], color="#1f77b4", alpha=0.7)
    plt.plot(history['val_loss'], color="#ff7f0e", linestyle="--")
    plt.plot(min_loss_x, min_loss, marker='o', markersize=3, color="#1f77b4", alpha=0.7, label='Inline label')
    plt.plot(min_val_loss_x, min_val_loss, marker='o', markersize=3, color="#ff7f0e", alpha=7, label='Inline label')
    plt.title('Model loss', fontsize=20)
    plt.ylabel('Loss', fontsize=16)
    plt.xlabel('Epoch', fontsize=16)
    plt.legend(['Train', 
                'Test', 
                ('%.3f' % min_loss), 
                ('%.3f' % min_val_loss)], 
                loc='upper right', 
                fancybox=True, 
                framealpha=0.9, 
                shadow=True, 
                borderpad=1)

    if (x_ticks_vertical):
        plt.xticks(np.arange(0, len(history['loss']), 5.0), rotation='vertical')
    else:
        plt.xticks(np.arange(0, len(history['loss']), 5.0))

    plt.show()

    # summarize history for accuracy, display max
    plt.figure(figsize=(16,6))
    plt.plot(history['accuracy'], alpha=0.7)
    plt.plot(history['val_accuracy'], linestyle="--")
    plt.plot(max_accuracy_x, max_accuracy, marker='o', markersize=3, color="#1f77b4", alpha=7)
    plt.plot(max_val_accuracy_x, max_val_accuracy, marker='o', markersize=3, color="orange", alpha=7)
    plt.title('Model accuracy', fontsize=20)
    plt.ylabel('Accuracy', fontsize=16)
    plt.xlabel('Epoch', fontsize=16)
    plt.legend(['Train', 
                'Test', 
                ('%.2f' % max_accuracy), 
                ('%.2f' % max_val_accuracy)], 
                loc='upper left', 
                fancybox=True, 
                framealpha=0.9, 
                shadow=True, 
                borderpad=1)
    plt.figure(num=1, figsize=(10, 6))

    if (x_ticks_vertical):
        plt.xticks(np.arange(0, len(history['accuracy']), 5.0), rotation='vertical')
    else:
        plt.xticks(np.arange(0, len(history['accuracy']), 5.0))

    plt.show()
    
def plot_confusion_matrix(cm,
                          classes, 
                          normalized=False, 
                          title=None, 
                          cmap=plt.cm.Blues,
                          size=(10,10)):
    fig, ax = plt.subplots(figsize=size)
    im = ax.imshow(cm, interpolation='nearest', cmap=cmap)
    ax.figure.colorbar(im, ax=ax)

    # We want to show all ticks...
    ax.set(xticks=np.arange(cm.shape[1]),
           yticks=np.arange(cm.shape[0]),
           # ... and label them with the respective list entries
           xticklabels=classes, yticklabels=classes,
           title=title,
           ylabel='True label',
           xlabel='Predicted label')

    # Rotate the tick labels and set their alignment.
    plt.setp(ax.get_xticklabels(), rotation=45, ha="right",
             rotation_mode="anchor")

    # Loop over data dimensions and create text annotations.
    fmt = '.2f' if normalized else 'd'
    thresh = cm.max() / 2.
    for i in range(cm.shape[0]):
        for j in range(cm.shape[1]):
            ax.text(j, i, format(cm[i, j], fmt),
                    ha="center", va="center",
                    color="white" if cm[i, j] > thresh else "black")

    fig.tight_layout()
    plt.show()
       
def acc_per_class(np_probs_array):    
    accs = []
    for idx in range(0, np_probs_array.shape[0]):
        correct = np_probs_array[idx][idx].astype(int)
        total = np_probs_array[idx].sum().astype(int)
        acc = (correct / total) * 100
        accs.append(acc)
    return accs
# Load best saved model
model = load_model(model_path)
model_evaluation_report(model, X_train, y_train_encoded, X_test, y_test_encoded)
plot_train_history(history, x_ticks_vertical=True)
# Predict probabilities for test set
y_probs = model.predict(X_test, verbose=0)

# Get predicted labels
yhat_probs = np.argmax(y_probs, axis=1)
y_trues    = np.argmax(y_test_encoded, axis=1)
y_probs.shape
objects      = ['Caramel', 'Corn_Flakes', 'No-Sounds', 'Rice','Vitamines']
# Sets decimal precision (for printing output only)
np.set_printoptions(precision=2)

# Compute confusion matrix data
cm = confusion_matrix(y_trues, yhat_probs)

plot_confusion_matrix(cm,
                      objects, 
                      normalized=False, 
                      title="Model Performance", 
                      cmap=plt.cm.Blues,
                      size=(12,12))
# Find per-class accuracy from the confusion matrix data
accuracies = helpers.acc_per_class(cm)

pd.DataFrame({
    'CLASS': objects,
    'ACCURACY': accuracies
}).sort_values(by="ACCURACY", ascending=False)
# Build classification report
re = classification_report(y_trues, yhat_probs, labels=[0,1,2,3,4], target_names=objects)
print(re)<jupyter_output>precision    recall  f1-score   support

     Caramel       1.00      1.00      1.00      1039
 Corn_Flakes       1.00      1.00      1.00      1017
   No-Sounds       1.00      1.00      1.00       949
        Rice       1.00      1.00      1.00      1005
   Vitamines       1.00      1.00      1.00      1031

    accuracy                           1.00      5041
   macro avg       1.00      1.00      1.00      5041
weighted avg       1.00      1.00      1.00      5041<jupyter_text>PART 3. TSNE Representation<jupyter_code>objects      = ['Caramel', 'Corn_Flakes','Rice' , 'No-Sounds','Vitamines']
model2 = load_model(model_path)
config = model.layers[-1].get_config()
weights = [x.numpy() for x in model.layers[-1].weights]
config['activation'] = tf.keras.activations.linear
config['name'] = 'logits'
new_layer = tf.keras.layers.Dense(**config)(model.layers[-2].output)
new_model = tf.keras.Model(inputs=[model.input], outputs=[new_layer])
new_model.layers[-1].set_weights(weights)
test_embeddings    = tf.zeros((0, 5), dtype=tf.float32)
test_predictions   = []
for x in X_test : 
    x = x.reshape(1, num_rows, num_columns, num_channels)
    embeddings = new_model(x)
    prediction = model(x)
    
    test_predictions.append(tf.math.argmax(prediction,1))
    test_embeddings = tf.concat((test_embeddings, embeddings), 0)
embeddings_array = np.array(test_embeddings)
predictions_array = np.array(test_predictions).flatten()
""" 2D Plot"""
tsne = TSNE(2, verbose=1)
tsne_proj = tsne.fit_transform(embeddings_array)

cmap = cm.get_cmap('tab20')
fig, ax = plt.subplots(figsize=(16,10))
num_categories = 5
for lab in range(num_categories):
    indices = predictions_array == lab
    ax.scatter(tsne_proj[indices,0],
               tsne_proj[indices,1],
               c=np.array(cmap(lab)).reshape(1,4), 
               label = objects[lab] ,
               alpha=0.5)

plt.title('TSNE : Embedded representation of the Sounds', fontsize=20)  
ax.legend(fontsize='large', markerscale=2)
plt.show()
""" 3D Plot"""
tsne = TSNE(3, verbose=1)
tsne_proj = tsne.fit_transform(embeddings_array)
cmap = cm.get_cmap('tab20')
ax = plt.figure(figsize=(16,10)).gca(projection='3d')
num_categories = 5
for lab in range(num_categories):
    indices = predictions_array == lab
    ax.scatter(tsne_proj[indices, 0],
               tsne_proj[indices, 1],
               tsne_proj[indices, 2],
               c=np.array(cmap(lab)).reshape(1, 4),
               label=objects[lab],
               alpha=0.5)
plt.title('TSNE : Embedded representation of the Sounds', fontsize=20)  
ax.legend(fontsize='large', markerscale=2) 
ax.grid(False)
plt.show()<jupyter_output>[t-SNE] Computing 91 nearest neighbors...
[t-SNE] Indexed 5041 samples in 0.002s...
[t-SNE] Computed neighbors for 5041 samples in 0.120s...
[t-SNE] Computed conditional probabilities for sample 1000 / 5041
[t-SNE] Computed conditional probabilities for sample 2000 / 5041
[t-SNE] Computed conditional probabilities for sample 3000 / 5041
[t-SNE] Computed conditional probabilities for sample 4000 / 5041
[t-SNE] Computed conditional probabilities for sample 5000 / 5041
[t-SNE] Computed conditional probabilities for sample 5041 / 5041
[t-SNE] Mean sigma: 1.129091
[t-SNE] KL divergence after 250 iterations with early exaggeration: 67.755959
[t-SNE] KL divergence after 1000 iterations: 0.892822
{'fraction_non_alphanumeric': 0.11376092430751, 'fraction_numerical': 0.027847726262775884, 'mean_word_length': 2.8461830611469807, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 26, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '4', 'dir': 'jupyter-structured-clean-dedup', 'id': '369860', 'license': 'MIT', 'n_tokens_mistral': 7163, 'n_tokens_neox': 6606, 'n_words': 1422, 'path': 'Audio_record_and_classifier_framework/Features_extraction_classification.ipynb', 'repo_name': 'stanFurrer/Multimodal-solution-for-grasp-stability-prediction', 'size': '20253'}