cosg 1.2.1

pip install -U cosg

Fixes and new options. Read "Re-run if affected" first: four changes move results that may have been reported. The dependencies are unchanged.

Re-run if affected

• use_raw=True reported the wrong gene names, in COSG 1.0.4 to 1.2.0. Affected: cosg.cosg(adata, use_raw=True) on an AnnData whose adata.raw holds more genes than adata, which is scanpy's usual layout after selecting variable genes. COSG scored adata.raw.X but named its columns by adata.var_names, so a marker at column j of raw was reported as the gene at position j of adata (a group's markers g2, g3 came back as g32, g33), or raised IndexError past its end. The names, fold changes and p-values now describe raw's genes. In 1.0.1 to 1.0.3 use_raw had no effect (adata.X was scored, and named correctly). Re-run any use_raw=True analysis made with 1.0.4 or later.
• Stratified p-values with a group alone in a batch. Affected: cosg.cosg(..., batch_key=..., calculate_pvalues=True) where some batch holds cells of only one group (sorted populations per sample, one type per donor). That batch's sum entered the statistic but not its mean, so on null data the group's p-values were near 1e-80 and every gene it expresses was called. A batch that holds none of a group's cells, or nothing but them, carries no information about the group and is now left out of its test (statistic, mean, variance and saddlepoint together); the p-values match a stratified permutation reference. Strata of one or two cells now count when they inform the group, as in the permutation null. Runs without batch_key, or whose every batch mixes groups, are unchanged; re-run the others.
• The normalized scores of a group that scores above zero at fewer than 5 % of the features. Affected: IQR-normalized scores, that is iqr_normalize on a cytome (on by default for the dict, long and dense outputs) and cosg.iqrLogNormalize, mostly on peaks and tiles, where a rare or shallowly sequenced cell type can score above zero at a few percent of them. The spread of such a group, between its q_upper and q_lower quantiles (the 95th and 75th percentiles by default), is zero: a cytome path then divided its scores by 1e-6, which put them about log(scale / 1e-6) above other groups', and iqrLogNormalize gave every column without a spread the smallest spread of another column. A group without a spread now keeps a scale of its own: its q_upper quantile when that is positive, else the 5th percentile of its positive scores (a group positive at 5 to 5.3 % of the features has its spread floored there, so the scale does not jump as the share crosses 5 %). The normalized scores of other groups are unchanged, and no ranking within a group changes. Re-run if you compared normalized scores across groups and one of them was this sparse.
• No expression filter on the GPU batch path, in 1.1.0 to 1.2.0. cosg.cosg(adata, device='gpu', batch_key=...) applied no expression filter, so remove_lowly_expressed and expressed_pct had no effect there and a gene detected in a handful of a group's cells could be one of its markers. The filter is now applied after that path, as after the CPU's batch path, and removes the same genes. Re-run such an analysis; its markers change wherever a gene was below the floor.

Other changes to results

