SVCFit

SVCFit is a fast and scalable computational tool designed to estimate the Structural Variant Cellular Fraction (SVCF) of inversions, deletions, tandem duplications, and translocations. Developed for the R environment, SVCFit integrates structural variant (SV) calls with Copy Number Variation (CNV) and Single Nucleotide Polymorphism (SNP) data to provide accurate cellular fraction estimates.

SVCFit also handles hemizygous (single-copy) chromosomes — for example chrX and chrY in a male subject — where heterozygous germline SNPs do not exist and the standard allele-copy-ratio correction cannot be estimated. On these chromosomes SVCFit switches to a ploidy-aware estimator that uses a read-depth–derived mean copy number (cn_bar) in place of SNP-based phasing. See Hemizygous chromosomes (chrX / chrY).

Resources

• Open access data: It is available on mendeley (doi: 10.17632/2nhhdjx225.6)

• Protected Data: Available via European Genome-phenome Archive (EGAD00001001343).

• Prostate mixture scripts: GitHub Repository

Installation

SVCFit is hosted on GitHub. You can install it directly within R using the remotes package.

Note: Installation requires a GitHub Personal Access Token (PAT) because the repository is hosted on GitHub.\

if (!requireNamespace("remotes", quietly = TRUE))
    install.packages("remotes")

# 1. Setup GitHub Credentials (if not already configured)
if (!requireNamespace("usethis", quietly = TRUE))
    install.packages("usethis")

# Create a token in your browser
usethis::create_github_token() 

# Store the token (paste when prompted)
credentials::set_github_pat()

# 2. Install SVCFit
remotes::install_github("KarchinLab/SVCFit", build_vignettes = TRUE, dependencies = TRUE)

Input Requirements

1. Structural variants

SVCFit accepts standard Variant Call Format (VCF) files. By default, the parser is optimized for VCFs produced by the SVTyper package [1].

CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	normal	tumor
chr1	1000	INV:6:0:1:0:0:0	T		100	PASS	END=1500;SVTYPE=INV;SVLEN=500;…	GT:PR:SR:…	0/1:76,0:70,0:…	0/1:76,0:70,0:…
chr2	5000	DEL:7:0:1:0:0:0	G		100	PASS	END=5300;SVTYPE=DEL;SVLEN=300;…	GT:PR:SR:…	0/1:76,0:70,0:…	0/1:76,0:70,0:…

Required INFO fields:

• SVTYPE (e.g. INV, DEL, DUP, BND)
• END

2. Copy number variants

SVCFit currently utilizes copy number calls from the FACETS package [2]. The tool uses allele-specific copy number, total copy number, and the cellular fraction (cncf) to annotate SVs with overlapping CNVs. No modification is required for standard FACETS output.

3. Heterozygous SNP near SV

SVCFit requires heterozygous SNP calls to phase overlapping CNVs. These can be generated using GATK4 [3] HaplotypeCaller and filtered using bcftools [4].

For computational efficiency, VCF can be filtered to include only heterozygous SNPs within 500bp of SV breakpoints.

The following code is not included in SVCFit and should be run separately.

# 1. Call SNPs using GATK
gatk --java-options "-Xmx4g" HaplotypeCaller \
    -R $ref \
    -I $normal_bam \
    -O $snp_dir/SNP.vcf.gz

# 2. Filter for heterozygous SNPs using bcftools
bcftools view -v snps -g het -Oz -o $snp_dir/het_snp.vcf.gz $snp_dir/SNP.vcf.gz
tabix -p vcf $snp_dir/het_snp.vcf.gz

4. Heterozygous SNP on SV supporting reads

To infer SV phasing, SVCFit specifically examines heterozygous SNPs found on reads that support the structural variant. This is done following the steps: 1) Extract SV-supporting reads from the BAM file using samtools [4]. 2) Generate a read pileup using bcftools restricted to the heterozygous SNP positions identified in step 3.

The following code is not included in SVCFit and should be run separately.

# Extract SV supporting reads
samtools view -f 1 -F 2 -b $tumor_bam > $snp_dir/sup_$samp_name.bam
samtools index $snp_dir/sup_$samp_name.bam

