Release notes

v1.2.6

pip install -U piaso-tools

piaso.tl.leiden is PIASO’s own: parallel, and the same clusters on any number of threads

Leiden now runs in PIASO’s Rust extension by default. It is the same algorithm (Traag, Waltman and van Eck 2019), with the same objective and the same randomised refinement as igraph’s, and it reaches the same modularity. What differs is how it is run:

10 iterations	cells	PIASO	igraph
mouse developing visual cortex	200,061	3.4 s, 0.11 GB	21.8 s, 0.32 GB
whole mouse brain	2,341,350	26.3 s, 1.3 GB	553.7 s, 3.7 GB

Measured on a 10-core workstation (20 threads); memory is what the call adds to the process.

One graph and one random_state give one partition. The number of threads, the cores of the machine and the layout of the graph in memory do not enter the result: the labels are identical at 1, 2, 4, 10 and 20 threads. Parallel clustering usually gives that up, because threads see each other’s moves in an order that changes from run to run. Here a round of moves reads the state the round began with, every random number is a hash of the seed and the cell, and edge weights are summed as integers, so nothing depends on which thread did what first.

piaso.tl.leiden(adata, resolution=1.0)                  # all the cores this process may use
piaso.tl.leiden(adata, resolution=1.0, n_threads=4)     # the same labels
piaso.tl.leiden(adata, resolution=1.0, backend="igraph")  # the partitions of 1.2.5

What to know when upgrading:

• Clusters change. The partitions differ from those of 1.2.5 by about as much as another seed would change them (ARI 0.82 to 0.88 between the two on the datasets above, where two seeds of igraph agree at 0.80 to 0.89). Pass backend="igraph" to reproduce an earlier result.
• So does everything that clusters. leiden_local and runGDR call piaso.tl.leiden, so their results follow. leiden_local takes backend= and n_threads= and passes them on.
• The graph must be symmetric. A directed graph is refused, with the pair of nodes whose two directions differ. igraph’s path took the upper triangle and clustered it without a word.
• n_iterations stays at 10. Two runs give 96 to 99 % of the partition of ten in under half the time, which is worth knowing at tens of millions of cells; ten make two seeds agree a little more.
• Ctrl-C stops a run. A graph of fewer than 32,768 nodes is clustered on one thread, which is faster there.

The Leiden at scale tutorial runs it on 200,000 cells and shows the clusters beside igraph’s.

Reading a store you cannot write

The plotting functions open a .cytome read-only, so plotting from a downloaded or shared store leaves the file exactly as it was. Statistics that are normally cached in the store (cell depth, INFOG and TF-IDF parameters) are computed for the call instead when it cannot be written. Functions that compute and then write (leiden, leiden_local, neighbors, umap, runSVD, runGDR, infog, run_TFIDF) check first, and a read-only store fails in a second with a message saying why, not at the end of the run.

The motif scanner: twice as fast, and safe after fork

piaso.pp.scan_motifs(backend="rust") scores a window by walking it and the motif’s columns together, with no index arithmetic in the inner loop: 300 motifs against 100,000 sequences of 300 bp take 5.7 s instead of 10.4 s, and peak memory is 1.03 GB instead of 1.26 GB. The scores are identical to the bit.

It also hung for ever in a process forked after it had run, which is Python’s multiprocessing on Linux. It now runs on a thread pool of its own, releases Python’s lock while it works (other Python threads keep running) and stops on Ctrl-C. Every parallel function of the extension now builds its pool the same way, and n_threads=0 means the same in each: RAYON_NUM_THREADS if set, else the cores the process may use.

runSVD on a cytome no longer opens the file a second time

The SVD engine used to open the store with its own copy of SQLite while Python’s connection held it. The two could not see each other’s locks, and the engine could reset the shared index of the write-ahead log that Python had mapped; the next large write then ended the process with a bus error (SIGBUS). The engine now computes on the compressed chunks the open dataset hands it (ds.chunk_source), so no second SQLite touches the file. The result is identical to the bit, and cache_chunks=True now works with the engine, keeping the compressed chunks within half the memory available.

Also in this release

