init

LICENSE.txt

@@ -0,0 +1,16 @@
+# =======
+# License BSD-3
+# =======
+# https://opensource.org/licenses/BSD-3-Clause
+#
+# Copyright 2018 Josh L. Espinoza
+#
+# Redistribution and use in source and binary forms, with or without modification, are permitted provided that the following conditions are met:
+#
+# 1. Redistributions of source code must retain the above copyright notice, this list of conditions and the following disclaimer.
+#
+# 2. Redistributions in binary form must reproduce the above copyright notice, this list of conditions and the following disclaimer in the documentation and/or other materials provided with the distribution.
+#
+# 3. Neither the name of the copyright holder nor the names of its contributors may be used to endorse or promote products derived from this software without specific prior written permission.
+#
+# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

MANIFEST.in

@@ -0,0 +1,10 @@
+graft clairvoyance
+include setup.py
+include LICENSE.txt
+include README.md
+global-exclude *.py[cod]
+global-exclude Icon*
+global-exclude .DS_Store
+exclude *.egg-info
+exclude dist
+exclude conda_builds

README.md

@@ -0,0 +1,243 @@
+
+```
+
+ _______        _______ _____  ______ _    _  _____  __   __ _______ __   _ _______ _______
+ |       |      |_____|   |   |_____/  \  /  |     |   \_/   |_____| | \  | |       |______
+ |_____  |_____ |     | __|__ |    \_   \/   |_____|    |    |     | |  \_| |_____  |______
+                                                                                           
+                                        
+
+```
+#### Description
+
+Reimplementation for `Clairvoyance` from [Espinoza & Dupont et al. 2021](https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1008857).  The updated version includes regression support, support for all linear/tree-based models, and improved visualizations.  `Clairvoyance` is currently in active development.   
+
+
+#### Details:
+`import clairvoyance as cy`
+
+`__version__ = "2022.12.27"`
+
+#### Installation
+
+```
+pip install clairvoyance_feature_selection
+conda install -c jolespin clairvoyance
+```
+
+#### Citation
+
+Espinoza JL, Dupont CL, O’Rourke A, Beyhan S, Morales P, Spoering A, et al. (2021) Predicting antimicrobial mechanism-of-action from transcriptomes: A generalizable explainable artificial intelligence approach. PLoS Comput Biol 17(3): e1008857. https://doi.org/10.1371/journal.pcbi.1008857
+
+#### Usage
+
+
+
+##### Feature selection based on classification tasks
+Here's a basic classifcation using a `LogisticRegression` model and a grid search for different `C` and `penalty` parameters. We add 996 noise variables within the range of values as the original Iris features. After that we normalize them so their scale is standardized.  In this case, we are optimizing for `accuracy`.
+
+```python
+import clairvoyance as cy
+import numpy as np
+import pandas as pd
+from sklearn.datasets import load_iris
+from sklearn.linear_model import LogisticRegression
+
+# Load iris dataset
+X, y = load_iris(return_X_y=True, as_frame=True)
+X.columns = X.columns.map(lambda j: j.split(" (cm")[0].replace(" ","_"))
+
+# Relabel targets
+target_names = load_iris().target_names
+y = y.map(lambda i: target_names[i])
+
+# Add 996 noise features (total = 1000 features) in the same range of values as the original features
+number_of_noise_features = 996
+vmin = X.values.ravel().min()
+vmax = X.values.ravel().max()
+X_noise = pd.DataFrame(
+    data=np.random.RandomState(0).randint(low=int(vmin*10), high=int(vmax*10), size=(150, number_of_noise_features))/10,
+    columns=map(lambda j:"noise_{}".format(j+1), range(number_of_noise_features)),
+)
+
+X_iris_with_noise = pd.concat([X, X_noise], axis=1)
+X_normalized = X_iris_with_noise - X_iris_with_noise.mean(axis=0).values
+X_normalized = X_normalized/X_normalized.std(axis=0).values
+
+# Specify model algorithm and parameter grid
+estimator=LogisticRegression(max_iter=1000, solver="liblinear", multi_class="ovr")
+param_grid={
+    "C":[1e-10] + (np.arange(1,11)/10).tolist(),
+    "penalty":["l1", "l2"],
+}
+
+# Instantiate model
+clf = cy.ClairvoyanceClassification(
+    n_jobs=-1, 
+    scorer="accuracy", 
+    n_draws=10, 
+    estimator=estimator, 
+    param_grid=param_grid, 
+    verbose=1,
+)
+clf.fit(X_normalized, y)#, sort_hyperparameters_by=["C", "penalty"], ascending=[True, False])
+history = clf.recursive_feature_inclusion(early_stopping=10)
+history.head()
+
+```
+![](images/1.png)
+
+```python
+clf.plot_scores(title="Iris", xtick_rotation=90)
+clf.plot_weights()
+clf.plot_weights(weight_type="cross_validation")
+```
+![](images/2.png)
+
+There are still a few noise variables, though with much lower weight, suggesting our classifier is modeling noise.  We can add an additional penalty where a change in score must exceed a threshold to add a new feature during the recursive feature inclusion algorithm. 
+
+```python
+history = clf.recursive_feature_inclusion(early_stopping=10, minimum_improvement_in_score=0.05)
+clf.plot_scores(title="Iris", xtick_rotation=90)
+clf.plot_weights()
+clf.plot_weights(weight_type="cross_validation")
+```
+![](images/3.png)
+
+Now let's do a binary classification but optimize `fbeta` score instead of `accuracy`.  Instead of a fixed penalty, we are going to use a custom penalty that scales with the number of features included. 
+
+```python
+from sklearn.metrics import fbeta_score
+
+# Let's do a binary classification
+y_notsetosa = y.map(lambda x: {True:"not_setosa", False:x}[x != "setosa"]) 
+
+# Let's also use a FBeta scorer
+scorer = make_scorer(fbeta_score, average="binary", beta=0.5, pos_label="setosa")
+
+# Instantiate model
+clf_binary = cy.ClairvoyanceClassification(
+    n_jobs=-1, 
+    scorer="accuracy", 
+    n_draws=10, 
+    estimator=estimator, 
+    param_grid=param_grid, 
+    verbose=1,
+)
+
+# Let's also prefer lower C values and l1 over l2 (i.e., stronger regularization and sparsity)
+clf_binary.fit(X_normalized, y, sort_hyperparameters_by=["C", "penalty"], ascending=[True, True]) 
+
+# Instead of adding a fixed penalty for adding new features, let's add a function that scales with the number of features
+history = clf_binary.recursive_feature_inclusion(early_stopping=10, additional_feature_penalty=lambda n: 1e-3*n**2)
+history.head()
+```
+![](images/4.png)
+
+```python
+clf_binary.plot_scores(title="Iris (Binary)", xtick_rotation=90)
+clf_binary.plot_weights()
+clf_binary.plot_weights(weight_type="cross_validation")
+```
+![](images/5.png)
+
+##### Feature selection based on regression tasks
+Here's a basic regression using a `DecisionTreeRegressor` model and a grid search for different `min_samples_leaf` and `min_samples_split` parameters. We add 87 noise variables and normalize all of the features so their scale is standardized.  In this case, we are optimizing for `neg_root_mean_squared_error`.  We are using a validation set of ~16% of the data during our recursive feature inclusion. 
+
+```python
+from sklearn.datasets import load_boston
+from sklearn.tree import DecisionTreeRegressor
+from sklearn.model_selection import train_test_split
+
+# Load Boston data
+boston = load_boston()
+X = pd.DataFrame(boston.data, columns=boston.feature_names)
+y = pd.Series(boston.target)
+
+number_of_noise_features = 100 - X.shape[1]
+X_noise = pd.DataFrame(np.random.RandomState(0).normal(size=(X.shape[0], number_of_noise_features)),  columns=map(lambda j: f"noise_{j}", range(number_of_noise_features)))
+X_boston_with_noise = pd.concat([X, X_noise], axis=1)
+X_normalized = X_boston_with_noise - X_boston_with_noise.mean(axis=0).values
+X_normalized = X_normalized/X_normalized.std(axis=0).values
+
+# Let's fit the model but leave a held out validation set
+X_training, X_validation, y_training, y_validation = train_test_split(X_normalized, y, random_state=0, test_size=0.1618)
+
+# Get parameters
+estimator = DecisionTreeRegressor(random_state=0)
+param_grid = {"min_samples_leaf":[1,2,3,5,8],"min_samples_split":{ 0.1618, 0.382, 0.5, 0.618}}
+
+# Fit model
+reg = cy.ClairvoyanceRegression(name="Boston", n_jobs=-1, n_draws=10, estimator=estimator, param_grid=param_grid, verbose=1)
+reg.fit(X_training, y_training)
+history = reg.recursive_feature_inclusion(early_stopping=10, X=X_validation, y=y_validation)
+history.head()
+```
+![](images/6.png)
+
+```python
+reg.plot_scores(title="Boston", xtick_rotation=90)
+reg.plot_weights()
+reg.plot_weights(weight_type="cross_validation")
+```
+![](images/7.png)
+
+Let's use the weighted fitted with a `DecisionTreeRegressor` but use an ensemble `RandomForestRegressor` for the actual feature inclusion algorithm. 
+
+```python
+from sklearn.ensemble import RandomForestRegressor
+history = reg.recursive_feature_inclusion(early_stopping=10, estimator=RandomForestRegressor(random_state=0), X=X_validation, y=y_validation)
+reg.plot_scores(title="Boston", xtick_rotation=90)
+reg.plot_weights()
+reg.plot_weights(weight_type="cross_validation")
+```
+![](images/8.png)
+
+##### Recursive feature selection based on classification tasks
+Here we are running `Clairvoyance` recursively identifying several feature sets that work with different hyperparameters to get a range of feature sets to select from in the end.  We will iterate through all of the hyperparamater configurations, recursively feed in the data using different percentiles of the weights, and use different score thresholds from the random draws.  The recursive usage is similar to the legacy implementation used in [Espinoza & Dupont et al. 2021](https://journals.plos.org/ploscompbiol/article?id=10.1371/journal.pcbi.1008857) (which is still provided as an executable).
+
+```python
+# Get the iris data again
+X_normalized = X_iris_with_noise - X_iris_with_noise.mean(axis=0).values
+X_normalized = X_normalized/X_normalized.std(axis=0).values
+target_names = load_iris().target_names
+y = pd.Series(load_iris().target)
+y = y.map(lambda i: target_names[i])
+
+# Specify model algorithm and parameter grid
+estimator=LogisticRegression(max_iter=1000, solver="liblinear", multi_class="ovr")
+param_grid={
+    "C":[1e-10] + (np.arange(1,11)/10).tolist(),
+    "penalty":["l1", "l2"],
+}
+
+# Instantiate model
+rci = cy.ClairvoyanceRecursive(
+    n_jobs=-1, 
+    scorer="accuracy", 
+    n_draws=10, 
+    estimator=estimator, 
+    param_grid=param_grid, 
+    percentiles=[0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.925, 0.95, 0.975, 0.99],
+    minimum_scores=[-np.inf, 0.382, 0.5],
+    verbose=0,
+)
+rci.fit(X_normalized, y, sort_hyperparameters_by=["C", "penalty"], ascending=[True, True])
+rci.plot_recursive_feature_selection()
+
+```
+
+We observe a nice separate around 10 features, so let's use that as a maximum. 
+
+![](images/9.png)
+
+```python
+# Plot the features with a maximum of 10
+rci.plot_recursive_feature_selection(max_features=10)
+
+# Filter out all the results that have more than 10 features
+rci.results_.query("number_of_features <= 10").sort_values("score", ascending=False).head()
+```
+
+![](images/10.png)
+