diff --git a/examples/methods/plot_fairness_da.py b/examples/methods/plot_fairness_da.py new file mode 100644 index 00000000..31d75b57 --- /dev/null +++ b/examples/methods/plot_fairness_da.py @@ -0,0 +1,195 @@ +# %% +import matplotlib.pyplot as plt +import numpy as np +from fairlearn.datasets import fetch_acs_income +from sklearn.ensemble import RandomForestRegressor +from sklearn.metrics import mean_absolute_error +from sklearn.model_selection import train_test_split +from sklearn.preprocessing import StandardScaler + +from skada import LinearOTMapping, source_target_split + +# %% +# Load the dataset +df = fetch_acs_income(as_frame=True) +X_df = df["frame"] +y = df["target"].to_numpy() +sample_domain = X_df["SEX"].to_numpy() # 1: male, 2: female + +# Drop target and sensitive attribute to get feature matrix +X = X_df.drop(columns=["PINCP", "SEX"]).to_numpy() + +# Take 10% of data while preserving the distribution +X, _, y, _, sample_domain, _ = train_test_split( + X, y, sample_domain, test_size=0.9, stratify=sample_domain, random_state=42 +) + +# Normalize features +X = StandardScaler().fit_transform(X) + +# Re-label domains: source=1 (e.g. male), target=2 (e.g. female) +sample_domain = np.where(sample_domain == 1, 1, -1) + + +# %% + +X_source, X_target, y_source, y_target = source_target_split( + X, y, sample_domain=sample_domain +) +print(f"Source domain size: {X_source.shape[0]}") +print(f"Target domain size: {X_target.shape[0]}") +print(f"Source domain income mean: {y_source.mean()}") +print(f"Target domain income mean: {y_target.mean()}") +plt.figure(figsize=(10, 4)) +plt.subplot(1, 2, 1) +plt.hist(y_source, bins=50, alpha=0.5, label="Male (source)") +plt.hist(y_target, bins=50, alpha=0.5, label="Female (target)") +plt.legend() +plt.title("Income distribution before adaptation") +plt.xlabel("PINCP") + +# %% + +# Train regressor on source + +clf = RandomForestRegressor(n_estimators=5, random_state=31415) +clf.fit(X_source, y_source) +# Evaluate performance using MAE +mae_source = mean_absolute_error(y_source, clf.predict(X_source)) / y_source.mean() +mae_target = mean_absolute_error(y_target, clf.predict(X_target)) / y_target.mean() + +# Scatterplot of predicted vs true (no decision boundary in regression) +y_pred_source = clf.predict(X_source) +y_pred_target = clf.predict(X_target) + +plt.figure(2, figsize=(10, 4)) +plt.subplot(1, 2, 1) +plt.scatter(y_source, y_pred_source, alpha=0.5, label="Source") +plt.plot([y_source.min(), y_source.max()], [y_source.min(), y_source.max()], "k--") +plt.xlabel("True income") +plt.ylabel("Predicted income") +plt.title(f"Source (MAE={mae_source:.2f})") +plt.legend() + +plt.subplot(1, 2, 2) +plt.scatter(y_target, y_pred_target, alpha=0.5, label="Target") +plt.plot([y_target.min(), y_target.max()], [y_target.min(), y_target.max()], "k--") +plt.xlabel("True income") +plt.ylabel("Predicted income") +plt.title(f"Target (MAE={mae_target:.2f})") +plt.legend() + +plt.tight_layout() +plt.show() + +# %% +# ---------------------------------- +# Build OTDA pipeline for regression +clf_otda = LinearOTMapping( + RandomForestRegressor(n_estimators=100, random_state=31415, alpha=0.5) +) + +# modify y such that for the target domain there are only nan +y_for_fit = np.where(sample_domain == 1, y, np.nan) +clf_otda.fit(X, y_for_fit, sample_domain=sample_domain) + + +# Evaluate Mean Absolute Error (MAE) scores +mae_source = mean_absolute_error(y_source, clf_otda.predict(X_source)) / y_source.mean() +mae_target = mean_absolute_error(y_target, clf_otda.predict(X_target)) / y_target.mean() + +print(f"Mean Absolute Error (MAE) - Source: {mae_source:.2f}") +print(f"Mean Absolute Error (MAE) - Target: {mae_target:.2f}") + +# Predict +y_pred_source_ot = clf_otda.predict(X_source) +y_pred_target_ot = clf_otda.predict(X_target) + +# Plot predictions +plt.figure(3, figsize=(10, 4)) + +plt.subplot(1, 2, 1) +plt.scatter(y_source, y_pred_source_ot, alpha=0.5, label="Source") +plt.plot([y_source.min(), y_source.max()], [y_source.min(), y_source.max()], "k--") +plt.xlabel("True income") +plt.ylabel("Predicted income") +plt.title(f"OTDA Source (MAE={mae_source:.2f})") +plt.legend() + +plt.subplot(1, 2, 2) +plt.scatter(y_target, y_pred_target_ot, alpha=0.5, label="Target") +plt.plot([y_target.min(), y_target.max()], [y_target.min(), y_target.max()], "k--") +plt.xlabel("True income") +plt.ylabel("Predicted income") +plt.title(f"OTDA Target (MAE={mae_target:.2f})") +plt.legend() + +plt.tight_layout() +plt.show() + +# %% + + +def compute_demographic_parity_difference(y_pred, sensitive_attr): + """Compute the demographic parity difference between two groups.""" + group_1 = y_pred[sensitive_attr == 1] + group_2 = y_pred[sensitive_attr == -1] + + p1 = np.mean(group_1) + p2 = np.mean(group_2) + + # Scale the difference by the overall mean prediction + return np.abs(p1 - p2) + + +def compute_error_difference(y_pred, y_true, sensitive_attr): + group_1 = y_pred[sensitive_attr == 1] + group_2 = y_pred[sensitive_attr == -1] + + group_1_true = y_true[sensitive_attr == 1] + group_2_true = y_true[sensitive_attr == -1] + + error_1 = np.linalg.norm(group_1 - group_1_true) + error_2 = np.linalg.norm(group_2 - group_2_true) + + return np.abs(error_1 - error_2) + + +# %% +# Concatenate source predictions (before OTDA) with target predictions (before OTDA) +y_stacked_no_ot = np.concatenate([y_pred_source, y_pred_target]) + +# Concatenate source predictions (before OTDA) with target predictions (after OTDA) +y_stacked_ot = np.concatenate([y_pred_source, y_pred_target_ot]) + +# %% + +print( + "Demographic parity difference before OTDA:", + compute_demographic_parity_difference( + y_stacked_no_ot, sensitive_attr=sample_domain + ), +) +print( + "Demographic parity difference after OTDA:", + compute_demographic_parity_difference(y_stacked_ot, sensitive_attr=sample_domain), +) + +# compute percrentage improvement in demographic parity difference +improvement_dp = ( + compute_demographic_parity_difference(y_stacked_no_ot, sensitive_attr=sample_domain) + - compute_demographic_parity_difference(y_stacked_ot, sensitive_attr=sample_domain) +) / compute_demographic_parity_difference(y_stacked_no_ot, sensitive_attr=sample_domain) + +print(f"Percentage improvement in demographic parity difference: {improvement_dp:.2%}") +# %% + +print( + "Error difference before OTDA:", + compute_error_difference(y_stacked_no_ot, y, sensitive_attr=sample_domain), +) +print( + "Error difference after OTDA:", + compute_error_difference(y_stacked_ot, y, sensitive_attr=sample_domain), +) +# %% diff --git a/notes_fairness.md b/notes_fairness.md new file mode 100644 index 00000000..fc8a1293 --- /dev/null +++ b/notes_fairness.md @@ -0,0 +1,9 @@ +# Notes Fairness + +- Idée 1: utiliser les méthodes de DA pour faire de la fairness. + - Importer les datasets depuis fairlearn + - Regarder si les méthodes de DA marchent + - Marie: add datasets + - Roman: format datasets + - Linus: metrics +- Idée 2: regarder le setup un attribut sensible + domaine \ No newline at end of file diff --git a/pyproject.toml b/pyproject.toml index b6960294..33652377 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -69,6 +69,7 @@ doc = [ "torch", "torchvision", "skorch", + "fairlearn", ] [project.urls] diff --git a/skada/_mapping.py b/skada/_mapping.py index cb1caf90..6b96fa34 100644 --- a/skada/_mapping.py +++ b/skada/_mapping.py @@ -82,7 +82,13 @@ def fit_transform(self, X, y=None, *, sample_domain=None, **params): return self.transform(X, sample_domain=sample_domain, allow_source=True) def transform( - self, X, y=None, *, sample_domain=None, allow_source=False, **params + self, + X, + y=None, + *, + sample_domain=None, + allow_source=False, + **params, ) -> np.ndarray: # xxx(okachaiev): implement auto-infer for sample_domain X, sample_domain = check_X_domain( @@ -100,6 +106,7 @@ def transform( X_adapt, _ = source_target_merge( X_source, X_target, sample_domain=sample_domain ) + X_adapt = self.alpha * X_adapt + (1 - self.alpha) * X return X_adapt @abstractmethod @@ -141,11 +148,13 @@ def __init__( metric="sqeuclidean", norm=None, max_iter=100_000, + alpha=1.0, ): super().__init__() self.metric = metric self.norm = norm self.max_iter = max_iter + self.alpha = alpha def _create_transport_estimator(self): return da.EMDTransport( @@ -155,7 +164,9 @@ def _create_transport_estimator(self): ) -def OTMapping(base_estimator=None, metric="sqeuclidean", norm=None, max_iter=100000): +def OTMapping( + base_estimator=None, metric="sqeuclidean", norm=None, max_iter=100000, alpha=1.0 +): """OTmapping pipeline with adapter and estimator. See [6]_ for details. @@ -188,7 +199,7 @@ def OTMapping(base_estimator=None, metric="sqeuclidean", norm=None, max_iter=100 base_estimator = SVC(kernel="rbf") return make_da_pipeline( - OTMappingAdapter(metric=metric, norm=norm, max_iter=max_iter), + OTMappingAdapter(metric=metric, norm=norm, max_iter=max_iter, alpha=alpha), base_estimator, ) @@ -234,6 +245,7 @@ def __init__( norm=None, max_iter=1000, tol=10e-9, + alpha=1.0, ): super().__init__() self.reg_e = reg_e @@ -241,6 +253,7 @@ def __init__( self.norm = norm self.max_iter = max_iter self.tol = tol + self.alpha = alpha def _create_transport_estimator(self): return da.SinkhornTransport( @@ -259,6 +272,7 @@ def EntropicOTMapping( max_iter=1000, reg_e=1.0, tol=1e-8, + alpha=1.0, ): """EntropicOTMapping pipeline with adapter and estimator. @@ -298,7 +312,12 @@ def EntropicOTMapping( return make_da_pipeline( EntropicOTMappingAdapter( - metric=metric, norm=norm, max_iter=max_iter, reg_e=reg_e, tol=tol + metric=metric, + norm=norm, + max_iter=max_iter, + reg_e=reg_e, + tol=tol, + alpha=alpha, ), base_estimator, ) @@ -352,6 +371,7 @@ def __init__( max_iter=10, max_inner_iter=200, tol=10e-9, + alpha=1.0, ): super().__init__() self.reg_e = reg_e @@ -361,6 +381,7 @@ def __init__( self.max_iter = max_iter self.max_inner_iter = max_inner_iter self.tol = tol + self.alpha = alpha def _create_transport_estimator(self): assert self.norm in ["lpl1", "l1l2"], "Unknown norm" @@ -388,6 +409,7 @@ def ClassRegularizerOTMapping( reg_e=1.0, reg_cl=0.1, tol=1e-8, + alpha=1.0, ): """ClassRegularizedOTMapping pipeline with adapter and estimator. @@ -435,6 +457,7 @@ def ClassRegularizerOTMapping( reg_e=reg_e, reg_cl=reg_cl, tol=tol, + alpha=alpha, ), base_estimator, ) @@ -468,10 +491,11 @@ class LinearOTMappingAdapter(BaseOTMappingAdapter): adaptation. arXiv preprint arXiv:1905.10155. """ - def __init__(self, reg=1e-08, bias=True): + def __init__(self, reg=1e-08, bias=True, alpha=1.0): super().__init__() self.reg = reg self.bias = bias + self.alpha = alpha def _create_transport_estimator(self): return da.LinearTransport(reg=self.reg, bias=self.bias) @@ -481,6 +505,7 @@ def LinearOTMapping( base_estimator=None, reg=1.0, bias=True, + alpha=1.0, ): """Returns a the linear OT mapping method with adapter and estimator. @@ -514,6 +539,7 @@ def LinearOTMapping( LinearOTMappingAdapter( reg=reg, bias=bias, + alpha=alpha, ), base_estimator, ) @@ -595,11 +621,12 @@ class MultiLinearMongeAlignmentAdapter(BaseAdapter): """ - def __init__(self, reg=1e-08, bias=True, test_time=False): + def __init__(self, reg=1e-08, bias=True, test_time=False, alpha=1.0): super().__init__() self.reg = reg self.bias = bias self.test_time = test_time + self.alpha = alpha def fit(self, X, y=None, *, sample_domain=None): """Fit adaptation parameters. @@ -685,7 +712,13 @@ def fit_transform(self, X, y=None, sample_domain=None, **params): return self.transform(X, sample_domain=sample_domain, allow_source=True) def transform( - self, X, y=None, *, sample_domain=None, allow_source=False, **params + self, + X, + y=None, + *, + sample_domain=None, + allow_source=False, + **params, ) -> np.ndarray: X, sample_domain = check_X_domain( X, sample_domain, allow_multi_source=True, allow_multi_target=True @@ -697,11 +730,13 @@ def transform( A, b = self.mappings_[domain] X_adapt[sel] = X[sel].dot(A) + b + X_adapt = self.alpha * X_adapt + (1 - self.alpha) * X + return X_adapt def MultiLinearMongeAlignment( - base_estimator=None, reg=1e-08, bias=True, test_time=False + base_estimator=None, reg=1e-08, bias=True, test_time=False, alpha=1.0 ): """MultiLinearMongeAlignment pipeline with adapter and estimator. @@ -747,7 +782,9 @@ def MultiLinearMongeAlignment( base_estimator = LogisticRegression() return make_da_pipeline( - MultiLinearMongeAlignmentAdapter(reg=reg, bias=bias, test_time=test_time), + MultiLinearMongeAlignmentAdapter( + reg=reg, bias=bias, test_time=test_time, alpha=alpha + ), base_estimator, ) @@ -837,10 +874,11 @@ class CORALAdapter(BaseAdapter): In Advances in Computer Vision and Pattern Recognition, 2017. """ - def __init__(self, reg="auto", assume_centered=False): + def __init__(self, reg="auto", assume_centered=False, alpha=1.0): super().__init__() self.reg = reg self.assume_centered = assume_centered + self.alpha = alpha def fit(self, X, y=None, sample_domain=None): """Fit adaptation parameters. @@ -897,7 +935,13 @@ def fit_transform(self, X, y=None, *, sample_domain=None, **params): return self.transform(X, sample_domain=sample_domain, allow_source=True) def transform( - self, X, y=None, *, sample_domain=None, allow_source=False, **params + self, + X, + y=None, + *, + sample_domain=None, + allow_source=False, + **params, ) -> np.ndarray: X, sample_domain = check_X_domain( X, @@ -927,6 +971,8 @@ def transform( X_adapt, _ = source_target_merge( X_source_adapt, X_target_adapt, sample_domain=sample_domain ) + X_adapt = self.alpha * X_adapt + (1 - self.alpha) * X + return X_adapt @@ -934,6 +980,7 @@ def CORAL( base_estimator=None, reg="auto", assume_centered=False, + alpha=1.0, ): """CORAL pipeline with adapter and estimator. @@ -968,7 +1015,7 @@ def CORAL( base_estimator = SVC(kernel="rbf") return make_da_pipeline( - CORALAdapter(reg=reg, assume_centered=assume_centered), + CORALAdapter(reg=reg, assume_centered=assume_centered, alpha=alpha), base_estimator, ) @@ -1017,7 +1064,9 @@ class MMDLSConSMappingAdapter(BaseAdapter): In ICML, 2013. """ - def __init__(self, gamma, reg_k=1e-10, reg_m=1e-10, tol=1e-5, max_iter=100): + def __init__( + self, gamma, reg_k=1e-10, reg_m=1e-10, tol=1e-5, max_iter=100, alpha=1.0 + ): super().__init__() self.gamma = gamma self.reg_k = reg_k @@ -1026,6 +1075,7 @@ def __init__(self, gamma, reg_k=1e-10, reg_m=1e-10, tol=1e-5, max_iter=100): self.max_iter = max_iter self.W_ = None self.B_ = None + self.alpha = alpha def _mapping_optimization(self, X_source, X_target, y_source): """Mapping optimization""" @@ -1149,7 +1199,13 @@ def fit_transform(self, X, y=None, sample_domain=None, **params): return self.transform(X, sample_domain=sample_domain, allow_source=True) def transform( - self, X, y=None, *, sample_domain=None, allow_source=False, **params + self, + X, + y=None, + *, + sample_domain=None, + allow_source=False, + **params, ) -> np.ndarray: X, sample_domain = check_X_domain(X, sample_domain, allow_source=allow_source) @@ -1180,11 +1236,18 @@ def transform( X_adapt, _ = source_target_merge( X_source_adapt, X_target, sample_domain=sample_domain ) + X_adapt = self.alpha * X_adapt + (1 - self.alpha) * X return X_adapt def MMDLSConSMapping( - base_estimator=None, gamma=1.0, reg_k=1e-10, reg_m=1e-10, tol=1e-5, max_iter=100 + base_estimator=None, + gamma=1.0, + reg_k=1e-10, + reg_m=1e-10, + tol=1e-5, + max_iter=100, + alpha=1.0, ): """MMDLSConSMapping pipeline with adapter and estimator. @@ -1220,7 +1283,12 @@ def MMDLSConSMapping( return make_da_pipeline( MMDLSConSMappingAdapter( - gamma=gamma, reg_k=reg_k, reg_m=reg_m, tol=tol, max_iter=max_iter + gamma=gamma, + reg_k=reg_k, + reg_m=reg_m, + tol=tol, + max_iter=max_iter, + alpha=alpha, ), base_estimator, ) diff --git a/skada/tests/test_mapping.py b/skada/tests/test_mapping.py index 455e1530..7c0c3549 100644 --- a/skada/tests/test_mapping.py +++ b/skada/tests/test_mapping.py @@ -43,20 +43,26 @@ [ make_da_pipeline(OTMappingAdapter(), LogisticRegression()), OTMapping(), + OTMapping(alpha=0.5), make_da_pipeline(EntropicOTMappingAdapter(), LogisticRegression()), EntropicOTMapping(), + EntropicOTMapping(alpha=0.5), make_da_pipeline( ClassRegularizerOTMappingAdapter(norm="lpl1"), LogisticRegression() ), ClassRegularizerOTMapping(), + ClassRegularizerOTMapping(alpha=0.5), make_da_pipeline( ClassRegularizerOTMappingAdapter(norm="l1l2"), LogisticRegression() ), ClassRegularizerOTMapping(norm="l1l2"), + ClassRegularizerOTMapping(norm="l1l2", alpha=0.5), make_da_pipeline(LinearOTMappingAdapter(), LogisticRegression()), LinearOTMapping(), + LinearOTMapping(alpha=0.5), make_da_pipeline(MultiLinearMongeAlignmentAdapter(), LogisticRegression()), MultiLinearMongeAlignment(), + MultiLinearMongeAlignment(alpha=0.5), make_da_pipeline(CORALAdapter(), LogisticRegression()), pytest.param( CORALAdapter(reg=None), @@ -64,6 +70,7 @@ ), make_da_pipeline(CORALAdapter(reg=0.1), LogisticRegression()), CORAL(), + CORAL(alpha=0.5), pytest.param( make_da_pipeline(MMDLSConSMappingAdapter(gamma=1e-3), SVC()), marks=pytest.mark.skipif(not torch, reason="PyTorch not installed"),