In the following notebook we will show how you can use the CARLA library.

How to use CARLA¶

[1]:

from IPython.display import display

import warnings
warnings.filterwarnings("ignore")

%load_ext autoreload
%autoreload 2

Data¶

Before we can do anything we need some data. Using CARLA, you have several options to handle data.

You could import one of the datasets from our OnlineCatalog.
However, you may want to use your own data instead. This can easily be done by using the CsvCatalog.

Using the OnlineCatalog¶

Using the OnlineCatalog is very easy. Currently, we support four data sets: “heloc”, “adult”, “compas”, and “give_me_credit”. In the examples below, we will use the adult data set. Below, we demonstrate how you can use the OnlineCatalog.

[2]:

from carla.data.catalog import OnlineCatalog

# load catalog dataset
data_name = "adult"
dataset = OnlineCatalog(data_name)

Using TensorFlow backend.

[INFO] Using Python-MIP package version 1.12.0 [model.py <module>]

Below, we take a look at how you can add your own data to CARLA.

Using the CsvCatalog¶

For the CsvCatalog there are 5 attributes. The file_path should be the path of the csv file you want to use. Then we have two different types of features, continuous and categorical, of which some can be immutable. Finally, the target attribute is the column which contains the targets/labels. For the adult income data set, this will be “income”, i.e., whether an individual earned more or less than $50.000.

Note that when using the CsvCatalog the data should already be cleaned; e.g., your .csv file should not contain any NaNs. Moreover, also make sure that the categorical variables are binary encoded, i.e., $x_j \in \{0,1\}$, if feature $j$ is a categorical variable (e.g., “workclass_private”). We are currently working on extensions to this.

[3]:

from carla.data.catalog import CsvCatalog

continuous = ["age", "fnlwgt", "education-num", "capital-gain", "hours-per-week", "capital-loss"]
categorical = ["marital-status", "native-country", "occupation", "race", "relationship", "sex", "workclass"]
immutable = ["age", "sex"]

dataset = CsvCatalog(file_path="adult.csv",
                     continuous=continuous,
                     categorical=categorical,
                     immutables=immutable,
                     target='income')

display(dataset.df.head())

	age	fnlwgt	education-num	capital-gain	...	occupation_Other	race_White	relationship_Non-Husband	sex_Male	workclass_Private
0	0.301370	0.044131	0.800000	0.02174	...	0.0	1.0	1.0	1.0	0.0
1	0.452055	0.048052	0.800000	0.00000	...	0.0	1.0	0.0	1.0	0.0
2	0.287671	0.137581	0.533333	0.00000	...	1.0	1.0	1.0	1.0	1.0
3	0.493151	0.150486	0.400000	0.00000	...	1.0	0.0	0.0	1.0	1.0
4	0.150685	0.220635	0.800000	0.00000	...	0.0	0.0	1.0	0.0	1.0

5 rows × 14 columns

ML Classifier¶

Now that we have the data loaded we also need a classification model. Again, you have two options:

You can easily define your own model. In our model documentation we describe how you can do that.
Here we will show how you can train one of our catalog models. Depending on your data and your use-case you might need to tweak the training hyperparameters. Note that if you are using the OnlineCatalog for your data, you can just set load_online to True and you don’t have to train anything!

For example, for the ANN used here we need to define the following hyperparameters: - learning rate - number of epochs - batch size - sizes of the hidden layers.

After defining the model using the MLModelCatalog, just call the train method with those parameters and you are good to go!

[4]:

from carla.models.catalog import MLModelCatalog

training_params = {"lr": 0.002, "epochs": 10, "batch_size": 1024, "hidden_size": [18, 9, 3]}

ml_model = MLModelCatalog(
    dataset,
    model_type="ann",
    load_online=False,
    backend="pytorch"
)

ml_model.train(
    learning_rate=training_params["lr"],
    epochs=training_params["epochs"],
    batch_size=training_params["batch_size"],
    hidden_size=training_params["hidden_size"]
)

Loaded model from /home/johan/carla/models/custom/ann_layers_18_9_3.pt
test accuracy for model: 0.8382208387942333

Counterfactual Explanations and Algorithmic Recourse¶

Now that we have both the data, and a model we can start using CARLA to generate counterfactuals. Again, you have two options:

You can pick a recourse method from the catalog.
Or you can implement one yourself using our recourse interface. If you would like to add a new method to the library, just submit a pull-request. :)

In the following example, we are getting negatively labeled samples for which we would like to find counterfactuals.

[5]:

from carla.models.negative_instances import predict_negative_instances
import carla.recourse_methods.catalog as recourse_catalog

factuals = predict_negative_instances(ml_model, dataset.df)
test_factual = factuals.iloc[:5]