• P-values deep in the tail. The saddlepoint stopped at a normal p of 1e-30, which left the sparse markers, exactly the ones whose normal p is wrong by orders of magnitude, at their normal value: a gene in 12 cells, all of a 100-cell group of 5,000, had log10 p -129.7 against an exact -20.7. Sparse entries (fewer than 30 expressing cells expected in the group) are now refined however deep (-20.2 here; within 0.1 to 0.8 log10 of the exact tail); dense genes that deep keep the normal approximation's value, as before. No p is reported below the permutation floor, 1 / C(n, n_c) (per batch with batch_key), the probability of the observed labeling itself. neg_log10_pvals of sparse and of the very strongest markers fall; decisions at any usual alpha do not change. pvalue_method='exact' no longer reports 0 for a strong marker whose tail lies past its table's cut (the saddlepoint takes those entries).
• The groups in the result are the groups scored that have cells. groups= with a categorical groupby returned every category, the unselected ones with arbitrary genes, and calculate_logfoldchanges then raised IndexError; an unused category did the same with groups='all'. Such categories are now left out, as on the cytome paths. Fold changes are those of the groups scored (reference may be a label outside groups) and no longer need a categorical groupby. PIASO's supervised runGDR with batch_key loses, for a batch that lacks a type, the gene set of arbitrary genes that type used to contribute.
• One rule for a group's top features on every path: by score, ties by feature index (the paths chose among tied scores differently, e.g. 10 tied peaks gave {p1, p4, p5} or {p0, p1, p2}). n_genes_user above the number of features returns them all (it raised). The ndarray output of the feature-batched cytome paths lists a group's next features by score, zeros and then -1, instead of "" and 0.0.
• "Expressed" means a nonzero value on every path. The streaming, GPU and cytome paths counted stored entries, explicit zeros included, and dense streaming input counted values above 0. Results change only for matrices that store zeros or hold negative values.
• expressed_pct=0 (or a group's 0) has one meaning on every path: the expressed_min_num_cells_in_target_group minimum alone, 3 cells by default. It used to switch the filter off altogether, the minimum included, on the default chunked CPU path, the GPU's chunked path and every path of run_cosg_cytome, and to keep the minimum on the others. Results change only for a 0 on those paths, where a feature seen in one or two cells of a group is no longer among its markers; remove_lowly_expressed=False switches the filter off.
• Float64 sums on a cytome's unbatched path (the RNA default), as on the AnnData paths: its scores move at the level of float32 rounding of the sums (1.5e-6 over a store of many chunks). The feature-batched paths keep float32 sums. log1p normalization is vectorized (80 times faster on the column chunks of a column-major matrix), with the same numbers to the bit.
• A default call writes nothing into a cytome. A path is opened read-only (cytome 0.3.5 or later), and log1p takes each cell's depth from the data, the row sums of the chunk it normalizes, instead of the RNA_cell_depth cache, which was trusted whenever its length matched (a matrix replaced by one of the same shape kept the old depths). The depths are the ones the cache held when it was fresh, so results do not move, and the first call on a store is faster. PIASO's infog and tfidf statistics are recomputed by each call on a path rather than cached in the store, which costs a store without them one more pass over the matrix per call (about 25 % more time in our test); pass a Dataset opened read-write to keep caching them. write_to_cytome=True writes only its markers, and a store it cannot write fails before the computation.
• plotMarkerDotplot draws the plot itself. scanpy's sc.pl.dotplot, which it called, cut the group labels beside the brackets to twice the number of genes per group (with swap_axes=True and top_n_genes=2, 'L6 CT', 'L6 IT' and 'L6 IT Car3' all read 'L6 .') and drew the dendrogram on a scale of its own, out of register with the rows. The default renderer (backend='native') needs no scanpy, keeps whole labels and draws the dendrogram on the rows. save is now a file path (scanpy pasted it after figures/dotplot_, so an absolute path raised FileNotFoundError). The function returns None with show=True (the default), as before, and the dot grid's axes with show=False. New arguments: show, dpi, expression_cutoff and backend.

New

• A detection floor per group. expressed_pct also takes one fraction per group: a dict {group: fraction}, a pandas Series indexed by group, or a sequence in the order of the result's groups. The share of cells that suits an abundant, deeply sequenced type need not suit a rare or shallow one, and on peaks and tiles a floor read from the data can be an order of magnitude below 0.1. Every path honours it: cosg.cosg on CSR, CSC or dense input, chunked or not, with batch_key or p-values, on the CPU or the GPU; run_cosg_cytome unbatched, batched by chromosome or by disk, on row- or column-major matrices, and on the GPU; and score_from_sums.
  • A group's floor removes features from that group's markers only, after the cosine and the penalty, which do not involve it, so the same value for every group gives the result of that value as a number, to the bit. One thing couples groups: under scale_rule='positive_median' a group with fewer than 20 positive scores takes the median of the other groups' scales, which their floors move.
  • With a dict or a Series every group with cells needs a fraction: a group left out raises rather than silently taking 0.1. A key that names no label of groupby raises; a key for a label not scored in the call (left out by groups or cell_mask) is ignored, so one table of every type serves any mask. A table of per-type cell counts is passed as fractions (count - 0.5) / n_cells, which make the threshold "count cells or more" whatever the rounding.
  • The floors used are recorded: in adata.uns[key_added]['params']['expressed_pct'] as a one-record array with a field per group, and in run_cosg_cytome's result under 'expressed_pct', as {group: fraction}, which write_to_cytome stores with the markers. A run with a single number returns exactly what it did before.
  • expressed_pct is checked on the cytome paths too: a value outside [0, 1] raises there as it did on the AnnData path. expressed_pct=None means 0.1 on every path.
• Column-major cytome matrices. A matrix stored as column (CSC) chunks only, which earlier versions could not read, is scored a block of features at a time: a column chunk holds every entry of its features, so memory follows the block, not the number of features times the number of groups. On counts the scores are those of the same matrix stored by rows.
• cosg.score_from_sums(dot, sq, group_sizes, nnz) returns COSG's raw scores from the sums every streaming path accumulates: per feature and group, the sum of the feature's values; per feature, the sum of its squares; the group sizes; and, for the expression filter, the nonzero counts. Sums of parts add up, so the sums of count-split halves, of shards or of added cells can be scored without the matrix, and per-batch sums (batch_group_sizes) give the batch-corrected scores of batch_key.
• scale_rule='positive_median' (run_cosg_cytome, cosg.cosg on a cytome, iqrLogNormalize) scales each group by the median of its positive scores, for sparse features (peaks, tiles), where most groups score above zero at 5 to 25 % of the features and the spread over all features is set by how many of them score at all (depth, the expression filter, the size of the feature set) more than by the group's own scores. A group with fewer than 20 positive scores takes the median of the other groups' scales. The default stays 'iqr'. A scale divides a group's scores, so no rule changes a ranking within a group.
• max_accumulator_bytes (run_cosg_cytome, 2 GiB by default): on a row-major matrix, chromosome batching holds every chromosome's per-feature sums in one pass only while they fit; beyond it the chromosomes are packed into passes that fit, so memory stops growing with the number of features times the number of groups (a whole-brain atlas's tiles by its subclasses would need about 13 GB). The scores are the same to the bit.
• neg_log10_pvals_adj, beside neg_log10_pvals: -log10 of the Benjamini-Hochberg adjusted p, computed in log space, so it stays finite where pvals_adj underflows to 0. neg_log10_pvals is of the raw p, so a line at 1.3 on it is raw 0.05, not FDR 0.05; this one is safe to threshold.

