Add files via upload

Author: makramchahine

Dockerfile

@@ -0,0 +1,16 @@
+FROM dolphonie1/causal_skipper:0.01
+
+# add custom shap lib. This dockerfile should be built from the root of this repo
+# with the shap dir inside
+WORKDIR /src
+RUN mkdir /src/drone_causality && mkdir /src/shap
+COPY ./shap /src/shap
+RUN pip install /src/shap
+
+# install opencv deps for shap
+RUN apt-get update && export DEBIAN_FRONTEND=noninteractive && apt-get install ffmpeg libsm6 libxext6  -y
+
+RUN pip install seaborn
+
+# add current repo contents
+COPY . /src/drone_causality/

analysis/analyze_study.py

@@ -0,0 +1,49 @@
+import argparse
+
+import joblib
+import optuna
+from optuna import Study
+from optuna.visualization import plot_optimization_history, plot_intermediate_values, plot_parallel_coordinate, \
+    plot_contour, plot_param_importances, plot_slice
+
+
+def analyze_study(study_name: str, storage_name: str):
+    study = optuna.create_study(study_name=study_name, storage=storage_name, load_if_exists=True,
+                                direction="minimize")
+    print(study.best_trial.params)
+    graph_study(study)
+
+
+def analyze_local(file_path: str):
+    study = joblib.load(file_path)
+    print(study.best_trial.params)
+    graph_study(study)
+
+
+def graph_study(study: Study):
+    fig = plot_optimization_history(study)
+    fig.show()
+    fig2 = plot_intermediate_values(study)
+    fig2.show()
+    fig3 = plot_parallel_coordinate(study)
+    fig3.show()
+    fig4 = plot_contour(study)
+    fig4.show()
+    fig5 = plot_param_importances(study)
+    fig5.show()
+    fig6 = plot_slice(study)
+    fig6.show()
+    best_params = study.best_params
+    print(best_params)
+    print(study.best_value)
+
+
+if __name__ == "__main__":
+    parser = argparse.ArgumentParser()
+    parser.add_argument("storage_path", type=str)
+    parser.add_argument("--study_name", type=str, default=None)
+    args = parser.parse_args()
+    if args.study_name is not None:
+        analyze_study(args.study_name, args.storage_path)
+    else:
+        analyze_local(args.storage_path)

analysis/calculate_output_noise.py

@@ -0,0 +1,39 @@
+# Created by Patrick Kao at 4/4/22
+import argparse
+
+import numpy as np
+from PIL import Image
+from tensorflow import keras
+
+from utils.vis_utils import parse_params_json
+from utils.model_utils import load_model_from_weights, get_readable_name, generate_hidden_list
+
+
+def calculate_output_noise(params_path: str, input_img_path: str, noise: float, n_trials=50):
+    input_img = np.expand_dims(np.array(Image.open(input_img_path), dtype=float), axis=0)
+    model_variances = {}
+    var_layer = keras.layers.GaussianNoise(stddev=noise)
+    for local_path, model_path, model_params in parse_params_json(params_path):
+        model = load_model_from_weights(model_params, model_path)
+        imgs = []
+        for _ in range(n_trials):
+            # might want to try with nonzero hiddens at some point
+            hiddens = generate_hidden_list(model=model, return_numpy=True)
+            noise_img = var_layer(input_img, training=True)
+            output = model.predict([noise_img, *hiddens])
+            imgs.append(output[0])  # don't save hidden output
+
+        channel_variances = np.var(imgs, axis=0)
+        avg_variance = np.mean(channel_variances)
+        model_variances[get_readable_name(model_params)] = avg_variance
+
+    return model_variances
+
+
+if __name__ == "__main__":
+    parser = argparse.ArgumentParser()
+    parser.add_argument("params_path")
+    parser.add_argument("input_img")
+    parser.add_argument("noise", type=float)
+    args = parser.parse_args()
+    print(calculate_output_noise(args.params_path, args.input_img, args.noise))

analysis/grad_cam.py

