(feat): match pbmc3k tutorial to seurat's

[ilan-gold]

TODO:

• Figure out best way for ranking genes so that we recover meaningful results. At the moment, running the method on the full gene list yields a lot of RP genes for the CD4 cluster ,which I would guess is basically noise. But it seems that seurat uses the full list. Separately I can't seem to figure out how to get the scores from seurat - they are claimed to be present but I don't seem them. I don't think they are just lfc.
• Marker gene documentation CD8A and CD8B are not present in the ranked genes either here or in seurat but are noted as marker genes. So I think we should just change that table and note that some genes are not present the ranked genes (maybe explain why? talk to Rahul again?)
• PCA Rahul's PCA is quite similar but not exact. It would be nice maybe for them to have arpack since an R implementation exists: https://search.r-project.org/CRAN/refmans/igraph/html/arpack.html and
• Clustering Same as above, especially since igraph is available in R: https://igraph.org/r/doc/cluster_leiden.html

rendered

[review-notebook-app]

Check out this pull request on 

See visual diffs & provide feedback on Jupyter Notebooks.

---

Powered by ReviewNB

[ilan-gold]

@flying-sheep Still very rough, but looking for some feedback given the above "outstanding" issues, especially on framing the reproducibility aspect

[flying-sheep]

Blast from the past: I also use an ARPACKy PCA in destiny, using RSpectra:

• Non-deterministic behaviour with set.seed() theislab/destiny#5
• theislab/destiny@b4e8dad

If nothing changed in the space, that might be their way forward as well, but of course I don’t know if RSpectra’s PCA is 100% identical to ARPACK.

[review-notebook-app]

View / edit / reply to this conversation on ReviewNB

flying-sheep commented on 2025-05-15T15:26:26Z
----------------------------------------------------------------

there seems to be no output from print_header

[review-notebook-app]

View / edit / reply to this conversation on ReviewNB

flying-sheep commented on 2025-05-15T15:26:27Z
----------------------------------------------------------------

I don’t really get what “up to ties” means

[review-notebook-app]

View / edit / reply to this conversation on ReviewNB

flying-sheep commented on 2025-05-15T15:26:28Z
----------------------------------------------------------------

Line #2.    adata_subset_hvg = adata[:, adata.var["highly_variable"]].copy()

hmm, maybe explain that you’re using that subset for a while until you go back to the non-subset one?

I think it’s maybe a bit confusing that there are two adata objects being used interspersedly. I think modifying a notebook like that can easily result in copy-pasting the wrong name.

ilan-gold commented on 2025-05-27T11:45:30Z
----------------------------------------------------------------

A couple of things about this:

1. If you don't do the subset, the marker genes found in 0 vs. rest are ribosomal proteins, which is not immediately clear from the Seurat equivalent but becomes immediately clear if you score the genes. Hence my mention to Rahul of this issue about providing some (transparent) way to rank. And I don't think ribosomal proteins are particularly helpful, just a guess. Th2 output:

2. Plotting the marker_genes at the end includes non-HVG genes.

3. The "How can I remove unwanted sources of variation" part of the seurat tutorial uses HVG for scaling but the part above it does not i.e., scaling without regressing out. But we can't regress out on a subset of the feature space.

ilan-gold commented on 2025-05-27T11:46:23Z
----------------------------------------------------------------

Genuinely unsure how to proceed here, we could stop doing the regressing out (and just do scaling), and then report the ribosomal protein genes. But I'd be curious to hear what Rahul has to say.

[review-notebook-app]

View / edit / reply to this conversation on ReviewNB

flying-sheep commented on 2025-05-15T15:26:28Z
----------------------------------------------------------------

We should just switch to sc.tl.marker_gene_overlap instead of changing these around everytime this file is touched.

[flying-sheep]

see above

[ilan-gold]

A couple of things about this:

1. If you don't do the subset, the marker genes found in 0 vs. rest are ribosomal proteins, which is not immediately clear from the Seurat equivalent but becomes immediately clear if you score the genes. Hence my mention to Rahul of this issue about providing some (transparent) way to rank. And I don't think ribosomal proteins are particularly helpful, just a guess. Th2 output:

2. Plotting the marker_genes at the end includes non-HVG genes.

3. The "How can I remove unwanted sources of variation" part of the seurat tutorial uses HVG for scaling but the part above it does not i.e., scaling without regressing out. But we can't regress out on a subset of the feature space.

---

View entire conversation on ReviewNB

[ilan-gold]

Genuinely unsure how to proceed here, we could stop doing the regressing out (and just do scaling), and then report the ribosomal protein genes. But I'd be curious to hear what Rahul has to say.

---

View entire conversation on ReviewNB