Skip to content

TEST_DATA_SPEC.md

Latest commit

 

History

History
332 lines (243 loc) · 7.64 KB

TEST_DATA_SPEC.md

File metadata and controls

332 lines (243 loc) · 7.64 KB
Raw

Test Data Specification

This document specifies the synthetic test data used for development and validation.

Overview

The test dataset simulates sequencing data from closely related bacterial strains with known variants, enabling validation of the detection algorithms.

Reference Genome

Properties

Property Value
Length 100,000 bp
GC content 50%
Structure Single circular chromosome (represented as linear)

Generation Method

import random

def generate_reference(length: int, gc_content: float, seed: int) -> str:
    random.seed(seed)
    gc_bases = ['G', 'C']
    at_bases = ['A', 'T']
    
    sequence = []
    for _ in range(length):
        if random.random() < gc_content:
            sequence.append(random.choice(gc_bases))
        else:
            sequence.append(random.choice(at_bases))
    
    return ''.join(sequence)

Use seed=42 for reproducibility.

Core Genes

Properties

Property Value
Number of genes 20
Size range 500-1000 bp each
Total coverage ~15,000 bp (15% of genome)
Distribution Evenly spaced across genome

Positions

Genes are placed at regular intervals to avoid overlap and ensure even distribution:

Gene  Start   End     Length
----  -----   ---     ------
1     1000    1750    750
2     5000    5600    600
3     10000   10800   800
4     15000   15500   500
5     20000   20900   900
...
20    95000   95700   700

Exact positions will be generated programmatically with random lengths in 500-1000 range.

Requirements

  • No variants (SNPs, indels) should be placed within core genes
  • Core genes must be single-copy in all samples
  • Used for robust coverage normalization

Sample A (Baseline)

Sample A contains reads from the unmodified reference genome.

Properties

Property Value
Source Reference genome (no variants)
Coverage 100x
Read length 150 bp
Insert size 300-500 bp (mean 400, sd 50)

Sample B (Variant Sample)

Sample B contains reads from a modified genome with known variants.

Variants

SNPs (10 total)

ID Position Ref Alt Notes
SNP1 3500 A G Outside core genes
SNP2 8200 C T Outside core genes
SNP3 12500 G A Outside core genes
SNP4 18000 T C Outside core genes
SNP5 25000 A T Outside core genes
SNP6 35000 G C Outside core genes
SNP7 48000 C A Outside core genes
SNP8 62000 T G Outside core genes
SNP9 78000 A C Outside core genes
SNP10 92000 G T Outside core genes

Exact positions will be chosen to avoid core genes.

Duplication (1 region)

Property Value
Position 40000-42000
Length 2000 bp
Copy number in A 1
Copy number in B 2
Expected detection ~2x k-mer counts in region

The duplicated region is tandem (inserted immediately after original).

Deletion (1 region)

Property Value
Position 70000-71000
Length 1000 bp
Copy number in A 1
Copy number in B 0
Expected detection Missing k-mers (unique to A)

Sample C (Second Variant Sample)

Sample C provides a third sample for testing all-vs-all comparison.

Variants

SNPs (8 total, partially overlapping with B)

ID Position Ref Alt Shared with B?
SNP1 3500 A G Yes (same as B)
SNP2 8200 C T Yes (same as B)
SNP3 22000 A G No (unique to C)
SNP4 30000 T A No (unique to C)
SNP5 45000 G C No (unique to C)
SNP6 55000 C T No (unique to C)
SNP7 85000 A G No (unique to C)
SNP8 98000 T C No (unique to C)

Structural Variants

Type Position Length Notes
Duplication 50000-51500 1500 bp Different region than B
Deletion 82000-82500 500 bp Different region than B

Read Simulation

Paired-End Parameters

Parameter Value
Read length 150 bp
Insert size mean 400 bp
Insert size std dev 50 bp
Insert size min 300 bp
Insert size max 500 bp

Error Model

Illumina-like errors with position-dependent rate:

def error_rate(position: int, read_length: int) -> float:
    """
    Error rate increases toward end of read.
    
    Position 0: ~0.1% error rate
    Position 149: ~0.5% error rate
    """
    base_rate = 0.001
    end_rate = 0.005
    fraction = position / read_length
    return base_rate + (end_rate - base_rate) * fraction

Error Types

Type Proportion
Substitution 95%
Insertion 2.5%
Deletion 2.5%

For simplicity, the initial implementation can use substitutions only.

Quality Scores

Quality scores follow error probability:

def quality_score(error_prob: float) -> int:
    """Convert error probability to Phred score."""
    import math
    if error_prob <= 0:
        return 40  # Max quality
    return min(40, int(-10 * math.log10(error_prob)))

Position-dependent quality:

  • Start of read: Q30-Q40 (0.1% - 0.01% error)
  • End of read: Q23-Q30 (0.5% - 0.1% error)

Coverage Calculation

For 100x coverage of 100,000 bp genome with 150 bp reads:

Total bases needed = 100,000 * 100 = 10,000,000
Bases per read pair = 300 (150 * 2)
Read pairs needed = 10,000,000 / 300 = 33,333

Generate ~33,500 read pairs per sample.

Output Files

Directory Structure

test_data/
├── reference.fasta
├── core_genes.fasta
├── sample_a_R1.fastq.gz
├── sample_a_R2.fastq.gz
├── sample_b_R1.fastq.gz
├── sample_b_R2.fastq.gz
├── sample_c_R1.fastq.gz
├── sample_c_R2.fastq.gz
└── variants.tsv

variants.tsv Format

Tab-separated file with ground truth:

sample	type	chrom	start	end	ref	alt	copy_number_a	copy_number_b
B	SNP	chr1	3500	3500	A	G	1	1
B	SNP	chr1	8200	8200	C	T	1	1
...
B	DUP	chr1	40000	42000	.	.	1	2
B	DEL	chr1	70000	71000	.	.	1	0
C	SNP	chr1	3500	3500	A	G	1	1
...

FASTA Format

Standard FASTA with 80-character line width:

>reference length=100000
ACGTACGTACGTACGT...
>core_gene_01 start=1000 end=1750
ACGTACGT...

FASTQ Format

Standard FASTQ with gzip compression:

@read_00001/1
ACGTACGTACGTACGT...
+
IIIIIIIIIIIIIII...

Read names encode:

  • Read number
  • /1 or /2 for paired end
  • Optional: source position for debugging

Validation Criteria

The test data should satisfy these properties:

Coverage Distribution

For single-copy regions:

  • Mean coverage: 100x ± 5x
  • Coverage should follow approximately Poisson distribution

Variant Detection

After running the full pipeline:

  1. SNPs: Each SNP should produce 2k-1 differential k-mers (k-mers overlapping the SNP position)

  2. Duplication:

    • ~2000 k-mers at ~2x copy number in B vs A
    • Should merge into single region of ~2000 bp

  3. Deletion:

    • ~1000 k-mers unique to A (absent in B)
    • Should merge into single region of ~1000 bp

Distance Matrix

Expected relative distances:

  • A-B: moderate (SNPs + dup + del)
  • A-C: moderate (SNPs + dup + del)
  • B-C: smaller (share some SNPs, different structural variants)

Reproducibility

All random processes must use fixed seeds:

Process Seed
Reference generation 42
SNP position selection 43
Read pair generation (A) 100
Read pair generation (B) 101
Read pair generation (C) 102
Error introduction Per-read: read_number + 1000

This ensures identical output across runs.