---
title: "Fair Interpolated Transport for Group-Fair Clustering"
author: "Jesse S. Ghashti, Warren Hare, John R. J. Thompson"
date: "`r format(Sys.Date())`"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Fair Interpolated Transport for Group-Fair Clustering}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.align = "center",
  fig.width = 6,
  fig.height = 5,
  message = FALSE,
  warning = FALSE
)
library(FITclust)
library(ggplot2)
```

# Introduction

The **FITclust** package implements **Fair Interpolated Transport** (FIT), an
algorithm-agnostic preprocessing framework for group-fair clustering. FIT moves
each group-conditional empirical distribution along its Wasserstein-2 geodesic
toward a shared barycenter at a tunable transport intensity `t` in the interval
from zero to one, then selects the smallest intensity at which a soft-fairness
violation falls within a tolerance `deltaFair`. This package accompanies the
paper _Algorithm-Agnostic Group-Fair Clustering via Fair Interpolated Transport_ 
by Ghashti, Hare, and Thompson (2026).

Features of this package include

- a Wasserstein-2 barycenter supported on the full sample, computed by the
  Alvarez-Esteban et al. (2016) fixed-point iteration,

- a McCann interpolation that produces the transported data at any intensity
  without recomputing the transport plans,

- three soft clustering families that return fuzzy membership matrices rather
  than hard assignments alone, centroid based, graph based, and model based, and

- a minimum-intervention rule that reports the smallest intensity meeting the
  fairness tolerance together with the full swept diagnostics.

## Package Overview

There are three clustering functions, one per family. Each takes the
data, the protected-group labels, and the number of clusters, and each carries
out the full transport sweep internally.

1. Fit the centroid family with `fitSKM()`, built on fuzzy c-means.

2. Fit the graph family with `fitSSC()`, built on soft normalized spectral
   clustering.

3. Fit the model family with `fitSMM()`, built on a Gaussian mixture.

The transport is conducted with the functions `buildTransport()`,
`wassersteinBarycenter()`, and `computeTransportMaps()`, and the fairness
diagnostics through `softViolation()` and `hardViolation()`, so the individual
steps can be examined on their own.

# Installation

Install the latest release version of **FITclust** from
[GitHub](https://github.com/ghashti-j/FITclust) or with the following:

```{r, eval = FALSE}
library(devtools)
install_github("ghashti-j/FITclust")
library(FITclust)
```

# Sample Usage

We follow the step-by-step demonstration from the Appendix A of Ghashti, Hare, and
Thompson (2026). The example is a bivariate mixture with two well-separated
clusters and a protected attribute that is imbalanced within each cluster.

1. First we generate the demonstration data. There are two clusters centreed at
`(-3, -3)` and `(3, 3)`, and two protected groups offset in opposite directions
on the second coordinate. The group counts are asymmetric within each cluster,
one hundred against two hundred, so the overall group marginals are balanced
while the baseline clustering is not fair.

```{r}
set.seed(42)
demoData <- rbind(
  data.frame(x1 = rnorm(100, -3, 1), x2 = rnorm(100, -3 - 0.25, 1), cluster = 1L, group = 0L),
  data.frame(x1 = rnorm(200, -3, 1), x2 = rnorm(200, -3 + 0.25, 1), cluster = 1L, group = 1L),
  data.frame(x1 = rnorm(200,  3, 1), x2 = rnorm(200,  3 - 0.25, 1), cluster = 2L, group = 0L),
  data.frame(x1 = rnorm(100,  3, 1), x2 = rnorm(100,  3 + 0.25, 1), cluster = 2L, group = 1L)
)
dataMat <- as.matrix(demoData[, c("x1", "x2")])
groupVec <- demoData$group
trueCluster <- demoData$cluster
cat("n =", nrow(dataMat),
    " group counts =", paste(table(groupVec), collapse = "/"),
    " cluster counts =", paste(table(trueCluster), collapse = "/"), "\n")
