7
$\begingroup$

What tool could I use to check the ethnicity of a human bam or vcf file? I would like to use the results as a QC check to know whether a given sample or set of samples match the ethnicity information annotated in the meta data or not.

Seen so far:

The examples work, but starting from a GRCh38 .vcf file (from GATK UG) gives me an error:

Read    mytest.vcf   :1 Individual(s)Exception in thread "main" java.lang.NumberFormatException: For input string: ""
        at java.lang.NumberFormatException.forInputString(NumberFormatException.java:65)
        at java.lang.Integer.parseInt(Integer.java:592)
        at java.lang.Integer.valueOf(Integer.java:766)
        at Individual.get_Individual(Individual.java:73)
        at similarity.main(similarity.java:80)

It won't work with GATK UG .vcf files straightaway, it needs to reprocess .bam files using the ANGSD pipeline:

This tool does not give ethnicity checks but only genotype checks. Running all-vs-all vcfs vs bams can create a semi-matrix of similarities.

$\endgroup$

1 Answer 1

8
$\begingroup$

The main difficulty here is the use of GRCh38. Unfortunately, despite the fact that it's more than four years old, the major ethnicity-labeled public datasets (1000 Genomes, gnomAD when allele frequencies are enough) still aren't available for that reference. It is necessary to perform a liftover operation, or just use overlapping rsIDs and hope for the best.

Let's suppose you go with overlapping rsIDs, and a list of those rsIDs, one per line is in 'rsids.txt'. (Due to strand flips between reference builds, you may want to restrict this list so that all allele codes match, and it's also reasonable to throw out A/T and C/G SNPs.) Then, the following process would work:

  1. Download ADMIXTURE (https://www.genetics.ucla.edu/software/admixture/download.html ), plink (https://www.cog-genomics.org/plink/1.9/ ), and/or chr1-chr22 of 1000 Genomes phase 3 (ftp://ftp.1000genomes.ebi.ac.uk/vol1/ftp/release/20130502/ ) if you don't already have them.

  2. Use plink to extract the overlapping higher-MAF rsIDs from each 1000 Genomes phase 3 VCF ("plink --vcf ... --extract rsids.txt --maf 0.05 --make-bed --out ..."), and then merge the resulting per-chromosome filesets ("plink --merge-list ... --out merged_phase3_subset").

  3. ADMIXTURE prefers a dataset with around 100k variants. Unless you have a small list of overlapping rsIDs (in which case you had better use another method), you should have more than that. plink's LD-pruning function is a good way to select a subset for ADMIXTURE's use: "plink --bfile merged_phase3_subset --indep-pairwise 500kb 1 0.2; plink --bfile merged_phase3_subset --extract plink.prune.in --make-bed --out admixture_data"

Adjust the 0.2 threshold as necessary to keep the right number of variants.

  1. Run ADMIXTURE in unsupervised mode ("admixture admixture_data.bed 5 -j8"; adjust the -j parameter depending on the number of processor cores). This generates an admixture_data.5.P file with population allele frequencies, and an admixture_data.5.Q file with sample population assignments. Verify that the .Q file corresponds to the 1000 Genomes phase 3 'superpopulations'.

  2. Convert your data to plink-format if necessary, keeping just the overlapping rsIDs and 'downgrading' to GRCh37 coordinates. Make sure your variants are sorted in GRCh37-coordinate order, and the allele order also matches that in your reference dataset ("--a2-allele admixture_data.bim 6 2" during your final --make-bed operation will do the trick); otherwise ADMIXTURE won't do the right thing. Then run ADMIXTURE in projection mode ("cp admixture_data.5.P my_data.5.P.in; admixture -P my_data.bed 5"). my_data.5.Q will then have the ethnicity estimates you're looking for.

$\endgroup$

Your Answer

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

Not the answer you're looking for? Browse other questions tagged or ask your own question.