Source

python-ml-explain / explainer.py

Full commit
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.patches as mpatches

def explain_instance(model, data, instance, iterations = 300):
  """
  Explains a dataset instance using the specified model. The returned
  explanation is a dictionary containing contribution of each individual
  feature and the model's prediction.

  NOTE: It currently only works for classification models!

  :param model: Scikit-learn model instance with enabled probability
    estimation
  :param data: Dataset used for estimating feature range
  :param instance: Instance to explain
  :param iterations: Number of iterations in Monte Carlo simulation
  """
  cls_typ = instance.dtypes[-1].type
  features = list(data.columns[:-1])
  actual = np.asarray(instance[data.columns[-1]])[0]
  instance = instance[features]
  prediction = cls_typ(model.predict(np.asarray(instance))[0])
  p_index = sorted(data[data.columns[-1]].unique()).index(prediction)
  data = data[features]

  nFeatures = len(data.columns)
  explanation = []
  for feature in data.columns:
    contribution = 0.0
    for j in xrange(iterations):
      perm = np.random.choice(data.columns, nFeatures, replace = False)
      tmp = instance.copy()
      idx = 0
      while perm[idx] != feature:
        tmp[perm[idx]] = select_random_value(data[perm[idx]])
        idx += 1

      tmp2 = tmp.copy()
      tmp2[perm[idx]] = select_random_value(data[perm[idx]])
      contribution += model.predict_proba(tmp)[0][p_index] - model.predict_proba(tmp2)[0][p_index]

    contribution /= iterations
    explanation.append((feature, contribution, np.asarray(instance[feature])[0]))

  return dict(explanation = explanation, prediction = prediction, actual = actual,
    model_name = model.__class__.__name__)

def explain_value(model, data, feature, value, iterations = 300):
  """
  Explains a single value of a single feature. The returned explanation
  is a dictionary containing the contribution mean and standard deviation.
  
  NOTE: It currently only works for discrete features!
  
  :param model: Scikit-learn model instance with enabled probability
    estimation
  :param data: Dataset used for estimating feature range
  :param feature: Feature to explain
  :param value: Feature value to explain
  :param iterations: Number of iterations in Monte Carlo simulation
  """
  cls_typ = data.dtypes[-1].type
  orig_data = data
  data = data[data.columns[:-1]]
  
  contribs = []
  for j in xrange(iterations):
    # Use first instance to get attribute format
    instance1 = data[0:1]
    # Replace all attributes with random values
    for ifeature in data.columns:
      instance1[ifeature] = select_random_value(data[ifeature])
    # Make another instance and replace the chosen feature with a
    # pre-selected value
    instance2 = instance1.copy()
    instance2[feature] = value
    # Append contribution
    contribs.append(model.predict_proba(instance2)[0][0] - \
      model.predict_proba(instance1)[0][0])
  
  return dict(mean = np.mean(contribs), std = np.std(contribs))

def explain_model(model, data, resolution = 50, iterations = 300):
  """
  Explains the complete model (all values of all features).
  
  NOTE: It currently only works for classification models!
  
  :param model: Scikit-learn model instance with enabled probability
    estimation
  :param data: Dataset used for estimating feature range
  :param resolution: Resolution for continous attributes
  :param iterations: Number of iterations in Monte Carlo simulation
  """
  explanation = []
  for feature in data.columns[:-1]:
    means, stds = [], []
    
    if data[feature].dtype.kind == 'i':
      values = sorted(data[feature].unique())
    elif data[feature].dtype.kind == 'f':
      values = np.linspace(0, 1, resolution)
    
    for value in values:
      e = explain_value(model, data, feature, value, iterations)
      means.append(e['mean'])
      stds.append(e['std'])
    
    explanation.append((feature, dict(values = values, means = means, stds = stds)))
  
  return explanation

def select_random_value(feature):
  """
  Selects a random value from a feature's range. It assumes that
  integer-typed features are discrete and float-typed features
  are numeric.

  :param feature: DataFrame describing feature's values
  """
  if feature.dtype.kind == 'i':
    return np.random.choice(feature.unique(), 1)
  elif feature.dtype.kind == 'f':
    return np.random.rand() * (feature.max() - feature.min()) + feature.min()
 
  raise TypeError, "Unsupported feature type!"

