Repository: owkin/PyDESeq2
Pull request #26 (2022-12-28 13:51:43)

CI: implement tests for design factors
This PR implements tests on the design factor in the clinical data. More precisely:

- The `DeseqDataset` class now tests whether the design factor contains NaN values and checks the design factors takes exactly 2 values. 
- Pytests now test that an error is thrown in case a `DeseqDataset` is initialised with NaNs, or not exactly two values in the design factor.

DO NOT MERGE: there are a few tests that remain to be implemented 


Patch:
diff --git a/pydeseq2/DeseqDataSet.py b/pydeseq2/DeseqDataSet.py
index 72cf35c..17bf63f 100644
--- a/pydeseq2/DeseqDataSet.py
+++ b/pydeseq2/DeseqDataSet.py
@@ -44,7 +44,13 @@ class DeseqDataSet:
         Must be indexed by sample barcodes.
 
     design_factor : str
-        Name of the column of clinical to be used as a design variable. (default: 'high_grade').
+        Name of the column of clinical to be used as a design variable.
+        (default: 'high_grade').
+
+    reference_level : str
+        The factor to use as a reference. Must be one of the values taken by the design.
+        If None, the reference will be chosen alphabetically (last in order).
+        (default: None).
 
     min_mu : float
         Threshold for mean estimates. (default: 0.5).
@@ -54,7 +60,8 @@ class DeseqDataSet:
 
     max_disp : float
         Upper threshold for dispersion parameters.
-        NB: The threshold that is actually enforced is max(max_disp, len(counts)). (default: 10).
+        NB: The threshold that is actually enforced is max(max_disp, len(counts)).
+        (default: 10).
 
     refit_cooks : bool
         Whether to refit cooks outliers. (default: True).
@@ -149,6 +156,7 @@ class DeseqDataSet:
         counts,
         clinical,
         design_factor="high_grade",
+        reference_level=None,
         min_mu=0.5,
         min_disp=1e-8,
         max_disp=10,
@@ -165,13 +173,26 @@ class DeseqDataSet:
 
         # Test counts before going further
         test_valid_counts(counts)
-
         self.counts = counts
+
+        # Import clinical data and convert design_column to string
         self.clinical = clinical
+        if not self.counts.index.identical(self.clinical.index):
+            raise ValueError(
+                "The count matrix and clinical data should have the same index."
+            )
         self.design_factor = design_factor
+        if self.clinical[self.design_factor].isna().any():
+            raise ValueError("NaNs are not allowed in the design factor.")
+        self.clinical[self.design_factor] = self.clinical[self.design_factor].astype(str)
+
         # Build the design matrix (splits the dataset in two cohorts)
         self.design_matrix = build_design_matrix(
-            self.clinical, design_factor, expanded=False, intercept=True
+            self.clinical,
+            design_factor,
+            ref=reference_level,
+            expanded=False,
+            intercept=True,
         )
         self.min_mu = min_mu
         self.min_disp = min_disp
diff --git a/pydeseq2/utils.py b/pydeseq2/utils.py
index 60edf03..98013d7 100644
--- a/pydeseq2/utils.py
+++ b/pydeseq2/utils.py
@@ -173,14 +173,33 @@ def build_design_matrix(
         A DataFrame with experiment design information (to split cohorts).
         Indexed by sample barcodes.
     """
+
     design_matrix = pd.get_dummies(clinical_df[design])
+
+    if design_matrix.shape[-1] == 1:
+        raise ValueError(
+            f"The design factor takes only one value: "
+            f"{design_matrix.columns.values.tolist()}, but two are necessary."
+        )
+    elif design_matrix.shape[-1] > 2:
+        raise ValueError(
+            f"The design factor takes more than two values: "
+            f"{design_matrix.columns.values.tolist()}, but should have exactly two."
+        )
     # Sort columns alphabetically
     design_matrix = design_matrix.reindex(
         sorted(design_matrix.columns, reverse=True), axis=1
     )
-    if ref is not None:  # Put reference level first
-        ref_level = design_matrix.pop(ref)
-        design_matrix.insert(0, ref, ref_level)
+    if ref is not None:  # Put reference level last
+        try:
+            ref_level = design_matrix.pop(ref)
+        except KeyError as e:
+            print(
+                "The reference level must correspond to"
+                " one of the design factor values."
+            )
+            raise e
+        design_matrix.insert(1, ref, ref_level)
     if not expanded:  # drop last factor
         design_matrix.drop(columns=design_matrix.columns[-1], axis=1, inplace=True)
     if intercept:


Test patch:
diff --git a/tests/test_edge_cases.py b/tests/test_edge_cases.py
index 3878d66..febbb29 100644
--- a/tests/test_edge_cases.py
+++ b/tests/test_edge_cases.py
@@ -51,6 +51,7 @@ def test_zero_genes():
     assert res_df.loc[zero_genes].padj.isna().all()
 
 
+# Tests on the count matrix
 def test_nan_counts():
     """Test that a ValueError is thrown when the count matrix contains NaNs."""
     counts_df = pd.DataFrame({"gene1": [0, np.nan], "gene2": [4, 12]})
@@ -89,3 +90,60 @@ def test_non_negative_counts():
 
     with pytest.raises(ValueError):
         DeseqDataSet(counts_df, clinical_df, design_factor="condition")
+
+
+# Tests on the clinical data (design factors)
+def test_nan_factors():
+    """Test that a ValueError is thrown when the design factor contains NaNs."""
+    counts_df = pd.DataFrame({"gene1": [0, 1], "gene2": [4, 12]})
+    clinical_df = pd.DataFrame({"condition": [0, np.NaN]})
+
+    with pytest.raises(ValueError):
+        DeseqDataSet(counts_df, clinical_df, design_factor="condition")
+
+
+def test_one_factors():
+    """Test that a ValueError is thrown when the design factor takes only one value ."""
+    counts_df = pd.DataFrame({"gene1": [0, 1], "gene2": [4, 12]})
+    clinical_df = pd.DataFrame({"condition": [0, 0]})
+
+    with pytest.raises(ValueError):
+        DeseqDataSet(counts_df, clinical_df, design_factor="condition")
+
+
+def test_too_many_factors():
+    """Test that a ValueError is thrown when the design factor takes
+    more than two values."""
+    counts_df = pd.DataFrame({"gene1": [0, 1, 5], "gene2": [4, 12, 8]})
+    clinical_df = pd.DataFrame({"condition": [0, 1, 2]})
+
+    with pytest.raises(ValueError):
+        DeseqDataSet(counts_df, clinical_df, design_factor="condition")
+
+
+def test_reference_level():
+    """Test that a ValueError is thrown when the reference level is not one of the
+    design factor values."""
+    counts_df = pd.DataFrame({"gene1": [0, 1], "gene2": [4, 12]})
+    clinical_df = pd.DataFrame({"condition": [0, 1]})
+
+    with pytest.raises(KeyError):
+        DeseqDataSet(
+            counts_df, clinical_df, design_factor="condition", reference_level="control"
+        )
+
+
+def test_indexes():
+    """Test that a ValueError is thrown when the count matrix and the clinical data
+    don't have the same index."""
+    counts_df = pd.DataFrame(
+        {"gene1": [0, 1], "gene2": [4, 12]}, index=["sample1", "sample2"]
+    )
+    clinical_df = pd.DataFrame({"condition": [0, 1]}, index=["sample01", "sample02"])
+
+    with pytest.raises(ValueError):
+        DeseqDataSet(
+            counts_df,
+            clinical_df,
+            design_factor="condition",
+        )