Fixed

• P-values on a cytome: strong markers got exactly 0.0; p is now floored at the smallest positive float64, and neg_log10_pvals and neg_log10_pvals_adj are reported, as on the AnnData path (they agree to 1e-6 on the same counts).
• device='gpu' on a cytome (and device='auto' on a machine with CuPy) dropped calculate_pvalues without a word and did not apply iqr_normalize. The GPU path now shares the CPU's post-processing: every output_format and the scale options work there. What it cannot do (p-values, feature batching, a column-major matrix) raises for device='gpu' and runs on the CPU for device='auto'.
• Scores that were not finite: the chunked GPU path (chosen automatically when the matrix exceeds 70 % of free GPU memory) marked filtered features with -inf instead of -1, which iqrLogNormalize turned into NaN. With cpu_chunk_size=0, batch_key and a dense matrix, a group too small in every batch made every score NaN; it now scores 0, as on the other batch paths. gpu_chunk_size=0 (the documented monolithic transfer) raised ZeroDivisionError.
• Inputs are checked before they are scored: a NaN or an infinity in the matrix made its gene the top marker of every group on the default path; non-finite values are now refused, naming the features (on a cytome, at no extra cost). A cell without a label raised KeyError or TypeError; it is refused with a count. mu=NaN made every score NaN. A call that failed wiped the earlier result under the same key_added; it now leaves it in place.
• cell_mask on a cytome: a boolean mask shorter than the cells was accepted and the last cells silently left out, and integer indices with batch_key raised IndexError. Both forms are validated (length, range, duplicates) and give the same result in any order.
• A feature name held twice in a cytome's dense output: a 10x GRCh38 reference holds 10 gene symbols twice (HSPA14, TBCE, MATR3, ...). output_format='dense' mapped names to rows, so one row of each pair stayed 0 and the other took either feature's score; it is now filled by feature index, each row with its own score. The dict and long outputs, keyed by name, keep the better score of such a pair and warn.
• Chromosome batching assumed each chromosome's features contiguous; on a store not ordered by chromosome, features were scored twice and the IQR scale padded to more features than there are. The batches are now runs of one chromosome. A column-major matrix with the default layer='auto' crashed on the unbatched path. Column chunks are normalized with their own features' TF-IDF and INFOG factors.
• A groupby column missing from a cytome's cells table was read by SQLite as a string, so every cell took that name as its label and the run reported one group named after it, with a plausible marker table. A missing column now raises, naming the columns there and a match that differs only in case; such a name is accepted, as SQLite resolves it.
• min_total_count on the feature-batched path of run_cosg_cytome was dropped with expressed_pct=0 or remove_lowly_expressed=False. It now applies whenever it is above 0, as documented.
• A group label holding "/" (L2/3 IT) left an AnnData that adata.write could not save: h5ad and zarr keys cannot hold "/", and anndata 0.13 refuses to write them by default. The by-group frame adata.uns[key]['COSG'] now writes "/" as "%2F" in its column names (names::L2%2F3 IT); indexByGene, and with it the plots that read the frame, reads it back as "/", for frames from earlier versions too, and the record arrays (names, scores, ...) keep the labels exactly. Labels without "/" are written as before.
• run_cosg_cytome(output_format='dense') raised "assignment destination is read-only" with pandas 3, whose copy-on-write makes a frame's .values read-only.
• run_cosg_cytome checks mu and n_genes_user as cosg.cosg does, and expressed_min_num_cells_in_target_group=0 is accepted on both paths. COSGR, the R package, applies expressed_pct without the 3-cell minimum, so for groups under about 30 cells it can report other markers.

Reproducing 1.2.0

• Normalized scores of the sparse groups above: iqr_fallback='constant' on a cytome (run_cosg_cytome, or cosg.cosg on a cytome) and iqr_fallback='smallest' in iqrLogNormalize.
• The dot plot: plotMarkerDotplot(..., backend='scanpy'), which needs pip install 'cosg[dotplot]'.
• The other changes have no switch back.

Requirements

Unchanged. The cosg[cytome] extra still accepts cytome 0.2.3 or later; with cytome before 0.3.5 a path is opened for writing, as in earlier releases.

Documentation

• Marker genes and their significance
• COSG on a cytome
• COSG across batches
• COSG on the GPU
• COSG on spatial data