Merge pull request #7 from mckenziephagen/pr_model

Pr model

Author: arokem

.gitignore

@@ -0,0 +1,3 @@
+*ipynb* 
+*cache*
+*egg*  

src/model_fc/models.py

@@ -14,20 +14,20 @@
 def run_model(train_ts, test_ts, n_rois, model, **kwargs):
     """Calculate a model based functional connectivity matrix.
 
-
     train_ts: training timeseries
     test_ts: testing timeseries
     n_rois: number of rois in parcellation
     model: model object
 
-
     """
     assert train_ts.shape[1] == n_rois == test_ts.shape[1]
-    fc_mat = np.zeros((n_rois, n_rois))
+    fc_mat = np.empty((n_rois, n_rois))
 
-    inner_rsq_dict = {"train": [], "test": []}
+    results_dict = {}
 
     for target_idx in range(train_ts.shape[1]):
+        results_dict[f"node_{target_idx}"] = {}
+        print(f"*****ECHO {target_idx}***********")
         y_train = np.array(train_ts[:, target_idx])
         X_train = np.delete(train_ts, target_idx, axis=1)
 
@@ -36,13 +36,15 @@ def run_model(train_ts, test_ts, n_rois, model, **kwargs):
 
         model.fit(X=X_train, y=y_train)
 
-        fc_mat[target_idx, :] = np.insert(model.coef_, target_idx, 0)
+        fc_mat[target_idx, :] = np.insert(model.coef_, target_idx, 1)
         test_rsq, train_rsq = eval_metrics(X_train, y_train, X_test, y_test, model)
 
-        inner_rsq_dict["test"].append(test_rsq)
-        inner_rsq_dict["train"].append(train_rsq)
+        results_dict[f"node_{target_idx}"]["train_r2"] = train_rsq
+        results_dict[f"node_{target_idx}"]["test_r2"] = test_rsq
+
+    results_dict["fc_matrix"] = fc_mat
 
-    return (fc_mat, inner_rsq_dict, model)
+    return results_dict
 
 
 def eval_metrics(X_train, y_train, X_test, y_test, model):
@@ -58,25 +60,19 @@ def init_model(
     model_str, max_iter, random_state, stability_selection=16, selection_frac=0.7
 ):
     """Initialize model object for FC calculations."""
-    if model_str == "uoi-lasso":
+    if model_str == "uoiLasso":
         uoi_lasso = UoI_Lasso(estimation_score="BIC")
-        comm = MPI.COMM_WORLD
-
         uoi_lasso.selection_frac = selection_frac
         uoi_lasso.stability_selection = stability_selection
         uoi_lasso.copy_X = True
         uoi_lasso.estimation_target = None
         uoi_lasso.logger = None
-        uoi_lasso.warm_start = False
-        uoi_lasso.comm = comm
-        uoi_lasso.random_state = 1
+        uoi_lasso.warm_start = True
         uoi_lasso.n_lambdas = 100
         uoi_lasso.max_iter = max_iter
-        uoi_lasso.random_state = random_state
-
         model = uoi_lasso
 
-    elif model_str == "lasso-cv":
+    elif model_str == "lassoCV":
         lasso = LassoCV(
             fit_intercept=True,
             cv=5,
@@ -86,8 +82,12 @@ def init_model(
         )
 
         model = lasso
+    elif model_str == "ridgeCV":
+        ridge = RidgeCV(fit_intercept=True)
+
+        model = ridge
 
-    elif model_str == "lasso-bic":
+    elif model_str == "lassoBIC":
         lasso = LassoLarsIC(criterion="bic", fit_intercept=True, max_iter=max_iter)
 
         model = lasso
@@ -119,10 +119,10 @@ def fit(self, X, y):
 
         # Calculate the Pearson coefficient between y and every column of x
         corrmatrix = np.corrcoef(np.concatenate([X, y[:, None]], -1).T)
-        self.coeff_ = corrmatrix[-1, :-1]
+        self.coef_ = corrmatrix[-1, :-1]
         # Calculate the linear combination of columns of X with these
         # coefficients
-        naive_pred_y = X @ self.coeff_
+        naive_pred_y = X @ self.coef_
         # Calculate the scale parameter to match the variance of y
         self.scale_ = np.mean(y / naive_pred_y)
         self.is_fitted_ = True
@@ -132,4 +132,4 @@ def fit(self, X, y):
     def predict(self, X):
         check_is_fitted(self)
         X = validate_data(self, X, reset=False)
-        return (X @ self.coeff_) / self.scale_
+        return X @ self.coef_ * self.scale_

test/test_models.py

@@ -1,3 +1,5 @@
+import math
+
 import numpy as np
 from model_fc.models import PearsonRegressor
 from sklearn.utils.estimator_checks import parametrize_with_checks
@@ -16,3 +18,16 @@ def test_PearsonRegressor():
     pr.fit(X, y)
     pred = pr.predict(X)
     assert pred.shape == y.shape
+
+
+def test_PearsonRegressorScale():
+    y = np.random.randn(100)
+    X = np.zeros((100, 10))
+    for col in range(X.shape[1]):
+        # no noise so that our scaling == 1 exactly
+        X[:, col] = y
+
+    pr = PearsonRegressor()
+    pr.fit(X, y)
+    # assert scale is 1/ncols
+    assert math.isclose(pr.scale_, 0.1)