# Generate pileup at known heterozygous sites
# Note: pos_$samp_name.bed should contain the positions from het_snp.vcf.gz above
bcftools mpileup -f $ref -a DP,AD -A \
    -R $snp_dir/pos_$samp_name.bed \
    $snp_dir/sup_$samp_name.bam -Ov > $snp_dir/het_on_sv_$samp_name.vcf

Example Preprocessing Pipeline

The commands below provide a complete worked example for generating SVCFit inputs from raw BAM files. Each tumor-versus-normal biopsy pair is processed as a separate sample. Tool versions used in the manuscript are listed in the Tool Versions table.

The expected input layout:

sample/
├── tumor.bam        (+ .bai)
├── normal.bam       (+ .bai)
├── ref.fa           (+ .fai, .dict)
└── intervals.bed    (optional)

1. Trim and align reads (if starting from FASTQ)

trim_galore --paired tumor_R1.fq.gz tumor_R2.fq.gz -o trimmed/
bwa mem -t 8 ref.fa trimmed/tumor_R1_val_1.fq.gz trimmed/tumor_R2_val_2.fq.gz \
  | samtools sort -@ 4 -o tumor.bam -
samtools index tumor.bam
# repeat for the matched normal

2. Mark duplicates and recalibrate (GATK4[3])

gatk MarkDuplicates -I tumor.bam -O tumor.md.bam -M tumor.metrics.txt
gatk BaseRecalibrator -I tumor.md.bam -R ref.fa --known-sites known.vcf.gz -O tumor.bqsr.table
gatk ApplyBQSR        -I tumor.md.bam -R ref.fa --bqsr-recal-file tumor.bqsr.table -O tumor.recal.bam
# repeat for the matched normal

3. Somatic SV calling (Manta[5])

configManta.py --tumorBam tumor.recal.bam --normalBam normal.recal.bam \
               --referenceFasta ref.fa --runDir manta_run/
manta_run/runWorkflow.py
# output: manta_run/results/variants/somaticSV.vcf.gz

4. SV genotyping (SVtyper[1])

# Ensure CIPOS and CIEND INFO fields are present; SVtyper requires them.
svtyper -B tumor.recal.bam -i tumor.vcf -o tumor.gt.vcf

The output VCF has AO (SV-supporting reads) and RO (reference reads) per breakpoint, read by SVCFit as BPC and BEC after scaling by mean read depth.

For multi-caller consensus (recommended for clinical samples), merge with SURVIVOR:

# vcf_list.txt contains one column where each row is the path to VCF from a SV caller 
SURVIVOR merge vcf_list.txt 500 1 1 1 0 30 consensus.vcf

5. Germline heterozygous SNP detection (GATK4 + bcftools)

# 1. Call SNPs using GATK
gatk --java-options "-Xmx4g" HaplotypeCaller \
    -R $ref \
    -I $normal_bam \
    -O $snp_dir/SNP.vcf.gz

# 2. Filter for heterozygous SNPs using bcftools (can be filtered for computation efficiency)
bcftools view -v snps -g het -Oz -o $snp_dir/het_snp.vcf.gz $snp_dir/SNP.vcf.gz
tabix -p vcf $snp_dir/het_snp.vcf.gz

6. Allele-specific copy-number profile (FACETS[2])

snp-pileup -g -q15 -Q20 -P100 -r25,0 normal.het.vcf.gz tumor.snp.csv \
           normal.recal.bam tumor.recal.bam

library(facets)
rcmat = readSnpMatrix(SNP_file)
xx=preProcSample(rcmat,ndepth=20, gbuild="hg38", cval=50)
oo = procSample(xx, cval=150, dipLogR=NULL)
ooo= procSample(xx, cval=500, dipLogR=oo$dipLogR)
fit = emcncf(ooo)
write.table(fit$cncf, file = "tumor.facets.tsv", sep = "\t", quote = FALSE, row.names = FALSE)

7. Heterozygous SNPs pileup on SV supporting reads

# 1. Extract SV supporting reads
samtools view -f 1 -F 2 -b $tumor_bam > $snp_dir/sup_$samp_name.bam
samtools index $snp_dir/sup_$samp_name.bam