• leiden_local on a cytome makes every group’s store in one scan of the matrix (ds.split), instead of one subset per group. The labels are unchanged.
• piaso.data.chrom_sizes(genome) gives the chromosome sizes of hg38, mm10, Mmul_10, mCalJac1, hg19 and mm39 by name, or of any chrom.sizes or .fai file, reduced to the primary chromosomes.
• piaso.data.fetch_dataset checks a download before it takes its final name, and never replaces a file that is already there: a downloaded .cytome is a working file, and checking it on every call re-downloaded it over the results written into it.
• color= given a NumPy array of gene names draws one panel per name. It was read as one value per cell.
• Fixed: plot_embeddings_split without groups= could give categories each other’s colours. Each panel’s palette was built in alphabetical order but applied in the column’s category order, so wherever the two differed (a categorical order, or one stored with ds.set_categories) the colours were swapped, with no warning. Every panel now colours each category as plotEmbedding does. Figures made with 1.2.5 or earlier without groups= are worth a second look.

Requires cytome>=0.3.6.

v1.2.5

pip install -U piaso-tools

runSVD on a cytome runs its passes in Rust, and converges by default

The streaming SVD’s cost is the pass over the store: one sparse product per pass, twenty-odd passes per solve. Those passes now run in a Rust engine that reads the chunks straight out of the file, in whichever layout the matrix is stored (row-major or column-major, any integer width or float, and the zlib blobs older files carry under the zstd label), blocks each chunk against the dense operand so every random access stays in cache, and lets each thread own its slice of the output. Every entry of every product is one sum in a fixed order, so the result does not depend on the thread count and reruns are bit-identical. Measured on a 37,609-cell, 317,726-peak matrix, one round of the solver went from 27 s to 9.5 s; on a million-cell tile matrix, from 199 s to 135 s. n_threads=0 (the default) uses every core the process is allowed; PIASO_SVD_PROFILE=1 prints, per pass, how long was spent reading, bucketing and multiplying.

The solver changed with it. method="auto", now the default, runs block Lanczos with a stopping rule (tol=1e-3) and reports what it did:

Solver: block Lanczos (at most 12 rounds = 26 passes)
Solver done: block Lanczos, 9 rounds, 20 passes (converged)
Solver time: 43.0 s, of which 41.6 s in 20 passes over the store and 1.4 s in numpy between them

The old power iteration at its default of seven rounds was not converged on a peak matrix; it needs about twenty-two. method="power" with n_iter still runs exactly that many rounds when you want a fixed budget, and method="krylov" with n_iter does the same for Lanczos.

Two more defaults for peak and tile matrices on a cytome. TF-IDF is applied on the fly to raw ATAC or tile counts unless auto_tfidf=False; the pass engine weights each entry as it reads it, so no normalised layer is written. The guard that refuses to normalise an already-normalised layer runs first: inferred TF-IDF steps aside with a warning on a layer that does not look like raw counts, and an explicit auto_tfidf=True raises. And the reads take no file locks: once the write-ahead log has been checkpointed, the engine opens the store as an immutable snapshot, which is what lets it run at full speed from a network file system. On a cluster node reading from a network mount, a twenty-pass solve went from 5 min 13 s to 44.9 s; from the node’s local disk the same solve is 29 s. The installation page has a section on running on a cluster.

leiden_local re-embeds every group with runSVD

dr_method="X_svd" is the new default: each group is embedded by runSVD on its own selected features and clustered from there. On a cytome the per-group work streams from disk for any modality, with the TF-IDF runSVD applies to peaks or tiles and the store’s highly_variable selection, and each group’s subset carries only the one matrix it needs. "X_svd_full" uses every feature rather than the selected ones; "X_pca" keeps the INFOG route. min_cells, svd_method, n_iter and svd_tol pass through, and the function no longer touches scanpy. infog_svd likewise defaults to the infog layer, and where there is no highly_variable column it uses every gene and says so.

Feature selection by name, and two primate genomes

piaso.pp.selectFeatures(ds, features) flags the named features (peak ids, gene symbols, Ensembl ids, whatever the cytome stores) in one call, writing the boolean column runSVD and leiden_local read by default. match="all" and min_matched= turn a partial match into an error.

piaso.data.fetch_genome knows rhesus macaque (Mmul_10) and common marmoset (mCalJac1), with chromosome names as the Cell Ranger ARC references carry them; piaso.data.makeGenomeFiles builds the same files for any assembly from a GTF and a chromosome-sizes file. A genome without a GTF source no longer announces one.

Every reader takes a column-major matrix

Peak and tile matrices stored column-major (the layout the SVD engine prefers; cytome 0.3.3 reads and subsets it as such) are read that way throughout: score, calculateCellMetrics, calculateFeatureMetrics, normalize_log1p, the ligand-receptor readers, the TF-IDF raw-count guard and the feature reads behind plotting all take the columns when that is what the store holds, instead of asking for row chunks a column-major matrix does not have.

