I am trying to reproduce the results from this paper "Human genomic regions with exceptionally high levels of population differentiation identified from 911 whole-genome sequences".

Specifically I want to extract high differentiated ( or less) SNPs from 1k Genome data set representing high proportion of the population present.

The paper mentioned using Vcffixup from Vcflib (https://github.com/vcflib/vcflib) but there might still be other things I am missing.

  • Update 2 (MWE)

So this is my code snippet.

## to dowload datafile, from linux terminal, use curl or wget
curl -O https://1000genomes.s3.amazonaws.com/release/20110521/ALL.chr22.phase1_release_v3.20101123.snps_indels_svs.genotypes.vcf.gz
### unzip the file
gunzip dowloaded file
### take the first few lines from the vcf file.
### since actual file is huge, the name of sfile for testing is 
### myfile.vcf
head -200 downloadedfile > myfile.vcf 

Now the next section was done in python environment (jupyter notebook).

## using the library pysam
from pysam import VariantFile
import pandas as pd 
with VariantFile(test_datafile) as vcf_reader:
for record in vcf_reader:
    sample_file = record.samples.keys() 
    alleles = [record.samples[x].allele_indices for x in record.samples]
    break ## just to make sure my data is as intended

Now to the VCF tools part - vcffixup was used in the paper.

## Install vcflib
   sudo apt-get install libvcflib-tools libvcflib-dev
### export path
   export PATH=$PATH:/usr/lib/vcflib/bin 
### refresh enviroment variables
   source ~/.bashrc 
## run vcffixup
  vcffixup myfile.vcf

A sample of vcffixup output is hown below: enter image description here

So the questions questions I have is:

  1. how to compute derived allele frequency differences at all variants for pairs of continent and pairs of population within continent


  • 1
    $\begingroup$ Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking. $\endgroup$
    – Community Bot
    Aug 24, 2022 at 21:27
  • $\begingroup$ Thanks for the edit! But how do we choose the snps? And what do you mean by present? Present in which sample of the multi sample vcf file? And do you care about zygosity? I mean, the vcf is already what you want:it has 0/0 when a variant is absent. $\endgroup$
    – terdon
    Aug 25, 2022 at 23:06


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.