```

We plot the data by protected group. The two groups overlap heavily in each
cluster, so the imbalance is not visually obvious, yet it is enough to make an
unconstrained clustering unfair.

```{r, fig.align='center'}
ggplot(demoData, aes(x1, x2, shape = factor(group), fill = factor(group))) +
  geom_point(size = 2, colour = "black", stroke = 0.3, alpha = 0.7) +
  scale_shape_manual("Group", values = c("0" = 21, "1" = 24)) +
  scale_fill_manual("Group", values = c("0" = "#8ABF69", "1" = "#D08890")) +
  labs(x = expression(x[1]), y = expression(x[2])) +
  coord_fixed() + theme_bw() +
  theme(panel.grid = element_blank(), legend.position = "bottom")
```

2. We construct the transport interpolation. `buildTransport()`
computes the Wasserstein-2 barycenter on the full sample and the group transport
maps, and returns a closure `fn` that produces the transported data at any
intensity.

```{r}
alphaVec <- resolveAlpha("uniform", groupVec, sort(unique(groupVec)))
transport <- buildTransport(dataMat, groupVec, alphaVec, verbose = FALSE)
cat("barycenter atoms =", nrow(transport$barycenter),
    " converged =", transport$baryConverged,
    " iterations =", transport$baryIter, "\n")
```

We can compare the soft-fairness violation of a baseline fuzzy c-means fit on the
original data against a fit on the fully transported data at intensity one.

```{r}
baseFit <- fcm(dataMat, numClusters = 2, numStart = 5)
fullFit <- fcm(transport$fn(1), numClusters = 2, numStart = 5)
cat("Delta_soft at t = 0:", round(softViolation(baseFit$membership, groupVec), 3), "\n")
cat("Delta_soft at t = 1:", round(softViolation(fullFit$membership, groupVec), 3), "\n")
```

3. We now run the full centroid-family procedure with `fitSKM()`, sweeping the
transport intensity on a grid at the tolerance `deltaFair = 0.05`. The
returned object reports the selected intensity `tOptimal`, the baseline and final
violations, and the swept diagnostics in `history`.

```{r}
set.seed(1)
fitCentroid <- fitSKM(dataMat, groupVec, numClusters = 2, deltaFair = 0.05,
                      tSeq = seq(0, 1, by = 0.02), verbose = FALSE)
cat("t* =", fitCentroid$tOptimal,
    " Delta_soft:", round(fitCentroid$violationBaseline, 3),
    "->", round(fitCentroid$violation, 3), "\n")
```

4. The `history` data frame records the soft and hard violations at every grid
point. Plotting the soft-violation trace against the intensity shows where it
first crosses the tolerance, which is the selected `tOptimal`.

```{r, fig.align='center'}
hist <- fitCentroid$history
ggplot(hist, aes(t, violationSoft)) +
  geom_line() + geom_point(size = 1) +
  geom_hline(yintercept = 0.05, linetype = "dashed", colour = "#E31A1C") +
  geom_vline(xintercept = fitCentroid$tOptimal, linetype = "dotted", colour = "grey30") +
  labs(x = expression(t), y = expression(Delta[soft](t))) +
  theme_bw() + theme(panel.grid = element_blank())
```

5. We can count how many observations are reassigned between the baseline and
the fair solution, after aligning the two labelings so that cluster identities
match. This is the practical footprint of the fairness intervention.

```{r}
alignLabels <- function(current, reference) {
  overlap <- table(current, reference)
  mapping <- apply(overlap, 1, which.max)
  as.integer(mapping[as.character(current)])
}
baseLabels <- fitCentroid$clustersBaseline
fairLabels <- alignLabels(fitCentroid$clusters, baseLabels)
cat("reassigned:", sum(fairLabels != baseLabels),
    "of", length(baseLabels),
    sprintf("(%.1f%%)", 100 * mean(fairLabels != baseLabels)), "\n")
```

6. Finally we visualize the fair clustering in the original coordinates, with
observations colored by their fair cluster label and shaped by protected group.

```{r, fig.align='center'}
plotDF <- data.frame(x1 = dataMat[, 1], x2 = dataMat[, 2],
                     cluster = factor(fairLabels), group = factor(groupVec))
