# Best Practices for Querying Large Genomic Variant Datasets in TileDB-VCF: A Parallel Query Guide

> Optimize genomic variant queries with TileDB-VCF. Learn best practices for parallel querying of large datasets achieving sub-second query times with partitioning and streaming.

- Repository: [K-Dense/scientific-agent-skills](https://github.com/K-Dense-AI/scientific-agent-skills)
- Tags: best-practices
- Published: 2026-05-14

---

**TileDB-VCF achieves sub-second queries on terabyte-scale variant datasets by partitioning genomic regions and sample cohorts across parallel workers with explicit memory budgets and streaming iterators.**

When working with the `K-Dense-AI/scientific-agent-skills` repository, the TileDB-VCF skill documented in [`scientific-skills/tiledbvcf/SKILL.md`](https://github.com/K-Dense-AI/scientific-agent-skills/blob/main/scientific-skills/tiledbvcf/SKILL.md) provides the canonical reference for implementing high-performance variant queries. This guide covers the exact patterns used to parallelize queries across billions of genomic records while maintaining predictable memory consumption.

## Understanding the TileDB-VCF Sparse Array Model

TileDB-VCF stores VCF/BCF data as **sparse multidimensional arrays**, where genomic coordinates (`chrom`, `pos_start`, `pos_end`) serve as dimensions and each sample is stored as an attribute. This layout enables TileDB to fetch only the sub-array covering requested regions and samples—a design that naturally supports parallelism.

The data model follows **1-based coordinates** according to the VCF specification (first base = 1). When defining query regions, always use 1-based positioning as implemented in the core library referenced in [`scientific-skills/tiledbvcf/SKILL.md`](https://github.com/K-Dense-AI/scientific-agent-skills/blob/main/scientific-skills/tiledbvcf/SKILL.md).

## Configuring Memory and Tile Partitions

Before executing parallel queries, instantiate a `ReadConfig` object to cap resource usage and define partition boundaries. According to the source documentation, three critical parameters govern performance:

- **`memory_budget`**: Caps RAM per query worker (typically 1–2 GB per worker for moderate loads, increased for high-throughput jobs)
- **`region_partition`**: Splits the genome into independent chunks (e.g., `(0, 3095677412)` for whole-genome runs)
- **`sample_partition`**: Divides the sample dimension across workers (e.g., `(0, 10000)` for up to 10,000 samples per partition)

```python
import tiledbvcf

config = tiledbvcf.ReadConfig(
    memory_budget=2048,
    region_partition=(0, 3095677412),
    sample_partition=(0, 10000),
)
ds = tiledbvcf.Dataset(uri="my_dataset", mode="r", cfg=config)

```

Setting these explicitly prevents out-of-memory crashes on large datasets and allows the TileDB engine to optimize tile traversal for concurrent access patterns.

## Parallel Query Strategies

### Region-Based Parallelism

TileDB-VCF accepts lists of **disjoint genomic regions**, with each region processed independently. For optimal throughput, break the genome into approximately **10 Mb chunks** to maximize parallelism while keeping individual result sets manageable.

```python
regions = [
    "chr1:1-10000000",
    "chr1:10000001-20000000",
    "chr2:1-10000000",
    # ...

]

def query_region(region):
    return ds.read(
        attrs=["sample_name", "pos_start", "pos_end", "fmt_GT"],
        regions=[region],
        samples=["sample1", "sample2"]
    )

# Execute across 8 workers

from multiprocessing import Pool
with Pool(processes=8) as pool:
    results = pool.map(query_region, regions)

```

### Sample-Based Parallelism

When querying across thousands of samples, partition the sample list into cohorts of **1,000 to 5,000 samples** per worker. Retrieve the full sample list via `ds.samples()`, then chunk accordingly:

```python
all_samples = ds.samples()
sample_chunks = [all_samples[i:i+1000] for i in range(0, len(all_samples), 1000)]

def query_sample_chunk(sample_subset):
    return ds.read(
        attrs=["sample_name", "fmt_GT"],
        regions=["chr1:1-5000000"],
        samples=sample_subset
    )

with Pool(processes=4) as pool:
    results = pool.map(query_sample_chunk, sample_chunks)

```

### Grid-Based Parallelism

For the largest workloads, create a **grid of (region × sample) tiles** and assign each tile to a separate process. This approach fully saturates CPU and I/O bandwidth by combining both partitioning strategies.

```python
regions = [f"chr1:{i}-{i+9999999}" for i in range(1, 250000001, 10000000)]
sample_chunks = [all_samples[i:i+1000] for i in range(0, len(all_samples), 1000)]
tasks = [(r, s) for r in regions for s in sample_chunks]

def query_tile(region, sample_subset):
    return ds.read(
        attrs=["sample_name", "fmt_GT"],
        regions=[region],
        samples=sample_subset
    )

with Pool(processes=12) as pool:
    for df in pool.starmap(query_tile, tasks):
        process(df)  # e.g., aggregate allele frequencies

```

### Streaming Large Results

When queries return millions of rows, materializing the entire result set in memory causes instability. Use **`read_iter()`** with an explicit `batch_size` to stream results:

```python
iterator = ds.read_iter(
    attrs=["sample_name", "pos_start", "pos_end", "fmt_GT"],
    regions=regions,
    samples=all_samples,
    batch_size=100_000,  # rows per batch

)

for batch in iterator:
    write_to_parquet(batch)  # or compute incremental statistics

```

Streaming is essential for result sets exceeding **10⁶ rows** and works seamlessly with any parallel partitioning strategy.

## Optimizing Ingestion for Parallel Reads

Parallel query performance depends on how the dataset was created. Ensure the array was built using **parallel ingestion** so the data layout is optimally tiled for concurrent reads. The `tiledbvcf store` command automatically parallelizes over input VCF files:

```bash
tiledbvcf store --uri my_dataset --samples sample1.vcf.gz,sample2.vcf.gz

```

For maximum efficiency, split large VCF collections into batches (e.g., 100 files per batch) and run multiple `store` commands concurrently. The [`scientific-skills/tiledbvcf/SKILL.md`](https://github.com/K-Dense-AI/scientific-agent-skills/blob/main/scientific-skills/tiledbvcf/SKILL.md) file notes that ingestion parallelism directly determines read scalability.

## Cloud Storage Configuration

TileDB-VCF reads directly from S3, Azure Blob, and GCS URIs. Pass VFS configuration via `tiledb_config` to optimize remote access:

```python
cfg = {"vfs.s3.region": "us-east-1", "vfs.s3.use_virtual_addressing": "true"}
ds = tiledbvcf.Dataset("s3://my-bucket/dataset", mode="r", tiledb_config=cfg)

```

Ensure the compute environment has proper credentials (IAM roles, `AWS_ACCESS_KEY_ID`, etc.) before executing parallel queries across cloud-hosted datasets.

## Performance Verification Checklist

Before running production workloads, verify these settings:

- **Memory budget** matches available RAM per worker (`config.memory_budget`)
- **Region chunks** average ~10 Mb (adjust if execution time varies significantly between chunks)
- **Sample partitions** contain 1,000–5,000 samples each
- **Streaming enabled** for expected outputs >1,000,000 rows via `read_iter()`
- **Cloud credentials** accessible (validate with `ds.samples()`)
- **Parallel ingestion** used during dataset creation (check ingestion logs for parallel execution confirmation)

## Summary

- **Partition the query space** using 10 Mb genomic regions and sample cohorts of 1,000–5,000 to maximize parallel efficiency.
- **Set explicit memory budgets** via `tiledbvcf.ReadConfig` to prevent OOM errors in multi-worker environments.
- **Stream large results** using `read_iter()` with defined batch sizes instead of `read()` when expecting millions of rows.
- **Ensure parallel ingestion** was used during dataset creation to optimize the underlying array layout for concurrent access.
- **Configure cloud VFS settings** explicitly for S3/Azure/GCS workloads to ensure credential and region compatibility.

## Frequently Asked Questions

### What is the optimal region chunk size for parallel TileDB-VCF queries?

The [`scientific-skills/tiledbvcf/SKILL.md`](https://github.com/K-Dense-AI/scientific-agent-skills/blob/main/scientific-skills/tiledbvcf/SKILL.md) documentation recommends approximately **10 Mb chunks** (e.g., `chr1:1-10000000`) as the optimal balance between parallelism and memory usage. Smaller chunks increase scheduling overhead, while larger chunks risk memory exhaustion when querying high-density variant regions across many samples.

### How do I prevent out-of-memory errors when querying billions of variants?

Set a conservative `memory_budget` (typically 1–2 GB per worker) in `tiledbvcf.ReadConfig`, partition samples into groups of 1,000–5,000, and use `read_iter()` with a manageable `batch_size` (e.g., 100,000 rows) instead of `read()`. This combination caps resident memory regardless of total dataset size.

### Can TileDB-VCF execute parallel queries directly against S3 storage?

Yes. Pass a `tiledb_config` dictionary with VFS parameters (e.g., `vfs.s3.region`) to the `Dataset` constructor, then execute the same parallel query patterns used for local files. TileDB’s internal engine handles concurrent reads from cloud storage, though you must ensure the compute environment has valid AWS/GCP/Azure credentials configured.

### What is the difference between `read()` and `read_iter()` in TileDB-VCF?

`read()` materializes the entire query result into a single DataFrame in memory, suitable for small result sets. `read_iter()` returns a Python iterator that yields DataFrame batches of size `batch_size`, enabling processing of multi-gigabyte results with constant memory usage. For queries returning more than one million rows, `read_iter()` is the required approach to maintain system stability.