---
title: "xaiHydro Methods Reference"
subtitle: "SHAP, LIME, PDP and ALE — Mathematical Background and Implementation"
author:
  - name: "Sadikul Islam"
    affiliation: "Department of Hydrology, [Your Institution], India"
    email: "sadikul.islam@institution.ac.in"
date: "`r Sys.Date()`"
output:
  rmarkdown::html_vignette:
    toc: true
    toc_depth: 3
    number_sections: true
    fig_width: 7
    fig_height: 4
vignette: >
  %\VignetteIndexEntry{xaiHydro Methods Reference}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
bibliography: references.bib
link-citations: true
---

```{r setup, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment  = "#>",
  warning  = FALSE,
  message  = FALSE
)
```

---

> **Author**: Sadikul Islam  
> **Package**: xaiHydro v0.1.0  
> **Reference**: Islam, S. (2026). Explainable AI for Hydro-Climate Models.
> In: *Hydro-Climate Analytics: Remote Sensing, AI and Geospatial Modelling*
> (Springer).

---

# Overview

This vignette documents the mathematical foundations of each XAI method
implemented in `xaiHydro`, their computational implementation, and
practical guidance for hydro-climate applications. It serves as the
**Methods** section reference for the accompanying book chapter.

---

# Method 1: SHAP (SHapley Additive exPlanations)

## Mathematical definition

For a prediction $\hat{f}(\mathbf{x})$ and background dataset
$\mathcal{D}$, the SHAP value for feature $j$ and observation
$\mathbf{x}$ is defined as the weighted average of marginal contributions
across all feature coalitions [@lundberg2017unified]:

$$
\phi_j(\hat{f}, \mathbf{x}) = \sum_{S \subseteq \mathcal{F} \setminus \{j\}}
\frac{|S|!\,(|\mathcal{F}| - |S| - 1)!}{|\mathcal{F}|!}
\left[\hat{f}_{S \cup \{j\}}(\mathbf{x}_{S \cup \{j\}}) -
      \hat{f}_{S}(\mathbf{x}_{S})\right]
$$

where $\mathcal{F}$ is the full set of features and $S$ ranges over all
subsets excluding feature $j$.

## Key properties

| Property | Definition | Hydrology relevance |
|----------|-----------|---------------------|
| **Efficiency** | $\sum_j \phi_j = \hat{f}(\mathbf{x}) - \mathbb{E}[\hat{f}]$ | Contributions explain full gap from mean discharge |
| **Symmetry** | Symmetric features receive equal attribution | Ensures fair credit between correlated rain indices |
| **Dummy** | Non-contributing features receive $\phi_j = 0$ | Zero-importance variables identified automatically |
| **Linearity** | Additive under model combination | Supports ensemble model analysis |

## Computational implementation

`xaiHydro` uses a built-in permutation-based Monte Carlo SHAP estimator (Strumbelj and Kononenko, 2014) requiring no additional packages:

$$
\hat{\phi}_j = \frac{1}{B} \sum_{b=1}^{B}
\left[\hat{f}(\mathbf{x}^{(b)}_{+j}) - \hat{f}(\mathbf{x}^{(b)}_{-j})\right]
$$

where $\mathbf{x}^{(b)}_{+j}$ and $\mathbf{x}^{(b)}_{-j}$ are randomly
constructed versions of $\mathbf{x}$ with and without feature $j$
replaced by background values.

```r
# B = nsim Monte Carlo replicates
shap_vals <- hydro_shap(exp, nsim = 50, seed = 42)
```

**Practical guidance for hydrology:**

- Use `nsim ≥ 50` for stable importance rankings in final figures.
- Background dataset should span the full observed range of predictors
  (all seasons, wet and dry years).
- For daily streamflow models with 8–15 predictors, `nsim = 100` runs
  in approximately 2–3 minutes on a standard laptop.

---

# Method 2: LIME

## Mathematical definition

LIME fits a locally weighted linear surrogate $g \in G$ that approximates
the complex model $\hat{f}$ in the neighbourhood of observation
$\mathbf{x}'$ [@ribeiro2016should]:

$$
\xi(\mathbf{x}') = \arg\min_{g \in G}
\mathcal{L}(\hat{f},\, g,\, \pi_{\mathbf{x}'}) + \Omega(g)
$$

where $\pi_{\mathbf{x}'}(\mathbf{z}) = \exp(-d(\mathbf{x}', \mathbf{z})^2 / \sigma^2)$
is an exponential kernel weighting neighbourhood samples by proximity, and
$\Omega(g)$ is a complexity penalty (LASSO for linear $g$).

## Implementation

```r
lime_result <- hydro_lime(
  explainer    = exp,
  new_obs      = X[storm_day, , drop = FALSE],
  n_features   = 6,       # number of features in local model
  kernel_width = 0.75     # σ — controls neighbourhood size
)
```

## Selecting `kernel_width`

| `kernel_width` | Neighbourhood | Best for |
|---------------|---------------|----------|
| 0.25–0.5 | Very local | Highly non-linear models; single event explanation |
| 0.75 (default) | Moderate | General use |
| 1.0–2.0 | Broad | Smooth models; more stable but less local |

**Note:** Check local R² in the plot subtitle. An R² < 0.5 indicates
the surrogate does not capture local model behaviour adequately; try
increasing `kernel_width` or reducing `n_features`.

---

# Method 3: Partial Dependence Plots (PDP)

## Mathematical definition

The partial dependence function for feature set $S$ is
[@friedman2001greedy]:

$$
\hat{f}_S(\mathbf{x}_S) =
\mathbb{E}_{\mathbf{x}_C}\left[\hat{f}(\mathbf{x}_S,\, \mathbf{x}_C)\right]
= \int \hat{f}(\mathbf{x}_S,\, \mathbf{x}_C)\, d\mathbb{P}(\mathbf{x}_C)
$$

estimated by averaging over the observed marginal distribution:

$$
\hat{f}_S(\mathbf{x}_S) \approx \frac{1}{n} \sum_{i=1}^{n}
\hat{f}(\mathbf{x}_S,\, \mathbf{x}^{(i)}_C)
$$

## Limitation: correlated predictors

When features in $S$ and $C$ are correlated, the marginal average
integrates over covariate combinations that never occur in practice
(e.g., very high precipitation with very low soil moisture in monsoon
hydrology). ALE avoids this.

---

# Method 4: Accumulated Local Effects (ALE)

## Mathematical definition

ALE resolves the PDP extrapolation problem by integrating local
conditional effects over observed data slices [@apley2020visualizing]:

$$
\hat{f}^{\text{ALE}}(x_j) =
\int_{x_{j,\min}}^{x_j}
\mathbb{E}\left[
  \frac{\partial \hat{f}(X)}{\partial X_j}
  \Big| X_j = z_j
\right] dz_j - c
$$

where $c$ is a centering constant ensuring $\mathbb{E}[\hat{f}^{\text{ALE}}] = 0$.

## When to use ALE vs PDP

```r
# ALE — recommended when predictors are correlated
ale <- hydro_pdp(exp, variable = c("precipitation", "soil_moisture"),
                 type = "accumulated")

# PDP — acceptable when predictors are approximately independent
pdp <- hydro_pdp(exp, variable = "ndvi",
                 type = "partial")
```

**Decision rule**: Compute the Pearson correlation matrix of your
predictors. If any |r| > 0.5 between the target predictor and others,
prefer ALE.

---

# Method 5: Permutation variable importance

## Definition

Permutation importance for feature $j$ is [@breiman2001random]:

$$
\text{VI}_j = \frac{1}{B} \sum_{b=1}^{B}
\left[L\!\left(\hat{f}, \tilde{\mathbf{X}}^{(b)}_j\right) - L(\hat{f}, \mathbf{X})\right]
$$

where $\tilde{\mathbf{X}}^{(b)}_j$ is the dataset with feature $j$
randomly permuted in the $b$-th repetition, and $L$ is the loss function
(RMSE by default).

```r
imp <- hydro_importance(exp, loss_function = "rmse", B = 20)
```

---

# Choosing the right XAI method

| Question | Recommended method | Function |
|---|---|---|
| Which features matter globally? | Permutation importance + SHAP summary | `hydro_importance()`, `hydro_shap()` |
| How does feature X affect predictions on average? | PDP or ALE | `hydro_pdp()` |
| Why was *this specific day* predicted as a flood? | SHAP waterfall or breakdown | `plot_shap_waterfall()`, `hydro_breakdown()` |
| What simple rule explains this event locally? | LIME | `hydro_lime()` |
| Are there unexpected residual patterns? | Residual diagnostics | `hydro_residuals()` |

---

# Reporting checklist for book chapters

When using `xaiHydro` results in a publication:

- [ ] Report `nsim` used for SHAP computation (affects reproducibility).
- [ ] State whether PDP or ALE was used and justify choice (correlation structure).
- [ ] Report local R² for LIME explanations.
- [ ] Report `B` (repetitions) for permutation importance.
- [ ] Include a residual diagnostic panel to validate the base model.
- [ ] Cite the underlying method papers (see References below).
- [ ] Set and report `seed` for all stochastic computations.

---

# References