def plot_instance_explanation(result, value_maps = None, filename = "output.png"):
  """
  Plots an instance explanation generated by `explain_instance`.

  :param result: Explanation generated by `explain_instance`
  :param filename: File where the output visualization should be
    saved to
  """
  plt.clf()
  fig = plt.figure(figsize = (8, 6))
  bar_ax = plt.axes((0.1, 0.1, 0.85, 0.7))
  result['explanation'] = result['explanation'][::-1]
  N = len(result['explanation'])
  contributions = [contribution for _, contribution, _ in result['explanation']]
  maxc = np.max(np.abs(contributions)) + 0.1

  bar_ax.barh(np.arange(N) + 0.3, contributions, height = 0.55)
  
  bar_ax.hlines(np.arange(N) + 0.05, -maxc, maxc, linestyles = 'dashed')
  bar_ax.axvline(x = 0.0, color = 'black', linewidth = 2)

  bar_ax.set_yticks(np.arange(len(result['explanation'])) + 0.5)
  bar_ax.set_yticklabels([feature for feature, _, _ in result['explanation']])
  bar_ax.set_xlim(left = -maxc, right = maxc)
  for t in bar_ax.get_xticklines(): t.set_marker(None)
  for t in bar_ax.get_yticklines(): t.set_marker(None)
  bar_ax.set_frame_on(False)

  info_ax = plt.axes((0.1, 0.84, 0.85, 0.1))
  info_ax.set_frame_on(False)
  info_ax.set_axis_off()

  p = mpatches.Rectangle((0.0, 0.01), 0.99, 0.99, facecolor = "lightblue",
    edgecolor = "black", alpha = 0.5)
  info_ax.add_patch(p)
  
  info_ax.text(0.01, 0.6, "Prediction: R = %s" % result['prediction'],
    fontsize = 13, fontweight = 'bold')
  info_ax.text(0.01, 0.2, "Actual value: R = %s" % result['actual'],
    fontsize = 13, fontweight = 'bold')
  info_ax.text(0.98, 0.6, result['model_name'], fontsize = 13,
    horizontalalignment = 'right')
  
  for idx, (feature, contribution, value) in enumerate(result['explanation']):
    bar_ax.text(maxc, 0.5 + idx, "%.2f" % contribution, horizontalalignment = 'right',
      bbox = dict(facecolor='red', alpha = 0.5, pad = 10.0))
  
    value = round(value, 2) if isinstance(value, float) else value
    if value_maps is not None and value in value_maps.get(feature, {}):
      value = value_maps[feature][value]
    
    bar_ax.text(-maxc + maxc/15., 0.5 + idx, "%s" % value, horizontalalignment = 'left',
      bbox = dict(facecolor='green', alpha = 0.5, pad = 15.0))

  plt.savefig(filename)
  plt.close(fig)

def plot_model_explanation(result, filename = "output.png"):
  """
  Plots an instance explanation generated by `explain_model`.

  :param result: Explanation generated by `explain_model`
  :param filename: File where the output visualization should be
    saved to
  """
  plt.clf()
  fig = plt.figure(figsize = (8, 12))
  height = 0.97 / len(result)
  result = result[::-1]
  for idx, (feature, data) in enumerate(result):
    maxc = np.max(data['means']) + 0.1
    minv, maxv = np.min(data['values']), np.max(data['values'])
    subplot = plt.axes((0.1, 0.05 + idx*height, 0.85, height - 0.05))
    subplot.set_title(feature)
    subplot.plot(data['values'], data['means'], 'bo')
    subplot.axhline(y = 0.0, color = 'black', linewidth = 1, linestyle = '--')
    subplot.set_ylim(bottom = -maxc, top = maxc)
    subplot.set_xlim(left = minv-0.02, right = maxv+0.02)
    
    subplot = subplot.twinx()
    subplot.set_axis_off()
    subplot.plot(data['values'], data['stds'], 'go')
    subplot.set_ylim(bottom = -1.0, top = 1.0)
    subplot.set_xlim(left = minv-0.02, right = maxv+0.02)
  
  plt.savefig(filename)
  plt.close(fig)