# 2. Generate pileup at known heterozygous sites that's on SV supporting reads
# Note: pos_$samp_name.bed should contain the positions from het_snp.vcf.gz above
bcftools mpileup -f $ref -a DP,AD -A \
    -R $snp_dir/pos_$samp_name.bed \
    $snp_dir/sup_$samp_name.bam -Ov > $snp_dir/het_on_sv_$samp_name.vcf

Caller-specific quirks

• Manta[5]: writes INFO/CIPOS and INFO/CIEND natively — no reformatting needed.
• Delly[6]: writes INFO/CIPOS but not INFO/CIEND — append CIEND=-50,50 (or the per-call confidence interval if available).
• GRIDSS[7]: uses paired-end / split-read counts in dedicated INFO fields; convert to LUMPY-style INFO/MATEID, INFO/CIPOS, and INFO/CIEND before SVtyper.
• Manta + Delly + GRIDSS consensus via SURVIVOR[8]: SURVIVOR drops CIPOS/CIEND from some merged records — backfill them before genotyping.

Usage Workflow

The SVCFit pipeline consists of three main steps: Extraction, Characterization, and Calculation.

1. Extract SV and SNP Information — extract_info()

Load and preprocess input VCF and CNV files. This step parses metadata, processes breakends (BND), and handles heterozygous SNPs. extract_info() internally performs:

1. Load input data — load_data()
2. Process BND events — proc_bnd()
3. Parse SV metadata — parse_sv_info()
4. Parse heterozygous SNPs — parse_het_snps(), parse_snp_on_sv()

info <- extract_info(
  p_het = "path/to/het_snps.vcf",
  p_onsv = "path/to/snps_on_sv.vcf",
  p_sv = "path/to/structural_variants.vcf",
  p_cnv = "path/to/cnv_file.txt",
  chr_lst = NULL,
  flank_del = 50, 
  QUAL_thresh = 100, 
  min_alt = 2, 
  tumor_only = FALSE
)

Function Arguments

extract_info()

Argument	Type	Default	Description
p_het	Character	—	Path to VCF of heterozygous SNPs.
p_onsv	Character	—	Path to VCF of SNPs overlapping SV-supporting reads.
p_sv	Character	—	Path to SV VCF.
p_cnv	Character	—	Path to CNV file.
chr_lst	Character	NULL	Chromosomes to include.
flank_del	numeric	50	Max distance to consider deletion overlapping a BND.
QUAL_thresh	numeric	100	Minimum QUAL score.
min_alt	numeric	2	Minimum alternative reads.
tum_only	Logical	—	Whether SVs come from tumor-only calling.

Output: A list of data frames containing parsed SV + SNP information.

2. Annotate SVs Using CNV and SNP Information — characterize_sv()

This step integrates CNV and heterozygous SNPs to infer phasing, zygosity, and overlapping CNV. characterize_sv() internally performs:

1. Assign SV IDs to SNPs — assign_svids()
2. Summarizes phasing + zygosity — sum_sv_info()
3. Assign CNV to SV — assign_cnv()
4. Annotate overlapping CNV — annotate_cnv(), parse_snp_on_sv()

sv_char <- characterize_sv(
  sv_phase = info$sv_phase, 
  sv_info = info$sv_info, 
  cnv = info$cnv,
  flank_snp = 500,
  flank_cnv = 1000
)

Function Arguments

characterize_sv()

Argument	Type	Default	Description
sv_phase	data.frame	—	Phasing/zygosity from SNPs.
sv_info	data.frame	—	Parsed SV metadata.
cnv	data.frame	—	CNV data.
flank_snp	numeric	500	Max assignment distance for SNPs.
flank_cnv	numeric	1000	Max assignment distance for CNVs.

3. Calculate SVCF for Structural Variants — calc_svcf()

This step computes the Structural Variant Cellular Fraction (SVCF). and returns an annotated VCF file in data.frame format.

svcf_out <- calc_svcf(
  anno_sv_cnv = sv_char$anno_sv_cnv,
  sv_info     = sv_char$sv_info,
  thresh      = 0.1,
  samp        = "SampleID",
  exper       = "ExperimentID"
)

Function Arguments