ggplot(plotDF, aes(x1, x2, shape = group, fill = cluster)) +
  geom_point(size = 2, colour = "black", stroke = 0.3) +
  scale_shape_manual("Group", values = c("0" = 21, "1" = 24)) +
  scale_fill_manual("Cluster", values = c("1" = "#4E9BC7", "2" = "#F4A460")) +
  labs(x = expression(x[1]), y = expression(x[2])) +
  coord_fixed() + theme_bw() +
  theme(panel.grid = element_blank(), legend.position = "bottom") +
  guides(fill = guide_legend(override.aes = list(shape = 22)),
         shape = guide_legend(override.aes = list(fill = "grey60")))
```

# The Three Clustering Families

We fit all three on the demonstration data and compare the selected intensities and
the violation reductions.

```{r}
set.seed(1)
fitGraph <- fitSSC(dataMat, groupVec, numClusters = 2, deltaFair = 0.05,
                   tSeq = seq(0, 1, by = 0.02), verbose = FALSE)
fitModel <- fitSMM(dataMat, groupVec, numClusters = 2, deltaFair = 0.05,
                   tSeq = seq(0, 1, by = 0.02), verbose = FALSE)
summaryTab <- data.frame(
  Family = c("Centroid (fitSKM)", "Graph (fitSSC)", "Model (fitSMM)"),
  tOptimal = c(fitCentroid$tOptimal, fitGraph$tOptimal, fitModel$tOptimal),
  DeltaSoftBaseline = round(c(fitCentroid$violationBaseline,
                              fitGraph$violationBaseline,
                              fitModel$violationBaseline), 3),
  DeltaSoftFair = round(c(fitCentroid$violation,
                          fitGraph$violation,
                          fitModel$violation), 3)
)
summaryTab
```

# Functions

The transport and fairness building blocks can be called on their own.

* Transport construction
    + `buildTransport()`: barycenter, transport maps, and the interpolation closure
    + `wassersteinBarycenter()`: the barycenter alone, with its transport plans
    + `computeTransportMaps()`: group-to-barycenter maps from existing plans
    + `resolveAlpha()`: barycenter weights, `"uniform"`, `"proportional"`, or `"inverseProportional"`
* Fairness diagnostics
    + `softViolation()`: worst-case soft group-share deviation on a membership matrix
    + `hardViolation()`: worst-case hard group-share deviation on a partition
* Base clustering routines
    + `fcm()`: fuzzy c-means
    + `spectralEmbed()`: self-tuning spectral embedding
    + `emGMM()`: Gaussian mixture by expectation maximization

**References**

* Alvarez-Esteban, P. C., del Barrio, E., Cuesta-Albertos, J. A., and C. Matran
  (2016). A fixed-point approach to barycenters in Wasserstein space. _Journal
  of Mathematical Analysis and Applications, 441_(2), 744-762.

* J. C. Bezdek (1981). _Pattern Recognition with Fuzzy Objective Function
  Algorithms_. Plenum Press, New York.

* Feldman, M., Friedler, S. A., Moeller, J., Scheidegger, C., and S.
  Venkatasubramanian (2015). Certifying and removing disparate impact.
  _Proceedings of the 21st ACM SIGKDD International Conference on Knowledge
  Discovery and Data Mining_, 259-268.

* Ghashti, J. S., Hare, W., and J. R. J. Thompson (2026). Algorithm-Agnostic Group-Fair Clustering via Fair Interpolated Transport. Submitted.

* Kleindessner, M., Samadi, S., Awasthi, P., and J. Morgenstern (2019).
  Guarantees for spectral clustering with fairness constraints. _Proceedings of
  the 36th International Conference on Machine Learning_, 3458-3467.

* McCann, R. J. (1997). A convexity principle for interacting gases. _Advances
  in Mathematics, 128_(1), 153-179.

* McLachlan, G. and T. Krishnan (2008). _The EM Algorithm and Extensions_,
  Second Edition. John Wiley & Sons.

* Ng, A., Jordan, M., and Y. Weiss (2001). On spectral clustering: Analysis and
  an algorithm. _Advances in Neural Information Processing Systems, 14_.

* Zelnik-Manor, L., and P. Perona (2004). Self-tuning spectral clustering.
  _Advances in Neural Information Processing Systems, 17_.
