resolve nan issue with merge priority

pinellolab · Oct 3, 2023 · bbfc8d2 · bbfc8d2
1 parent 7446cbe
commit bbfc8d2
Show file tree

Hide file tree

Showing 5 changed files with 21 additions and 14 deletions.
diff --git a/bean/annotate/filter_alleles.py b/bean/annotate/filter_alleles.py
@@ -491,7 +491,8 @@ def _map_alleles_to_filtered(
             # prioritize allele with most mean counts
             merge_priority = guide_filtered_allele_counts.mean(
                 axis=1, numeric_only=True
-            ).values
+            )
+            merge_priority = merge_priority.fillna(0)
             if is_cn_allele:
                 guide_raw_counts["allele_mapped"] = guide_raw_counts[allele_col].map(
                     lambda allele: allele.map_to_closest(

diff --git a/bean/annotate/translate_allele.py b/bean/annotate/translate_allele.py
@@ -419,7 +419,7 @@ def translate_allele(
                 )
             return get_allele_aa_change_multi_genes(
                 allele,
-                gene_names=fasta_file_dict,
+                gene_names=gene_names,
                 include_synonymous=include_synonymous,
             )
         elif fasta_file:
@@ -432,11 +432,19 @@ def translate_allele(
                 fasta_file=fasta_file,
                 include_synonymous=include_synonymous,
             )
-        return get_allele_aa_change_multi_genes(
-            allele,
-            fasta_file_dict=fasta_file_dict,
-            include_synonymous=include_synonymous,
-        )
+        elif fasta_file_dict:
+            return get_allele_aa_change_multi_genes(
+                allele,
+                fasta_file_dict=fasta_file_dict,
+                include_synonymous=include_synonymous,
+            )
+        else:
+            fasta_file_name = os.path.dirname(be.__file__) + "/annotate/ldlr_exons.fa"
+            return get_allele_aa_change_single_gene(
+                allele,
+                fasta_file=fasta_file_name,
+                include_synonymous=include_synonymous,
+            )
     except RefBaseMismatchException as e:
         if not allow_ref_mismatch:
             raise e
@@ -446,7 +454,7 @@ def translate_allele(
 
 def translate_allele_df(
     allele_df,
-    include_synonymouns=True,
+    include_synonymous=True,
     allow_ref_mismatch=True,
     fasta_file=None,
     fasta_file_dict: Dict[str, str] = None,
@@ -457,7 +465,7 @@ def translate_allele_df(
     translated_alleles = allele_df.allele.map(
         lambda a: be.translate_allele(
             a,
-            include_synonymouns=include_synonymouns,
+            include_synonymous=include_synonymous,
             allow_ref_mismatch=allow_ref_mismatch,
             fasta_file=fasta_file,
             fasta_file_dict=fasta_file_dict,

diff --git a/bean/annotate/utils.py b/bean/annotate/utils.py
@@ -267,7 +267,7 @@ def check_args(args):
         + int(args.translate_fastas_csv is not None)
         + int(args.translate_gene is not None)
         + int(args.translate_genes_list is not None)
-        != 1
+        > 1
     ):
         raise ValueError(
             "Invalid arguments: You should specify exactly one of --translate-fasta, --translate-fastas-csv, --translate-gene, translate-genes-list to translate alleles."

diff --git a/bean/framework/Edit.py b/bean/framework/Edit.py
@@ -235,9 +235,7 @@ def map_to_closest(
                         raise ValueError(
                             f"merge_priority length {len(merge_priority)} is not the same as allele_list length {len(allele_list)}"
                         )
-                    nt_max_idx = nt_max_idx[
-                        np.argmax(merge_priority[nt_max_idx], skipna=True)
-                    ]
+                    nt_max_idx = nt_max_idx[np.nanargmax(merge_priority[nt_max_idx])]
                 else:
                     nt_max_idx = nt_max_idx[0]
                 if nt_jaccards[nt_max_idx] > jaccard_threshold:

diff --git a/bean/mapping/GuideEditCounter.py b/bean/mapping/GuideEditCounter.py
@@ -743,7 +743,7 @@ def _check_names_filter_fastq(self, filter_by_qual=False):
             )
 
         R1, R2 = self._get_seq_records()
-
+        info("Done loading reads for quality filtering")
         _check_readname_match(R1, R2)
         if filter_by_qual:
             self.n_reads_after_filtering = self._filter_read_quality(R1, R2)