A breaking change: getMarkers(as_dict=True) returns the dictionary

It used to return (DataFrame, dict). Almost every use of it wanted the dict — it is the form predictCellTypeByMarker takes — and paid for that by unpacking a pair and discarding the first element. Asking for a dict now gets a dict:

marker_db = piaso.tl.getMarkers(study="AllenWholeMouseBrain_isocortex",
                                as_dict=True)
piaso.tl.predictCellTypeByMarker(adata, marker_gene_set=marker_db)

as_dict='both' returns the old pair. It is deprecated and will be removed in the next release. as_dict=False, the default, is unchanged.

The change announces itself once per session on the first as_dict=True, because a caller that keeps unpacking does not always fail loudly: with exactly two cell types, df, d = getMarkers(as_dict=True) succeeds and binds two cell-type names. Search your code for , ... = piaso.tl.getMarkers( — the tutorials are updated.

Dot plots: the same numbers, two orders of magnitude faster

piaso.pl.dotplot computed each group’s mean and expressed fraction with one Python comparison per (feature, group) pair, each scanning every cell. It now builds a one-hot group indicator once and takes two sparse products. The numbers are identical to 5e-16. A 100,000-cell, 60-feature, 30-group plot goes from 68 s to under half a second; at 20,000 cells it is 9.1 s to 0.05 s.

A missing groupby column is an error, not a wrong answer

SQLite has a compatibility misfeature: a double-quoted token that does not resolve to a column is treated as a string literal rather than raising. So SELECT "Leiden" FROM cells on a cytome with no Leiden column returns the string 'Leiden' for every row — and the caller sees one group, named after the typo, with a marker table that looks entirely plausible. That is how a COSG run reported 1 groups for a column that was never there.

Every place a user-supplied column name reaches such a query now looks it up first: cosg.run_cosg_cytome, piaso.tl.runSCALAR, and the new fragment-length functions. The error lists the columns that do exist. A case variant is still accepted, because SQLite column names are case-insensitive and leiden really does name Leiden. piaso.utils.require_cells_column is the shared check.

groups= draws a legend of what is on the plot

plotEmbedding(color='CellTypes', groups=['A','B']) draws A and B in their palette colours and everything else as one grey cloud. The legend was built from the full category list, so on a forty-cluster column it listed forty entries in forty colours, thirty-eight of which appeared nowhere on the figure — while the grey covering most of the points had no entry at all.

It now lists the highlighted groups and one final entry for the grey, labelled by the new na_label (default 'other') and shown only when something was actually left out. Highlighted groups keep the colour they have on an unfiltered plot.

color= can be the values, not only a column name

Anything computed on the fly — a regulon activity vector, a per-cell score never written back, a mask — had to be assigned into obs first. Passing the array reached color in adata.obs.columns and came back as TypeError: unhashable type: 'numpy.ndarray'.

plotEmbedding and plot_embeddings_split now accept a numpy array, a pandas Series or a Categorical as one value per cell. Numeric values get a colorbar, everything else a category legend; a Categorical keeps its own order, and a named Series titles the panel. A wrong length names both numbers. A plain list/tuple still means one panel per colour.

A Series is aligned by name when its index resolves to the cell names — anything that came out of a groupby, a sort or a .loc colours the right cells instead of being silently permuted. A RangeIndex, or an index sharing no labels with the cells, is used positionally; an index that resolves some of its labels and not others is refused, naming the count, because neither reading is safe. Duplicate cell names make alignment ambiguous, so those fall back to position with a warning.

predictCellTypeByMarker can weight its markers

A marker list treats a textbook marker and a weak, broadly expressed one identically. marker_gene_weights= passes per-gene weights through to piaso.tl.score(gene_weights=): a dict keyed by cell type, a DataFrame read column by column, or a sequence in set order.

The weights are used as given, with no normalisation per cell type. score() divides by median(weights) * n_genes and builds its control sets from the same weights, so rescaling one cell type’s weights cannot change which cell type wins — measured at 9e-16 across a 1000× rescale. Normalising would be a step that provably changes nothing.

What a rescale cannot survive is a median of zero: the divisor becomes zero and every score becomes NaN without a word. score() now rejects that, and non-finite weights, naming the gene set. Mismatched lengths, a missing cell type and an unknown one are each refused with both numbers.

Smaller things

• predictCellTypeByGDR fits on the shared label column. With a query that lacks the reference’s label column, the classifier was fitted on a column that only the reference rows carried and failed with a KeyError; it now fits on the combined column both sides share.

• Twenty-two public functions have complete docstrings, checked by a test against their signatures: every parameter present, every default the one the code uses, the cytome-only parameters and the deprecated aliases documented where they live. The remaining functions are on the same test’s backlog.

• leiden_local stops warning about its own work. It called neighbors without passing the modality on, so an ATAC run warned on every call that it had found the embedding PIASO itself had just written under another modality. The argument is forwarded; graphs and labels are unchanged, and a genuine cross-modality guess still warns.

• A per-cell column will not land on someone else’s. SQLite column names are case-insensitive, and cytome resolves a new spelling to an existing column with only a warning — so writing tss_score into a cytome already carrying a TSS_score column from another tool replaced it. piaso.utils.write_cells_column is the shared guard: the same spelling overwrites (re-running a QC step replaces its own column), a spelling that differs only in case is refused unless overwrite=True.

• The TF-IDF cache no longer grows with the dataset. compute_tfidf_stats kept its per-cell depth and per-feature IDF vectors in ds.metadata, which is a JSON document rewritten whole on every save — large enough on a real ATAC cytome for cytome to warn about its size. The depth is now a cells column, the IDF stays the feature column it was already written to, and metadata keeps the scale factor and the two column names. At 20K cells by 30K peaks that is 98 bytes instead of 1.0 MB. Cytomes written the old way are read unchanged.

v1.2.4

pip install -U piaso-tools

SCALAR: a detection floor, applied in the test

piaso.tl.runSCALAR gains expressed_pct (default 0.1) and groupby. An interaction is tested only if the ligand is detected in at least that fraction of the sender’s cells and the receptor in the receiver’s, and each gene’s matched control genes are drawn from genes that clear the same floor in the same cell type. The null therefore never contains a control that is not expressed where it is used, and the k × k support stays full. Rows that fail the floor are returned with testable = False, p_value = 1 and are excluded from FDR; two new columns, ligand_expressed_pct and receptor_expressed_pct, carry the fractions the floor was applied to. groupby names the cell-type column and is required when the floor is on; expressed_pct=None reproduces the previous behaviour exactly.

Why in the test and not in the matrix: zeroing lowly-detected genes in the specificity matrix leaves the observed pair conditioned on passing the floor while its controls are zeroed wherever they fall below it, so most of the null becomes zeros and the p-value drifts toward “detected in more cells than its peers”. Measured on the SCALAR tutorial data, a quarter of the calls that approach produces are significant only because of the zero padding.

runSCALAR now validates its input. NaN specificity scores raise (they are what cosg.indexByGene leaves for genes outside a top-N list unless set_nan_to_zero=True, and a NaN observation would otherwise score at the p-value floor). Negative scores raise (COSG marks genes failing its own expressed_pct filter with −1, and two such sentinels multiply to +1, the largest product in the table). Existing callers that pass a clean matrix and no groupby will now see an error asking for one; add groupby= or expressed_pct=None.

The SCALAR tutorial’s recipe changes with it: COSG at mu=10 over all genes, raw scores, no COSG-side filter, and the floor applied by runSCALAR. Any per-cell-type linear rescaling of the specificity matrix (for instance dividing by an IQR) leaves every p-value unchanged, because the observed product and its null scale together; cosg.iqrLogNormalize applies a log1p on top and is not recommended as SCALAR input.

SCALAR: a readable score, and a magnitude that compares across pairs

The interaction score is a product of two specificities, each below 1, so a pair of two 0.02 specificities scored 0.0004 — a number no one can hold, and one that is not comparable across sender–receiver pairs because specificity is calibrated within a cell type. Four columns fix both without touching a single p-value:

• ligand_specificity and receptor_specificity, the two inputs.
• specificity_geomean, their geometric mean: a monotone function of the product, so it ranks exactly as the score does, on the specificity scale. The tutorial’s plots are sized by it.
• null_mean and log2_enrichment: the mean of the interaction’s own matched null and the log2 ratio of the observation to it (a 10⁻¹² guard against division by zero, nothing more), so every pair is referenced to what expressed, expression-matched gene pairs achieve between the same two cell types. This is the number to compare across pairs. On the tutorial data the largest raw score in the table (PTPRC→CD22) is ten times its null; a somatostatin pair scoring 0.0008 is six thousand times its null.

statistic="min" is a second test statistic: the smaller of the two specificities, so both sides must beat the controls’ weaker sides — the “are both specific” question. Its exact null factorises into two one-sided counts. The default remains the product.

piaso.tl.umap and piaso.tl.leiden no longer return a value

Breaking. Both wrote their result and returned it on AnnData, while returning None on a cytome — one function with two contracts, and a 161,027-row array echoed into every notebook cell that called it. They now return None on both backends; read the result from adata.obsm[key_added] or adata.obs[key_added], as with every other piaso.tl function.

The two calls that have nowhere to write still return their result: data=None with a knn_result (in-memory), and leiden(..., cell_mask=...), whose labels cover the masked cells only and would mis-align if written full-length.

If you have emb = piaso.tl.umap(adata) in a script, it now binds None.

TF lists: one source that works, and a file that is checked

piaso.data.fetch_tf_list(species, source=...) replaces fetch_animaltfdb_tf_list as the entry point, defaulting to the cisTarget/SCENIC+ TF universe (~1,890 human, ~1,860 mouse).

AnimalTFDB was removed as a source. Its only host (guolab.wchscu.cn) sits behind a WAF that answers every non-browser client — urllib, requests and curl alike, under any User-Agent — with a captcha page or 405, and the old bioinfo.life.hust.edu.cn mirror no longer resolves. source="animaltfdb" now raises and says so; fetch_animaltfdb_tf_list, which shipped in 1.2.3, raises without making a request and will be deleted in 1.3. An AnimalTFDB table saved from a browser still works: pass it to load_tf_list(path=...), whose TSV parser reads its Symbol column.

Any TF-list file is now validated. The captcha page arrived as HTTP 200 and used to parse into a dozen plausible-looking “symbols” that silently shrank the TF universe. A file that starts with markup, or that yields tokens no gene symbol could be, is refused; a downloaded catalogue must also yield at least 50 symbols, and one that does not is deleted rather than cached. No minimum applies to a file you pass yourself — a ten-TF list is a legitimate thing to hand in.

For regulon inference you need neither: cytorete.tl.inferRegulon with tf_list=None already restricts to TFs that have a motif and are expressed in the data, which is the only set that can produce an edge.

Fixed: the source distribution was too large for PyPI

maturin sdist builds from the cargo package, which swept in docs/ — the rendered tutorials and their images. The archive came to 407 MB against 2.1 MB of actual source, and PyPI refuses any file over 100 MiB.

The three wheels upload first and succeed, so the effect was a red publish step at the end of every release since 1.2.0 and no sdist on PyPI for any version. Cargo.toml now excludes the documentation and the other non-source trees; the archive is 590 kB, installs from source, and passes twine check.

Fixed: a dict palette set the colours and lost the order

palette={"Neuron": ..., "Glia": ...} coloured the categories correctly and then re-sorted them alphabetically, so the sequence written at the call site was dropped without a word. A dict states an order — Python has preserved insertion order since 3.7 — and it now becomes the category order, ahead of a stored categorical order, because the call site is the more specific of the two. A partial dict leads and the rest follow alphabetically; a list palette is unaffected, since a list cannot state an order.

plot_embeddings_split had a second half of the same bug: it resolved one order across all panels and then did not pass it to the renderer, which re-sorted each panel on its own.

Fixed: legend_loc="none" left the colour bar behind

For a continuous colour the colour bar is the legend, and legend_loc never reached it — plotEmbedding drew one anyway, and plot_embeddings_split drew one per panel. Both now honour "none". Categorical legends were already suppressed correctly.

umap-learn must be 0.5.8 or newer

Earlier versions call scikit-learn’s check_array(force_all_finite=...), which scikit-learn renamed in 1.6 and removed in 1.8. The combination raises a TypeError from three frames inside scikit-learn with no mention of umap, which reads like a PIASO bug; it is an upstream one, fixed upstream in umap-learn 0.5.8. The requirement is now umap-learn>=0.5.8, and piaso.tl.umap catches that specific failure and names the package to upgrade, for environments that already have the bad pair.

pandas 3: supported, not excluded

pandas 3 turns on future.infer_string, so a string column is StringDtype rather than object. Two checks had been asking numpy about that dtype:

• np.issubdtype(StringDtype, np.number) raises rather than returning False, so the numeric-versus-categorical dispatch crashed and took every categorical embedding with it. It now asks pandas, which understands both its own extension dtypes and numpy’s. Booleans stay categorical — the one place the two libraries disagree, pandas calling bool numeric.
• A string cluster label could return from a .cytome as an integer. That was in cytome’s entity writer, which reached TEXT only through is_object_dtype and otherwise defaulted to INTEGER; fixed in cytome 0.3.2, which this release now requires. Nothing raised at the time, because SQLite stores what it is handed and only the read is typed — the damage appeared later, wherever a label was compared to a string or looked up in a stored category order.

There is no upper bound on pandas. An earlier build of this release capped it at <3; that was wrong. The problem was ours, it is fixed, and the suite passes on 3.x. A cap would have made PIASO uninstallable beside anything that requires pandas 3, which in this ecosystem is increasingly everything.

Both behaviours are covered by tests that run on pandas 2 as well, through the future.infer_string option that has existed since 2.1 — so this stays tested wherever it is installed rather than only where pandas 3 happens to be.

Fixed: piaso.tl.umap(data=None, knn_result=...) crashed

The in-memory path — documented, and advertised in this release as one of the two calls that still return their result — asked None for .obsm and died with AttributeError. It had been invisible because the test covering it skips wherever the installed umap-learn and scikit-learn disagree, which is most machines. Calling it with no knn_result now explains itself instead of raising from three frames down.

Fixed: importCellRanger reported a broken install

It imported the ATAC fragment leg unconditionally, so in the released package every call failed with ModuleNotFoundError, including modality='rna', which needs nothing from it. RNA-only now works; asking for fragments says which method is unavailable and why.

Fixed: embedding dots were drawn as rings, at 3x the size requested

plotEmbedding, plot_embeddings_split and scatter never passed linewidths to ax.scatter, so matplotlib stroked every marker at rcParams['patch.linewidth'] — 1.0 pt, in edgecolors='face'.

On a cell-sized point that stroke is most of the dot. point_size=0.2 is a disc 0.45 pt across, so it drew 1.45 pt wide: 3.2x the size asked for. With alpha < 1 the stroke and the fill composite twice around the rim, so each point rendered as an annulus with a washed-out core whose colour was neither the palette colour nor the background. Measured on a 600 dpi render, ink peaked two pixels off-centre.

All three scatter calls now pass linewidths=0, edgecolors='none', and both are exposed as parameters for anyone who wants a deliberate outline. The same fix is applied in piaso.pl.scatter and the plotByCluster jitter points. No visible change on large markers; the workaround plt.rc_context({"patch.linewidth": 0}) is no longer needed.

Fixed: plot_embeddings_split silently ignored extra keywords

It accepted **kwargs and never forwarded them, so styling a split panel had no effect and no error. Unrecognised keywords now raise TypeError, as plotEmbedding does. show is accepted as an alias for show_figure — it was among the dropped keywords, which meant show=False displayed the figure anyway.

Codebook motifs: about twice the TF universe, opt-in

piaso.data.fetch_codebook() and load_codebook() add the representative PWM set of Jolma, Laverty, Fathi et al., Nature 657, 275-283 (2026) — 1,421 human TFs against JASPAR 2024 CORE vertebrates’ 754, covering all but eleven of them. On SEA-AD middle temporal gyrus that lifts the count of TFs with a motif and expression from 483 to 1,023. The archive is ~1 MB, downloaded from the publisher under the article’s CC-BY-4.0 licence on first use and cached; nothing is redistributed with PIASO.

cytorete.tl.inferRegulon(..., motif_db="codebook") uses it. JASPAR remains the default, and that is now measured rather than assumed. Both were run end-to-end on the same 20,000 SEA-AD nuclei:

	JASPAR 2024	Codebook
TFs with motifs in the data	740	1,408
regulons	415	982 (412 of JASPAR’s 415 among them)
zinc-finger regulons	65 (16%)	385 (39%)
agreement with the TF’s own specificity	26.7%	12.7%
runtime	13 min	56 min

The result is mixed rather than bad. Codebook finds more agreeing regulons in absolute terms (125 against 111) and fixes cases JASPAR gets wrong — OLIG2 now tops OPC rather than oligodendrocytes, and SOX10, absent from JASPAR, tops oligodendrocytes. But it dilutes: agreement halves, the zinc-finger share more than doubles, and neuronal factors such as RORB and LHX6 fall in rank. Reach for it when the question is glial, zinc-finger, or coverage-limited; keep JASPAR when it is a survey.

The loader repairs what the archive carries: construct names (CASZ1.FL, ZNF729.DBD3) are filed under the gene symbol, three clone identifiers are dropped, and five CIS-BP conversions with all-zero positions are cleaned and renormalised.

Removed: cytorete.tl.markerMotifEnrichment

Never released, and removed after measurement rather than for tidiness. Two cutoff-free designs were tried on SEA-AD: a Fisher test on top-N markers against a GC- and length-matched background returned 0, 1 and 5 significant pairs at 50, 100 and 200 markers, with Jaccard 0.00 between the first two; a PIASOscore of each TF’s motif gene set against matched controls agreed with the TF’s own expression in 3-4% of cases.

The same agreement measure on inferRegulon’s 415 regulons, same cells and same motifs, is 26.7%. A PWM is a family address, not a TF address, so a motif-only statistic cannot say which cell type a TF acts in — the TF’s expression profile is what separates SPI1 from ELF1. §10 of the motif analysis tutorial now reports both measurements and hands over to inferRegulon, which adds exactly that step. cytorete.pp.build_cistrome still returns the TF × gene matrix for anyone who wants to run their own test.

Sankey: node heights, and one size key instead of two

piaso.pl.sankey gains min_node_size (height of the shortest category block, in points) and node_scale ('linear', 'sqrt', 'log'), the node-level parallel of min_flow_width and flow_scale. Both are solved the same affine way, so blocks stay monotone in group size and nothing is tied.

This fixes a distortion the default layout has always had: a node is exactly as tall as the ribbons landing on it, so under min_flow_width its height is m * size + b * fan_out, and two equally large categories are drawn at different heights when one splits more ways than the other. Measured on a two-category example at a 20 pt baseline, the four-way split draws twice as tall as the equal one-way split.

What it costs is stated in the docstring and enforced in the key: each node’s ribbons are renormalised to fill its height, so the count-to-thickness map becomes per-node and min_flow_width stops being exact.

show_width_legend / width_legend_loc are now show_size_legend / size_legend_loc, and draw one key rather than two. Only one map can be global at a time, so only one key can be true: the key describes flows when the flow map is non-proportional, and switches to describing nodes — suppressing the flow key — as soon as the nodes are sized.

Confusion matrices: hierarchical ordering

order_query / order_reference accept 'hclust' (average linkage on that axis’s own profiles) and 'hclust_symmetric' (one clustered order on both axes, for a square table over one label set). 'auto' now prefers the symmetric order when the two axes carry the same labels — two independent orders moved a real diagonal off the diagonal — and is unchanged otherwise.

Fixed

• The GDR SettingWithCopyWarning. Each batch was scored by slicing an anndata view, whose constructor rewrites categorical obs columns through _remove_unused_categories — an anndata-internal pandas idiom, and pure overhead here since the scorer reads only the matrix and var_names. The batch object is now built directly. This also removes two warnings.catch_warnings blocks, one of them inside a ThreadPoolExecutor worker, where mutating the process-global filter list races between threads.

• plotLigandReceptorLollipop tick labels were formatted to one decimal, so an axis spanning a few 10⁻⁴ — what neuronal ligand–receptor scores look like — read 0.0 at every tick. The precision now follows the axis range.

v1.2.3

pip install -U piaso-tools

Fixed

• piaso.tl.predictCellTypeByGDR discarded its own result. The function copies the query internally to avoid mutating .X, and the copy rebound the local name, so the prediction was written to the copy. It printed “All finished. The predicted cell types are saved as CellTypes_gdr in adata.obs” and the caller’s AnnData came back without that column. The cytome path was unaffected — it writes through a separate handle, which is why this survived. If you called it on an AnnData and found no prediction, this was why.

• plot_embeddings_split(vmin_pct=…, vmax_pct=…) was accepted and ignored. The parameters arrived through **kwargs and were dropped, so a request for 10/90 silently produced the full range. Nothing raised and the figure still rendered — only the contrast was wrong. Now implemented to match plotEmbedding, with explicit vmin/vmax still taking precedence.

Added

• piaso.pp.alignSpatialCoordinates. Sections are placed on their chips independently, so a split plot renders each panel at its own offset; this centres each sample on its own centroid, and with_std=True also equalises apparent section size. Works on AnnData, an open cytome, or a path, and takes groupby= or batch_key= for the same column.

  This was described in the v1.2.2 notes but did not make the 1.2.2 wheel. It ships here.

• palette= accepts a colormap name for a categorical column, sampled across the categories in their resolved order. On ordered categories — ages, timepoints, stages — that turns the legend into a ramp instead of a set of unrelated hues.

v1.2.2

pip install -U piaso-tools

Spatial

• Tissue images and regions of interest. piaso.pl.plotEmbedding and plot_embeddings_split take image=True to draw a registered morphology image under a spatial embedding, with orientation and units handled: the image is placed in coordinate space, so nothing needs flipping or scaling. Needs cytome ≥ 0.2.6, which stores the image in the file.
• cell_mask= plots a subset — a boolean mask or integer indices — without writing a subset file. It pairs with cytome’s cells_in_region, so one rectangle selects its cells and, via spatial_images.crop, its pixels.
• piaso.pp.rotateSpatialCoordinates works directly on a cytome, not only on an AnnData, and rebuilds the coordinate index in the same call so region queries keep agreeing with what you plot.

Plotting

• Continuous colour now depends on where the value came from. A feature keeps PIASO’s sequential color_1; a numeric cell column — regulon activity, a QC metric, pseudotime — gets Spectral_r. An explicit cmap= still wins. Both are continuous, but drawing metadata on the expression ramp made a metadata panel read as expression.
• Fixed: plot_embeddings_split(fix_coordinate_ratio=True) is documented as the default but only the False branch acted, so panels were stretched to their subplot box and spatial tissue rendered squashed.
• Fixed: basis='spatial' resolved to the last matching embedding, so writing an aligned or backup copy beside it silently redirected existing calls. Candidates are ranked now — the name that is the key beats one that merely contains it.

v1.2.1

pip install -U piaso-tools

• piaso.tl.cospecificity_trans is public; it is the trans co-specificity used by cytorete to decide regulon edges.
• The GRN entry points that used to live in piaso.tl remain as lazy forwarders to cytorete, raising an actionable ImportError when it is not installed.
• Motif scanning (piaso.pp.scan_motifs), the PWM loaders and the SCREEN cCRE registry are part of the public API.

v1.2.0

pip install -U piaso-tools

cytome datasets and streaming

PIASO works directly on cytome files. Anywhere a function took an AnnData, it also takes a path to a .cytome file or an open dataset:

import piaso, cytome

ds = cytome.open("atlas.cytome")
piaso.tl.runGDR(ds, groupby="cell_type", layer="infog")

The matrix is read in chunks rather than loaded, so peak memory is set by the batch size instead of by the number of cells. cytome is installed automatically — there is no extra to remember.

A self-contained analysis workflow

From raw UMI counts to clusters, an embedding and marker genes:

import piaso, cosg

# Input must be RAW UMI counts. This reads adata.X by default; if your raw
# counts live in a layer, pass it explicitly:
#     piaso.tl.infog(adata, layer="counts", n_top_genes=3000)
piaso.tl.infog(adata, n_top_genes=3000)

# Pass layer="infog" — runSVD defaults to adata.X, which would silently
# ignore the normalization above.
piaso.tl.runSVD(adata, layer="infog", n_components=50, key_added="X_svd")

piaso.tl.neighbors(adata, use_rep="X_svd", n_neighbors=15)
piaso.tl.leiden(adata, resolution=1.0, key_added="leiden")
piaso.tl.umap(adata, use_rep="X_svd")

cosg.cosg(adata, groupby="leiden", key_added="cosg")
piaso.pl.embedding(adata, basis="X_umap", color="leiden")

Every step runs on a plain pip install piaso-tools — no scanpy required. scanpy remains an optional extra (pip install 'piaso-tools[scanpy]') for interoperability.

One key to watch when moving older code over: piaso.tl.leiden writes its labels to adata.obs['leiden'] by default, whereas the scanpy-based tutorials that preceded this release passed key_added='Leiden' explicitly. Downstream calls that name the column need the same spelling.

Motif scanning

piaso.pp.scan_motifs is a pure-numpy PWM scanner; scan_motifs_rust is a rayon-parallel implementation with the same contract. The motif-database loaders and .2bit sequence access live alongside it in piaso.data, so the whole workflow — fetch a genome, extract sequences, load PWMs, scan — is available from one install.

Reference data

piaso.data fetches and caches what analyses need: genome sequence and annotation, example datasets, motif databases (JASPAR, CIS-BP, cisTarget), and the SCREEN cCRE registry.

Requirements

Requires cosg>=1.1.0. pip install -U piaso-tools pulls it in.

Security

pyo3 upgraded to 0.29.2, closing GHSA-36hh-v3qg-5jq4 and GHSA-chgr-c6px-7xpp. PIASO called neither affected API, so exposure was nil, but the crate linked the code.

v1.1.0

Marker-gene-guided dimensionality reduction (GDR), INFOG normalization, gene-set scoring, cell-type annotation and label transfer, and the plotting suite. See the GitHub releases for the full history.