From f4e979f3c72230e66891364e0039fced17cf6446 Mon Sep 17 00:00:00 2001 From: Tyler Reddy Date: Tue, 23 Jun 2026 11:46:24 -0600 Subject: [PATCH 1/2] WIP, ENH: exact solve with array backends * Related to gh-68, but for the opposite case of the exact solve. * This is a (very) early draft/prototype of supporting the Python array API standard (https://data-apis.org/array-api/latest/) for `GFDLClassifier`. Some of the sklearn API conformance tests will currently fail, but all numerical correctness tests with the NumPy backend should continue to pass. * The basic idea is to take an early peak at our potential for speedups when using JIT compiled/GPU backend array types, since i.e., hyperopt can be slow on CPU alone, and also because we'd just like to deliver the most performant estimators that we can, which also helps us conduct numerical experiments for the mathematicians with faster turnaround times. --- .github/constraints/deps.txt | 1 + pyproject.toml | 3 ++- src/gfdl/model.py | 12 ++++++++---- 3 files changed, 11 insertions(+), 5 deletions(-) diff --git a/.github/constraints/deps.txt b/.github/constraints/deps.txt index b586d0c..4f42bd8 100644 --- a/.github/constraints/deps.txt +++ b/.github/constraints/deps.txt @@ -1,3 +1,4 @@ +array-api-compat numpy>=2.0.0 scipy>=1.13.0 scikit-learn>=1.5.0 diff --git a/pyproject.toml b/pyproject.toml index d566d9c..e091a9e 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -28,7 +28,8 @@ classifiers = [ dependencies = ["numpy>=2.0.0", "scikit-learn>=1.5.0", "scipy>=1.13.0", - "packaging>=24.0"] + "packaging>=24.0", + "array-api-compat"] [project.urls] source = "https://github.com/lanl/GFDL" diff --git a/src/gfdl/model.py b/src/gfdl/model.py index 56dd28a..9bae1dc 100644 --- a/src/gfdl/model.py +++ b/src/gfdl/model.py @@ -18,6 +18,7 @@ from sklearn.utils.metaestimators import available_if from sklearn.utils.multiclass import check_classification_targets, unique_labels from sklearn.utils.validation import check_is_fitted, validate_data +import array_api_compat from gfdl.activations import resolve_activation from gfdl.weights import resolve_weight @@ -47,9 +48,10 @@ def fit(self, X, Y): # Assumption : X, Y have been pre-processed. # X shape: (n_samples, n_features) # Y shape: (n_samples, n_classes-1) + xp = array_api_compat.array_namespace(X, Y) if self.reg_alpha is not None and self.reg_alpha < 0.0: raise ValueError("Negative reg_alpha. Expected range : None or [0.0, inf).") - hidden_layer_sizes = np.asarray(self.hidden_layer_sizes) + hidden_layer_sizes = xp.asarray(self.hidden_layer_sizes) if hidden_layer_sizes.min() < 1: raise ValueError("hidden_layer_sizes must be > 0, " f"got {hidden_layer_sizes}") @@ -95,7 +97,7 @@ def fit(self, X, Y): # or (n_samples, sum_hidden) if self.direct_links: Hs.append(X) - D = np.hstack(Hs) + D = xp.concat(Hs, axis=1) # beta shape: (sum_hidden+n_features, n_classes-1) # or (sum_hidden, n_classes-1) @@ -103,7 +105,7 @@ def fit(self, X, Y): # If reg_alpha is None, use direct solve using # MoorePenrose Pseudo-Inverse, otherwise use ridge regularized form. if self.reg_alpha is None: - self.coeff_ = np.linalg.pinv(D, rtol=self.rtol) @ Y + self.coeff_ = xp.linalg.pinv(D, rtol=self.rtol) @ Y else: ridge = Ridge(alpha=self.reg_alpha, fit_intercept=False) ridge.fit(D, Y) @@ -410,6 +412,7 @@ def fit(self, X, y): object Fitted estimator. """ + xp = array_api_compat.array_namespace(X, y) # shape: (n_samples, n_features) X, Y = validate_data(self, X, y) self.classes_ = unique_labels(Y) @@ -418,7 +421,8 @@ def fit(self, X, y): # (this is necessary for everything beyond binary classification) self.enc_ = OneHotEncoder(handle_unknown="ignore", sparse_output=False) # shape: (n_samples, n_classes-1) - Y = self.enc_.fit_transform(Y.reshape(-1, 1)) + Y = self.enc_.fit_transform(np.from_dlpack(Y).reshape(-1, 1)) + Y = xp.asarray(Y) # call base fit method super().fit(X, Y) From 34d44817521814c0335676c48a383599ab8ee900 Mon Sep 17 00:00:00 2001 From: Tyler Reddy Date: Wed, 24 Jun 2026 16:34:25 -0600 Subject: [PATCH 2/2] WIP, ENH, MAINT: PR 115 revisions * GPU/device-related shims for classifier `fit()` in above PR. Seems to allow basic fitting on GPU with torch but is probably still hacky/inefficient, and of course there are still not tests for xp/device support. * A few ruff-related fixes while in the neighborhood. --- src/gfdl/model.py | 10 +++++++--- 1 file changed, 7 insertions(+), 3 deletions(-) diff --git a/src/gfdl/model.py b/src/gfdl/model.py index 9bae1dc..d7faea9 100644 --- a/src/gfdl/model.py +++ b/src/gfdl/model.py @@ -2,6 +2,7 @@ Estimators for gradient free deep learning. """ +import array_api_compat import numpy as np import scipy from scipy.special import logsumexp @@ -18,7 +19,6 @@ from sklearn.utils.metaestimators import available_if from sklearn.utils.multiclass import check_classification_targets, unique_labels from sklearn.utils.validation import check_is_fitted, validate_data -import array_api_compat from gfdl.activations import resolve_activation from gfdl.weights import resolve_weight @@ -89,6 +89,8 @@ def fit(self, X, Y): Hs = [] H_prev = X for w, b in zip(self.W_, self.b_, strict=False): + w = xp.asarray(w, device=H_prev.device) + b = xp.asarray(b, device=H_prev.device) Z = H_prev @ w.T + b # (n_samples, n_hidden) H_prev = self._activation_fn(Z) Hs.append(H_prev) @@ -413,6 +415,7 @@ def fit(self, X, y): Fitted estimator. """ xp = array_api_compat.array_namespace(X, y) + y_device = y.device # shape: (n_samples, n_features) X, Y = validate_data(self, X, y) self.classes_ = unique_labels(Y) @@ -421,8 +424,9 @@ def fit(self, X, y): # (this is necessary for everything beyond binary classification) self.enc_ = OneHotEncoder(handle_unknown="ignore", sparse_output=False) # shape: (n_samples, n_classes-1) - Y = self.enc_.fit_transform(np.from_dlpack(Y).reshape(-1, 1)) - Y = xp.asarray(Y) + Y = self.enc_.fit_transform(np.from_dlpack(xp.asarray(Y, + device="cpu")).reshape(-1, 1)) + Y = xp.asarray(Y, device=y_device) # call base fit method super().fit(X, Y)