fixing and adding validation

Author: Jimmy-INL

pysensors/optimizers/_gqr.py

@@ -16,10 +16,10 @@ class GQR(QR):
     """
     General QR optimizer for sensor selection.
     Ranks sensors in descending order of "importance" based on
-    reconstruction performance. This is an extension that requires a more intrusive
+    reconstruction accuracy. This is an extension that requires a more intrusive
     access to the QR optimizer to facilitate a more adaptive optimization. This is a generalized version of cost constraints
-    in the sense that users can allow n constrained sensors in the constrained area.
-    if n = 0 this converges to the CCQR results. If no constraints it converges to QR results.
+    in the sense that users can allow `n_const_sensors` in the constrained area.
+    if n = 0 this converges to the CCQR results. and if no constrained region it should converge to the results from QR optimizer.
 
     See the following reference for more information
         Manohar, Krithika, et al.
@@ -46,8 +46,6 @@ def __init__(self):
         constraint_option : string,
             max_n_const_sensors : The number of sensors in the constrained region should be less than or equal to n_const_sensors.
             exact_n_const_sensors : The number of sensors in the constrained region should be exactly equal to n_const_sensors.
-        nx, ny : integer,
-            X, Y dimensions of the grid.
         """
         self.pivots_ = None
         self.idx_constrained = []
@@ -59,7 +57,7 @@ def __init__(self):
         self.ny = None
         self.r = 1
 
-    def fit(self,basis_matrix=None,**optimizer_kws):
+    def fit(self,basis_matrix,**optimizer_kws):
         """
         Parameters
         ----------
@@ -78,6 +76,13 @@ def fit(self,basis_matrix=None,**optimizer_kws):
         n_features, n_samples = basis_matrix.shape  # We transpose basis_matrix below
         max_const_sensors = len(self.idx_constrained) # Maximum number of sensors allowed in the constrained region
 
+        ## Assertions and checks:
+        # if self.n_sensors > n_features - max_const_sensors + self.nConstrainedSensors:
+        #     raise IOError ("n_sensors cannot be larger than n_features - all possible locations in the constrained area + allowed constrained sensors")
+        # if self.n_sensors > n_samples + self.nConstrainedSensors: ## Handling zero constraint?
+        #     raise IOError ("Currently n_sensors should be less than min(number of samples, number of modes) + number of constrained sensors,\
+        #                    got: n_sensors = {}, n_samples + const_sensors = {} + {} = {}".format(self.n_sensors,n_samples,self.nConstrainedSensors,n_samples+self.nConstrainedSensors))
+
         # Initialize helper variables
         R = basis_matrix.conj().T.copy()
         p = np.arange(n_features)

pysensors/utils/__init__.py

@@ -6,6 +6,8 @@
 from ._norm_calc import exact_n
 from ._norm_calc import max_n
 from ._norm_calc import predetermined
+from ._validation import determinant
+from ._validation import relative_reconstruction_error
 
 __all__ = [
     "constrained_binary_solve",
@@ -17,5 +19,7 @@
     "functional_constraints",
     "exact_n",
     "max_n",
-    "predetermined"
+    "predetermined",
+    "determinant",
+    "relative_reconstruction_error"
 ]

pysensors/utils/_norm_calc.py

@@ -30,9 +30,23 @@ def exact_n(lin_idx, dlens, piv, j, n_const_sensors, **kwargs): ##Will first for
     -------
     dlens : np.darray, shape [Variable based on j] with constraints mapped into it.
     """
-    didx = np.isin(piv[j:],lin_idx,invert=j<n_const_sensors)
-    dlens[didx] = 0
-    return dlens
+    if 'all_sensors' in kwargs.keys():
+        all_sensors = kwargs['all_sensors']
+    else:
+        all_sensors = []
+    if 'n_sensors' in kwargs.keys():
+        n_sensors = kwargs['n_sensors']
+    else:
+        n_sensors = len(all_sensors)
+    for i in range(n_sensors):
+        if np.isin(all_sensors[:n_sensors],lin_idx,invert=False).sum() < n_const_sensors:
+            if n_sensors >= j > (n_sensors - (n_const_sensors-1)):
+                didx = np.isin(piv[j:],lin_idx,invert=True)
+                dlens[didx] = 0
+        else:
+            max_n(lin_idx, dlens, piv, j, n_const_sensors, **kwargs)
+    return(dlens)
+
 
 def max_n(lin_idx, dlens, piv, j, n_const_sensors, **kwargs):
     """
@@ -109,9 +123,9 @@ def predetermined(lin_idx, dlens, piv, j, n_const_sensors, **kwargs):
 
 __norm_calc_type = {}
 __norm_calc_type[''] = unconstrained
-__norm_calc_type['exact_n_const_sensors'] = exact_n
-__norm_calc_type['max_n_const_sensors'] = max_n
-__norm_calc_type['predetermined_norm_calc'] = predetermined
+__norm_calc_type['exact_n'] = exact_n
+__norm_calc_type['max_n'] = max_n
+__norm_calc_type['predetermined'] = predetermined
 
 def returnInstance(cls, name):
   """

pysensors/utils/_validation.py

@@ -0,0 +1,53 @@
+"""
+Various utility functions for validation and computing reconstruction scores and errors.
+"""
+import numpy as np
+from scipy.sparse import csr_matrix
+
+def determinant(top_sensors, n_features, basis_matrix):
+    """
+    Function for calculating |C.T phi.T C phi|.
+
+    Parameters
+    ----------
+    top_sensors: np.darray,
+        Column indices of choosen sensor locations
+    n_features : int,
+        No. of features of dataset
+    basis_matrix : np.darray,
+        The basis matrix calculated by model.basis_matrix_
+    Returns
+    -------
+    optimality : Float,
+        The dterminant value obtained.
+    """
+
+    p = len(top_sensors) # Number of sensors
+    n,r = np.shape(basis_matrix) # state dimension X Number of modes
+    c = csr_matrix((p,n),dtype=np.int8)
+
+    for i in range(p):
+        c[i,top_sensors[i]] = 1
+    phi = basis_matrix
+    # optimality = np.linalg.det(( c @ phi).T @ (c@phi)) #np.log(np.linalg.det(phi.T @ c.T)) np.log(np.linalg.det((c@phi).T @ (c@phi)))
+    optimality = abs(np.linalg.det(c @ phi)) if p==r else abs(np.linalg.det(( c @ phi).T @ (c @ phi)))
+    # optimality = abs(np.linalg.det(c @ phi))
+    return optimality
+
+def relative_reconstruction_error(data, prediction):
+    """
+    Function for calculating relative error between actual data and the reconstruction
+
+    Parameters
+        ----------
+        data: np.darray,
+            The actual data from the dataset evaluated
+        prediction : np.darray,
+            The predicted values from model.predict(X[:,top_sensors])
+        Returns
+        -------
+        error_val : Float,
+            The relative error calculated.
+    """
+    error_val = (np.linalg.norm((data - prediction)/np.linalg.norm(data)))*100
+    return (error_val)