@@ -0,0 +1,99 @@
+# Created by Patrick Kao at 3/9/22
+from math import ceil
+from typing import Optional, Sequence, Union
+
+import tensorflow as tf
+from numpy import ndarray
+from tensorflow import Tensor
+from tensorflow.python.keras import Model
+from tensorflow.python.keras.layers import Conv2D
+from tensorflow.python.keras.models import Functional
+
+from utils.vis_utils import image_grid
+from utils.model_utils import load_model_from_weights, load_model_no_params, ModelParams
+
+
+def compute_gradcam(img: Union[Tensor, ndarray], grad_model: Functional, hiddens: Sequence[Tensor],
+                    pred_index: Optional[Sequence[Tensor]] = None):
+    heatmaps, hiddens = _compute_gradcam(img=img, grad_model=grad_model, hiddens=hiddens, pred_index=pred_index)
+    avg_heat = tf.math.add_n(heatmaps)
+    avg_heat = tf.expand_dims(avg_heat, axis=-1)
+    return avg_heat, hiddens
+
+
+def compute_gradcam_tile(img: Union[Tensor, ndarray], grad_model: Functional, hiddens: Sequence[Tensor],
+                         pred_index: Optional[Sequence[Tensor]] = None):
+    heatmaps, hiddens = _compute_gradcam(img=img, grad_model=grad_model, hiddens=hiddens, pred_index=pred_index)
+    num_rows = ceil(len(heatmaps) / 2)
+    return image_grid(imgs=heatmaps, rows=num_rows, cols=2), hiddens
+
+
+def _compute_gradcam(img: Union[Tensor, ndarray], grad_model: Functional, hiddens: Sequence[Tensor],
+                     pred_index: Optional[Sequence[Tensor]] = None):
+    """
+    Adaptation of grad-cam code at https://keras.io/examples/vision/grad_cam/ with
+    the following adjustments:
+
+    - because we want the impact of pixels not on a class decision but on any of the 4 axes, sum heatmaps for all 4 outputs
+    - because we don't care about positive or negative impact, drop the ReLU (wasn't in this implementation to begin with)
+    - Before adding heatmaps together, take absolute value of each heatmap because don't care if positive or negative contribution to direction
+    - idea (not implemented): instead of just upscaling, multiply grad-cam heatmap against visual backprop heatmap to model actual network weights
+    :return tuple of image tensor (height, width, 1) with shape of final conv layer, list of hidden vectors each of shape
+    (batch, hidden_dim)
+    """
+    if pred_index is None:
+        pred_index = range(grad_model.output_shape[1][-1])
+
+    # Then, we compute the gradient of the top predicted class for our input image
+    # with respect to the activations of the last conv layer
+    with tf.GradientTape() as tape:
+        out = grad_model([img, *hiddens])
+        last_conv_layer_output = out[0]
+        preds = out[1]
+        hiddens = out[2:]
+
+    heatmaps = []
+    # for each element of preds, compute gradient of last_conv_out wrt this element of pred, abs and sum these gradients
+    # strip batch dim
+    # jacobian shape 4x last_conv_layer_output.shape where each element is gradient, preds[:,i] wrt last_conv_layer_out
+    grads = tape.jacobian(preds, last_conv_layer_output)[0]
+    last_conv_layer_output = last_conv_layer_output[0]
+    for pred in pred_index:
+        # This is the gradient of the output neuron (top pred1icted or chosen)
+        # with regard to the output feature map of the last conv layer
+        grad = grads[pred]
+
+        # This is a vector where each entry is the mean intensity of the gradient
+        # over a specific feature map channel
+        pooled_grads = tf.reduce_mean(grad, axis=(0, 1, 2))
+
+        # We multiply each channel in the feature map array
+        # by "how important this channel is" with regard to the top predicted class
+        # then sum all the channels to obtain the heatmap class activation
+        heatmap = last_conv_layer_output @ pooled_grads[..., tf.newaxis]
+        heatmap = tf.squeeze(heatmap)
+
+        # patrick edit: absolute value heatmaps to not discount/cancel negative and positive contributions
+        heatmap = tf.math.abs(heatmap)
+
+        heatmaps.append(heatmap)
+
+    return heatmaps, hiddens
+
+
+def get_last_conv(model_path: str, model_params: Optional[ModelParams] = None) -> Model:
+    if model_params is not None:
+        model_params.single_step = True
+        vis_model = load_model_from_weights(model_params, checkpoint_path=model_path)
+    else:
+        vis_model = load_model_no_params(model_path, single_step=True)
+
+    # get last conv layer
+    # cleave off only convolutional head
+    conv_layers = [layer for layer in vis_model.layers if isinstance(layer, Conv2D)]
+
+    # First, we create a model that maps the input image to the activations
+    # of the last conv layer as well as the output predictions
+    return tf.keras.models.Model(
+        [vis_model.inputs], [conv_layers[-1].output, *vis_model.output]
+    )

