Support non-perturbed scRNA-seq and mouse data

.claude/skills/h5mu-structure/scripts/inspect_h5mu.py

@@ -45,7 +45,11 @@ def get_unique_examples(data, n=3):
             return ", ".join(_fmt_value(x) for x in uniques)
         elif sp.issparse(data):
             flat = data.toarray().flatten()
-            uniques = np.unique(flat[~np.isnan(flat)] if np.issubdtype(flat.dtype, np.floating) else flat)[:n]
+            uniques = np.unique(
+                flat[~np.isnan(flat)]
+                if np.issubdtype(flat.dtype, np.floating)
+                else flat
+            )[:n]
             return ", ".join(_fmt_value(x) for x in uniques)
         elif isinstance(data, np.ndarray):
             flat = data.flatten()
@@ -182,7 +186,7 @@ def format_mapping_section(name, mapping, indent="", is_last_section=False):
         lines.append(tree_line(is_last_section, f"{name}/", indent))
         child_indent = tree_child_indent(is_last_section, indent)
         for i, key in enumerate(items):
-            is_last_item = (i == len(items) - 1)
+            is_last_item = i == len(items) - 1
             val = mapping[key]
             if name in ("uns",):
                 desc = describe_uns_value(val)
@@ -207,9 +211,11 @@ def format_adata(adata, mod_name=None):
     obs_names_examples = ", ".join(str(x) for x in adata.obs_names[:3])
     all_obs_items = ["obs_names"] + obs_cols
     for i, col in enumerate(all_obs_items):
-        is_last = (i == len(all_obs_items) - 1)
+        is_last = i == len(all_obs_items) - 1
         if col == "obs_names":
-            lines.append(tree_line(is_last, f"obs_names: e.g. {obs_names_examples}", obs_indent))
+            lines.append(
+                tree_line(is_last, f"obs_names: e.g. {obs_names_examples}", obs_indent)
+            )
         else:
             desc = describe_column(adata.obs[col])
             lines.append(tree_line(is_last, f"{col} — {desc}", obs_indent))
@@ -221,9 +227,11 @@ def format_adata(adata, mod_name=None):
     var_names_examples = ", ".join(str(x) for x in adata.var_names[:3])
     all_var_items = ["var_names"] + var_cols
     for i, col in enumerate(all_var_items):
-        is_last = (i == len(all_var_items) - 1)
+        is_last = i == len(all_var_items) - 1
         if col == "var_names":
-            lines.append(tree_line(is_last, f"var_names: e.g. {var_names_examples}", var_indent))
+            lines.append(
+                tree_line(is_last, f"var_names: e.g. {var_names_examples}", var_indent)
+            )
         else:
             desc = describe_column(adata.var[col])
             lines.append(tree_line(is_last, f"{col} — {desc}", var_indent))
@@ -238,8 +246,12 @@ def format_adata(adata, mod_name=None):
         ("varp", adata.varp),
     ]
     for i, (sec_name, sec_data) in enumerate(sections):
-        is_last = (i == len(sections) - 1)
-        lines.extend(format_mapping_section(sec_name, sec_data, indent="", is_last_section=is_last))
+        is_last = i == len(sections) - 1
+        lines.extend(
+            format_mapping_section(
+                sec_name, sec_data, indent="", is_last_section=is_last
+            )
+        )
 
     return lines
 
@@ -278,7 +290,9 @@ def format_mudata(mdata, filename):
         obs_indent = tree_child_indent(False)
         for i, col in enumerate(obs_cols):
             desc = describe_column(mdata.obs[col])
-            lines.append(tree_line(i == len(obs_cols) - 1, f"{col} — {desc}", obs_indent))
+            lines.append(
+                tree_line(i == len(obs_cols) - 1, f"{col} — {desc}", obs_indent)
+            )
     else:
         lines.append(tree_line(False, "obs columns: (none beyond obs_names)"))
 
@@ -289,7 +303,9 @@ def format_mudata(mdata, filename):
         var_indent = tree_child_indent(False)
         for i, col in enumerate(var_cols):
             desc = describe_column(mdata.var[col])
-            lines.append(tree_line(i == len(var_cols) - 1, f"{col} — {desc}", var_indent))
+            lines.append(
+                tree_line(i == len(var_cols) - 1, f"{col} — {desc}", var_indent)
+            )
     else:
         lines.append(tree_line(False, "var columns: (none beyond var_names)"))
 
@@ -302,8 +318,12 @@ def format_mudata(mdata, filename):
         ("varp", mdata.varp),
     ]
     for i, (sec_name, sec_data) in enumerate(sections):
-        is_last = (i == len(sections) - 1)
-        lines.extend(format_mapping_section(sec_name, sec_data, indent="", is_last_section=is_last))
+        is_last = i == len(sections) - 1
+        lines.extend(
+            format_mapping_section(
+                sec_name, sec_data, indent="", is_last_section=is_last
+            )
+        )
 
     return "\n".join(lines) + "\n"
 

.claude/skills/perturbNMF-runner/references/02-evaluation.md

@@ -37,7 +37,8 @@ Reconstruction error and stability are computed automatically.
 | `--guide_annotation_path` | None | TSV with a `targeting` column to identify non-targeting controls (alternative to `--guide_annotation_key`) |
 | `--gene_names_key` | `symbol` | Column in data_guide["rna"].var with gene names |
 | `--FDR_method` | `StoreyQ` | FDR correction: `StoreyQ` or `BH` |
-| `--organism` | `human` | Species for enrichment |
+| `--organism` | `human` | Species for enrichment (`human` or `mouse`; mouse genes are mapped to human orthologs) |
+| `--ortholog_table` | None | Mouse-human ortholog TSV; default MGI + Ensembl, cached in `Resources/` |
 | `--n_top` | `300` | Top genes for enrichment tests |
 | `--guide_annotation_key` | `["non-targeting"]` | Non-targeting guide identifiers (accepts multiple values) |
 | `--use_cache` | flag | Load enrichr gene set libraries from cached JSON in `Resources/` instead of downloading; falls back to download + cache on miss |

.claude/skills/perturbNMF-runner/references/04-visualization.md

@@ -29,6 +29,8 @@
 | `--selected_k` | None | K value to highlight with a red dashed line |
 | `--stability_file` | None | Pre-computed stability/error TSV/NPZ (bypasses `cnmf.consensus()`; for torch-cNMF runs) |
 | `--run_program_dotplot` | off | Enable per-(K, sel_thresh) program dotplots (default off; requires the inference h5mu artifact per K/threshold) |
+| `--no_trait` | off | Skip trait (GWAS) enrichment inputs and plots (non-human data) |
+| `--no_perturbation` | off | Skip perturbation inputs and regulator plots (non-perturbed data, no guides) |
 
 Enrichment file patterns use `{k}` placeholder (and `{sample}` for perturbation). See `references/parameter-catalog.md` Section 4 for all optional params.
 

.claude/skills/perturbNMF-runner/references/05-annotation-summary.md

@@ -120,6 +120,7 @@ other stage (no notebook required). Submit one job per K to cover multiple K val
 | `--perturbation_file_name` | `perturbation_association_results` | Perturbation file stem (e.g. `CRT`) |
 | `--effect_size` | `log2FC` | Effect-size column (e.g. `approx_log2FC`) |
 | `--gene_names_key` | `symbol` | var column with gene symbols |
+| `--organism` | `human` | Species of the data; `mouse` annotates genes with their human orthologs' summaries |
 | `--non_targeting_key` | `non-targeting` | Negative-control target label(s) |
 
 See `references/parameter-catalog.md` (Section 12) for the full flag list, including

.claude/skills/perturbNMF-runner/references/parameter-catalog.md

@@ -219,7 +219,8 @@ Same as sk-cNMF except:
 | `--gwas_data_path` | str | None | Path to GWAS data (required only for trait enrichment) |
 | `--X_normalized_path` | str | None | Normalized counts h5ad (needed for explained variance). Path pattern: `<out_dir>/<run_name>/Inference/cnmf_tmp/Inference.norm_counts.h5ad` (prefix is always `Inference`, not the run_name) |
 | `--guide_annotation_path` | str | None | Guide annotation TSV |
-| `--organism` | str | `human` | Species for enrichment |
+| `--organism` | str | `human` | Species for enrichment (`human` or `mouse`; mouse genes are mapped to human orthologs) |
+| `--ortholog_table` | str | None | Mouse-human ortholog TSV; default MGI + Ensembl, cached in `Resources/` |
 | `--FDR_method` | str | `StoreyQ` | FDR method: `StoreyQ` or `BH` |
 | `--n_top` | int | `300` | Top genes for enrichment tests |
 | `--gene_names_key` | str | `symbol` | Column in data_guide["rna"].var with gene names |
@@ -271,6 +272,8 @@ Same as sk-cNMF except:
 | `--variance_col` | str | `Total` | No | Variance column name |
 | `--stability_file` | str | None | No | Pre-computed stability/error file (TSV or NPZ) |
 | `--run_program_dotplot` | flag | False | No | Enable per-(K, sel_thresh) program dotplots; requires `{output_dir}/{run_name}/Inference/adata/cNMF_{K}_{thresh}.h5mu` to exist for each pair |
+| `--no_trait` | flag | False | No | Skip trait (GWAS) enrichment files and plots (non-human data) |
+| `--no_perturbation` | flag | False | No | Skip perturbation files and regulator plots (non-perturbed data, no guides) |
 
 ---
 
@@ -616,6 +619,8 @@ cover multiple K values.
 | `--save_path` | str | `{out_dir}/{run_name}/Interpretation/Summary_table/{K}_{thresh}` | Output dir for the `.xlsx`, TSVs, and specificity/correlation/KD sidecars |
 | `--mdata_path` | str | `{out_dir}/{run_name}/Inference/adata/cNMF_{K}_{thresh}.h5mu` | Input MuData |
 | `--num_gene` | int | `300` | Top genes per program / per enriched term to keep |
+| `--organism` | str | `human` | Species of the data; with `mouse`, each gene gets the MyGene summaries of its human orthologs (listed in a `Human ortholog` column) |
+| `--ortholog_table` | str | None | Mouse-human ortholog TSV; default MGI + Ensembl, cached in `Resources/` |
 | `--perturbation_file_name` | str | `perturbation_association_results` | Perturbation result file stem (between `{K}_` and `_{Condition}.txt`); e.g. `CRT` |
 | `--Conditions` | str (nargs=\*) | all labels found in the data | Condition labels (values of obs[categorical_key]); `--Sample` is a deprecated alias |
 | `--effect_size` | str | `log2FC` | Effect-size column in perturbation files (e.g. `approx_log2FC`) |

.claude/skills/perturbNMF-runner/scripts/add_loadings_to_h5mu.py

@@ -36,83 +36,108 @@
 
 
 def parse_args():
-    p = argparse.ArgumentParser(description=__doc__,
-                                formatter_class=argparse.RawDescriptionHelpFormatter)
-    p.add_argument('--h5mu', required=True, help='Existing .h5mu to patch in place')
-    p.add_argument('--gene_spectra_score', required=True,
-                   help='cNMF gene_spectra_score TSV (K x genes)')
-    p.add_argument('--prog_key', default='cNMF', help='Program modality to patch')
-    p.add_argument('--data_key', default='rna',
-                   help='Gene modality to validate the spectra width against '
-                        '(skipped if absent)')
-    p.add_argument('--overwrite', action='store_true',
-                   help="Replace varm['loadings'] / uns['var_names'] if already present")
-    p.add_argument('--dry_run', action='store_true',
-                   help='Validate and report, but write nothing')
+    p = argparse.ArgumentParser(
+        description=__doc__, formatter_class=argparse.RawDescriptionHelpFormatter
+    )
+    p.add_argument("--h5mu", required=True, help="Existing .h5mu to patch in place")
+    p.add_argument(
+        "--gene_spectra_score",
+        required=True,
+        help="cNMF gene_spectra_score TSV (K x genes)",
+    )
+    p.add_argument("--prog_key", default="cNMF", help="Program modality to patch")
+    p.add_argument(
+        "--data_key",
+        default="rna",
+        help="Gene modality to validate the spectra width against (skipped if absent)",
+    )
+    p.add_argument(
+        "--overwrite",
+        action="store_true",
+        help="Replace varm['loadings'] / uns['var_names'] if already present",
+    )
+    p.add_argument(
+        "--dry_run", action="store_true", help="Validate and report, but write nothing"
+    )
     return p.parse_args()
 
 
 def main():
     args = parse_args()
 
-    print(f'Reading gene spectra score: {args.gene_spectra_score}', flush=True)
-    loadings = pd.read_csv(args.gene_spectra_score, sep='\t', index_col=0)
+    print(f"Reading gene spectra score: {args.gene_spectra_score}", flush=True)
+    loadings = pd.read_csv(args.gene_spectra_score, sep="\t", index_col=0)
     loadings.index = loadings.index.astype(str)
-    print(f'  loadings: {loadings.shape[0]} programs x {loadings.shape[1]} genes',
-          flush=True)
+    print(
+        f"  loadings: {loadings.shape[0]} programs x {loadings.shape[1]} genes",
+        flush=True,
+    )
 
-    mode = 'r' if args.dry_run else 'r+'
+    mode = "r" if args.dry_run else "r+"
     with h5py.File(args.h5mu, mode) as f:
-        mod = f'mod/{args.prog_key}'
+        mod = f"mod/{args.prog_key}"
         if mod not in f:
-            raise ValueError(f"modality '{args.prog_key}' not found in {args.h5mu} "
-                             f"(have: {list(f['mod'].keys())})")
+            raise ValueError(
+                f"modality '{args.prog_key}' not found in {args.h5mu} "
+                f"(have: {list(f['mod'].keys())})"
+            )
 
-        prog_ids = read_elem(f[f'{mod}/var']).index.astype(str)
+        prog_ids = read_elem(f[f"{mod}/var"]).index.astype(str)
         if set(prog_ids) != set(loadings.index):
             raise ValueError(
-                f'program ids differ between the h5mu ({len(prog_ids)}) and the '
-                f'spectra score ({loadings.shape[0]}): e.g. missing from spectra '
-                f'{sorted(set(prog_ids) - set(loadings.index))[:5]}')
+                f"program ids differ between the h5mu ({len(prog_ids)}) and the "
+                f"spectra score ({loadings.shape[0]}): e.g. missing from spectra "
+                f"{sorted(set(prog_ids) - set(loadings.index))[:5]}"
+            )
         # Row i of varm must be program var_names[i]; reindex rather than trust order.
         loadings = loadings.loc[prog_ids]
 
-        data_var = f'mod/{args.data_key}/var'
+        data_var = f"mod/{args.data_key}/var"
         if data_var in f:
             n_genes = read_elem(f[data_var]).shape[0]
             if n_genes != loadings.shape[1]:
                 raise ValueError(
                     f"the spectra score has {loadings.shape[1]} genes but modality "
-                    f"'{args.data_key}' has {n_genes} -- Stage 2a asserts these match")
+                    f"'{args.data_key}' has {n_genes} -- Stage 2a asserts these match"
+                )
             print(f"  validated against '{args.data_key}': {n_genes} genes", flush=True)
         else:
-            print(f"  note: no '{args.data_key}' modality, skipping the gene-count check",
-                  flush=True)
+            print(
+                f"  note: no '{args.data_key}' modality, skipping the gene-count check",
+                flush=True,
+            )
 
-        targets = [(f'{mod}/varm', 'loadings'), (f'{mod}/uns', 'var_names')]
-        existing = [f'{g}/{k}' for g, k in targets if f'{g}/{k}' in f]
+        targets = [(f"{mod}/varm", "loadings"), (f"{mod}/uns", "var_names")]
+        existing = [f"{g}/{k}" for g, k in targets if f"{g}/{k}" in f]
         if existing and not args.overwrite:
-            raise ValueError(f'already present: {existing}. Pass --overwrite to replace.')
+            raise ValueError(
+                f"already present: {existing}. Pass --overwrite to replace."
+            )
 
         if args.dry_run:
-            print('Dry run: validation passed, nothing written.', flush=True)
+            print("Dry run: validation passed, nothing written.", flush=True)
             return 0
 
         for group, key in targets:
-            if f'{group}/{key}' in f:
-                print(f'  overwriting {group}/{key}', flush=True)
-                del f[f'{group}/{key}']
-
-        print(f'Writing {mod}/varm/loadings and {mod}/uns/var_names', flush=True)
-        write_elem(f[f'{mod}/varm'], 'loadings',
-                   np.asarray(loadings.values, dtype=np.float64))
-        write_elem(f[f'{mod}/uns'], 'var_names',
-                   loadings.columns.astype(str).to_numpy())
-
-    print(f'Done. {args.h5mu} now carries loadings for {loadings.shape[0]} programs '
-          f'x {loadings.shape[1]} genes.', flush=True)
+            if f"{group}/{key}" in f:
+                print(f"  overwriting {group}/{key}", flush=True)
+                del f[f"{group}/{key}"]
+
+        print(f"Writing {mod}/varm/loadings and {mod}/uns/var_names", flush=True)
+        write_elem(
+            f[f"{mod}/varm"], "loadings", np.asarray(loadings.values, dtype=np.float64)
+        )
+        write_elem(
+            f[f"{mod}/uns"], "var_names", loadings.columns.astype(str).to_numpy()
+        )
+
+    print(
+        f"Done. {args.h5mu} now carries loadings for {loadings.shape[0]} programs "
+        f"x {loadings.shape[1]} genes.",
+        flush=True,
+    )
     return 0
 
 
-if __name__ == '__main__':
+if __name__ == "__main__":
     sys.exit(main())