display(test_factual)

	age	fnlwgt	education-num	capital-gain	...	occupation_Other	race_White	relationship_Non-Husband	sex_Male	workclass_Private
0	0.301370	0.044131	0.800000	0.02174	...	0.0	1.0	1.0	1.0	0.0
2	0.287671	0.137581	0.533333	0.00000	...	1.0	1.0	1.0	1.0	1.0
3	0.493151	0.150486	0.400000	0.00000	...	1.0	0.0	0.0	1.0	1.0
4	0.150685	0.220635	0.800000	0.00000	...	0.0	0.0	1.0	0.0	1.0
6	0.438356	0.100061	0.266667	0.00000	...	1.0	0.0	1.0	0.0	1.0

5 rows × 14 columns

Wachter et al (2018) (gradient method)¶

The recourse objective function looks as follows:

\begin{align} \delta_x^* & = argmin_{\delta_x, x+ \delta_x \in \mathcal{A}} \, \ell \big(h(x + \delta_x), 0.5)\big) + \lambda \cdot \, d(x + \delta_x, x), %\\ \end{align}

where $\lambda \geq 0$ is a trade-off parameter, $0.5$ is the probabilistic target, $\mathcal{A}$ is the feasible set of actions, and $\ell(\cdot, \cdot)$ is the binary-cross-entropy loss. The first term on the right-hand-side ensures that the model prediction corresponding to the counterfactual i.e., $h(x + \delta_x )$ is close to the favorable outcome with classification prediction $1$. The second term encourages low-cost recourses; for example, Wachter et al (2018) propose $\ell_1$ or $\ell_2$ distances to ensure that the distance between the factual instance $x$ and the counterfactual $\check{x} = x + \delta_x^*$ is small.

[6]:

hyperparams = {"loss_type": "BCE", "binary_cat_features": True}
recourse_method = recourse_catalog.Wachter(ml_model, hyperparams)
df_cfs = recourse_method.get_counterfactuals(test_factual)

display(df_cfs)

[INFO] Counterfactual Explanation Found [wachter.py wachter_recourse]
[INFO] Counterfactual Explanation Found [wachter.py wachter_recourse]
[INFO] Counterfactual Explanation Found [wachter.py wachter_recourse]
[INFO] Counterfactual Explanation Found [wachter.py wachter_recourse]
[INFO] Counterfactual Explanation Found [wachter.py wachter_recourse]

	age	capital-gain	capital-loss	education-num	fnlwgt	...	race_White	relationship_Non-Husband	sex_Male	workclass_Private	income
0	0.350025	0.070708	0.048941	0.849036	0.000182	...	1.0	1.0	0.0	0.042775	1.0
2	0.414486	0.125881	0.127741	0.658535	0.033042	...	0.0	1.0	0.0	0.873237	1.0
3	0.550512	0.058995	0.058867	0.459005	0.210247	...	0.0	0.0	1.0	0.950902	1.0
4	0.257907	0.106133	0.107396	0.905325	0.150003	...	0.0	1.0	0.0	0.894212	1.0
6	0.605591	0.167044	0.167439	0.433625	-0.056802	...	0.0	1.0	0.0	0.834089	1.0

5 rows × 14 columns

CCHVAE by Pawelczyk et al (2020) (manifold method)¶

Let $g: \mathcal{Z} \to \mathcal{X}$ be the decoder of a generative model (e.g., VAE). Let $e: \mathcal{X} \to \mathcal{Z}$ be the correpsonding encoder. We then encode the factual input x, for which we wish to find a counterfactual explanation, as follows: $e(x)=z$, and conduct the search in the latent space. Manifold-based methods solve an objective function that looks as follows:

\begin{align} \delta_z^* & = argmin_{\delta_z, g(z + \delta_z) \in \mathcal{A}} \, \ell \big(h(g(z + \delta_z), 0.5)\big) + \cdot \, d(z + \delta_z, z), %\\ \end{align}

where $\lambda \geq 0$ is a trade-off parameter, $0.5$ is the probabilistic target, and $\ell(\cdot, \cdot)$ is the binary-cross-entropy loss, and $\mathcal{A}$ is the feasible set of actions. The first term on the right-hand-side ensures that the model prediction corresponding to the counterfactual i.e., $h(g(z + \delta_z))$ is close to the favorable outcome with classification label $1$. The second term encourages low-cost recourses;

For example, Pawelczyk et al (2020) use random search in the latent space to approximate the above objective function, while Joshi et al (2019) use a gradient-based algorithm on a variant of the above objective function. We refer to the respective papers for more details

[7]:

hyperparams = {
    "data_name": dataset.name,
    "n_search_samples": 100,
    "p_norm": 1,
    "step": 0.1,
    "max_iter": 1000,
    "clamp": True,
    "binary_cat_features": True,
    "vae_params": {
        "layers": [len(ml_model.feature_input_order), 512, 256, 8],
        "train": True,
        "lambda_reg": 1e-6,
        "epochs": 5,
        "lr": 1e-3,
        "batch_size": 32,
    },
}

