text: #Ref: Sreenivas Sarwar Anik
  
  """
  Skin cancer lesion classification using the HAM10000 dataset
  
  Dataset link:
  https://www.kaggle.com/kmader/skin-cancer-mnist-ham10000
  Data description: 
  https://arxiv.org/ftp/arxiv/papers/1803/1803.10417.pdf
  
  The 7 classes of skin cancer lesions included in this dataset are:
  Melanocytic nevi (nv)
  Melanoma (mel)
  Benign keratosis-like lesions (bkl)
  Basal cell carcinoma (bcc) 
  Actinic keratoses (akiec)
  Vascular lesions (vas)
  Dermatofibroma (df)
  """
  
  
  import matplotlib.pyplot as plt
  import numpy as np
  import pandas as pd
  import os
  from glob import glob
  import seaborn as sns
  from PIL import Image
  
  np.random.seed(42)
  from sklearn.metrics import confusion_matrix
  
  import keras
  from keras.utils.np_utils import to_categorical # used for converting labels to one-hot-encoding
  from keras.models import Sequential
  from keras.layers import Dense, Dropout, Flatten, Conv2D, MaxPool2D, BatchNormalization
  from sklearn.model_selection import train_test_split
  from scipy import stats
  from sklearn.preprocessing import LabelEncoder
  
  skin_df = pd.read_csv('data/HAM10000/HAM10000_metadata.csv')
  
  SIZE=32
  
  # label encoding to numeric values from text
  le = LabelEncoder()
  le.fit(skin_df['dx'])
  LabelEncoder()
  print(list(le.classes_))
   
  skin_df['label'] = le.transform(skin_df["dx"]) 
  print(skin_df.sample(10))
  
  
  # Data distribution visualization
  fig = plt.figure(figsize=(12,8))
  
  ax1 = fig.add_subplot(221)
  skin_df['dx'].value_counts().plot(kind='bar', ax=ax1)
  ax1.set_ylabel('Count')
  ax1.set_title('Cell Type');
  
  ax2 = fig.add_subplot(222)
  skin_df['sex'].value_counts().plot(kind='bar', ax=ax2)
  ax2.set_ylabel('Count', size=15)
  ax2.set_title('Sex');
  
  ax3 = fig.add_subplot(223)
  skin_df['localization'].value_counts().plot(kind='bar')
  ax3.set_ylabel('Count',size=12)
  ax3.set_title('Localization')
  
  ax4 = fig.add_subplot(224)
  sample_age = skin_df[pd.notnull(skin_df['age'])]
  sns.distplot(sample_age['age'], fit=stats.norm, color='red');
  ax4.set_title('Age')
  
  plt.tight_layout()
  plt.show()
  
  
  # Distribution of data into various classes 
  from sklearn.utils import resample
  print(skin_df['label'].value_counts())
  
  #Balance data.
  # Many ways to balance data... you can also try assigning weights during model.fit
  #Separate each classes, resample, and combine back into single dataframe
  
  df_0 = skin_df[skin_df['label'] == 0]
  df_1 = skin_df[skin_df['label'] == 1]
  df_2 = skin_df[skin_df['label'] == 2]
  df_3 = skin_df[skin_df['label'] == 3]
  df_4 = skin_df[skin_df['label'] == 4]
  df_5 = skin_df[skin_df['label'] == 5]
  df_6 = skin_df[skin_df['label'] == 6]
  
  n_samples=500 
  df_0_balanced = resample(df_0, replace=True, n_samples=n_samples, random_state=42) 
  df_1_balanced = resample(df_1, replace=True, n_samples=n_samples, random_state=42) 
  df_2_balanced = resample(df_2, replace=True, n_samples=n_samples, random_state=42)
  df_3_balanced = resample(df_3, replace=True, n_samples=n_samples, random_state=42)
  df_4_balanced = resample(df_4, replace=True, n_samples=n_samples, random_state=42)
  df_5_balanced = resample(df_5, replace=True, n_samples=n_samples, random_state=42)
  df_6_balanced = resample(df_6, replace=True, n_samples=n_samples, random_state=42)
  
  #Combined back to a single dataframe
  skin_df_balanced = pd.concat([df_0_balanced, df_1_balanced, 
                                df_2_balanced, df_3_balanced, 
                                df_4_balanced, df_5_balanced, df_6_balanced])
  
  #Check the distribution. All classes should be balanced now.
  print(skin_df_balanced['label'].value_counts())
  
  
  #Now time to read images based on image ID from the CSV file
  #This is the safest way to read images as it ensures the right image is read for the right ID
  image_path = {os.path.splitext(os.path.basename(x))[0]: x
                       for x in glob(os.path.join('data/HAM10000/', '*', '*.jpg'))}
  
  #Define the path and add as a new column
  skin_df_balanced['path'] = skin_df['image_id'].map(image_path.get)
  #Use the path to read images.
  skin_df_balanced['image'] = skin_df_balanced['path'].map(lambda x: np.asarray(Image.open(x).resize((SIZE,SIZE))))
  
  
  n_samples = 5  # number of samples for plotting
  # Plotting
  fig, m_axs = plt.subplots(7, n_samples, figsize = (4*n_samples, 3*7))
  for n_axs, (type_name, type_rows) in zip(m_axs, 
                                           skin_df_balanced.sort_values(['dx']).groupby('dx')):
      n_axs[0].set_title(type_name)
      for c_ax, (_, c_row) in zip(n_axs, type_rows.sample(n_samples, random_state=1234).iterrows()):
          c_ax.imshow(c_row['image'])
          c_ax.axis('off')
  
  #Convert dataframe column of images into numpy array
  X = np.asarray(skin_df_balanced['image'].tolist())
  X = X/255.  # Scale values to 0-1. You can also used standardscaler or other scaling methods.
  Y=skin_df_balanced['label']  #Assign label values to Y
  Y_cat = to_categorical(Y, num_classes=7) #Convert to categorical as this is a multiclass classification problem
  #Split to training and testing
  x_train, x_test, y_train, y_test = train_test_split(X, Y_cat, test_size=0.25, random_state=42)
  
  #Define the model.
  #I've used autokeras to find out the best model for this problem.
  #You can also load pretrained networks such as mobilenet or VGG16
  
  num_classes = 7
  
  model = Sequential()
  model.add(Conv2D(256, (3, 3), activation="relu", input_shape=(SIZE, SIZE, 3)))
  #model.add(BatchNormalization())
  model.add(MaxPool2D(pool_size=(2, 2)))  
  model.add(Dropout(0.3))
  
  model.add(Conv2D(128, (3, 3),activation='relu'))
  #model.add(BatchNormalization())
  model.add(MaxPool2D(pool_size=(2, 2)))  
  model.add(Dropout(0.3))
  
  model.add(Conv2D(64, (3, 3),activation='relu'))
  #model.add(BatchNormalization())
  model.add(MaxPool2D(pool_size=(2, 2)))  
  model.add(Dropout(0.3))
  model.add(Flatten())
  
  model.add(Dense(32))
  model.add(Dense(7, activation='softmax'))
  model.summary()
  
  model.compile(loss='categorical_crossentropy', optimizer='Adam', metrics=['acc'])
  
  
  # Train
  #You can also use generator to use augmentation during training.
  
  batch_size = 16 
  epochs = 50
  
  history = model.fit(
      x_train, y_train,
      epochs=epochs,
      batch_size = batch_size,
      validation_data=(x_test, y_test),
      verbose=2)
  
  score = model.evaluate(x_test, y_test)
  print('Test accuracy:', score[1])
  
  
  #plot the training and validation accuracy and loss at each epoch
  loss = history.history['loss']
  val_loss = history.history['val_loss']
  epochs = range(1, len(loss) + 1)
  plt.plot(epochs, loss, 'y', label='Training loss')
  plt.plot(epochs, val_loss, 'r', label='Validation loss')
  plt.title('Training and validation loss')
  plt.xlabel('Epochs')
  plt.ylabel('Loss')
  plt.legend()
  plt.show()
  
  
  acc = history.history['acc']
  val_acc = history.history['val_acc']
  plt.plot(epochs, acc, 'y', label='Training acc')
  plt.plot(epochs, val_acc, 'r', label='Validation acc')
  plt.title('Training and validation accuracy')
  plt.xlabel('Epochs')
  plt.ylabel('Accuracy')
  plt.legend()
  plt.show()
  
  
  # Prediction on test data
  y_pred = model.predict(x_test)
  # Convert predictions classes to one hot vectors 
  y_pred_classes = np.argmax(y_pred, axis = 1) 
  # Convert test data to one hot vectors
  y_true = np.argmax(y_test, axis = 1) 
  
  #Print confusion matrix
  cm = confusion_matrix(y_true, y_pred_classes)
  
  fig, ax = plt.subplots(figsize=(6,6))
  sns.set(font_scale=1.6)
  sns.heatmap(cm, annot=True, linewidths=.5, ax=ax)
  
  
  #PLot fractional incorrect misclassifications
  incorr_fraction = 1 - np.diag(cm) / np.sum(cm, axis=1)
  plt.bar(np.arange(7), incorr_fraction)
  plt.xlabel('True Label')
  plt.ylabel('Fraction of incorrect predictions')
  
  
          
          