Argument	Type	Default	Description
anno_sv_cnv	data.frame	—	CNV-annotated SVs.
sv_info	data.frame	—	Parsed SV info.
thresh	numeric	0.1	Threshold for SV-before-CNV inference.
samp	character	—	Sample name.
exper	character	—	Experiment name.
hemizygous_chr	character	NULL	Chromosomes single-copy in the germline (e.g. c("chrX","chrY")). NULL = diploid-only behavior.
hemi_cn_bar	data.frame / numeric	NULL	Read-depth mean copy number (cn_bar) for hemizygous rows — a data.frame(CHROM, POS, cn_bar) or numeric vector. Without it (or for an unmatched row) hemizygous duplications are left unresolved; deletions, inversions, and translocations fall back to SVCF = VAF.
hemi_bg_cn	data.frame / numeric	NULL	Flanking (background) copy number for the CNV-first deletion form; optional.
hemi_dup_r	numeric	2	Copies in carrier cells for a hemizygous tandem duplication; SVCFs reported as upper bounds.
zero_ref_allowlist	data.frame	NULL	Hemizygous rows with sv_ref = 0 that BAM evidence confirms are genuine clonal losses (SVCF = VAF = 1).

Output: An annotated VCF-like data frame with additional fields for VAF, Rbar, r, and SVCF. On hemizygous chromosomes it also carries pl (local normal ploidy), svcf_status, sv_cnv_order, and svcf_is_bound.

1. VAF: variant allele frequency
2. Rbar: average break interval count in a sample
3. r: inferred integer copy number of break intervals
4. SVCF: structural variant cellular fraction.

4. Build tumor evolution tree — build_tree()

This step build the tumor evolutionary tree based on SV clusters obtained from Dirichlet process Gaussian Mixture Model (DP-GMM).Currently, this step is optimized for two sample longitudinal data.

output <- cluster_data(
  pair_path,
  pur_path,
  data_dir,
  pair_num = 1)
clone2=output[[3]]

build_tree(
  clones,
  lineage_precedence_thresh=0.2, 
  sum_filter_thresh=0.2)

Function Arguments

cluster_data()

Argument	Type	Default	Description
pair_path	character	—	Path to a tab-separated file with columns for ‘pre_BAT sample’ and ‘on_BAT sample’.
pur_path	character	—	Path to a tab-separated file with columns for ‘sample’ and ‘purity’.
data_dir	character	—	Path to the directory containing per-sample SVCF output files.
Kmax	numeric	10	Maximum number of clusters for DP-GMM.
n_steps	numeric	100	Number of DP-GMM iterations.
thr_min_w	numeric	0.01	Minimum cluster weight threshold.
random_state	integer	0	Random seed for reproducibility.
concentration	numeric	1	Dirichlet concentration parameter.
min_n	numeric	5	Minimum cluster size to retain.
min_dist	numeric	0.2	Minimum distance for merging nearby clusters.
pair_num	numeric	1	The identifier (index or ID) for the specific sample pair (patient) being analyzed.
pairs	integer vector	NULL	Subset of pair IDs to process; defaults to all pairs.
exclude_pairs	integer vector	integer(0)	Pair IDs to exclude from analysis.
deduplicate	Logical	TRUE	Whether to deduplicate events before clustering.
ccf_floor	numeric	0.1	Minimum CCF value before flooring.

build_tree()

Argument	Type	Default	Description
clones	data.frame	—	SV clustering result.
lineage_precedence_thresh	numeric	0.2	Maximum violation of lineage precedence rule.
sum_filter_thresh	numeric	0.2	Maximum violation of sum condition rule.
linear_penalty	numeric	0	Penalty applied to linear (chain) topologies during tree scoring.

Output: A tumor evolutionary tree rooted at the germline (G). Node numbers correspond to SV cluster numbers. The branching depicts the chronological occurrence of SV clusters.

5. Simulation & Benchmarking

SVCFit includes utility functions for processing simulation data from VISOR and attaching “ground truth” labels to structural variants for benchmarking.

5.1 read clonal assignment

truth <- load_truth(
  truth_path = "path/to/truth_beds", 
  overlap = FALSE
  )

This function has 1 arguments:

Argument	Type	Default	Description
truth_path	Character	N/A	Path to BED files storing true structural variant information with clonal assignment. Each BED file should be named like "c1.bed, c2.bed", etc for non-overlapping simulations and "c11.bed, c22.bed", etc for overlapping simulations. Structural variants should be saved in separate BED files if they belong to different (sub)clones.
overlap	Logical	FALSE	Whether the simulation has SV-CNV overlap.