cchvae = recourse_catalog.CCHVAE(ml_model, hyperparams)
df_cfs = cchvae.get_counterfactuals(test_factual)

display(df_cfs)

[INFO] Start training of Variational Autoencoder... [models.py fit]
[INFO] [Epoch: 0/5] [objective: 0.381] [models.py fit]
[INFO] [ELBO train: 0.38] [models.py fit]
[INFO] [ELBO train: 0.14] [models.py fit]
[INFO] [ELBO train: 0.12] [models.py fit]
[INFO] [ELBO train: 0.12] [models.py fit]
[INFO] [ELBO train: 0.12] [models.py fit]
[INFO] ... finished training of Variational Autoencoder. [models.py fit]

	age	capital-gain	capital-loss	education-num	fnlwgt	...	race_White	sex_Male	workclass_Private	income
0	0.296011	0.036202	0.039719	0.601064	0.12027	...	1.0	1.0	0.737570	1.0
2	0.296000	0.036201	0.039718	0.601064	0.12027	...	1.0	1.0	0.737570	1.0
3	0.296004	0.036203	0.039718	0.601062	0.12027	...	1.0	1.0	0.737568	1.0
4	0.296001	0.036203	0.039718	0.601063	0.12027	...	1.0	1.0	0.737570	1.0
6	0.296014	0.036203	0.039718	0.601063	0.12027	...	1.0	1.0	0.737570	1.0

5 rows × 14 columns

FOCUS by Lucic et al (2021) (tree method)¶

Our library also supports sklearn and xgboost tree-based classifiers such as Random Forests, Decision Trees or Gradient Boosted Decision Trees. Those classifiers are needed for methods, which explicitly require the use of tree models (e.g., FeatureTweak and FOCUS).

A lot of recourse methods assume that the prediction model is differentiable. These methods do therefore not work on tree-based models. The FOCUS method avoids this problem by using a differentiable approximation to the original tree model. For a normal decision tree $t$ a node $j$ is activated if its parent node $p_j$ is activated and feature $x_{f_j}$ is above the threshold $\theta_j$.

\[\begin{split}t_j(x) = \begin{cases} 1 & \text{ if j is the root,} \\ t_{p_j}(x) \cdot 𝟙 [x_{f_j} > \theta_j] & \text{ if j is a left child,} \\ t_{p_j}(x) \cdot 𝟙 [x_{f_j} \leq \theta_j] & \text{ if j is a right child} \end{cases}\end{split}\]

This can be replaced by a differentiable approximation, where $\text{sig}$ is the sigmoid function.

\[\begin{split}\tilde{t}_j(x) = \begin{cases} 1 & \text{ if j is the root,} \\ \tilde{t}_{p_j}(x) \cdot \text{sig}(\theta_j - x_{f_j}) & \text{ if j is a left child,} \\ \tilde{t}_{p_j}(x) \cdot \text{sig}(x_{f_j} - \theta_j) & \text{ if j is a right child} \end{cases}\end{split}\]

Using that an approximation to the original tree model can be made, with the same structure and thresholds. However, now the model is differentiable, and we can do recourse as e.g. shown in the Wachter method.

Note that in the original work the authors removed the categorical features from the data, and we do the same.

You can use one of our catalog forest models if you want. However, we will take this as an opportunity to show how to easily implement your own model.

[8]:

ml_model = MLModelCatalog(dataset, "forest", backend="sklearn", load_online=False)
ml_model.train(max_depth=2, n_estimators=5, force_train=True)

factuals = predict_negative_instances(ml_model, dataset.df)
test_factual = factuals.iloc[:5]

display(test_factual)

balance on test set 0.24115334207077327, balance on test set 0.2338630406290957
model fitted with training score 0.7983016601135867 and test score 0.8052096985583224

	age	fnlwgt	education-num	...	occupation_Other	race_White	relationship_Non-Husband	sex_Male	workclass_Private
4	0.150685	0.220635	0.800000	...	0.0	0.0	1.0	0.0	1.0
5	0.273973	0.184219	0.866667	...	0.0	1.0	1.0	0.0	1.0
11	0.178082	0.087281	0.800000	...	0.0	0.0	0.0	1.0	0.0
12	0.082192	0.074410	0.800000	...	0.0	1.0	1.0	0.0	1.0
16	0.109589	0.111271	0.533333	...	1.0	1.0	1.0	1.0	0.0

5 rows × 14 columns

You just need to import our MLModel api, train a model, and define some methods. These are the methods shown in our Black-Box-Model example.

[9]:

