add counterfactual code

Pilot1/NT3/make_csv.py

@@ -0,0 +1,61 @@
+import pandas as pd
+import pickle
+import argparse
+import glob, os
+from pathlib import Path
+import matplotlib.pyplot as plt
+
+def get_args():
+    parser = argparse.ArgumentParser()
+    parser.add_argument("-f",type=str, help="Run folder")
+    parser.add_argument("-c1", type=str, help="cluster 1 name")
+    parser.add_argument("-c2", type=str, help="cluster 2 name")
+    args = parser.parse_args()
+    return args
+
+def main():
+    args = get_args()
+    l1 = []
+    l2 = []
+    runs = glob.glob(args.f+"/EXP000/*/")
+    print(runs)
+    for r in runs:
+        print(r)
+        global_data = pd.read_csv(r+"training.log")
+        val_abs = global_data['val_abstention'].iloc[-1]
+        val_abs_acc = global_data['val_abstention_acc'].iloc[-1]
+        if os.path.exists(r+"cluster_trace.pkl"):
+            cluster_data = pickle.load(open(r+"cluster_trace.pkl", "rb"))
+        else: 
+            continue
+        polluted_abs = cluster_data['Abs polluted']
+        val_abs_cluster = cluster_data['Abs val cluster']
+        val_abs_acc_cluster = cluster_data['Abs val acc']
+        ratio = float(r[-8:-5])
+        if args.c1 in r:
+            l1.append([ratio, val_abs, val_abs_acc, val_abs_cluster, val_abs_acc_cluster, polluted_abs])
+        elif args.c2 in r:
+            l2.append([ratio, val_abs, val_abs_acc, val_abs_cluster, val_abs_acc_cluster, polluted_abs])
+
+    df1 = pd.DataFrame(l1, columns=['Noise Fraction', 'Val Abs', 'Val Abs Acc', 'Val Abs Cluster', 'Val Abs Acc Cluster', 'Polluted Abs'])
+    df2 = pd.DataFrame(l2, columns=['Noise Fraction', 'Val Abs', 'Val Abs Acc', 'Val Abs Cluster', 'Val Abs Acc Cluster', 'Polluted Abs'])
+    print(df1)    
+    df1.to_csv("cluster_1.csv")
+    df2.to_csv("cluster_2.csv")
+    plt.plot(df1['Noise Fraction'], df1['Val Abs'], marker='o', label='Val Abs')
+    plt.plot(df1['Noise Fraction'], df1['Val Abs Acc'],  marker='o',label='Val Abs Acc')
+    plt.plot(df1['Noise Fraction'], df1['Val Abs Cluster'],  marker='o',label='Val Abs Cluster')
+    plt.plot(df1['Noise Fraction'], df1['Val Abs Acc Cluster'],  marker='o',label='Val Abs Acc Cluster')
+    plt.xlabel("Noise fraction")
+    plt.legend()
+    plt.savefig('c1.png')
+
+    plt.plot(df2['Noise Fraction'], df2['Val Abs'],  marker='o',label='Val Abs')
+    plt.plot(df2['Noise Fraction'], df2['Val Abs Acc'],  marker='o',label='Val Abs Acc')
+    plt.plot(df2['Noise Fraction'], df2['Val Abs Cluster'],  marker='o',label='Val Abs Cluster')
+    plt.plot(df2['Noise Fraction'], df2['Val Abs Acc Cluster'],  marker='o',label='Val Abs Acc Cluster')
+    plt.xlabel("Noise Fraction")
+    plt.legend()
+    plt.savefig('c2.png')
+if __name__ == "__main__":
+    main()

Pilot1/NT3/nt3_abstention_keras2.py

@@ -16,6 +16,7 @@
 
 import nt3 as bmk
 import candle
+import pickle
 
 additional_definitions = abs_definitions
 
@@ -51,7 +52,13 @@ def initialize_parameters(default_model='nt3_noise_model.txt'):
     gParameters = candle.finalize_parameters(nt3Bmk)
 
     return gParameters
-
+
+def load_data_cf(cf_path):
+    # Pickle file holds the test train split and cf index info                                                       
+    print("Loading data...")
+    X_train, X_test, Y_train, Y_test, polluted_inds, cluster_inds = pickle.load(open(cf_path, 'rb'))
+    print('done')
+    return X_train, Y_train, X_test, Y_test, polluted_inds, cluster_inds
 
 def load_data(train_path, test_path, gParameters):
 
@@ -86,6 +93,38 @@ def load_data(train_path, test_path, gParameters):
 
     return X_train, Y_train, X_test, Y_test
 