The file path should follow this structure:

root/
├── true_clone/
│   ├── c1.bed/
│   ├── c2.bed/
│   ├── c3.bed/
│   └── .../

Parent nodes should always have lower number in name than its children (i.e. c1.bed instead of c3.bed) and all child node bed file should conatin its ancestors mutations.

5.2 attach clonal assignment to output

svcf_truth <- attach_truth(svcf_out, truth)

This function has 2 arguments:

Variable	Type	Default	Description
svcf_out	DataFrame	N/A	The output from calc_svcf
truth	DataFrame	N/A	Stores the clone assignment for each structural variant designed in a simulation.

This appends the known clonal assignment to the calculated SVCF output for performance evaluation.

Hemizygous chromosomes (chrX / chrY)

On a diploid autosome SVCFit uses heterozygous germline SNPs to estimate the allele copy ratio and phase each SV against overlapping CNVs. A hemizygous chromosome — a male X or Y, or any chromosome that is single-copy in the germline — has no heterozygous SNPs, so that route is undefined and the diploid conversion factor of 2 no longer applies. Left uncorrected, a clonal hemizygous SV is estimated at roughly twice its true cellular fraction and is then silently dropped downstream.

SVCFit handles these chromosomes with a ploidy-aware estimator. The key quantity is cn_bar, the mean copy number of the locus per cell, measured from read depth (cn_bar = R * psi_sample / 2, where R is the tumor/normal depth ratio and psi_sample the autosomal ploidy scaling). Total alleles per cell is always cn_bar, and with VAF = BPC / (BPC + BEC) the two orderings of an SV relative to an overlapping CNV are:

H1  (SV precedes CNV):  SVCF = cn_bar * VAF - (cn_bar - 1)
H2  (CNV precedes SV):  SVCF = cn_bar * VAF

Both reduce to SVCF = VAF at cn_bar = 1 (copy-neutral). The ordering is decided by the sign of the H1 form — no integer copy number or CNV cellular fraction is needed. A tandem duplication is its own copy-number change and uses SVCF = (cn_bar - 1) / (r - 1); because r is not identifiable from a single locus, these are reported as upper bounds.

Enabling it

Pass the single-copy chromosomes to run_svcfit() (or calc_svcf()), together with a per-SV cn_bar table from a depth segmentation of the chromosome:

result <- run_svcfit(
  p_het = p_het, p_onsv = p_onsv, p_sv = p_sv, p_cnv = p_cnv,
  samp = "SampleID", exper = "ExperimentID",
  hemizygous_chr = c("chrX", "chrY"),   # single-copy in this subject's germline
  hemi_cn_bar    = cn_bar_table         # data.frame(CHROM, POS, cn_bar), from read depth
)

With hemi_cn_bar = NULL, hemizygous SVs are treated as copy-neutral and resolve to SVCF = VAF — except tandem duplications, which need cn_bar and are left unresolved (final_svcf = NA, flagged hemizygous_dup_needs_cn_bar in svcf_status). A deletion or inversion on a genuinely copy-altered segment is not detectable without cn_bar on a hemizygous chromosome (there are no heterozygous SNPs to flag it), so it too falls back to SVCF = VAF rather than being singled out. Autosomes are unaffected: with hemizygous_chr = NULL every result is identical to the diploid path.

Obtaining cn_bar (chrX/chrY depth segmentation with DNAcopy)

SVCFit does not compute cn_bar itself — it consumes it. On the autosomes copy number comes from FACETS, but FACETS cannot fit a hemizygous chromosome: a male X has no heterozygous germline SNPs to supply the allelic-imbalance signal FACETS relies on, so it returns a single whole-chromosome segment with an inflated, unusable total copy number (tcn.em of 3–9 in our cohort, contradicted by read depth). The containing SV’s own SVCFit estimate cannot be used either — it is hemizygous and copy-altered, so that would be circular.

cn_bar is therefore measured directly from read depth, segmented with the Bioconductor DNAcopy package (circular binary segmentation, CBS). The recipe:

1. Bin chrX read depth (e.g. 100 kb bins) in both the tumor and the matched-normal BAM (samtools depth).
2. Normalize each BAM by its own autosomal 2-copy baseline (the median depth of a few known-diploid autosomal regions), then take the tumor/normal ratio. Because the matched normal is also single-copy on chrX, this ratio is the tumor’s mean copy number per cell relative to the germline.
3. Scale by psi_sample / 2 to recover absolute copy number, where psi_sample = purity * psi_tumor + (1 - purity) * 2 is the sample mean autosomal ploidy (from the FACETS autosomal fit, which is trustworthy). This step is not optional: across our cohort psi_sample / 2 ran from 0.87 to 1.72.
4. Segment log2(cn_bar) per bin with DNAcopy and report each segment’s mean as cn_bar = 2^seg.mean.

library(DNAcopy)
# cn_bar_bin : per-bin (tumor/normal depth ratio) * (psi_sample / 2)
# pos        : bin start positions on the hemizygous contig
cna <- CNA(log2(cn_bar_bin), rep("chrX", length(pos)), pos,
           data.type = "logratio", sampleid = "chrX")
seg <- segment(smooth.CNA(cna), alpha = 0.01, min.width = 2)$output
seg_cn_bar <- 2^seg$seg.mean          # mean copy number per cell, per segment

cn_bar is the final product and is never rounded to an integer — the corrected hemizygous forms take it directly (there is no integer copy number c or CNV cellular fraction f_CNV to solve for; read depth alone cannot separate them on a hemizygous locus). A segment at true germline copy number reads near 1.0; if none does, either the chromosome is wholly altered or psi_sample is wrong.

DNAcopy’s own defaults split a flat chrX into dozens of noise segments on real data; the alpha, min.width, and undo.SD knobs control that and should be tuned against a known copy-neutral control before use.

The format SVCFit actually needs — DNAcopy is replaceable

SVCFit never sees DNAcopy or its segment objects. It consumes a plain per-SV table, so any method that yields a mean copy number per cell — a different segmenter (CBS, PSCBS, HMMcopy, GATK ModelSegments, CNVkit, Battenberg/ASCAT on the autosome-analogous signal), or even a single tumor/normal depth ratio per breakpoint — works as long as the output is coerced into this shape:

hemi_cn_bar — data.frame with exactly these columns:

Column	Type	Meaning
CHROM	character	Contig name, matching the SV’s CHROM exactly (e.g. "chrX", or "X").
POS	integer	The SV breakpoint position, matching the SV’s POS exactly.
cn_bar	numeric	Mean copy number of that locus per cell (not rounded to an integer).

The lookup is an exact (CHROM, POS) match against the SV rows — one row per SV, keyed on the breakpoint, not per segment. That is the one manual step: if your source produces segment intervals (chrom, start, end, cn_bar), assign each SV the cn_bar of the segment spanning its POS before passing the table. An SV whose breakpoint has no matching row gets no depth correction — a tandem duplication is then left unresolved (final_svcf = NA, flagged hemizygous_dup_needs_cn_bar), while a deletion, inversion, or translocation falls back to the copy-neutral SVCF = VAF. A bare numeric vector is also accepted, but it is recycled across all SV rows in internal order (not just the hemizygous ones), so in practice only a single value — the same cn_bar for every hemizygous SV — is reliable; use the keyed data frame otherwise.

# minimal hemi_cn_bar — however you produced cn_bar, this is all SVCFit reads
hemi_cn_bar <- data.frame(
  CHROM  = c("chrX",     "chrX"),
  POS    = c(31200000L,  67500000L),   # exact SV breakpoint positions
  cn_bar = c(1.02,       1.54)         # mean copies/cell at each locus
)

The optional hemi_bg_cn (flanking/background copy number for the CNV-first deletion form) uses the same convention with columns CHROM, POS, bg_cn. DNAcopy is only needed for this hemizygous depth step — it is not a dependency of the autosomal pipeline.

Estimator helpers

These building blocks are exported so a cn_bar and read counts can be scored directly, without the full VCF pipeline:

