Changes between Initial Version and Version 1 of ChipBasedQcPipeline

Timestamp:

Sep 12, 2010 11:28:23 PM (14 years ago)

Author:

Morris Swertz

Comment:

--

ChipBasedQcPipeline

@@ +1 @@
+[[TOC()]]
+
+= Protocol for comparison between sequencing (VCF) and chip data =
+
+Created: Yurii Aulchenko, 2010.09.12
+
+Modified: (please fill in your name and date)
+
+Aims of application of the protocol:
+ * Establish custom pipeline for Chip-based QC.
+ * Check quality of sequence data. 
+ * Identify factors affecting quality of sequencing (e.g. batch effects). 
+ * Establish (preliminary) thresholds of quality metrics maximizing sensitivity and specificity of genotype calling.
+ * Using above thresholds, establish the false-positive and false-negative rates for variants discovered in our study (if we do not take trio structure into account). 
+ * Check if these rates are in agreement with theoretically expected (thus we do not miss any important experimental factor). 
+
+Below a workflow is provided. This document assumes VCF v.4 format (http://1000genomes.org/wiki/doku.php?id=1000_genomes:analysis:vcf4.0) is used; “+” strand is used in VCF. It is assumed that chip data come in ??? format (). This protocol does not specify exact tools to be used; as most of operations are data manipulations, it will be up to the involved analysts to decide what tool may be more convenient for them. The actual implementation of the workflow should allow reproducing the results or apply the same workflow to other data. This is not only important from good practice point of view, but also keeping in mind that more data will come to the same pipeline in the future. 
+
+== CHIP-VCF BUILD AND DBSNP MATCHING TABLE ==
+
+Before starting, check what is the reference build used by VCF files. Next, what is dbSNP version used to name the SNPs. If any of these (and/or strand is not always ‘+’ in chip data), arrange tab-separated file containing conversion table (name: chip_data_conversion_table_yyyy.mm.dd.txt):
+
+Header line:
+
+CHRC    POSC    SNPC    STRANDC A1C     A2C     CHRV    POSV    SNPV    A1V     A2V
+
+Next lines should all contain 12 tab-delimited values. There should be no missing values. 
+ * CHR[C/V]: chromosome SNP is located at, according to Chip [C] and VCF [V] build versions <integer from 1 to 22 for autosomes, “X” for X-chromosome, and “Y” for Y-chromosome>
+ * POS[C/V]: chromosomal position, according to Chip [C] and VCF [V] build versions <integer>
+ * SNP[C/V]: SNP rs-name, according to Chip [C] and VCF [V] dbSNP versions <alphanumeric>
+ * STRANDC: strand in chip annotation <single character, either “+” or “-“>
+ * A1C: first allele in chip annotation <single character, either “A”, “C”, “G” or “T”>
+ * A2C: second allele in chip annotation <single character, either “A”, “C”, “G” or “T”>
+ * A1V: (translated) first chip allele (A1C) according to “+” strand on VCF build <single character, either “A”, “C”, “G” or “T”>
+ * A2V: (translated) second chip allele (A2C) according to “+” strand on VCF build <single character, either “A”, “C”, “G” or “T”>
+
+Questions:
+ * Q: what is the build and dbSNP version used by chip and VCF?
+ * Q: how many SNPs changed the name in VCF build?
+ * Q: how many SNPs changed the strand in VCF build?
+ * Q: please provide a 2x2 table (name change/not) x (strand change/not)
+
+Note that in future, samples typed using a number of different Chip platforms will be coming in. Therefore above step should not assume a particular chip is used!
+
+== UPDATED CHIP GENOTYPES ==
+
+Using above described translation table, generate updated chip genotypes file (name: chip_genotypes_yyyy.mm.dd.txt)
+
+This is a tab-delimited text file containing a table. The header line is 
+
+ID      SNPV    QUALCHIP        A1VCHIP A2VCHIP GTCHIP
+
+Next lines should all contain 5 tab-delimited values. Use “.” (dot) for missing. 
+ * ID: sample ID (genotyped individual’s code) <alphanumeric>
+ * SNPV: SNP rs-name, according to VCF dbSNP version <alphanumeric>
+ * QUALCHIP: calling quality for the individual genotype
+ * A1VCHIP: first allele the personal genotype, translated according to “+” strand on VCF build <single character, either “A”, “C”, “G” or “T”>
+ * A2VCHIP: second allele the personal genotype, translated according to “+” strand on VCF build <single character, either “A”, “C”, “G” or “T”>
+GTCHIP: genotype with alleles in alphabetic order, <two characters, each either “A”, “C”, “G” or “T”>
+
+Questions:
+ * Q: do all SNPs in chip data have rs-number?
+ * Q: what alleles are observed in chip data? Only A/T/G/C?
+ * Q: are all SNPs bi-allelic?
+
+== EXTRACTION OF CHIP SNPS FROM VCF FILE ==
+
+From VCF, extract only lines containing SNPs also observed in the chip (see SNPV column of “chip_data_conversion_table_yyyy.mm.dd.txt”)
+
+Parse extracted lines, and arrange “Annotation” and “Genotypic” tables
+
+Annotation table (name: VCF_annotation_yyyy.mm.dd.txt). Tab-delimited file with header lines (and consequently extracting following columns from VCF):
+
+CHROM   POS     ID      REF     ALT     QUAL    FILTER  INFO
+
+At the beginning of the file, add meta-info from VCF file
+
+Genotypic table (name: VCF_genotypes_yyyy.mm.dd.txt). Tab-delimited file containing following information. Header line:
+
+ID      SNPV    GTVCF   GQ      DP      BATCH   ????
+
+Next lines should all contain XXX tab-delimited values. Use “.” (dot) for missing. 
+ * ID: sample ID (genotyped individual’s code) <alphanumeric>
+ * SNPV: SNP rs-name, according to VCF dbSNP version <alphanumeric>
+ * GTVCF: genotype with alleles in alphabetic order, <two characters, each either “A”, “C”, “G” or “T”>. This can be done by mapping the numbers provided in VCF GT field to REF and ALT and then ordering.
+ * GQ, DP: directly from VCF file
+BATCH …
+
+Merge chip and VCF genotypic tables (“chip_genotypes_yyyy.mm.dd.txt” and “VCF_genotypes_yyyy.mm.dd.txt”) using ID and SNPV as key variables. Keep all chip genotypes, substituting missing (“.”) when no information is available from VCF. Name the table “merged_chip_and_VCF_genotypes_yyy.mm.dd.txt”.
+
+Questions:
+ * Q: What is count and proportion of genotypes that do not match between GTCHIP and GTVCF? How much these counts/proportions changes if dropping rows with QUALCHIP < X (vary X)? How much these counts/proportions changes if dropping rows with GQ (DP) < X (vary X)?
+ * Q: What is proportion of false-positive and false-negative findings in our study, if we do not take trio structure into account?
+ * Q: Find out QC metrics thresholds maximizing specificity and sensitivity.
+
+Update the table with variable “CHIPVCFMISMATCH” (1 if mismatch, 0 for match, missing (“.”) if any is missing). 
+
+ * Q: Explore, which variables are significant predictors of mismatch using multiple logistic regression. 
+
+
+== CHIP SNPS MISSING FROM VCF ==
+
+Write the list of the chip SNPs not in VCF into the file “list_of_chip_snps_missing_in_VCF_yyyy.mm.dd.txt” (single column containing SNPV name). This should be done when matching chip SNPs with VCF SNPs (see section “EXTRACT CHIP SNPs FROM VCF”)
+
+ * Q: How many variants do we miss in VCF (how many SNPs in file list_of_chip_snps_missing_in_VCF_yyyy.mm.dd.txt)?
+
+For each SNP in list_of_chip_snps_missing_in_VCF_yyyy.mm.dd.txt, based on chip_genotypes_yyyy.mm.dd.txt derive frequency from chip data and arrange the following table (name: annot_chip_snps_missing_in_VCF_yyyy.mm.dd.txt). The header line should contain 
+
+SNPV     A1V     A2V     FREQA1V
+
+Each next line should contain 4 values delimited by tab; SNPV, A1V, and A2V explained above (the same as in “chip_data_conversion_table_yyyy.mm.dd.txt” file). FREQA1V is a floating-point frequency of allele “A1V”. 
+
+ * Q: Does the distribution of frequency of missed variants match the expected under the assumption that we miss at random because of limited #chromosomes and coverage (for each trio we read two chromosomes at 12x and 2 chromosomes at 24x) 
+