Skip to content

Windowed fst_wc disagrees with the scalar fst_weir_cockerham #133

Description

@nspope

Came up during #132

windowed_analysis(..., statistics=['fst_wc'], populations=[p1, p2]) routes to the fused CUDA kernel (_fused_windowed_twopop_kernel), which computes a haploid Weir & Cockerham FST (h_bar = 0, each haplotype treated as its own unit). The scalar divergence.fst_weir_cockerham is a diploid estimator: it pairs consecutive haplotypes into individuals and uses per-allele observed heterozygosity (the c component). These are different estimators, so windowed fst_wc does not match the scalar fst_weir_cockerham.

window: fused fst_wc (haploid)   scalar fst_weir_cockerham (diploid)
     0:  -0.000442               -0.000534      (mismatch)
 25000:   0.012111                0.014604      (mismatch)
 50000:   0.015716                0.017199      (mismatch)

Reproduce: request fst_wc per window via windowed_statistics_fused(hm, bp_bins=..., statistics=('fst_wc',), pop1=..., pop2=...) and compare each window to divergence.fst_weir_cockerham on the window subset.

If this is intentional, then the function needs to be corrected for multiallelic data; otherwise one of the scalar/windowed paths should be rewritten so they match.

Metadata

Metadata

Assignees

No one assigned

    Labels

    bugSomething isn't working

    Type

    No type

    Fields

    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions