Fix slide_cells edge cases in region_extract (#60)

scFates/plot/milestones.py

@@ -100,7 +100,6 @@ def milestones(
         X = proj[list(adata.uns["graph"]["milestones"].values()), :]
         adata_m = sc.AnnData(
             X,
-            dtype=X.dtype,
             obs=dict(mil=list(adata.uns["graph"]["milestones"].keys())),
             obsm={basis: X},
         )

scFates/tools/bifurcation_tools.py

@@ -671,7 +671,9 @@ def gamfit(sdf):
     fitted = (fitted - fitted.min()) / (fitted.max() - fitted.min())
     fitted = fitted[df_t.index]
     changes = wins.apply(
-        lambda x: (fitted.loc[df_t[(df_t > x[0]) & (df_t < x[1])].index]).diff().sum(),
+        lambda x: (fitted.loc[df_t[(df_t > x["start"]) & (df_t < x["end"])].index])
+        .diff()
+        .sum(),
         axis=1,
     )
 

scFates/tools/dendrogram.py

@@ -91,7 +91,7 @@ def dendrogram(adata: anndata.AnnData, crowdedness: float = 1):
 
     # Assign positions to cells based on sorted order
     for seg_id, pos_idx in seg_to_pos.items():
-        df.loc[df.seg == seg_id, "seg_pos"] = pos_idx
+        df.loc[df.seg == seg_id, "seg_pos"] = str(pos_idx)
 
     df.seg_pos = df.seg_pos.astype(int).astype("category")
 

scFates/tools/pseudotime.py

@@ -227,7 +227,7 @@ def pseudotime(
                     )
                     adata.obs.loc[cells_front.index, "seg"] = cells_front.values
 
-    milestones = pd.Series(index=adata.obs_names, dtype=str)
+    milestones = pd.Series(index=adata.obs_names, dtype=object)
     for seg in pp_seg.n:
         cell_seg = adata.obs.loc[adata.obs["seg"] == seg, "t"]
         if len(cell_seg) > 0:

scFates/tools/slide_cors.py

@@ -158,9 +158,12 @@ def region_extract(pt_cur, segs_cur, nbranch):
                 freq = freq + [cell_probs]
             return freq
         else:
+            # if the first index is 0, the very first principal point already
+            # exceeds the window; keep it so the slice is not empty (see #60)
+            end = inds[0] + 1 if inds[0] == 0 else inds[0]
             pps_region = pp_next[
                 np.argsort(graph["pp_info"].loc[pp_next, "time"].values)
-            ][: inds[0]]
+            ][:end]
             if mapping:
                 cell_probs = adata.obsm["X_R"][:, pps_region].sum(axis=1)
             else:
@@ -238,9 +241,17 @@ def region_extract(pt_cur, segs_cur, nbranch):
                         .max()
                     )
                     segs_cur2 = seg_branch2
-                    res1 = region_extract(pt_cur1, segs_cur1, nbranch)
-                    res2 = region_extract(pt_cur2, segs_cur2, nbranch)
-                    return freq + res1 + res2
+                    # a branch may have no remaining cells for the current
+                    # pseudotime, yielding a NaN pt_cur; only recurse into
+                    # branches that still have cells to extract (see #60)
+                    res_final = freq
+                    if not np.isnan(pt_cur1):
+                        res1 = region_extract(pt_cur1, segs_cur1, nbranch)
+                        res_final = res_final + res1
+                    if not np.isnan(pt_cur2):
+                        res2 = region_extract(pt_cur2, segs_cur2, nbranch)
+                        res_final = res_final + res2
+                    return res_final
 
     pt_cur = graph["pp_info"].loc[pps, "time"].min()
 
@@ -363,12 +374,9 @@ def slide_cors(
     def gather_cor(i, geneset):
         freq = freqs[i][adata.obs_names]
         with np.errstate(divide="ignore", invalid="ignore"):
-            cormat = pd.DataFrame(
-                DescrStatsW(X_r.values, weights=freq).corrcoef,
-                index=genesets,
-                columns=genesets,
-            )
-        np.fill_diagonal(cormat.values, np.nan)
+            cormat_arr = np.array(DescrStatsW(X_r.values, weights=freq).corrcoef)
+        np.fill_diagonal(cormat_arr, np.nan)
+        cormat = pd.DataFrame(cormat_arr, index=genesets, columns=genesets)
         return cormat.loc[:, geneset].mean(axis=1)
 
     gather = partial(gather_cor, geneset=genesetA)