+def evaluate_cf(model, nb_classes, output_dir, X_train, X_test, Y_train, Y_test, polluted_inds, cluster_inds, gParameters):
+    if len(polluted_inds) > 0:
+        y_pred = model.predict(X_test)
+        abstain_inds = []
+        for i in range(y_pred.shape[0]):
+            if np.argmax(y_pred[i]) == nb_classes:
+                abstain_inds.append(i)
+
+        # Cluster indices and polluted indices are wrt to entire train + test dataset                            
+        # whereas y_pred only contains test dataset so add offset for correct indexing                            
+        offset_testset = Y_train.shape[0]
+        abstain_inds=[i+offset_testset for i in abstain_inds]
+        polluted_percentage = np.sum([el in polluted_inds for el in abstain_inds])/np.max([len(abstain_inds),\
+1])
+        print("Percentage of abstained samples that were polluted {}".format(polluted_percentage))
+
+        cluster_inds_test = list(filter(lambda cluster_inds: cluster_inds >= offset_testset, cluster_inds))
+        cluster_inds_test_abstain = [el in abstain_inds for el in cluster_inds_test]
+        cluster_percentage = c = np.sum(cluster_inds_test_abstain)/len(cluster_inds_test)
+        print("Percentage of cluster (in test set) that was abstained {}".format(cluster_percentage))
+
+        unabstain_inds = []
+        for i in range(y_pred.shape[0]):
+            if np.argmax(y_pred[i]) != nb_classes and (i+offset_testset in cluster_inds_test):
+                unabstain_inds.append(i)
+        # Make sure number of unabstained indices in cluster test set plus number of abstainsed indices in cluster test set                                                                                                        
+        # equals number of indices in cluster in the test set                                                     
+        assert(len(unabstain_inds)+np.sum(cluster_inds_test_abstain) == len(cluster_inds_test))
+        score_cluster = 1 if len(unabstain_inds)==0 else model.evaluate(X_test[unabstain_inds], Y_test[unabstain_inds])[1]
+        print("Accuracy of unabastained cluster {}".format(score_cluster))
+        if gParameters['noise_save_cf']:
+            pickle.dump({'Abs polluted': polluted_percentage, 'Abs val cluster': cluster_percentage, 'Abs val acc': score_cluster}, open("{}/cluster_trace.pkl".format(output_dir), "wb"))
 
 def run(gParameters):
 
@@ -96,7 +135,10 @@ def run(gParameters):
     train_file = candle.get_file(file_train, url + file_train, cache_subdir='Pilot1')
     test_file = candle.get_file(file_test, url + file_test, cache_subdir='Pilot1')
 
-    X_train, Y_train, X_test, Y_test = load_data(train_file, test_file, gParameters)
+    if gParameters['noise_cf'] is not None:
+        X_train, Y_train, X_test, Y_test, polluted_inds, cluster_inds  = load_data_cf(gParameters['noise_cf'])
+    else:
+        X_train, Y_train, X_test, Y_test = load_data(train_file, test_file, gParameters)
 
     # only training set has noise
     X_train, Y_train = candle.add_noise(X_train, Y_train, gParameters)
@@ -274,6 +316,8 @@ def run(gParameters):
 
     score = model.evaluate(X_test, Y_test, verbose=0)
 
+    if gParameters['noise_cf'] is not None:
+        evaluate_cf(model, nb_classes, output_dir, X_train, X_test, Y_train, Y_test, polluted_inds, cluster_inds, gParameters)
     alpha_trace = open(output_dir + "/alpha_trace", "w+")
     for alpha in abstention_cbk.alphavalues:
         alpha_trace.write(str(alpha) + '\n')

Pilot1/NT3/nt3_baseline_keras2.py

@@ -15,9 +15,9 @@
 
 import nt3 as bmk
 import candle
+import pickle
 
