Integrating datasets by merging#

Dataset integration starts by placing compatible assays in one datastore. DataStoreMerge aligns their feature order, carries selected metadata, and records the source of each cell. It does not alter expression values or correct the joint representation. This guide builds that uncorrected baseline first.

1. Load compatible source stores#

The control and interferon beta stimulated Kang PBMC stores use the same RNA feature space. Their publication recipe physically removes cells without an imported cell-type label before running source-level quality control. The remaining I cell key records that quality-control selection.

import pandas as pd

import scarf

scarf.configure_output(level="ERROR", progress=False)

repository = scarf.cytebase.connect("scarf_docs")
ctrl_path = repository.download_dataset(
    name="kang_15K_pbmc_rnaseq",
    destination="scarf_datasets",
    zarr=True,
)
stim_path = repository.download_dataset(
    name="kang_14K_ifnb-pbmc_rnaseq",
    destination="scarf_datasets",
    zarr=True,
)

ds_ctrl = scarf.DataStore(f"{ctrl_path}/data.zarr", nthreads=4)
ds_stim = scarf.DataStore(f"{stim_path}/data.zarr", nthreads=4)

Confirm assay type, cell counts, and feature counts before merging. DataStoreMerge validates the feature axes; matching gene symbols alone do not establish compatible genome builds or quantification conventions.

pd.DataFrame(
    [
        {
            "source": label,
            "assay": type(store.RNA).__name__,
            "cells": store.cells.N,
            "active cells": int(store.cells.fetch_all("I").sum()),
            "features": store.RNA.feats.N,
        }
        for label, store in (("ctrl", ds_ctrl), ("stim", ds_stim))
    ]
)
source assay cells active cells features
0 ctrl RNAassay 8487 8486 35635
1 stim RNAassay 10111 10111 35635

2. Merge counts and metadata#

names supplies the source labels, source_column names their metadata column, and prepend_text keeps imported metadata names distinct from columns authored in the merged store. reset_cell_filter=False preserves the source quality-control selections.

merged_path = "scarf_datasets/kang_dataset_merging.zarr"
scarf.DataStoreMerge(
    datasets=[ds_ctrl, ds_stim],
    zarr_path=merged_path,
    names=["ctrl", "stim"],
    assays=["RNA"],
    prepend_text="orig",
    reset_cell_filter=False,
    source_column="sample_id",
    overwrite=True,
).dump()

merged = scarf.DataStore(merged_path, nthreads=4)

sample_id records the source label. Columns imported from the sources keep the orig_ prefix so their origin remains explicit.

The merged active population contains labelled cells from both sources.

merged.cells.to_pandas_dataframe(
    ["sample_id", "orig_cluster_labels"],
    key="I",
).groupby("sample_id")["orig_cluster_labels"].agg(
    cells="count",
    cell_types="nunique",
)
cells cell_types
sample_id
ctrl 8486 13
stim 10111 13

3. Open the rebuilt uncorrected baseline#

The catalog’s merged store is rebuilt with the merge recipe above and a labelled standard RNA run. Open that frozen run instead of repeating PCA, graph construction, clustering, and UMAP in this merge tutorial. Its graph uses 21 neighbours so the correction methods on the next page can branch from the same baseline.

prepared_path = repository.download_dataset(
    name="kang_29K_ctrl-ifnb_pbmc_rnaseq",
    destination="scarf_datasets",
    zarr=True,
)
ds = scarf.DataStore(f"{prepared_path}/data.zarr", nthreads=4)
baseline = ds.pipeline.open(label="docs_default")
sorted(baseline)
['analysis_cell_selection',
 'ann_index',
 'cluster_selection',
 'clusters',
 'connectivity_map',
 'embedding_initialization',
 'feature_universe',
 'highly_variable_features',
 'input_cell_selection',
 'leiden_1.0',
 'neighbors',
 'normalized',
 'pca',
 'umap']

The durable run maps each output name to its exact artifact. Requested metadata and results remain in its frozen view.

One plotting call compares source identity, imported cell types, and the exact clustering artifact on the same layout.

ds.plots.embedding(
    layout=baseline["umap"],
    color_by=["sample_id", "orig_cluster_labels", baseline["clusters"]],
    n_columns=3,
)
../_images/a87cbc686de1d80abe8fa54b9b0cc96c55e71f4b82f428755484102f50e16181.png

A proportional composition plot makes source dominance within the uncorrected Leiden clusters explicit.

pd.crosstab(
    baseline.cells.fetch("clusters"),
    baseline.cells.fetch("sample_id"),
    normalize="index",
)
col_0 ctrl stim
row_0
1 0.007531 0.992469
2 0.004461 0.995539
3 0.002390 0.997610
4 0.403670 0.596330
5 0.988064 0.011936
6 0.005602 0.994398
7 0.472826 0.527174
8 0.002604 0.997396
9 0.391771 0.608229
10 0.013255 0.986745
11 0.520000 0.480000
12 0.966443 0.033557
13 0.974338 0.025662
14 0.410256 0.589744
15 0.360360 0.639640
16 0.004975 0.995025
17 0.997170 0.002830
18 0.993827 0.006173
19 0.998350 0.001650
20 0.993289 0.006711

iLISI summarizes local source mixing on a zero-to-one scale. Zero means the median neighbourhood effectively contains cells from only one source. One is the maximum mixing score across the observed sources.

uncorrected_ilisi = ds.metric_ilisi(
    batch_colname="sample_id",
    neighbors=baseline["neighbors"],
    perplexity=7,
)
{"uncorrected iLISI": round(uncorrected_ilisi, 3)}
{'uncorrected iLISI': 0.0}

The stimulated sample received interferon beta, and PBMC cell types do not all respond identically to that treatment. Source-associated structure can therefore include biological response as well as technical variation. An interferon-response gene such as ISG15 makes that stim-enriched program visible on the same uncorrected layout.

Inspect treatment-linked expression separately before interpreting the source mixing as purely technical.

This page establishes the uncorrected observation; Correcting batch effects with Harmony compares how partial PCA and Harmony change it. Keep uncorrected counts for condition-level differential expression.