analysis/input_grad.py

@@ -0,0 +1,32 @@
+# Created by Patrick Kao at 3/10/22
+from typing import Sequence, Union
+
+import tensorflow as tf
+from numpy import ndarray
+from tensorflow import Tensor
+from tensorflow.python.keras.models import Functional
+
+
+def compute_input_grad(img: Union[Tensor, ndarray], model: Functional, hiddens: Sequence[Tensor]):
+    """
+    Computes gradients of model output with respect to img
+    :param img:
+    :param model:
+    :param hiddens:
+    :return: tuple of image tensor (height, width, 1) with shape of input img, list of hidden vectors each of shape
+    (batch, hidden_dim)
+    """
+    with tf.GradientTape() as tape:
+        tape.watch(img)
+        out = model([img, *hiddens])
+        preds = out[0]
+        hiddens = out[1:]
+
+    grads = tape.jacobian(preds, img)[0]  # shape: 4 x 1 x height x width x channels
+    grads = tf.math.abs(grads)  # take absolute value so + and - impacts don't cancel each other out
+
+    heatmap = tf.math.reduce_sum(grads, axis=0)  # shape 1 x height x width x channels
+    heatmap = tf.squeeze(heatmap, axis=0)
+    # convert heatmap to black and white by summing channels
+    heatmap = tf.math.reduce_sum(heatmap, axis=-1, keepdims=True)
+    return heatmap, hiddens, None

analysis/lipschitz_constant.py