-
-def initialize_parameters(default_model='nt3_default_model.txt'):
+def initialize_parameters(default_model='nt3_noise_model.txt'):
 
     # Build benchmark object
     nt3Bmk = bmk.BenchmarkNT3(
@@ -238,6 +238,10 @@ def run(gParameters):
 
         print("json %s: %.2f%%" % (loaded_model_json.metrics_names[1], score_json[1] * 100))
 
+
+    if gParameters['noise_save_cf']:
+        model.save('{}/{}.autosave.model'.format(output_dir, model_name))
+        pickle.dump([X_train, X_test, Y_train, Y_test], open('{}/{}.autosave.data.pkl'.format(output_dir, model_name), "wb"))
     return history
 
 

Pilot1/NT3/nt3_cf/README.md

@@ -0,0 +1,30 @@
+NT3 with counterfactuals:
+Code to generate counterfactual examples given an input model and dataset in pkl format. \
+Clusters and thresholds counterfactuals, injects noise into dataset \
+Workflow: 
+1) Generate counterfactuals using cf_nb.py
+```
+python cf_nb.py
+```
+
+2) Create threshold pickle files using threshold.py (provide a threshold value between 0 and 1, see --help) 
+```
+python threshold.py -d ../nt3.autosave.data.pkl -c cf_redo_all_reformat.pkl -t 0.9 -o threshold_0.9.pkl
+```
+
+3) Cluster threshold files using gen_clusters.py
+```
+python gen_clusters.py -t_value 0.9 -t threshold_0.9.pkl
+```
+
+4) Inject noise into dataset using inject_noise.py (provide a scale value to modify the amplitude of the noise, see --help)
+```
+python inject_noise.py -t threshold_0.9.pkl -c1 cf_class_0_cluster0.pkl -c2 cf_class_1_cluster0.pkl -scale 1.0 -r True -d ../nt3.autosave.data.pkl -f cf_failed_inds.pkl -o noise_data 
+```
+
+Abstention with counterfactuals:
+Code located in abstention/
+Workflow:
+1) Run abstention model with nt3_abstention_keras2_cf.py, pass in a pickle file with X (with noise), y (this is the output of 4) above)
+2) For a sweep use run_abstention_sweep.sh
+3) To collect metrics (abstention, cluster abstention) run make_csv.py

Pilot1/NT3/nt3_cf/abstention/make_csv.py

@@ -0,0 +1,41 @@
+import pandas as pd
+import pickle
+import argparse
+import glob, os
+from pathlib import Path
+
+def get_args():
+    parser = argparse.ArgumentParser()
+    parser.add_argument("-f",type=str, help="Run folder")
+    parser.add_argument("-c1", type=str, help="cluster 1 name")
+    parser.add_argument("-c2", type=str, help="cluster 2 name")
+    args = parser.parse_args()
+    return args
+
+def main():
+    args = get_args()
+    l1 = []
+    l2 = []
+    runs = glob.glob(args.f+"/EXP000/*/")
+    print(runs)
+    for r in runs:
+        global_data = pd.read_csv(r+"training.log")
+        val_abs = global_data['val_abstention'].iloc[-1]
+        val_abs_acc = global_data['val_abstention_acc'].iloc[-1]
+        cluster_data = pickle.load(open(r+"cluster_trace.pkl", "rb"))
+        polluted_abs = cluster_data['Abs polluted']
+        val_abs_cluster = cluster_data['Abs val cluster']
+        val_abs_acc_cluster = cluster_data['Abs val acc']
+        ratio = float(r[-4:-1])
+        if args.c1 in r:
+            l1.append([ratio, val_abs, val_abs_acc, val_abs_cluster, val_abs_acc_cluster, polluted_abs])
+        elif args.c2 in r:
+            l2.append([ratio, val_abs, val_abs_acc, val_abs_cluster, val_abs_acc_cluster, polluted_abs])
+
+    df1 = pd.DataFrame(l1, columns=['Noise Fraction', 'Val Abs', 'Val Abs Acc', 'Val Abs Cluster', 'Val Abs Acc Cluster', 'Polluted Abs'])
+    df2 = pd.DataFrame(l2, columns=['Noise Fraction', 'Val Abs', 'Val Abs Acc', 'Val Abs Cluster', 'Val Abs Acc Cluster', 'Polluted Abs'])
+    print(df1)    
+    df1.to_csv("cluster_1.csv")
+    df2.to_csv("cluster_2.csv")
+if __name__ == "__main__":
+    main()

Pilot1/NT3/nt3_cf/abstention/run_abstention_sweep.sh

@@ -0,0 +1,5 @@
+#!/bin/bash
+for filename in /vol/ml/shahashka/xai-geom/nt3/nt3.data*; do
+    python nt3_abstention_keras2_cf.py --cf_noise $filename --output_dir cf_sweep_0906 --run_id ${filename:40:21} --epochs 100
+    #cp cf_sweep_0902/EXP000/RUN000/training.log ${filename}_training_0902.log 
+done

Pilot1/NT3/nt3_cf/analyze.py

