adding plot_single_neuron script

plotting.py

@@ -12,6 +12,10 @@ import numpy as np
 from .keras_visualize_activations.read_activations import get_activations
 from .tfhelpers import get_grad_function, max_activation_wrt_input, create_random_event
 
+import logging
+logger = logging.getLogger(__name__)
+logger.addHandler(logging.NullHandler())
+
 """
 Some further plotting functions
 """
@@ -24,12 +28,12 @@ def plot_NN_vs_var_1D(plotname, means, scorefun, var_index, var_range, var_label
     "Plot the NN output vs one variable with the other variables set to the given mean values"
 
     # example: vary var1
-    print("Creating varied events (1d)")
+    logger.info("Creating varied events (1d)")
     sequence = np.arange(*var_range)
     events = np.tile(means, len(sequence)).reshape(-1, len(means))
     events[:,var_index] = sequence
 
-    print("Predicting scores")
+    logger.info("Predicting scores")
     scores = scorefun(events)
 
     fig, ax = plt.subplots()
@@ -43,20 +47,23 @@ def plot_NN_vs_var_1D(plotname, means, scorefun, var_index, var_range, var_label
 
 def plot_NN_vs_var_2D(plotname, means,
                       scorefun,
-                      var1_index, var1_range,
-                      var2_index, var2_range,
-                      var1_label=None,
-                      var2_label=None,
+                      varx_index,
+                      vary_index,
+                      nbinsx, xmin, xmax,
+                      nbinsy, ymin, ymax,
+                      varx_label=None,
+                      vary_label=None,
                       logscale=False,
                       ncontours=20,
                       only_pixels=False,
                       black_contourlines=False,
                       cmap="inferno"):
 
-    print("Creating varied events (2d)")
-    # example: vary var1 vs var2
-    sequence1 = np.arange(*var1_range)
-    sequence2 = np.arange(*var2_range)
+    logger.info("Creating varied events (2d)")
+
+    sequence1 = np.linspace(xmin, xmax, nbinsx)
+    sequence2 = np.linspace(ymin, ymax, nbinsy)
+
     # the following is a 2d array of events (so effectively 3D)
     events = np.tile(means, len(sequence1)*len(sequence2)).reshape(len(sequence2), len(sequence1), -1)
 
@@ -64,13 +71,13 @@ def plot_NN_vs_var_2D(plotname, means,
     # (probably there is a more clever way, but sufficient here)
     for i, y in enumerate(sequence2):
         for j, x in enumerate(sequence1):
-            events[i][j][var1_index] = x
-            events[i][j][var2_index] = y
+            events[i][j][varx_index] = x
+            events[i][j][vary_index] = y
 
     # convert back into 1d array
     events = events.reshape(-1, len(means))
 
-    print("Predicting scores")
+    logger.info("Predicting scores")
     scores = scorefun(events)
 
     # convert scores into 2d array
@@ -84,7 +91,7 @@ def plot_NN_vs_var_2D(plotname, means,
     if logscale:
         if zmin <= 0:
             zmin = 1e-5
-            print("Setting zmin to {}".format(zmin))
+            logger.info("Setting zmin to {}".format(zmin))
         lvls = np.logspace(math.log10(zmin), math.log10(zmax), ncontours)
         if only_pixels:
             pcm = ax.pcolormesh(sequence1, sequence2, scores, norm=matplotlib.colors.LogNorm(vmin=zmin, vmax=zmax), cmap=cmap, linewidth=0, rasterized=True)
@@ -104,10 +111,10 @@ def plot_NN_vs_var_2D(plotname, means,
         cbar = fig.colorbar(pcm, ax=ax, extend='max')
 
     cbar.set_label("NN output")
-    if var1_label is not None:
-        ax.set_xlabel(var1_label)
-    if var2_label is not None:
-        ax.set_ylabel(var2_label)
+    if varx_label is not None:
+        ax.set_xlabel(varx_label)
+    if vary_label is not None:
+        ax.set_ylabel(vary_label)
     fig.savefig(plotname)
     plt.close(fig)
 
@@ -166,15 +173,15 @@ def plot_NN_vs_var_2D_all(plotname, model, means,
     global_min = min([np.min(ar_layer) for ar_layer in acts[:-1]])
     global_max = max([np.max(ar_layer) for ar_layer in acts[:-1]])
 
-    print("global_min: {}".format(global_min))
-    print("global_max: {}".format(global_max))
+    logger.info("global_min: {}".format(global_min))
+    logger.info("global_max: {}".format(global_max))
 
     output_min_default = 0
     output_max_default = 1
 
     if global_min <= 0 and logz:
         global_min = log_default_ymin
-        print("Changing global_min to {}".format(log_default_ymin))
+        logger.info("Changing global_min to {}".format(log_default_ymin))
 
     ims = []
     for layer in range(layers):
@@ -231,10 +238,6 @@ def plot_hist_2D(plotname, xedges, yedges, hist, varx_label=None, vary_label=Non
 
 def plot_hist_2D_events(plotname, valsx, valsy, nbinsx, xmin, xmax, nbinsy, ymin, ymax, varx_label=None, vary_label=None, log=False):
 
-
-    print(valsx)
-    print(valsy)
-
     xedges = np.linspace(xmin, xmax, nbinsx)
     yedges = np.linspace(ymin, ymax, nbinsy)
 
@@ -312,6 +315,7 @@ if __name__ == "__main__":
     import logging
     logging.basicConfig()
     logging.getLogger("tfhelpers").setLevel(logging.DEBUG)
+    logging.getLogger(__name__).setLevel(logging.DEBUG)
 
     from .tfhelpers import get_single_neuron_function, get_max_activation_events
 
@@ -319,7 +323,7 @@ if __name__ == "__main__":
     # meme.setOptions(overrideCache="/scratch-local/nhartmann/meme_cache")
 
     if len(sys.argv) < 2:
-        c = ClassificationProject(os.path.expanduser("~/p/scripts/keras/008-allhighlevel/all_highlevel_985"))
+        c = ClassificationProject("/project/etp/nhartmann/p/keras/021-check-low-vs-high-fewvar/all_high")
     else:
         c = ClassificationProject(sys.argv[1])
 
@@ -345,9 +349,11 @@ if __name__ == "__main__":
 
         plot_NN_vs_var_2D("mt_vs_met.pdf", means=mean_signal,
                           scorefun=c.evaluate,
-                          var1_index=c.branches.index("met"), var1_range=(0, 1000, 10),
-                          var2_index=c.branches.index("mt"), var2_range=(0, 500, 10),
-                          var1_label="met [GeV]", var2_label="mt [GeV]")
+                          varx_index=c.branches.index("met"),
+                          vary_index=c.branches.index("mt"),
+                          nbinsx=100, xmin=0, xmax=1000,
+                          nbinsy=100, ymin=0, ymax=500,
+                          varx_label="met [GeV]", vary_label="mt [GeV]")
 
 
         plot_NN_vs_var_2D_all("mt_vs_met_all.pdf", means=mean_signal,
@@ -358,9 +364,11 @@ if __name__ == "__main__":
 
         plot_NN_vs_var_2D("mt_vs_met_crosscheck.pdf", means=mean_signal,
                           scorefun=get_single_neuron_function(c.model, layer=3, neuron=0, scaler=c.scaler),
-                          var1_index=c.branches.index("met"), var1_range=(0, 1000, 10),
-                          var2_index=c.branches.index("mt"), var2_range=(0, 500, 10),
-                          var1_label="met [GeV]", var2_label="mt [GeV]")
+                          varx_index=c.branches.index("met"),
+                          vary_index=c.branches.index("mt"),
+                          nbinsx=100, xmin=0, xmax=1000,
+                          nbinsy=100, ymin=0, ymax=500,
+                          varx_label="met [GeV]", vary_label="mt [GeV]")
 
 
     def test_max_act():
@@ -426,6 +434,29 @@ if __name__ == "__main__":
             log=True,
         )
 
+        plot_hist_2D_events(
+            "mt_vs_met_signal.pdf",
+            utrf_x_test[c.y_test==1][:,c.branches.index("met")],
+            utrf_x_test[c.y_test==1][:,c.branches.index("mt")],
+            100, 0, 1000,
+            100, 0, 500,
+            varx_label="met [GeV]",
+            vary_label="mt [GeV]",
+            log=True,
+        )
+
+        plot_hist_2D_events(
+            "mt_vs_met_backgound.pdf",
+            utrf_x_test[c.y_test==0][:,c.branches.index("met")],
+            utrf_x_test[c.y_test==0][:,c.branches.index("mt")],
+            100, 0, 1000,
+            100, 0, 500,
+            varx_label="met [GeV]",
+            vary_label="mt [GeV]",
+            log=True,
+        )
+
+
         # plot_hist_2D_events(
         #     "apl_vs_output_actmax.pdf",
         #     events[:,c.branches.index("LepAplanarity")],

scripts/plot_single_neuron.py

+#!/usr/bin/env python
+
+import sys
+import argparse
+
+import numpy as np
+
+from KerasROOTClassification import ClassificationProject
+from KerasROOTClassification.plotting import get_mean_event, plot_NN_vs_var_2D
+from KerasROOTClassification.tfhelpers import get_single_neuron_function
+
+parser = argparse.ArgumentParser(description='Create various 2D plots for a single neuron')
+parser.add_argument("project_dir")
+parser.add_argument("output_filename")
+parser.add_argument("varx")
+parser.add_argument("vary")
+parser.add_argument("-m", "--mode", choices=["mean_sig", "mean_bkg"], default="mean_sig")
+parser.add_argument("-l", "--layer", type=int, help="Layer index (takes last layer by default)")
+parser.add_argument("-n", "--neuron", type=int, default=0, help="Neuron index (takes first neuron by default)")
+parser.add_argument("--log", action="store_true", help="Plot in color in log scale")
+parser.add_argument("--contour", action="store_true", help="Interpolate with contours")
+parser.add_argument("-b", "--nbins", default=20, type=int, help="Number of bins in x and y direction")
+parser.add_argument("-x", "--xrange", type=float, nargs="+", help="xrange (low, high)")
+parser.add_argument("-y", "--yrange", type=float, nargs="+", help="yrange (low, high)")
+args = parser.parse_args()
+
+c = ClassificationProject(args.project_dir)
+
+layer = args.layer
+neuron = args.neuron
+
+if layer is None:
+    layer = c.layers
+
+varx_index = c.branches.index(args.varx)
+vary_index = c.branches.index(args.vary)
+
+x_test = c.x_test
+
+percentilesx = np.percentile(x_test[:,varx_index], [1,99])
+percentilesy = np.percentile(x_test[:,vary_index], [1,99])
+
+if args.xrange is not None:
+    if len(args.xrange) < 3:
+        args.xrange.append(args.nbins)
+    varx_range = args.xrange
+else:
+    varx_range = (percentilesx[0], percentilesx[1], args.nbins)
+
+if args.yrange is not None:
+    if len(args.yrange) < 3:
+        args.yrange.append(args.nbins)
+    vary_range = args.yrange
+else:
+    vary_range = (percentilesy[0], percentilesy[1], args.nbins)
+
+if args.mode == "mean_sig":
+    means = get_mean_event(c.x_test, c.y_test, 1)
+elif args.mode == "mean_bkg":
+    means = get_mean_event(c.x_test, c.y_test, 0)
+
+plot_NN_vs_var_2D(args.output_filename,
+                  means=means,
+                  varx_index=varx_index,
+                  vary_index=vary_index,
+                  scorefun=get_single_neuron_function(c.model, layer, neuron, scaler=c.scaler),
+                  xmin=varx_range[0], xmax=varx_range[1], nbinsx=varx_range[2],
+                  ymin=vary_range[0], ymax=vary_range[1], nbinsy=vary_range[2],
+                  varx_label=args.varx, vary_label=args.vary,
+                  logscale=args.log, only_pixels=(not args.contour))