@@ -0,0 +1,134 @@
+# Created by Patrick Kao at 4/4/22
+import argparse
+import json
+import os
+from pathlib import Path
+from typing import Sequence, Optional, Dict, Tuple
+
+import matplotlib.pyplot as plt
+import numpy as np
+from numpy import ndarray
+from tqdm import tqdm
+
+from keras_models import IMAGE_SHAPE
+from utils.data_utils import image_dir_generator
+from utils.model_utils import ModelParams, load_model_from_weights, generate_hidden_list, get_readable_name
+from utils.vis_utils import parse_params_json
+
+
+def calculate_lipschitz_constant(model_path: str, model_params: ModelParams, sequence_path: str,
+                                 reverse_channels: bool) -> Sequence[ndarray]:
+    model = load_model_from_weights(model_params, model_path)
+    hiddens = generate_hidden_list(model=model, return_numpy=True)
+    all_hiddens = []  # list of list of arrays with shape num_timesteps x num_hiddens x hidden_dim
+    for i, img in tqdm(enumerate(image_dir_generator(sequence_path, IMAGE_SHAPE, reverse_channels))):
+        all_hiddens.append(hiddens)
+        out = model.predict([img, *hiddens])
+        hiddens = out[1:]  # list num_hidden long, each el is hidden_dim,
+
+    # flatten batch dim
+    all_hiddens = [[np.squeeze(hid, axis=0) for hid in step_hid] for step_hid in all_hiddens]
+    # crete list with same shape as hidden vectors where contents are lipschitz values of each dimension
+    lip = [np.zeros_like(h) for h in all_hiddens[0]]
+    for i in range(len(all_hiddens) - 1):
+        current_hiddens = all_hiddens[i]
+        next_hiddens = all_hiddens[i + 1]
+        diff = [np.abs(n - c) for n, c in zip(next_hiddens, current_hiddens)]
+        lip = [np.maximum(l, d) for l, d in zip(lip, diff)]
+    return lip
+
+
+def graph_lipschitz_constant(lip_mean: Dict[str, ndarray], lip_std: Optional[Dict[str, ndarray]],
+                             display_result: bool = False,
+                             save_path: Optional[str] = None):
+    plt.clf()
+    plt.xlabel("Node Rank")
+    plt.ylabel("Lipschitz Constant of Hidden State Nodes")
+
+    # concat all hidden dims into 1d
+    for model_name, hiddens in lip_mean.items():
+        # hiddens shape: flattened_hidden_dim
+        sort_order = np.argsort(hiddens)
+
+        lip_sorted = hiddens[sort_order]
+        lip_x = np.linspace(0, 1, num=len(lip_sorted))
+        plt.plot(lip_x, lip_sorted, label=model_name)
+        if lip_std is not None:
+            std = lip_std[model_name]
+            std_sorted = std[sort_order]
+            plt.fill_between(lip_x, lip_sorted + std_sorted, lip_sorted - std_sorted, alpha=0.5)
+
+    plt.legend(loc="upper left")
+    # note that save needs to happen before plt.show() b/c plt.show() clears figures
+    if save_path is not None:
+        plt.savefig(save_path)
+
+    if display_result:
+        plt.show()
+
+
+def params_lipschitz_constant(datasets_json: str, params_path: str, display_result: bool = False,
+                              save_dir: Optional[str] = None) -> Tuple[Dict, Dict]:
+    with open(datasets_json, "r") as f:
+        datasets: Dict[str, Tuple[str, bool]] = json.load(f)
+
+    all_mean = {}
+    all_std = {}
+    if save_dir:
+        Path(save_dir).mkdir(parents=True, exist_ok=True)
+    for local_path, model_path, model_params in parse_params_json(params_path):
+        model_name = get_readable_name(model_params)
+        if model_name == "tcn":
+            # TCN has no hidden state
+            continue
+
+        lips = []
+        for dataset_name, (data_path, reverse_channels, csv_path) in datasets.items():
+            lip = calculate_lipschitz_constant(model_path, model_params, data_path, reverse_channels=reverse_channels)
+            lips.append(lip)
+
+        # lips shape: num_datasets x num_hiddens x hidden dim
+        # shape : num_dataset x flattened_hidden_dim
+        lips_flat = np.array([np.hstack(dataset) for dataset in lips])
+        lip_mean = np.mean(lips_flat, axis=0)
+        lip_std = np.std(lips_flat, axis=0)
+        model_lip_mean = {model_name: lip_mean}
+        model_lip_std = {model_name: lip_std}
+        if display_result or save_dir is not None:
+            save_path = os.path.join(save_dir, model_name) if save_dir is not None else None
+            graph_lipschitz_constant(lip_mean=model_lip_mean, lip_std=model_lip_std, display_result=display_result,
+                                     save_path=save_path)
+
+        all_mean.update(model_lip_mean)
+        all_std.update(model_lip_std)
+
+    # graph all lipschitz
+    if display_result or save_dir is not None:
+        graph_lipschitz_constant(lip_mean=all_mean, lip_std=all_std, display_result=display_result,
+                                 save_path=os.path.join(save_dir, "all_lipschitz"))
+
+    if save_dir is not None:
+        with open(os.path.join(save_dir, "lip_data.json"), "w") as f:
+            lip_data = {"means": convert_values_to_list(all_mean), "stds": convert_values_to_list(all_std)}
+            json.dump(lip_data, f)
+
+    return all_mean, all_std
+
+
+def convert_values_to_list(to_convert: Dict[str, ndarray]):
+    to_ret = {}
+    for key, np_arr in to_convert.items():
+        to_ret[key] = list(np_arr)
+
+    return to_ret
+
+
+if __name__ == "__main__":
+    parser = argparse.ArgumentParser()
+    parser.add_argument("datasets_json")
+    parser.add_argument("params_path")
+    parser.add_argument("--display_result", action="store_true")
+    parser.add_argument("--save_dir", type=str, default="lipschitz_out")
+    args = parser.parse_args()
+    params_lipschitz_constant(datasets_json=args.datasets_json, params_path=args.params_path,
+                              display_result=args.display_result, save_dir=args.save_dir)