@@ -0,0 +1,33 @@
+# Script to analyze perturbation by cluster
+# Plot the perturbations by cluster
+# Plot the pertubation centroids
+
+import os
+import pickle
+import matplotlib.pyplot as plt
+import numpy as np
+directory = 'clusters_0911_0.5/'
+orig_dataset = pickle.load(open("nt3.autosave.data.pkl", 'rb'))[0]
+cf_dataset = pickle.load(open("threshold_0905.pkl", 'rb'))['perturbation vector']
+for filename in os.listdir(directory):
+    if filename.startswith("cf_class_0") or filename.startswith("cf_class_1") :
+        data = pickle.load(open(os.path.join(directory, filename), 'rb'))
+        x_range = np.arange(len(data['centroid perturb vector']))
+        ind_in_cluster = data['sample indices in this cluster'][0:5]
+        fig,ax = plt.subplots(3, figsize=(20,15))
+        fig.suptitle("Perturbation Vectors for counterfactual class 1, cluster 1", fontsize=25)
+        for i,ax_i in zip(ind_in_cluster,ax):
+            d = cf_dataset[i]
+            ax_i.plot(x_range, d, label='perturbation vector')
+            ax_i.plot(x_range ,data['centroid perturb vector'], label='centroid')
+            #ax_i.axhline(y=0.5*np.max(np.abs(d)), color='r', linestyle='-')
+            #ax_i.axhline(y=-0.5*np.max(np.abs(d)), color='r', linestyle='-')
+            ax_i.axvline(x=9603, color='r', linestyle='-', linewidth=5, alpha=0.3)
+
+            ax_i.set_title("sample {}".format(i))
+            ax_i.legend()
+        fig.supxlabel("Feature index", fontsize=18)
+        plt.savefig("centroids_{}.png".format(filename))
+
+    else:
+        continue

Pilot1/NT3/nt3_cf/cf_nb.py

@@ -0,0 +1,56 @@
+import tensorflow as tf
+tf.get_logger().setLevel(40) # suppress deprecation messages
+tf.compat.v1.disable_v2_behavior() # disable TF2 behaviour as alibi code still relies on TF1 constructs
+from tensorflow.keras.models import Model, load_model
+import matplotlib.pyplot as plt
+import numpy as np
+import os
+os.environ["CUDA_VISIBLE_DEVICES"]="1"
+from time import time
+from alibi.explainers import CounterFactual, CounterFactualProto
+print('TF version: ', tf.__version__)
+print('Eager execution enabled: ', tf.executing_eagerly()) # False
+print(tf.test.is_gpu_available())
+import pickle
+model_nt3 = tf.keras.models.load_model('../nt3.autosave.model')
+with open('../nt3.autosave.data.pkl', 'rb') as pickle_file:
+        X_train,X_test,Y_train,Y_test = pickle.load(pickle_file)
+
+shape_cf = (1,) + X_train.shape[1:]
+print(shape_cf)
+target_proba = 0.9
+tol = 0.1 # want counterfactuals with p(class)>0.90
+target_class = 'other' # any class other than will do
+max_iter = 1000
+lam_init = 1e-1
+max_lam_steps = 20
+learning_rate_init = 0.1
+feature_range = (0,1)
+cf = CounterFactual(model_nt3, shape=shape_cf, target_proba=target_proba, tol=tol,
+                                        target_class=target_class, max_iter=max_iter, lam_init=lam_init,
+                                        max_lam_steps=max_lam_steps, learning_rate_init=learning_rate_init,
+                                        feature_range=feature_range)
+shape = X_train[0].shape[0]
+results=[]
+failed_inds = []
+X = np.concatenate([X_train,X_test])
+
+for i in np.arange(0,X.shape[0]):
+    print(i)
+    x_sample=X[i:i+1]
+    print(x_sample.shape)
+    start = time()
+    try:
+        explanation = cf.explain(x_sample)
+        print('Counterfactual prediction: {}, {}'.format(explanation.cf['class'], explanation.cf['proba']))
+        print("Actual prediction: {}".format(model_nt3.predict(x_sample)))
+        results.append([i, explanation.cf['X'],explanation.cf['class'], explanation.cf['proba']])
+        test = model_nt3.predict(explanation.cf['X'])
+        print(test, explanation.cf['proba'], explanation.cf['class'])
+    except:
+        print("Failed cf generation")
+        failed_inds.append(i)
+    if i%100 == 0 and i is not 0:
+        pickle.dump(results, open("cf_{}.pkl".format(i), "wb"))
+        results = []
+pickle.dump(failed_inds, open("cf_failed_inds.pkl", "wb"))