Function	Purpose
local_ploidy(chrom, hemizygous_chr)	Local normal ploidy (1 on a hemizygous chromosome, 2 otherwise).
classify_cn(cna, minor, pl)	Ploidy-aware DUP / norm / DEL classification (reduces to the diploid tests at pl = 2).
resolve_hemizygous_svcf(bpc, bec, cn_bar)	SVCF for an SV inside a CNV made by another event; picks H1/H2 by the sign rule.
hemizygous_dup_svcf(cn_bar, r = 2)	SVCF for a hemizygous tandem duplication (upper bound).
hemizygous_del_svcf(bpc, bec, cn_bar)	SVCF for a hemizygous deletion, with a depth-vs-read consistency check.
svcf_status(pl, cn_type, sv_ref)	Per-row status label so exclusions are countable, not silent.

resolve_hemizygous_svcf(bpc = 4, bec = 1, cn_bar = 2.0)
#>   svcf      ordering status  h1  h2
#> 1  0.6 sv_before_cnv     ok 0.6 1.6

hemizygous_dup_svcf(cn_bar = 1.5, r = 2)
#>   svcf is_upper_bound status
#> 1  0.5           TRUE     ok

Tutorial

library(SVCFit)
vignette("SVCFit_guide", package = "SVCFit")

Frequently Asked Questions

Can I use a different SV caller? Yes. Any caller that produces a per-SV breakpoint VCF compatible with SVtyper will work. Differences across callers in breakpoint detection sensitivity, split-read versus discordant-pair definitions, and quality filtering may produce different REF/ALT counts and therefore different SVCF estimates from the same data — see Discussion in the manuscript.

Can I skip FACETS and use Battenberg / ASCAT instead? Yes, as long as you provide per-segment total copy number and gain/loss classification in the same TSV format.

Do I need a matched normal? For the COMBAT analysis we used matched normals throughout. SVCFit can run on tumor-only when a matched normal is unavailable, but tumor purity must be supplied externally and ASCN inference becomes less reliable without germline-heterozygous SNP calls; this is treated as an unsupported configuration in the current release.

What about complex SVs (chromothripsis, BFB, chromoplexy)? The closed-form SVCF estimators cover deletions, tandem duplications, inversions, and the three classes of translocations. Multi-breakpoint complex SVs are not yet explicitly modeled — the per-breakpoint estimates are still produced, but their interpretation as a single cellular fraction is approximate. Future releases will add structure-aware handling.

Tool Versions

Tool versions used in the manuscript:

Step	Tool	Version
Read trimming	trim_galore	v0.6.1
Alignment	bwa mem	v0.7.19
MarkDuplicates / BQSR	GATK4	v4.6.2.0
Somatic SV (single-caller)	Manta	v1.6.0
Somatic SV (multi-caller)	Manta + Delly v1.5.0 + GRIDSS v2.13.2	—
Multi-caller merge	SURVIVOR	v1.0.7
SV genotyping	SVtyper	v0.7.1
Germline SNP calling	GATK4 HaplotypeCaller	v4.6.2.0
Het-SNP filter	bcftools	v1.20
Copy number	FACETS	v0.6.2
Copy number (hemizygous chrX/chrY)	DNAcopy (circular binary segmentation)	Bioconductor
Breakpoint read filter	samtools	v1.21
SNP pileup	bcftools mpileup	v1.20

Reference

1. Chiang, C. et al. SpeedSeq: ultra-fast personal genome analysis and interpretation. Nat Methods 12, 966–968 (2015).
2. Shen, R. & Seshan, V. E. FACETS: allele-specific copy number and clonal heterogeneity analysis tool for high-throughput DNA sequencing. Nucleic Acids Res 44, e131–e131 (2016).
3. Van Der Auwera, G. A. & O’Connor, B. D. Genomics in the Cloud. (O’Reilly Media, 2020).
4. Danecek, P. et al. Twelve years of SAMtools and BCFtools. GigaScience 10, giab008 (2021).
5. Chen, X. et al. Manta: rapid detection of structural variants and indels for germline and cancer sequencing applications. Bioinformatics 32, 1220–1222 (2016).
6. Rausch, T. et al. DELLY: structural variant discovery by integrated paired-end and split-read analysis. Bioinformatics 28, i333–i339 (2012).
7. Cameron, D. L. et al. GRIDSS2: comprehensive characterisation of somatic structural variation using single breakend variants and structural variant phasing. Genome Biol 22, 202 (2021).
8. Jeffares, D. C. et al. Transient structural variations have strong effects on quantitative traits and reproductive isolation in fission yeast. Nat. Commun. 8, 14061 (2017).