from carla import MLModel
import xgboost

class XGBoostModel(MLModel):
    """The default way of implementing XGBoost
    https://xgboost.readthedocs.io/en/latest/python/python_intro.html"""

    def __init__(self, data):
        super().__init__(data)

        # get preprocessed data
        df_train = self.data.df_train
        df_test = self.data.df_test

        x_train = df_train[self.data.continuous]
        y_train = df_train[self.data.target]
        x_test = df_test[self.data.continuous]
        y_test = df_test[self.data.target]

        self._feature_input_order = self.data.continuous

        param = {
            "max_depth": 2,  # determines how deep the tree can go
            "objective": "binary:logistic",  # determines the loss function
            "n_estimators": 5,
        }
        self._mymodel = xgboost.XGBClassifier(**param)
        self._mymodel.fit(
                x_train,
                y_train,
                eval_set=[(x_train, y_train), (x_test, y_test)],
                eval_metric="logloss",
                verbose=True,
            )

    @property
    def feature_input_order(self):
        # List of the feature order the ml model was trained on
        return self._feature_input_order

    @property
    def backend(self):
        # The ML framework the model was trained on
        return "xgboost"

    @property
    def raw_model(self):
        # The black-box model object
        return self._mymodel

    @property
    def tree_iterator(self):
        # make a copy of the trees, else feature names are not saved
        booster_it = [booster for booster in self.raw_model.get_booster()]
        # set the feature names
        for booster in booster_it:
            booster.feature_names = self.feature_input_order
        return booster_it

    # The predict function outputs
    # the continuous prediction of the model
    def predict(self, x):
        return self._mymodel.predict(self.get_ordered_features(x))

    # The predict_proba method outputs
    # the prediction as class probabilities
    def predict_proba(self, x):
        return self._mymodel.predict_proba(self.get_ordered_features(x))

[10]:

ml_model = XGBoostModel(dataset)

factuals = predict_negative_instances(ml_model, dataset.df)
test_factual = factuals.iloc[:5]

display(test_factual)

[0]     validation_0-logloss:0.58480    validation_1-logloss:0.58251
[1]     validation_0-logloss:0.52502    validation_1-logloss:0.52138
[2]     validation_0-logloss:0.48735    validation_1-logloss:0.48210
[3]     validation_0-logloss:0.46126    validation_1-logloss:0.45584
[4]     validation_0-logloss:0.44237    validation_1-logloss:0.43615

	age	fnlwgt	education-num	capital-gain	...	occupation_Other	race_White	relationship_Non-Husband	sex_Male	workclass_Private
0	0.301370	0.044131	0.800000	0.02174	...	0.0	1.0	1.0	1.0	0.0
1	0.452055	0.048052	0.800000	0.00000	...	0.0	1.0	0.0	1.0	0.0
2	0.287671	0.137581	0.533333	0.00000	...	1.0	1.0	1.0	1.0	1.0
3	0.493151	0.150486	0.400000	0.00000	...	1.0	0.0	0.0	1.0	1.0
4	0.150685	0.220635	0.800000	0.00000	...	0.0	0.0	1.0	0.0	1.0

5 rows × 14 columns

Below we start generating counterfactuals using FOCUS.

[11]:

hyperparams = {
    "optimizer": "adam",
    "lr": 0.001,
    "n_class": 2,
    "n_iter": 1000,
    "sigma": 1.0,
    "temperature": 1.0,
    "distance_weight": 0.01,
    "distance_func": "l1",
}

focus = recourse_catalog.FOCUS(ml_model, hyperparams)
df_cfs = focus.get_counterfactuals(test_factual)
display(df_cfs)

[WARNING] From /home/johan/Dropbox/Documents/Master/HiWi/CARLA/env/lib/python3.7/site-packages/tensorflow/python/ops/array_ops.py:1354: add_dispatch_support.<locals>.wrapper (from tensorflow.python.ops.array_ops) is deprecated and will be removed in a future version.
Instructions for updating:
Use tf.where in 2.0, which has the same broadcast rule as np.where [deprecation.py new_func]
[WARNING] Calling GradientTape.gradient on a persistent tape inside its context is significantly less efficient than calling it outside the context (it causes the gradient ops to be recorded on the tape, leading to increased CPU and memory usage). Only call GradientTape.gradient inside the context if you actually want to trace the gradient in order to compute higher order derivatives. [backprop.py gradient]

	age	fnlwgt	education-num	capital-gain	hours-per-week
0	0.301444	0.044131	0.800015	0.070616	0.397955
1	0.452064	0.048052	0.800001	0.070650	0.122511
2	0.287699	0.137581	0.533315	0.070581	0.397937
3	0.493142	0.150486	0.400060	0.070592	0.397903
4	0.150696	0.220635	0.799994	0.070879	0.397962