#21 add itemised trace to hERG QC

Author: kwabenantim

pcpostprocess/hergQC.py

@@ -123,38 +123,57 @@ def run_qc(self, voltage_steps, times,
         before = self.filter_capacitive_spikes(before, times, voltage_steps)
         after = self.filter_capacitive_spikes(after, times, voltage_steps)
 
+        QC = {label: [(False, None)] for label in self.qc_labels}
+
         if len(before) == 0 or len(after) == 0:
-            return False, [False for lab in self.qc_labels]
+            return QC
 
         if (None in qc_vals_before) or (None in qc_vals_after):
-            return False, False * self.no_QC
+            return QC
 
         qc1_1 = self.qc1(*qc_vals_before)
         qc1_2 = self.qc1(*qc_vals_after)
-        qc1 = [i and j for i, j in zip(qc1_1, qc1_2)]
 
-        qc2_1 = True
-        qc2_2 = True
+        QC['qc1.rseal'] = [qc1_1[0], qc1_2[0]]
+        QC['qc1.cm'] = [qc1_1[1], qc1_2[1]]
+        QC['qc1.rseries'] = [qc1_1[2], qc1_2[2]]
+
+        QC['qc2.raw'] = []
+        QC['qc2.subtracted'] = []
         for i in range(n_sweeps):
-            qc2_1 = qc2_1 and self.qc2(before[i])
-            qc2_2 = qc2_2 and self.qc2(before[i] - after[i])
+            qc2_1 = self.qc2(before[i])
+            qc2_2 = self.qc2(before[i] - after[i])
+
+            QC['qc2.raw'].append(qc2_1)
+            QC['qc2.subtracted'].append(qc2_2)
 
         qc3_1 = self.qc3(before[0, :], before[1, :])
         qc3_2 = self.qc3(after[0, :], after[1, :])
         qc3_3 = self.qc3(before[0, :] - after[0, :],
                          before[1, :] - after[1, :])
 
+        QC['qc3.raw'] = [qc3_1]
+        QC['qc3.E4031'] = [qc3_2]
+        QC['qc3.subtracted'] = [qc3_3]
+
         rseals = [qc_vals_before[0], qc_vals_after[0]]
         cms = [qc_vals_before[1], qc_vals_after[1]]
         rseriess = [qc_vals_before[2], qc_vals_after[2]]
         qc4 = self.qc4(rseals, cms, rseriess)
 
+        QC['qc4.rseal'] = [qc4[0]]
+        QC['qc4.cm'] = [qc4[1]]
+        QC['qc4.rseries'] = [qc4[2]]
+
         # indices where hERG peaks
         qc5 = self.qc5(before[0, :], after[0, :],
                        self.qc5_win)
 
         qc5_1 = self.qc5_1(before[0, :], after[0, :], label='1')
 
+        QC['qc5.staircase'] = [qc5]
+        QC['qc5.1.staircase'] = [qc5_1]
+
         # Ensure thats the windows are correct by checking the voltage trace
         assert np.all(
             np.abs(self.voltage[self.qc6_win[0]: self.qc6_win[1]] - 40.0))\
@@ -166,14 +185,17 @@ def run_qc(self, voltage_steps, times,
             np.abs(self.voltage[self.qc6_2_win[0]: self.qc6_2_win[1]] - 40.0))\
             < 1e-8
 
-        qc6, qc6_1, qc6_2 = True, True, True
+        QC['qc6.subtracted'] = []
+        QC['qc6.1.subtracted'] = []
+        QC['qc6.2.subtracted'] = []
         for i in range(before.shape[0]):
-            qc6 = qc6 and self.qc6((before[i, :] - after[i, :]),
-                                   self.qc6_win, label='0')
-            qc6_1 = qc6_1 and self.qc6((before[i, :] - after[i, :]),
-                                       self.qc6_1_win, label='1')
-            qc6_2 = qc6_2 and self.qc6((before[i, :] - after[i, :]),
-                                       self.qc6_2_win, label='2')
+            qc6 = self.qc6((before[i, :] - after[i, :]), self.qc6_win, label="0")
+            qc6_1 = self.qc6((before[i, :] - after[i, :]), self.qc6_1_win, label="1")
+            qc6_2 = self.qc6((before[i, :] - after[i, :]), self.qc6_2_win, label="2")
+
+            QC['qc6.subtracted'].append(qc6)
+            QC['qc6.1.subtracted'].append(qc6_1)
+            QC['qc6.2.subtracted'].append(qc6_2)
 
         if self._debug:
             fig = plt.figure(figsize=(8, 5))
@@ -199,37 +221,47 @@ def run_qc(self, voltage_steps, times,
             fig.savefig(os.path.join(self.plot_dir, 'qc_debug.png'))
             plt.close(fig)
 
-        # Make a flat list of all QC criteria (pass/fail bool)
-        QC = np.hstack([qc1, [qc2_1, qc2_2, qc3_1, qc3_2, qc3_3],
-                        qc4, [qc5, qc5_1, qc6, qc6_1, qc6_2]]).flatten()
+        # Check if all QC criteria passed
+        passed = all([x for qc in QC.values() for x, _ in qc])
 
-        passed = np.all(QC)
         return passed, QC
 
     def qc1(self, rseal, cm, rseries):
-
-        if any([x is None for x in (rseal, cm, rseries)]):
-            return False, False, False
-
         # Check R_seal, C_m, R_series within desired range
-        if rseal < self.rsealc[0] or rseal > self.rsealc[1] \
-                or not np.isfinite(rseal):
+        if (
+            rseal is None
+            or rseal < self.rsealc[0]
+            or rseal > self.rsealc[1]
+            or not np.isfinite(rseal)
+        ):
             self.logger.debug(f"rseal: {rseal}")
             qc11 = False
         else:
             qc11 = True
-        if cm < self.cmc[0] or cm > self.cmc[1] or not np.isfinite(cm):
+
+        if (
+            cm is None
+            or cm < self.cmc[0]
+            or cm > self.cmc[1]
+            or not np.isfinite(cm)
+        ):
             self.logger.debug(f"cm: {cm}")
             qc12 = False
         else:
             qc12 = True
-        if rseries < self.rseriesc[0] or rseries > self.rseriesc[1] \
-                or not np.isfinite(rseries):
+
+        if (
+            rseries is None
+            or rseries < self.rseriesc[0]
+            or rseries > self.rseriesc[1]
+            or not np.isfinite(rseries)
+        ):
             self.logger.debug(f"rseries: {rseries}")
             qc13 = False
         else:
             qc13 = True
-        return [qc11, qc12, qc13]
+
+        return [(qc11, rseal), (qc12, cm), (qc13, rseries)]
 
     def qc2(self, recording, method=3):
         # Check SNR is good
@@ -245,9 +277,11 @@ def qc2(self, recording, method=3):
 
         if snr < self.snrc or not np.isfinite(snr):
             self.logger.debug(f"snr: {snr}")
-            return False
+            result = False
+        else:
+            result = True
 
-        return True
+        return (result, snr)
 
     def qc3(self, recording1, recording2, method=3):
         # Check 2 sweeps similar
@@ -270,36 +304,51 @@ def qc3(self, recording1, recording2, method=3):
         rmsd = np.sqrt(np.mean((recording1 - recording2) ** 2))
         if rmsd > rmsdc or not (np.isfinite(rmsd) and np.isfinite(rmsdc)):
             self.logger.debug(f"rmsd: {rmsd}, rmsdc: {rmsdc}")
-            return False
-        return True
-
-    def qc4(self, rseals, cms, rseriess):
+            result = False
+        else:
+            result = True
 
-        if any([x is None for x in list(rseals) + list(cms) + list(rseriess)]):
-            return [False, False, False]
+        return (result, rmsd)
 
+    def qc4(self, rseals, cms, rseriess):
         # Check R_seal, C_m, R_series stability
         # Require std/mean < x%
         qc41 = True
         qc42 = True
         qc43 = True
-        if np.std(rseals) / np.mean(rseals) > self.rsealsc or not (
-                np.isfinite(np.mean(rseals)) and np.isfinite(np.std(rseals))):
-            self.logger.debug(f"d_rseal:  {np.std(rseals) / np.mean(rseals)}")
+
+        if None in list(rseals):
             qc41 = False
+            d_rseal = None
+        else:
+            d_rseal = np.std(rseals) / np.mean(rseals)
+            if d_rseal > self.rsealsc or not (
+                    np.isfinite(np.mean(rseals)) and np.isfinite(np.std(rseals))):
+                self.logger.debug(f"d_rseal:  {d_rseal}")
+                qc41 = False
 
-        if np.std(cms) / np.mean(cms) > self.cmsc or not (
-                np.isfinite(np.mean(cms)) and np.isfinite(np.std(cms))):
-            self.logger.debug(f"d_cm: {np.std(cms) / np.mean(cms)}")
+        if None in list(cms):
             qc42 = False
+            d_cm = None
+        else:
+            d_cm = np.std(cms) / np.mean(cms)
+            if d_cm > self.cmsc or not (
+                    np.isfinite(np.mean(cms)) and np.isfinite(np.std(cms))):
+                self.logger.debug(f"d_cm: {d_cm}")
+                qc42 = False
 
-        if np.std(rseriess) / np.mean(rseriess) > self.rseriessc or not (
-                np.isfinite(np.mean(rseriess))
-                and np.isfinite(np.std(rseriess))):
-            self.logger.debug(f"d_rseries: {np.std(rseriess) / np.mean(rseriess)}")
+        if None in list(rseriess):
             qc43 = False
+            d_rseries = None
+        else:
+            d_rseries = np.std(rseriess) / np.mean(rseriess)
+            if d_rseries > self.rseriessc or not (
+                    np.isfinite(np.mean(rseriess))
+                    and np.isfinite(np.std(rseriess))):
+                self.logger.debug(f"d_rseries: {d_rseries}")
+                qc43 = False
 
-        return [qc41, qc42, qc43]
+        return [(qc41, d_rseal), (qc42, d_cm), (qc43, d_rseries)]
 
     def qc5(self, recording1, recording2, win=None, label=''):
         # Check pharma peak value drops after E-4031 application
@@ -325,8 +374,11 @@ def qc5(self, recording1, recording2, win=None, label=''):
         if (max_diff < max_diffc) or not (np.isfinite(max_diff)
                                           and np.isfinite(max_diffc)):
             self.logger.debug(f"max_diff:  {max_diff}, max_diffc: {max_diffc}")
-            return False
-        return True
+            result = False
+        else:
+            result = True
+
+        return (result, max_diff)
 
     def qc5_1(self, recording1, recording2, win=None, label=''):
         # Check RMSD_0 drops after E-4031 application
@@ -355,8 +407,11 @@ def qc5_1(self, recording1, recording2, win=None, label=''):
         if (rmsd0_diff < rmsd0_diffc) or not (np.isfinite(rmsd0_diff)
                                               and np.isfinite(rmsd0_diffc)):
             self.logger.debug(f"rmsd0_diff: {rmsd0_diff}, rmsd0c: {rmsd0_diffc}")
-            return False
-        return True
+            result = False
+        else:
+            result = True
+
+        return (result, rmsd0_diff)
 
     def qc6(self, recording1, win=None, label=''):
         # Check subtracted staircase +40mV step up is non negative
@@ -377,8 +432,11 @@ def qc6(self, recording1, win=None, label=''):
         if (val < valc) or not (np.isfinite(val)
                                 and np.isfinite(valc)):
             self.logger.debug(f"qc6_{label} val:  {val}, valc: {valc}")
-            return False
-        return True
+            result = False
+        else:
+            result = True
+
+        return (result, val)
 
     def filter_capacitive_spikes(self, current, times, voltage_step_times):
         """

tests/test_herg_qc.py

@@ -89,4 +89,9 @@ def test_run_qc(self):
                                             qc_vals_after_well, n_sweeps=2)
 
                 logging.debug(well, passed)
-                self.assertTrue(passed)
+
+                trace = ""
+                for label, results in qcs.items():
+                    if any([x == False for x, _ in results]):
+                        trace += f"{label}: {results}\n"
+                self.assertTrue(passed, trace)