Genome Analysis Wiki - User contributions [en]

SeqShop: June 2014

2014-06-17T03:25:44Z

Gjun:

Sequence Analysis Workshop June 16-20, 2014

{| class="wikitable" cellpadding=5
|-
! style="background:#ABE" | Monday Lecture !! style="background:#ABE" | Practical
|-
| [[Media:2014.06_-_Abecasis_Workshop_Introduction.pdf|Design of Human Genetic Studies]]
| [[SeqShop: Design of Human Genetic Studies Practical]]
|-
| [[Media:Mapping and assembly.pdf|Sequence Mapping and Assembly]]
| [[SeqShop: Sequence Mapping and Assembly Practical]]
|-
! style="background:#ABE" | Tuesday Lecture !! style="background:#ABE" | Practical
|-
| [[Media:SNPcalling.pdf|Variant Calling and Filtering for SNPs]]
| [[SeqShop: Variant Calling and Filtering for SNPs Practical]]
|-
| Variant Calling and Filtering for INDELs
| [[SeqShop: Variant Calling and Filtering for INDELs Practical]]
|-
! style="background:#ABE" | Wednesday Lecture !! style="background:#ABE" | Practical
|-
| Analysis of Structural Variation
| [[SeqShop: Analysis of Structural Variation Practical]]
|-
! style="background:#ABE" | Thursday Lecture !! style="background:#ABE" | Practical
|-
| Estimates of Genetic Ancestry
| [[SeqShop: Estimates of Genetic Ancestry Practical]]
|-
| Genetic Association Analysis
| [[SeqShop: Genetic Association Analysis Practical]]
|}

File:SNPcalling.pdf

2014-06-17T03:24:37Z

Gjun:

VerifyIDintensity

2014-04-24T19:34:16Z

Gjun: /* Download verifyIDintensity */

[[Category:Software|VerifyIDintensity]]

'''verifyIDintensity''' is a software that detects and estimates sample contamination using intensity data from Illumina genotyping arrays using a mixture model.

== Download verifyIDintensity ==

* [https://github.com/gjun/verifyIDintensity GitHub repository]

== Build verifyID intensity ==

To build verifyIDintensity, run the following series of commands. You need [http://www.boost.org/ boost library] and [http://tclap.sourceforge.net/ tclap].

$ tar xzvf verifyIDintensity.tgz
$ make

== Basic Usage ==

verifyIDintensity [-t <float>] [-m <int>] -n <int> [-b <string>] [-s <string>] -i <string> [-v] [-p] [--] [--version] [-h]

== Options ==

-t <float>, --threshold <float>
Minimum allele frequency for likelihood estimation, default is 0.01

-m <int>, --marker <int>
(required) Number of markers

-n <int>, --number <int>
(required) Number of samples

-b <string>, --abf <string>
Allele frequency file (ABF), which is a plain text file with SNP_ID and Allele_B frequency.
SNP_IDs should be sorted in the same order as the intensity file

-s <string>, --stat <string>
Statistics file (created if not exist)

-i <string>, --in <string>
(required) Input pre-computed intensity (.adpc.bin) file

-v, --verbose
Turn on verbose mode

-p, --persample
Do per-sample analysis, default is per-marker analysis

--, --ignore_rest
Ignores the rest of the labeled arguments following this flag.

--version
Displays version information and exits.

-h, --help
Displays usage information and exits.

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== For sequence data ==

[[VerifyBamID]] software can estimate sample contamination from aligned sequence reads and population minor allele frequency

VerifyBamID

2013-02-26T05:55:02Z

Gjun: /* For array intensity data */

[[Category:Software|VerifyBamID]]

'''verifyBamID''' is a software that verifies whether the reads in particular file match previously known genotypes for an individual (or group of individuals), and checks whether the reads are contaminated as a mixture of two samples. '''verifyBamID''' can detect sample contamination and swaps when external genotypes are available. When external genotypes are not available, '''verifyBamID''' still robustly detects sample swaps.

== Download verifyBamID ==

To get a copy go to the [http://www.sph.umich.edu/csg/kang/verifyBamID/download VerifyBamID Download] download page.

== Join in verifyBamID mailing list ==

Please join in the [http://groups.google.com/group/verifybamid VerifyBamID Google Group] to ask / discuss / comment about verifyBamID.

== What's new ==

(2012/06/20)
* Fixed a bug of incorrect estimate of contamination when --chip-full option was used (Thanks to Richard Smith)
* Fixed a bug of incorrect per-readgroup output in --chip-* parameter

(2012/05/24)
* Fixed a bug of incorrect per-readgroup output (Thanks to Matthew Flickinger)
* '''(IMPORTANT)''' Add an option to remove either side of overlapping fragment. This option is turned on by default, and can be turned off usig --ignoreOverlapPair. If your sequence data has very short insert size, this update may increase the sensitivity of estimated contamination.
* Changes in the directory structure and Makefile

(2012/05/18) The new release of verifyBamID have undergone major change since the last version (as of 2011 April). Here are the highlights
* The genotype / allele frequency file is now based on VCF format rather than PLINK format.
* The reference sequence information is no longer required
* Uses Brent's method for precise estimation of contamination parameters
* Generate the depth distribution statistics.
* Estimated reference-bias parameters (useful mostly for ABI SOLiD sequence data)

== Build verifyBamID ==

The binary download of verifyBamID is available. You may use the version in Ubuntu 64-bit platform. To build verifyBamID, download the statgen library and run the following series of commands
tar xzvf verifyBamID.20120620.tar.gz
cd verifyBamID
make cloneLib (in the case ../libStatGen does not exist)
make
./bin/verifyBamID

Note that '''git clone''' command will create a directory ./libStatGen under your working directory, and '''make''' will create binary of verifyBamID under verifyBamID/bin/

verifyBamID is designed to be reasonably portable.

However, since development occurs only on Ubuntu 9.10 x86 and x64 platforms, and later, there are likely other portability issues.

Currently we support verifyBamID only on Ubuntu 9.10 and later on 64-bit processors.

== Basic Usage ==

A key step in any genetic analysis is to verify whether data being generated matches expectations. ''verifyBamID'' checks whether reads in a BAM file match previous genotypes for a specific sample. In addition, it detects possible sample mixture from population allele frequency only, which can be particularly useful when the genotype data is not available.

Using a mathematical model that relates observed sequence reads to an hypothetical true genotype, ''verifyBamID'' tries to decide whether sequence reads match a particular individual or are more likely to be contaminated (including a small proportion of foreign DNA), derived from a closely related individual, or derived from a completely different individual.

== Basic Usage Example ==

Here is a typical command line:

verifyBamID --vcf [input.vcf] --bam [input.bam] --out [output.prefix] --verbose --ignoreRG

where
[input.bam] is a BAM (Binary Alignment Map) file of a sequence reads
[input.vcf] is input VCF file containing individual genotypes or AF or AC/AN fields in the INFO field. gzipped VCF is also allowed.
[outPrefix] is output prefix of output files - [outPrefix].{selfRG,selfSM,bestRG,bestSM,depthRG,depthSM} will be created.

More detailed description of command line input is below

== Preparing input files ==

verifyBamID requires two input files - VCF file containing external genotypes or allele frequency information, and the BAM file.

=== VCF input genotype file ===

The input VCF file contains (1) external genotype information and/or (2) allele frequency information as AF entry or AC/AN entries in the INFO field. (See [http://www.1000genomes.org/wiki/Analysis/Variant%20Call%20Format/vcf-variant-call-format-version-41 | VCF specification] for further details). If neither information is provided, verifyBamID will not work properly.

If external genotype information is provided, sequence+array method will identify contamination and sample swaps by comparing the concordance between the external genotypes and the sequence reads. Additionally, sequence-only method will provide additional contamination estimates by modeling the sequence reads as mixture of two unknown samples based on the allele frequency information in the VCF file.

Input VCF file needs to meet several additional contraints need to meet in order to properly run verifyBamID.
* The VCF is assumed to be well-formed. For example, verifyBamID does not check whether REF allele actually matches with reference sequence.
* The VCF should only contain SNPs. Current version of verifyBamID does not accept INDELs, MNPs, Structural Variations, or other complex variants.
* The individual IDs in the VCF file, must be identical with the individual identifier in the BAM file. Otherwise, --smID option can override the sample ID information of the BAM file to the ID that matches to the individual IDs in the VCF file.
* IMPORTANT : For targeted sequencing data, it is important to subselect the markers to only include on-target markers in the genotype file. Off-target markers are not likely to have multiple non-duplicated reads at the marker position, and it may create artifacts in the analysis due to overlapping fragments.
* Currently, verifyBamID takes only autosomal chromosomes as input VCF.

An example input VCF file (without external genotype) is provided below. Note that AC and AC entries exists in the INFO field for the allele frequency information.

#CHROM POS ID REF ALT QUAL FILTER INFO
20 61651 SNP20-9651 C A . PASS CR=99.86851;GentrainScore=0.7055;HW=0.077647716;AN=2180;AC=11
20 63231 SNP20-11231 T G . PASS CR=99.93036;GentrainScore=0.7837;HW=0.035481825;AN=2182;AC=275
20 63244 rs6139074 A C . PASS CR=98.893394;GentrainScore=0.8001;HW=7.327299E-7;AN=2162;AC=501
20 63799 rs1418258 C T . PASS CR=99.75217;GentrainScore=0.8170;HW=0.6653377;AN=2182;AC=881

=== Input BAM file ===

verifyBamID requires a sorted, indexed, base quality recalibrated, and duplication-marked BAM file. It also requires to contain "@RG" header lines to annotation different readGroups (sequencing runs and lanes). The SM tag in the "@RG" header should match with one of the genotyped sample. Otherwise, verifyBamID may not be able to test whether the sequenced sample matches with genotyped sample, but will try to detect sample mixture from allele frequency, and will try to detect the best-matching sample among the genotyped sample.

== What the default option does ==

The default option of '''verifyBamID''' is the recommended setting for the most sequencing studies to provide a rapid and informative response. The default option provides the following features:
* --free-mix is turned on for estimating contamination using sequence-only method
* --chip-mix is turned on for estimating contamination or swap using sequence+array method, if the external genotype file is provided in the VCF
* --self is turnd on : The default option does not try to compare the sequence reads to identify the best matching individual (which is possible with --best option). It only compares with the external genotypes from the same individual to the sequenced individual.
* --maxDepth 20 is used without --precise option : The default option is intended for whole genome low coverage sequencing. For the targeted exome sequencing, --maxDepth 1000 and --precise is recommended.
* --ignoreRG is not a default option, but a recommended option, when you want to check the contamination for the entire BAM rather than examining each read group separately. This option will increase the computational efficiency especially in the case whether the sequence reads are multiplexed across many sequencing runs.

== Interpreting output files ==

=== Output files ===
When verifyBamID runs successfully, the following sets of files may be generated.
* [outPrefix].selfSM - Per-sample statistics describing how well the sample matches to the annotated sample.
* [outPrefix].depthSM - The depth distribution of the sequence reads per sample
* [outPrefix].selfRG - Per-readGroup statistics describing how well each lane matches to the annotated sample. (available only without --ignoreRG option)
* [outPrefix].depthRG - The depth distribution of the sequence reads per readGroup. (available only without --ignoreRG option)
* [outPrefix].bestSM - Per-sample best-match statistics with best-matching sample among the genotyped sample (available only with --best option)
* [outPrefix].bestRG - Per-readgroup best-match statistics with best-matching sample among the genotyped sample (available only with --best and without --ignoreRG option)

=== Column information in the output files ===
The .selfSM/.selfRG/.bestSM/.bestRG files have the following 19 columns per sample, or per readgroup (lane).

# SEQ_SM : Sample ID of the sequenced sample. Obtained from @RG header / SM tag in the BAM file
# RG : ReadGroup ID of sequenced lane. For [outPrefix].selfSM and [outPrefix].bestSM, these values are "ALL"
# CHIP_ID : Sample ID compared to in the genotype file. For [outPrefix].selfRG and [outPrefix].selfSM, these values should be identical to [SEQ_SM] or "NA" if the genotype of sequenced samples are unavailable. For [outPrefix].bestRG and [outPrefix].bestSM, these values should be the ID of best-matching sample among the genotype files compared to.
# # SNPs : # of SNPs passing the criteria from the VCF file
# # READS : Total # of reads loaded from the BAM file
# # AVG_DP : Average sequencing depth at the sites in the VCF file
# FREEMIX : Sequence-only estimate of contamination (0-1 scale)
# FREELK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence-only method
# FREELK0 : Log-ikelihood of the sequence reads given no contamination under sequence-only method
# FREE_RH : Estimated reference bias parameter Pr(refBase|HET) (when --free-refBias or --free-full is used)
# FREE_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --free-refBias or --free-full is used)
# CHIPMIX : Sequence+array estimate of contamination (NA if the external genotype is unavailable) (0-1 scale)
# CHIPLK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIPLK0 : Log-likelihood of the sequence reads given no contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIP_RH : Estimated reference bias parameter Pr(refBase|HET) (when --chip-refBias or --chip-full is used)
# CHIP_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --chip-refBias or --chip-full is used)
# DPREF : Depth (Coverage) of HomRef site (based on the genotypes of (SELF_SM/BEST_SM), passing mapQ, baseQual, maxDepth thresholds.
# RDPHET : DPHET/DPREF, Relative depth at Heterozygous site.
# RDPALT : DPHET/DPREF, Relative depth at HomAlt site.

=== A guideline to interpret output files ===

verifyBamID provides a series of information that is informative to determine whether the sample is possibly contaminated or swapped, but there is no single criteria that works for every circumstances. There are a few unmodeled factor in the estimation of [SELF-IBD]/[BEST-IBD] and [%MIX], so please note that the MLE estimation may not always exactly match to the true amount of contamination. Here we provide a guideline to flag potentially contaminated/swapped samples

* Each sample or lane can be checked in this way. When [CHIPMIX] >> 0.02 and/or [FREEMIX] >> 0.02, meaning 2% or more of non-reference bases are observed in reference sites, we recommend to examine the data more carefully for the possibility of contamination.
* We recommend to check each lane for the possibility of sample swaps. When [CHIPMIX] ~ 1 AND [FREEMIX] ~ 0, then it is possible that the sample is swapped with another sample. When [CHIPMIX] ~ 0 in .bestSM file, [CHIP_ID] might be actually the swapped sample. Otherwise, the swapped sample may not exist in the genotype data you have compared.
* When genotype data is not available but allele-frequency-based estimates of [FREEMIX] >= 0.03 and [FREELK1]-[FREELK0] is large, then it is possible that the sample is contaminated with other sample. We recommend to use per-sample data rather than per-lane data for checking this for low coverage data, because the inference will be more confident when there are large number of bases with depth 2 or higher.

== Command Line Options ==

The following parameters are available. Ones with "[]" are in effect:

Available Options
Input Files : --vcf [], --bam [], --subset [],
--smID []
VCF analysis options : --genoError [1.0e-03],
--minAF [0.01],
--minCallRate [0.50]
Individuals to compare with chip data : --site, --self, --best
Chip-free optimization options : --free-none, --free-mix [ON],
--free-refBias, --free-full
With-chip optimization options : --chip-none, --chip-mix [ON],
--chip-refBias, --chip-full
BAM analysis options : --ignoreRG, --ignoreOverlapPair,
--noEOF, --precise, --minMapQ [10],
--maxDepth [20], --minQ [13],
--maxQ [40], --grid [0.05]
Modeling Reference Bias : --refRef [1.00], --refHet [0.50],
--refAlt [0.00]
Output options : --out [], --verbose

Each option provides the following features:
* --vcf : specify required VCF file
* --bam : specify required BAM file (indexed with .bam.bai or .bai file)
* --subset : list of individual IDs to calculate the allele frequency. All individuals are used if unspecified
* --smID : If the individual ID in the BAM file and VCF file does not match, substitute the BAM file's ID into the specified argument
* --genoError : error rate of the external genotype file
* --minAF : minimum allele frequency of the markers to include
* --minAF : minimum call rate of the markers to include
* --site : If set, use only site information in the VCF and do not compare with the actual genotypes
* --self : Only compare the ID-matching individuals between the VCF and BAM file
* --best : Find the best matching individuals (.bestSM and .bestRG files will be produced). This option is substantially longer than the default option
* --free-none : Do not perform sequence-only method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence-only method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence-only method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence-only method
* --chip-none : Do not perform sequence+array method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence+array method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence+array method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence+array method
* --ignoreRG : ignore the read grouup level comparison and compare samples only (recommended for an expedited run)
* --ignoreOverlapPair : ignore overlapping pair end fragment covering the same base. Disabling this option may decrease the sensitivity of the method when the insert size is short (with slight gain in the computational speed)
* --noEOF : do not check the EOF marker of the BAM file (for earlier version of BAM)
* --precise : calculate the likelihood in log-scale for high-depth data (recommended when --maxDepth is greater than 20. Can be a little bit slower)
* --minMapQ : minimum mapping quality of the sequence reads to compare
* --minQ : minimum base quality to include
* --maxQ : maximum base quality to cap
* --grid : the grid interval to search the optimum before running Brent's algorithm.
* --refRef : Initial Pr(refBase|HOMREFGeno) parameter
* --refHet : Initial Pr(refBase|HETGeno) parameter
* --refAlt : Initial Pr(refBase|HOMALTGeno) parameter
* --out : output file prefix (required)
* --verbose : print the progress of the method on the screeen

== Principle of Operation ==

Each read group in a BAM file is evaluated independently. This means that in file with multiple read groups, problems will be flagged at the read group level (a plus). However, it also means that it might be hard to discern the correct assignment of read groups with very little data.

For each aligned base that overlaps a known genotype, we calculate the probability the probability that it was derived from a particular known genotype. This comparison considers only bases that overlap previously known genotypes and that meet the base quality and mapping quality thresholds.

Each individual in a pedigree has a different combination of genotypes, and bamGenotypeCheck will systematically search for the individual whose genotypes best match the observed read data.

For more about the technical details, see the page [[Verifying Sample Identities - Implementation]]

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== Contamination in Array Data ==

[[VerifyIDintensity]] software can estimate sample contamination from Illumina genotype array intensity data.

== Acknowledgements ==

VerifyBamID is a result from collaborative effort by Hyun Min Kang, Goo Jun, Matthew Flickinger, Mary Kate Wing, and Goncalo Abecasis. Please email to Hyun Min Kang [[mailto:hmkang@umich.edu| hmkang@umich.edu ]] for any questions.

VerifyBamID

2013-02-26T05:54:30Z

Gjun:

[[Category:Software|VerifyBamID]]

'''verifyBamID''' is a software that verifies whether the reads in particular file match previously known genotypes for an individual (or group of individuals), and checks whether the reads are contaminated as a mixture of two samples. '''verifyBamID''' can detect sample contamination and swaps when external genotypes are available. When external genotypes are not available, '''verifyBamID''' still robustly detects sample swaps.

== Download verifyBamID ==

To get a copy go to the [http://www.sph.umich.edu/csg/kang/verifyBamID/download VerifyBamID Download] download page.

== Join in verifyBamID mailing list ==

Please join in the [http://groups.google.com/group/verifybamid VerifyBamID Google Group] to ask / discuss / comment about verifyBamID.

== What's new ==

(2012/06/20)
* Fixed a bug of incorrect estimate of contamination when --chip-full option was used (Thanks to Richard Smith)
* Fixed a bug of incorrect per-readgroup output in --chip-* parameter

(2012/05/24)
* Fixed a bug of incorrect per-readgroup output (Thanks to Matthew Flickinger)
* '''(IMPORTANT)''' Add an option to remove either side of overlapping fragment. This option is turned on by default, and can be turned off usig --ignoreOverlapPair. If your sequence data has very short insert size, this update may increase the sensitivity of estimated contamination.
* Changes in the directory structure and Makefile

(2012/05/18) The new release of verifyBamID have undergone major change since the last version (as of 2011 April). Here are the highlights
* The genotype / allele frequency file is now based on VCF format rather than PLINK format.
* The reference sequence information is no longer required
* Uses Brent's method for precise estimation of contamination parameters
* Generate the depth distribution statistics.
* Estimated reference-bias parameters (useful mostly for ABI SOLiD sequence data)

== Build verifyBamID ==

The binary download of verifyBamID is available. You may use the version in Ubuntu 64-bit platform. To build verifyBamID, download the statgen library and run the following series of commands
tar xzvf verifyBamID.20120620.tar.gz
cd verifyBamID
make cloneLib (in the case ../libStatGen does not exist)
make
./bin/verifyBamID

Note that '''git clone''' command will create a directory ./libStatGen under your working directory, and '''make''' will create binary of verifyBamID under verifyBamID/bin/

verifyBamID is designed to be reasonably portable.

However, since development occurs only on Ubuntu 9.10 x86 and x64 platforms, and later, there are likely other portability issues.

Currently we support verifyBamID only on Ubuntu 9.10 and later on 64-bit processors.

== Basic Usage ==

A key step in any genetic analysis is to verify whether data being generated matches expectations. ''verifyBamID'' checks whether reads in a BAM file match previous genotypes for a specific sample. In addition, it detects possible sample mixture from population allele frequency only, which can be particularly useful when the genotype data is not available.

Using a mathematical model that relates observed sequence reads to an hypothetical true genotype, ''verifyBamID'' tries to decide whether sequence reads match a particular individual or are more likely to be contaminated (including a small proportion of foreign DNA), derived from a closely related individual, or derived from a completely different individual.

== Basic Usage Example ==

Here is a typical command line:

verifyBamID --vcf [input.vcf] --bam [input.bam] --out [output.prefix] --verbose --ignoreRG

where
[input.bam] is a BAM (Binary Alignment Map) file of a sequence reads
[input.vcf] is input VCF file containing individual genotypes or AF or AC/AN fields in the INFO field. gzipped VCF is also allowed.
[outPrefix] is output prefix of output files - [outPrefix].{selfRG,selfSM,bestRG,bestSM,depthRG,depthSM} will be created.

More detailed description of command line input is below

== Preparing input files ==

verifyBamID requires two input files - VCF file containing external genotypes or allele frequency information, and the BAM file.

=== VCF input genotype file ===

The input VCF file contains (1) external genotype information and/or (2) allele frequency information as AF entry or AC/AN entries in the INFO field. (See [http://www.1000genomes.org/wiki/Analysis/Variant%20Call%20Format/vcf-variant-call-format-version-41 | VCF specification] for further details). If neither information is provided, verifyBamID will not work properly.

If external genotype information is provided, sequence+array method will identify contamination and sample swaps by comparing the concordance between the external genotypes and the sequence reads. Additionally, sequence-only method will provide additional contamination estimates by modeling the sequence reads as mixture of two unknown samples based on the allele frequency information in the VCF file.

Input VCF file needs to meet several additional contraints need to meet in order to properly run verifyBamID.
* The VCF is assumed to be well-formed. For example, verifyBamID does not check whether REF allele actually matches with reference sequence.
* The VCF should only contain SNPs. Current version of verifyBamID does not accept INDELs, MNPs, Structural Variations, or other complex variants.
* The individual IDs in the VCF file, must be identical with the individual identifier in the BAM file. Otherwise, --smID option can override the sample ID information of the BAM file to the ID that matches to the individual IDs in the VCF file.
* IMPORTANT : For targeted sequencing data, it is important to subselect the markers to only include on-target markers in the genotype file. Off-target markers are not likely to have multiple non-duplicated reads at the marker position, and it may create artifacts in the analysis due to overlapping fragments.
* Currently, verifyBamID takes only autosomal chromosomes as input VCF.

An example input VCF file (without external genotype) is provided below. Note that AC and AC entries exists in the INFO field for the allele frequency information.

#CHROM POS ID REF ALT QUAL FILTER INFO
20 61651 SNP20-9651 C A . PASS CR=99.86851;GentrainScore=0.7055;HW=0.077647716;AN=2180;AC=11
20 63231 SNP20-11231 T G . PASS CR=99.93036;GentrainScore=0.7837;HW=0.035481825;AN=2182;AC=275
20 63244 rs6139074 A C . PASS CR=98.893394;GentrainScore=0.8001;HW=7.327299E-7;AN=2162;AC=501
20 63799 rs1418258 C T . PASS CR=99.75217;GentrainScore=0.8170;HW=0.6653377;AN=2182;AC=881

=== Input BAM file ===

verifyBamID requires a sorted, indexed, base quality recalibrated, and duplication-marked BAM file. It also requires to contain "@RG" header lines to annotation different readGroups (sequencing runs and lanes). The SM tag in the "@RG" header should match with one of the genotyped sample. Otherwise, verifyBamID may not be able to test whether the sequenced sample matches with genotyped sample, but will try to detect sample mixture from allele frequency, and will try to detect the best-matching sample among the genotyped sample.

== What the default option does ==

The default option of '''verifyBamID''' is the recommended setting for the most sequencing studies to provide a rapid and informative response. The default option provides the following features:
* --free-mix is turned on for estimating contamination using sequence-only method
* --chip-mix is turned on for estimating contamination or swap using sequence+array method, if the external genotype file is provided in the VCF
* --self is turnd on : The default option does not try to compare the sequence reads to identify the best matching individual (which is possible with --best option). It only compares with the external genotypes from the same individual to the sequenced individual.
* --maxDepth 20 is used without --precise option : The default option is intended for whole genome low coverage sequencing. For the targeted exome sequencing, --maxDepth 1000 and --precise is recommended.
* --ignoreRG is not a default option, but a recommended option, when you want to check the contamination for the entire BAM rather than examining each read group separately. This option will increase the computational efficiency especially in the case whether the sequence reads are multiplexed across many sequencing runs.

== Interpreting output files ==

=== Output files ===
When verifyBamID runs successfully, the following sets of files may be generated.
* [outPrefix].selfSM - Per-sample statistics describing how well the sample matches to the annotated sample.
* [outPrefix].depthSM - The depth distribution of the sequence reads per sample
* [outPrefix].selfRG - Per-readGroup statistics describing how well each lane matches to the annotated sample. (available only without --ignoreRG option)
* [outPrefix].depthRG - The depth distribution of the sequence reads per readGroup. (available only without --ignoreRG option)
* [outPrefix].bestSM - Per-sample best-match statistics with best-matching sample among the genotyped sample (available only with --best option)
* [outPrefix].bestRG - Per-readgroup best-match statistics with best-matching sample among the genotyped sample (available only with --best and without --ignoreRG option)

=== Column information in the output files ===
The .selfSM/.selfRG/.bestSM/.bestRG files have the following 19 columns per sample, or per readgroup (lane).

# SEQ_SM : Sample ID of the sequenced sample. Obtained from @RG header / SM tag in the BAM file
# RG : ReadGroup ID of sequenced lane. For [outPrefix].selfSM and [outPrefix].bestSM, these values are "ALL"
# CHIP_ID : Sample ID compared to in the genotype file. For [outPrefix].selfRG and [outPrefix].selfSM, these values should be identical to [SEQ_SM] or "NA" if the genotype of sequenced samples are unavailable. For [outPrefix].bestRG and [outPrefix].bestSM, these values should be the ID of best-matching sample among the genotype files compared to.
# # SNPs : # of SNPs passing the criteria from the VCF file
# # READS : Total # of reads loaded from the BAM file
# # AVG_DP : Average sequencing depth at the sites in the VCF file
# FREEMIX : Sequence-only estimate of contamination (0-1 scale)
# FREELK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence-only method
# FREELK0 : Log-ikelihood of the sequence reads given no contamination under sequence-only method
# FREE_RH : Estimated reference bias parameter Pr(refBase|HET) (when --free-refBias or --free-full is used)
# FREE_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --free-refBias or --free-full is used)
# CHIPMIX : Sequence+array estimate of contamination (NA if the external genotype is unavailable) (0-1 scale)
# CHIPLK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIPLK0 : Log-likelihood of the sequence reads given no contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIP_RH : Estimated reference bias parameter Pr(refBase|HET) (when --chip-refBias or --chip-full is used)
# CHIP_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --chip-refBias or --chip-full is used)
# DPREF : Depth (Coverage) of HomRef site (based on the genotypes of (SELF_SM/BEST_SM), passing mapQ, baseQual, maxDepth thresholds.
# RDPHET : DPHET/DPREF, Relative depth at Heterozygous site.
# RDPALT : DPHET/DPREF, Relative depth at HomAlt site.

=== A guideline to interpret output files ===

verifyBamID provides a series of information that is informative to determine whether the sample is possibly contaminated or swapped, but there is no single criteria that works for every circumstances. There are a few unmodeled factor in the estimation of [SELF-IBD]/[BEST-IBD] and [%MIX], so please note that the MLE estimation may not always exactly match to the true amount of contamination. Here we provide a guideline to flag potentially contaminated/swapped samples

* Each sample or lane can be checked in this way. When [CHIPMIX] >> 0.02 and/or [FREEMIX] >> 0.02, meaning 2% or more of non-reference bases are observed in reference sites, we recommend to examine the data more carefully for the possibility of contamination.
* We recommend to check each lane for the possibility of sample swaps. When [CHIPMIX] ~ 1 AND [FREEMIX] ~ 0, then it is possible that the sample is swapped with another sample. When [CHIPMIX] ~ 0 in .bestSM file, [CHIP_ID] might be actually the swapped sample. Otherwise, the swapped sample may not exist in the genotype data you have compared.
* When genotype data is not available but allele-frequency-based estimates of [FREEMIX] >= 0.03 and [FREELK1]-[FREELK0] is large, then it is possible that the sample is contaminated with other sample. We recommend to use per-sample data rather than per-lane data for checking this for low coverage data, because the inference will be more confident when there are large number of bases with depth 2 or higher.

== Command Line Options ==

The following parameters are available. Ones with "[]" are in effect:

Available Options
Input Files : --vcf [], --bam [], --subset [],
--smID []
VCF analysis options : --genoError [1.0e-03],
--minAF [0.01],
--minCallRate [0.50]
Individuals to compare with chip data : --site, --self, --best
Chip-free optimization options : --free-none, --free-mix [ON],
--free-refBias, --free-full
With-chip optimization options : --chip-none, --chip-mix [ON],
--chip-refBias, --chip-full
BAM analysis options : --ignoreRG, --ignoreOverlapPair,
--noEOF, --precise, --minMapQ [10],
--maxDepth [20], --minQ [13],
--maxQ [40], --grid [0.05]
Modeling Reference Bias : --refRef [1.00], --refHet [0.50],
--refAlt [0.00]
Output options : --out [], --verbose

Each option provides the following features:
* --vcf : specify required VCF file
* --bam : specify required BAM file (indexed with .bam.bai or .bai file)
* --subset : list of individual IDs to calculate the allele frequency. All individuals are used if unspecified
* --smID : If the individual ID in the BAM file and VCF file does not match, substitute the BAM file's ID into the specified argument
* --genoError : error rate of the external genotype file
* --minAF : minimum allele frequency of the markers to include
* --minAF : minimum call rate of the markers to include
* --site : If set, use only site information in the VCF and do not compare with the actual genotypes
* --self : Only compare the ID-matching individuals between the VCF and BAM file
* --best : Find the best matching individuals (.bestSM and .bestRG files will be produced). This option is substantially longer than the default option
* --free-none : Do not perform sequence-only method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence-only method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence-only method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence-only method
* --chip-none : Do not perform sequence+array method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence+array method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence+array method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence+array method
* --ignoreRG : ignore the read grouup level comparison and compare samples only (recommended for an expedited run)
* --ignoreOverlapPair : ignore overlapping pair end fragment covering the same base. Disabling this option may decrease the sensitivity of the method when the insert size is short (with slight gain in the computational speed)
* --noEOF : do not check the EOF marker of the BAM file (for earlier version of BAM)
* --precise : calculate the likelihood in log-scale for high-depth data (recommended when --maxDepth is greater than 20. Can be a little bit slower)
* --minMapQ : minimum mapping quality of the sequence reads to compare
* --minQ : minimum base quality to include
* --maxQ : maximum base quality to cap
* --grid : the grid interval to search the optimum before running Brent's algorithm.
* --refRef : Initial Pr(refBase|HOMREFGeno) parameter
* --refHet : Initial Pr(refBase|HETGeno) parameter
* --refAlt : Initial Pr(refBase|HOMALTGeno) parameter
* --out : output file prefix (required)
* --verbose : print the progress of the method on the screeen

== Principle of Operation ==

Each read group in a BAM file is evaluated independently. This means that in file with multiple read groups, problems will be flagged at the read group level (a plus). However, it also means that it might be hard to discern the correct assignment of read groups with very little data.

For each aligned base that overlaps a known genotype, we calculate the probability the probability that it was derived from a particular known genotype. This comparison considers only bases that overlap previously known genotypes and that meet the base quality and mapping quality thresholds.

Each individual in a pedigree has a different combination of genotypes, and bamGenotypeCheck will systematically search for the individual whose genotypes best match the observed read data.

For more about the technical details, see the page [[Verifying Sample Identities - Implementation]]

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== For array intensity data ==

[[VerifyIDintensity]] software can estimate sample contamination from Illumina genotype array intensity data.

== Acknowledgements ==

VerifyBamID is a result from collaborative effort by Hyun Min Kang, Goo Jun, Matthew Flickinger, Mary Kate Wing, and Goncalo Abecasis. Please email to Hyun Min Kang [[mailto:hmkang@umich.edu| hmkang@umich.edu ]] for any questions.

VerifyBamID

2013-02-26T05:54:03Z

Gjun:

[[Category:Software|VerifyBamID]]

'''verifyBamID''' is a software that verifies whether the reads in particular file match previously known genotypes for an individual (or group of individuals), and checks whether the reads are contaminated as a mixture of two samples. '''verifyBamID''' can detect sample contamination and swaps when external genotypes are available. When external genotypes are not available, '''verifyBamID''' still robustly detects sample swaps.

== Download verifyBamID ==

To get a copy go to the [http://www.sph.umich.edu/csg/kang/verifyBamID/download VerifyBamID Download] download page.

== Join in verifyBamID mailing list ==

Please join in the [http://groups.google.com/group/verifybamid VerifyBamID Google Group] to ask / discuss / comment about verifyBamID.

== What's new ==

(2012/06/20)
* Fixed a bug of incorrect estimate of contamination when --chip-full option was used (Thanks to Richard Smith)
* Fixed a bug of incorrect per-readgroup output in --chip-* parameter

(2012/05/24)
* Fixed a bug of incorrect per-readgroup output (Thanks to Matthew Flickinger)
* '''(IMPORTANT)''' Add an option to remove either side of overlapping fragment. This option is turned on by default, and can be turned off usig --ignoreOverlapPair. If your sequence data has very short insert size, this update may increase the sensitivity of estimated contamination.
* Changes in the directory structure and Makefile

(2012/05/18) The new release of verifyBamID have undergone major change since the last version (as of 2011 April). Here are the highlights
* The genotype / allele frequency file is now based on VCF format rather than PLINK format.
* The reference sequence information is no longer required
* Uses Brent's method for precise estimation of contamination parameters
* Generate the depth distribution statistics.
* Estimated reference-bias parameters (useful mostly for ABI SOLiD sequence data)

== Build verifyBamID ==

The binary download of verifyBamID is available. You may use the version in Ubuntu 64-bit platform. To build verifyBamID, download the statgen library and run the following series of commands
tar xzvf verifyBamID.20120620.tar.gz
cd verifyBamID
make cloneLib (in the case ../libStatGen does not exist)
make
./bin/verifyBamID

Note that '''git clone''' command will create a directory ./libStatGen under your working directory, and '''make''' will create binary of verifyBamID under verifyBamID/bin/

verifyBamID is designed to be reasonably portable.

However, since development occurs only on Ubuntu 9.10 x86 and x64 platforms, and later, there are likely other portability issues.

Currently we support verifyBamID only on Ubuntu 9.10 and later on 64-bit processors.

== Basic Usage ==

A key step in any genetic analysis is to verify whether data being generated matches expectations. ''verifyBamID'' checks whether reads in a BAM file match previous genotypes for a specific sample. In addition, it detects possible sample mixture from population allele frequency only, which can be particularly useful when the genotype data is not available.

Using a mathematical model that relates observed sequence reads to an hypothetical true genotype, ''verifyBamID'' tries to decide whether sequence reads match a particular individual or are more likely to be contaminated (including a small proportion of foreign DNA), derived from a closely related individual, or derived from a completely different individual.

== Basic Usage Example ==

Here is a typical command line:

verifyBamID --vcf [input.vcf] --bam [input.bam] --out [output.prefix] --verbose --ignoreRG

where
[input.bam] is a BAM (Binary Alignment Map) file of a sequence reads
[input.vcf] is input VCF file containing individual genotypes or AF or AC/AN fields in the INFO field. gzipped VCF is also allowed.
[outPrefix] is output prefix of output files - [outPrefix].{selfRG,selfSM,bestRG,bestSM,depthRG,depthSM} will be created.

More detailed description of command line input is below

== Preparing input files ==

verifyBamID requires two input files - VCF file containing external genotypes or allele frequency information, and the BAM file.

=== VCF input genotype file ===

The input VCF file contains (1) external genotype information and/or (2) allele frequency information as AF entry or AC/AN entries in the INFO field. (See [http://www.1000genomes.org/wiki/Analysis/Variant%20Call%20Format/vcf-variant-call-format-version-41 | VCF specification] for further details). If neither information is provided, verifyBamID will not work properly.

If external genotype information is provided, sequence+array method will identify contamination and sample swaps by comparing the concordance between the external genotypes and the sequence reads. Additionally, sequence-only method will provide additional contamination estimates by modeling the sequence reads as mixture of two unknown samples based on the allele frequency information in the VCF file.

Input VCF file needs to meet several additional contraints need to meet in order to properly run verifyBamID.
* The VCF is assumed to be well-formed. For example, verifyBamID does not check whether REF allele actually matches with reference sequence.
* The VCF should only contain SNPs. Current version of verifyBamID does not accept INDELs, MNPs, Structural Variations, or other complex variants.
* The individual IDs in the VCF file, must be identical with the individual identifier in the BAM file. Otherwise, --smID option can override the sample ID information of the BAM file to the ID that matches to the individual IDs in the VCF file.
* IMPORTANT : For targeted sequencing data, it is important to subselect the markers to only include on-target markers in the genotype file. Off-target markers are not likely to have multiple non-duplicated reads at the marker position, and it may create artifacts in the analysis due to overlapping fragments.
* Currently, verifyBamID takes only autosomal chromosomes as input VCF.

An example input VCF file (without external genotype) is provided below. Note that AC and AC entries exists in the INFO field for the allele frequency information.

#CHROM POS ID REF ALT QUAL FILTER INFO
20 61651 SNP20-9651 C A . PASS CR=99.86851;GentrainScore=0.7055;HW=0.077647716;AN=2180;AC=11
20 63231 SNP20-11231 T G . PASS CR=99.93036;GentrainScore=0.7837;HW=0.035481825;AN=2182;AC=275
20 63244 rs6139074 A C . PASS CR=98.893394;GentrainScore=0.8001;HW=7.327299E-7;AN=2162;AC=501
20 63799 rs1418258 C T . PASS CR=99.75217;GentrainScore=0.8170;HW=0.6653377;AN=2182;AC=881

=== Input BAM file ===

verifyBamID requires a sorted, indexed, base quality recalibrated, and duplication-marked BAM file. It also requires to contain "@RG" header lines to annotation different readGroups (sequencing runs and lanes). The SM tag in the "@RG" header should match with one of the genotyped sample. Otherwise, verifyBamID may not be able to test whether the sequenced sample matches with genotyped sample, but will try to detect sample mixture from allele frequency, and will try to detect the best-matching sample among the genotyped sample.

== What the default option does ==

The default option of '''verifyBamID''' is the recommended setting for the most sequencing studies to provide a rapid and informative response. The default option provides the following features:
* --free-mix is turned on for estimating contamination using sequence-only method
* --chip-mix is turned on for estimating contamination or swap using sequence+array method, if the external genotype file is provided in the VCF
* --self is turnd on : The default option does not try to compare the sequence reads to identify the best matching individual (which is possible with --best option). It only compares with the external genotypes from the same individual to the sequenced individual.
* --maxDepth 20 is used without --precise option : The default option is intended for whole genome low coverage sequencing. For the targeted exome sequencing, --maxDepth 1000 and --precise is recommended.
* --ignoreRG is not a default option, but a recommended option, when you want to check the contamination for the entire BAM rather than examining each read group separately. This option will increase the computational efficiency especially in the case whether the sequence reads are multiplexed across many sequencing runs.

== Interpreting output files ==

=== Output files ===
When verifyBamID runs successfully, the following sets of files may be generated.
* [outPrefix].selfSM - Per-sample statistics describing how well the sample matches to the annotated sample.
* [outPrefix].depthSM - The depth distribution of the sequence reads per sample
* [outPrefix].selfRG - Per-readGroup statistics describing how well each lane matches to the annotated sample. (available only without --ignoreRG option)
* [outPrefix].depthRG - The depth distribution of the sequence reads per readGroup. (available only without --ignoreRG option)
* [outPrefix].bestSM - Per-sample best-match statistics with best-matching sample among the genotyped sample (available only with --best option)
* [outPrefix].bestRG - Per-readgroup best-match statistics with best-matching sample among the genotyped sample (available only with --best and without --ignoreRG option)

=== Column information in the output files ===
The .selfSM/.selfRG/.bestSM/.bestRG files have the following 19 columns per sample, or per readgroup (lane).

# SEQ_SM : Sample ID of the sequenced sample. Obtained from @RG header / SM tag in the BAM file
# RG : ReadGroup ID of sequenced lane. For [outPrefix].selfSM and [outPrefix].bestSM, these values are "ALL"
# CHIP_ID : Sample ID compared to in the genotype file. For [outPrefix].selfRG and [outPrefix].selfSM, these values should be identical to [SEQ_SM] or "NA" if the genotype of sequenced samples are unavailable. For [outPrefix].bestRG and [outPrefix].bestSM, these values should be the ID of best-matching sample among the genotype files compared to.
# # SNPs : # of SNPs passing the criteria from the VCF file
# # READS : Total # of reads loaded from the BAM file
# # AVG_DP : Average sequencing depth at the sites in the VCF file
# FREEMIX : Sequence-only estimate of contamination (0-1 scale)
# FREELK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence-only method
# FREELK0 : Log-ikelihood of the sequence reads given no contamination under sequence-only method
# FREE_RH : Estimated reference bias parameter Pr(refBase|HET) (when --free-refBias or --free-full is used)
# FREE_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --free-refBias or --free-full is used)
# CHIPMIX : Sequence+array estimate of contamination (NA if the external genotype is unavailable) (0-1 scale)
# CHIPLK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIPLK0 : Log-likelihood of the sequence reads given no contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIP_RH : Estimated reference bias parameter Pr(refBase|HET) (when --chip-refBias or --chip-full is used)
# CHIP_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --chip-refBias or --chip-full is used)
# DPREF : Depth (Coverage) of HomRef site (based on the genotypes of (SELF_SM/BEST_SM), passing mapQ, baseQual, maxDepth thresholds.
# RDPHET : DPHET/DPREF, Relative depth at Heterozygous site.
# RDPALT : DPHET/DPREF, Relative depth at HomAlt site.

=== A guideline to interpret output files ===

verifyBamID provides a series of information that is informative to determine whether the sample is possibly contaminated or swapped, but there is no single criteria that works for every circumstances. There are a few unmodeled factor in the estimation of [SELF-IBD]/[BEST-IBD] and [%MIX], so please note that the MLE estimation may not always exactly match to the true amount of contamination. Here we provide a guideline to flag potentially contaminated/swapped samples

* Each sample or lane can be checked in this way. When [CHIPMIX] >> 0.02 and/or [FREEMIX] >> 0.02, meaning 2% or more of non-reference bases are observed in reference sites, we recommend to examine the data more carefully for the possibility of contamination.
* We recommend to check each lane for the possibility of sample swaps. When [CHIPMIX] ~ 1 AND [FREEMIX] ~ 0, then it is possible that the sample is swapped with another sample. When [CHIPMIX] ~ 0 in .bestSM file, [CHIP_ID] might be actually the swapped sample. Otherwise, the swapped sample may not exist in the genotype data you have compared.
* When genotype data is not available but allele-frequency-based estimates of [FREEMIX] >= 0.03 and [FREELK1]-[FREELK0] is large, then it is possible that the sample is contaminated with other sample. We recommend to use per-sample data rather than per-lane data for checking this for low coverage data, because the inference will be more confident when there are large number of bases with depth 2 or higher.

== Command Line Options ==

The following parameters are available. Ones with "[]" are in effect:

Available Options
Input Files : --vcf [], --bam [], --subset [],
--smID []
VCF analysis options : --genoError [1.0e-03],
--minAF [0.01],
--minCallRate [0.50]
Individuals to compare with chip data : --site, --self, --best
Chip-free optimization options : --free-none, --free-mix [ON],
--free-refBias, --free-full
With-chip optimization options : --chip-none, --chip-mix [ON],
--chip-refBias, --chip-full
BAM analysis options : --ignoreRG, --ignoreOverlapPair,
--noEOF, --precise, --minMapQ [10],
--maxDepth [20], --minQ [13],
--maxQ [40], --grid [0.05]
Modeling Reference Bias : --refRef [1.00], --refHet [0.50],
--refAlt [0.00]
Output options : --out [], --verbose

Each option provides the following features:
* --vcf : specify required VCF file
* --bam : specify required BAM file (indexed with .bam.bai or .bai file)
* --subset : list of individual IDs to calculate the allele frequency. All individuals are used if unspecified
* --smID : If the individual ID in the BAM file and VCF file does not match, substitute the BAM file's ID into the specified argument
* --genoError : error rate of the external genotype file
* --minAF : minimum allele frequency of the markers to include
* --minAF : minimum call rate of the markers to include
* --site : If set, use only site information in the VCF and do not compare with the actual genotypes
* --self : Only compare the ID-matching individuals between the VCF and BAM file
* --best : Find the best matching individuals (.bestSM and .bestRG files will be produced). This option is substantially longer than the default option
* --free-none : Do not perform sequence-only method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence-only method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence-only method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence-only method
* --chip-none : Do not perform sequence+array method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence+array method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence+array method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence+array method
* --ignoreRG : ignore the read grouup level comparison and compare samples only (recommended for an expedited run)
* --ignoreOverlapPair : ignore overlapping pair end fragment covering the same base. Disabling this option may decrease the sensitivity of the method when the insert size is short (with slight gain in the computational speed)
* --noEOF : do not check the EOF marker of the BAM file (for earlier version of BAM)
* --precise : calculate the likelihood in log-scale for high-depth data (recommended when --maxDepth is greater than 20. Can be a little bit slower)
* --minMapQ : minimum mapping quality of the sequence reads to compare
* --minQ : minimum base quality to include
* --maxQ : maximum base quality to cap
* --grid : the grid interval to search the optimum before running Brent's algorithm.
* --refRef : Initial Pr(refBase|HOMREFGeno) parameter
* --refHet : Initial Pr(refBase|HETGeno) parameter
* --refAlt : Initial Pr(refBase|HOMALTGeno) parameter
* --out : output file prefix (required)
* --verbose : print the progress of the method on the screeen

== Principle of Operation ==

Each read group in a BAM file is evaluated independently. This means that in file with multiple read groups, problems will be flagged at the read group level (a plus). However, it also means that it might be hard to discern the correct assignment of read groups with very little data.

For each aligned base that overlaps a known genotype, we calculate the probability the probability that it was derived from a particular known genotype. This comparison considers only bases that overlap previously known genotypes and that meet the base quality and mapping quality thresholds.

Each individual in a pedigree has a different combination of genotypes, and bamGenotypeCheck will systematically search for the individual whose genotypes best match the observed read data.

For more about the technical details, see the page [[Verifying Sample Identities - Implementation]]

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== For array intensity data ==

== For sequence data ==

[[VerifyIDintensity]] software can estimate sample contamination from Illumina genotype array intensity data.

== Acknowledgements ==

VerifyBamID is a result from collaborative effort by Hyun Min Kang, Goo Jun, Matthew Flickinger, Mary Kate Wing, and Goncalo Abecasis. Please email to Hyun Min Kang [[mailto:hmkang@umich.edu| hmkang@umich.edu ]] for any questions.

VerifyIDintensity

2013-02-26T05:53:04Z

Gjun:

[[Category:Software|VerifyIDintensity]]

'''verifyIDintensity''' is a software that detects and estimates sample contamination using intensity data from Illumina genotyping arrays using a mixture model.

== Download verifyIDintensity ==

* [[File:VerifyIDintensity.tgz]]

== Build verifyID intensity ==

To build verifyIDintensity, run the following series of commands. You need [http://www.boost.org/ boost library] and [http://tclap.sourceforge.net/ tclap].

$ tar xzvf verifyIDintensity.tgz
$ make

== Basic Usage ==

verifyIDintensity [-t <float>] [-m <int>] -n <int> [-b <string>] [-s <string>] -i <string> [-v] [-p] [--] [--version] [-h]

== Options ==

-t <float>, --threshold <float>
Minimum allele frequency for likelihood estimation, default is 0.01

-m <int>, --marker <int>
(required) Number of markers

-n <int>, --number <int>
(required) Number of samples

-b <string>, --abf <string>
Allele frequency file (ABF), which is a plain text file with SNP_ID and Allele_B frequency.
SNP_IDs should be sorted in the same order as the intensity file

-s <string>, --stat <string>
Statistics file (created if not exist)

-i <string>, --in <string>
(required) Input pre-computed intensity (.adpc.bin) file

-v, --verbose
Turn on verbose mode

-p, --persample
Do per-sample analysis, default is per-marker analysis

--, --ignore_rest
Ignores the rest of the labeled arguments following this flag.

--version
Displays version information and exits.

-h, --help
Displays usage information and exits.

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== For sequence data ==

[[VerifyBamID]] software can estimate sample contamination from aligned sequence reads and population minor allele frequency

VerifyBamID

2012-12-04T04:15:15Z

Gjun:

[[Category:Software|VerifyBamID]]

'''verifyBamID''' is a software that verifies whether the reads in particular file match previously known genotypes for an individual (or group of individuals), and checks whether the reads are contaminated as a mixture of two samples. '''verifyBamID''' can detect sample contamination and swaps when external genotypes are available. When external genotypes are not available, '''verifyBamID''' still robustly detects sample swaps.

== Download verifyBamID ==

To get a copy go to the [http://www.sph.umich.edu/csg/kang/verifyBamID/download VerifyBamID Download] download page.

== Join in verifyBamID mailing list ==

Please join in the [http://groups.google.com/group/verifybamid | VerifyBamID Google Group] to ask / discuss / comment about verifyBamID.

== What's new ==

(2012/06/20)
* Fixed a bug of incorrect estimate of contamination when --chip-full option was used (Thanks to Richard Smith)
* Fixed a bug of incorrect per-readgroup output in --chip-* parameter

(2012/05/24)
* Fixed a bug of incorrect per-readgroup output (Thanks to Matthew Flickinger)
* '''(IMPORTANT)''' Add an option to remove either side of overlapping fragment. This option is turned on by default, and can be turned off usig --ignoreOverlapPair. If your sequence data has very short insert size, this update may increase the sensitivity of estimated contamination.
* Changes in the directory structure and Makefile

(2012/05/18) The new release of verifyBamID have undergone major change since the last version (as of 2011 April). Here are the highlights
* The genotype / allele frequency file is now based on VCF format rather than PLINK format.
* The reference sequence information is no longer required
* Uses Brent's method for precise estimation of contamination parameters
* Generate the depth distribution statistics.
* Estimated reference-bias parameters (useful mostly for ABI SOLiD sequence data)

== Build verifyBamID ==

The binary download of verifyBamID is available. You may use the version in Ubuntu 64-bit platform. To build verifyBamID, download the statgen library and run the following series of commands
tar xzvf verifyBamID.20120620.tar.gz
cd verifyBamID
make cloneLib (in the case ../libStatGen does not exist)
make
./bin/verifyBamID

Note that '''git clone''' command will create a directory ./libStatGen under your working directory, and '''make''' will create binary of verifyBamID under verifyBamID/bin/

verifyBamID is designed to be reasonably portable.

However, since development occurs only on Ubuntu 9.10 x86 and x64 platforms, and later, there are likely other portability issues.

Currently we support verifyBamID only on Ubuntu 9.10 and later on 64-bit processors.

== Basic Usage ==

A key step in any genetic analysis is to verify whether data being generated matches expectations. ''verifyBamID'' checks whether reads in a BAM file match previous genotypes for a specific sample. In addition, it detects possible sample mixture from population allele frequency only, which can be particularly useful when the genotype data is not available.

Using a mathematical model that relates observed sequence reads to an hypothetical true genotype, ''verifyBamID'' tries to decide whether sequence reads match a particular individual or are more likely to be contaminated (including a small proportion of foreign DNA), derived from a closely related individual, or derived from a completely different individual.

== Basic Usage Example ==

Here is a typical command line:

verifyBamID --vcf [input.vcf] --bam [input.bam] --out [output.prefix] --verbose --ignoreRG

where
[input.bam] is a BAM (Binary Alignment Map) file of a sequence reads
[input.vcf] is input VCF file containing individual genotypes or AF or AC/AN fields in the INFO field. gzipped VCF is also allowed.
[outPrefix] is output prefix of output files - [outPrefix].{selfRG,selfSM,bestRG,bestSM,depthRG,depthSM} will be created.

More detailed description of command line input is below

== Preparing input files ==

verifyBamID requires two input files - VCF file containing external genotypes or allele frequency information, and the BAM file.

=== VCF input genotype file ===

The input VCF file contains (1) external genotype information and/or (2) allele frequency information as AF entry or AC/AN entries in the INFO field. (See [http://www.1000genomes.org/wiki/Analysis/Variant%20Call%20Format/vcf-variant-call-format-version-41 | VCF specification] for further details). If neither information is provided, verifyBamID will not work properly.

If external genotype information is provided, sequence+array method will identify contamination and sample swaps by comparing the concordance between the external genotypes and the sequence reads. Additionally, sequence-only method will provide additional contamination estimates by modeling the sequence reads as mixture of two unknown samples based on the allele frequency information in the VCF file.

Input VCF file needs to meet several additional contraints need to meet in order to properly run verifyBamID.
* The VCF is assumed to be well-formed. For example, verifyBamID does not check whether REF allele actually matches with reference sequence.
* The VCF should only contain SNPs. Current version of verifyBamID does not accept INDELs, MNPs, Structural Variations, or other complex variants.
* The individual IDs in the VCF file, must be identical with the individual identifier in the BAM file. Otherwise, --smID option can override the sample ID information of the BAM file to the ID that matches to the individual IDs in the VCF file.
* IMPORTANT : For targeted sequencing data, it is important to subselect the markers to only include on-target markers in the genotype file. Off-target markers are not likely to have multiple non-duplicated reads at the marker position, and it may create artifacts in the analysis due to overlapping fragments.
* Currently, verifyBamID takes only autosomal chromosomes as input VCF.

An example input VCF file (without external genotype) is provided below. Note that AC and AC entries exists in the INFO field for the allele frequency information.

#CHROM POS ID REF ALT QUAL FILTER INFO
20 61651 SNP20-9651 C A . PASS CR=99.86851;GentrainScore=0.7055;HW=0.077647716;AN=2180;AC=11
20 63231 SNP20-11231 T G . PASS CR=99.93036;GentrainScore=0.7837;HW=0.035481825;AN=2182;AC=275
20 63244 rs6139074 A C . PASS CR=98.893394;GentrainScore=0.8001;HW=7.327299E-7;AN=2162;AC=501
20 63799 rs1418258 C T . PASS CR=99.75217;GentrainScore=0.8170;HW=0.6653377;AN=2182;AC=881

=== Input BAM file ===

verifyBamID requires a sorted, indexed, base quality recalibrated, and duplication-marked BAM file. It also requires to contain "@RG" header lines to annotation different readGroups (sequencing runs and lanes). The SM tag in the "@RG" header should match with one of the genotyped sample. Otherwise, verifyBamID may not be able to test whether the sequenced sample matches with genotyped sample, but will try to detect sample mixture from allele frequency, and will try to detect the best-matching sample among the genotyped sample.

== What the default option does ==

The default option of '''verifyBamID''' is the recommended setting for the most sequencing studies to provide a rapid and informative response. The default option provides the following features:
* --free-mix is turned on for estimating contamination using sequence-only method
* --chip-mix is turned on for estimating contamination or swap using sequence+array method, if the external genotype file is provided in the VCF
* --self is turnd on : The default option does not try to compare the sequence reads to identify the best matching individual (which is possible with --best option). It only compares with the external genotypes from the same individual to the sequenced individual.
* --maxDepth 20 is used without --precise option : The default option is intended for whole genome low coverage sequencing. For the targeted exome sequencing, --maxDepth 1000 and --precise is recommended.
* --ignoreRG is not a default option, but a recommended option, when you want to check the contamination for the entire BAM rather than examining each read group separately. This option will increase the computational efficiency especially in the case whether the sequence reads are multiplexed across many sequencing runs.

== Interpreting output files ==

=== Output files ===
When verifyBamID runs successfully, the following sets of files may be generated.
* [outPrefix].selfSM - Per-sample statistics describing how well the sample matches to the annotated sample.
* [outPrefix].depthSM - The depth distribution of the sequence reads per sample
* [outPrefix].selfRG - Per-readGroup statistics describing how well each lane matches to the annotated sample. (available only without --ignoreRG option)
* [outPrefix].depthRG - The depth distribution of the sequence reads per readGroup. (available only without --ignoreRG option)
* [outPrefix].bestSM - Per-sample best-match statistics with best-matching sample among the genotyped sample (available only with --best option)
* [outPrefix].bestRG - Per-readgroup best-match statistics with best-matching sample among the genotyped sample (available only with --best and without --ignoreRG option)

=== Column information in the output files ===
The .selfSM/.selfRG/.bestSM/.bestRG files have the following 19 columns per sample, or per readgroup (lane).

# SEQ_SM : Sample ID of the sequenced sample. Obtained from @RG header / SM tag in the BAM file
# RG : ReadGroup ID of sequenced lane. For [outPrefix].selfSM and [outPrefix].bestSM, these values are "ALL"
# CHIP_ID : Sample ID compared to in the genotype file. For [outPrefix].selfRG and [outPrefix].selfSM, these values should be identical to [SEQ_SM] or "NA" if the genotype of sequenced samples are unavailable. For [outPrefix].bestRG and [outPrefix].bestSM, these values should be the ID of best-matching sample among the genotype files compared to.
# # SNPs : # of SNPs passing the criteria from the VCF file
# # READS : Total # of reads loaded from the BAM file
# # AVG_DP : Average sequencing depth at the sites in the VCF file
# FREEMIX : Sequence-only estimate of contamination
# FREELK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence-only method
# FREELK0 : Log-ikelihood of the sequence reads given no contamination under sequence-only method
# FREE_RH : Estimated reference bias parameter Pr(refBase|HET) (when --free-refBias or --free-full is used)
# FREE_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --free-refBias or --free-full is used)
# CHIPMIX : Sequence+array estimate of contamination (NA if the external genotype is unavailable)
# CHIPLK1 : Maximum log-likelihood of the sequence reads given estimated contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIPLK0 : Log-likelihood of the sequence reads given no contamination under sequence+array method (NA if the external genotypes are unavailable)
# CHIP_RH : Estimated reference bias parameter Pr(refBase|HET) (when --chip-refBias or --chip-full is used)
# CHIP_RA : Estimated reference bias parameter Pr(refBase|HOMALT) (when --chip-refBias or --chip-full is used)
# DPREF : Depth (Coverage) of HomRef site (based on the genotypes of (SELF_SM/BEST_SM), passing mapQ, baseQual, maxDepth thresholds.
# RDPHET : DPHET/DPREF, Relative depth at Heterozygous site.
# RDPALT : DPHET/DPREF, Relative depth at HomAlt site.

=== A guideline to interpret output files ===

verifyBamID provides a series of information that is informative to determine whether the sample is possibly contaminated or swapped, but there is no single criteria that works for every circumstances. There are a few unmodeled factor in the estimation of [SELF-IBD]/[BEST-IBD] and [%MIX], so please note that the MLE estimation may not always exactly match to the true amount of contamination. Here we provide a guideline to flag potentially contaminated/swapped samples

* Each sample or lane can be checked in this way. When [CHIPMIX] >> 0.02 and/or [FREEMIX] >> 0.02, meaning 2% or more of non-reference bases are observed in reference sites, we recommend to examine the data more carefully for the possibility of contamination.
* We recommend to check each lane for the possibility of sample swaps. When [CHIPMIX] ~ 1 AND [FREEMIX] ~ 0, then it is possible that the sample is swapped with another sample. When [CHIPMIX] ~ 0 in .bestSM file, [CHIP_ID] might be actually the swapped sample. Otherwise, the swapped sample may not exist in the genotype data you have compared.
* When genotype data is not available but allele-frequency-based estimates of [FREEMIX] >= 0.03 and [FREELK1]-[FREELK0] is large, then it is possible that the sample is contaminated with other sample. We recommend to use per-sample data rather than per-lane data for checking this for low coverage data, because the inference will be more confident when there are large number of bases with depth 2 or higher.

== Command Line Options ==

The following parameters are available. Ones with "[]" are in effect:

Available Options
Input Files : --vcf [], --bam [], --subset [],
--smID []
VCF analysis options : --genoError [1.0e-03],
--minAF [0.01],
--minCallRate [0.50]
Individuals to compare with chip data : --site, --self, --best
Chip-free optimization options : --free-none, --free-mix [ON],
--free-refBias, --free-full
With-chip optimization options : --chip-none, --chip-mix [ON],
--chip-refBias, --chip-full
BAM analysis options : --ignoreRG, --ignoreOverlapPair,
--noEOF, --precise, --minMapQ [10],
--maxDepth [20], --minQ [13],
--maxQ [40], --grid [0.05]
Modeling Reference Bias : --refRef [1.00], --refHet [0.50],
--refAlt [0.00]
Output options : --out [], --verbose

Each option provides the following features:
* --vcf : specify required VCF file
* --bam : specify required BAM file (indexed with .bam.bai or .bai file)
* --subset : list of individual IDs to calculate the allele frequency. All individuals are used if unspecified
* --smID : If the individual ID in the BAM file and VCF file does not match, substitute the BAM file's ID into the specified argument
* --genoError : error rate of the external genotype file
* --minAF : minimum allele frequency of the markers to include
* --minAF : minimum call rate of the markers to include
* --site : If set, use only site information in the VCF and do not compare with the actual genotypes
* --self : Only compare the ID-matching individuals between the VCF and BAM file
* --best : Find the best matching individuals (.bestSM and .bestRG files will be produced). This option is substantially longer than the default option
* --free-none : Do not perform sequence-only method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence-only method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence-only method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence-only method
* --chip-none : Do not perform sequence+array method to estimate parameters
* --free-mix : (default) Estimate contamination using sequence+array method with Brent's single dimensional optimization.
* --free-refBias : Estimate the refernece bias parameters using sequence+array method with Simplex method
* --free-full : Estimate both reference bias parameters and the contamination parameters using sequence+array method
* --ignoreRG : ignore the read grouup level comparison and compare samples only (recommended for an expedited run)
* --ignoreOverlapPair : ignore overlapping pair end fragment covering the same base. Disabling this option may decrease the sensitivity of the method when the insert size is short (with slight gain in the computational speed)
* --noEOF : do not check the EOF marker of the BAM file (for earlier version of BAM)
* --precise : calculate the likelihood in log-scale for high-depth data (recommended when --maxDepth is greater than 20. Can be a little bit slower)
* --minMapQ : minimum mapping quality of the sequence reads to compare
* --minQ : minimum base quality to include
* --maxQ : maximum base quality to cap
* --grid : the grid interval to search the optimum before running Brent's algorithm.
* --refRef : Initial Pr(refBase|HOMREFGeno) parameter
* --refHet : Initial Pr(refBase|HETGeno) parameter
* --refAlt : Initial Pr(refBase|HOMALTGeno) parameter
* --out : output file prefix (required)
* --verbose : print the progress of the method on the screeen

== Principle of Operation ==

Each read group in a BAM file is evaluated independently. This means that in file with multiple read groups, problems will be flagged at the read group level (a plus). However, it also means that it might be hard to discern the correct assignment of read groups with very little data.

For each aligned base that overlaps a known genotype, we calculate the probability the probability that it was derived from a particular known genotype. This comparison considers only bases that overlap previously known genotypes and that meet the base quality and mapping quality thresholds.

Each individual in a pedigree has a different combination of genotypes, and bamGenotypeCheck will systematically search for the individual whose genotypes best match the observed read data.

For more about the technical details, see the page [[Verifying Sample Identities - Implementation]]

== Reference ==

Please cite the following paper:

G. Jun, M. Flickinger, K. N. Hetrick, Kurt, J. M. Romm, K. F. Doheny, G. Abecasis, M. Boehnke,and H. M. Kang, ''Detecting and Estimating Contamination of Human DNA Samples in Sequencing and Array-Based Genotype Data'', American journal of human genetics doi:10.1016/j.ajhg.2012.09.004 (volume 91 issue 5 pp.839 - 848)

== Acknowledgements ==

VerifyBamID is a result from collaborative effort by Hyun Min Kang, Goo Jun, Matthew Flickinger, Mary Kate Trost, and Goncalo Abecasis. Please email to Hyun Min Kang [[mailto:hmkang@umich.edu| hmkang@umich.edu ]] for any questions.

File:VerifyIDintensity.tgz

2012-11-21T16:16:45Z

Gjun: uploaded a new version of "File:VerifyIDintensity.tgz"

VerifyIDintensity

2012-11-12T18:36:22Z

Gjun:

[[Category:Software|VerifyIDintensity]]

'''verifyIDintensity''' is a software that detects and estimates sample contamination using intensity data from Illumina genotyping arrays using a mixture model.

== Download verifyIDintensity ==

* [[File:VerifyIDintensity.tgz]]

== Build verifyID intensity ==

To build verifyIDintensity, run the following series of commands. You need [http://www.boost.org/ boost library] and [http://tclap.sourceforge.net/ tclap].

$ tar xzvf verifyIDintensity.tgz
$ make

== Basic Usage ==

verifyIDintensity [-t <float>] [-m <int>] -n <int> [-b <string>] [-s <string>] -i <string> [-v] [-p] [--] [--version] [-h]

== Options ==

-t <float>, --threshold <float>
Minimum allele frequency for likelihood estimation, default is 0.01

-m <int>, --marker <int>
(required) Number of markers

-n <int>, --number <int>
(required) Number of samples

-b <string>, --abf <string>
Allele frequency file (ABF), which is a plain text file with SNP_ID and Allele_B frequency.
SNP_IDs should be sorted in the same order as the intensity file

-s <string>, --stat <string>
Statistics file (created if not exist)

-i <string>, --in <string>
(required) Input pre-computed intensity (.adpc.bin) file

-v, --verbose
Turn on verbose mode

-p, --persample
Do per-sample analysis, default is per-marker analysis

--, --ignore_rest
Ignores the rest of the labeled arguments following this flag.

--version
Displays version information and exits.

-h, --help
Displays usage information and exits.

Polymutt

2012-10-05T16:21:31Z

Gjun: /* Creation of GLF files */

== Updates ==
The latest version of 0.04 is available for [[#Download | Download]].

== Introduction ==
* The program '''polymutt''' implemented a likelihood-based framework for calling '''single nucleotide variants''' and detecting '''''de novo''''' '''point mutation''' events in families for next-generation sequencing data. The program takes as input genotype likelihood format (GLF) files which can be generated following the [[#Creation of GLF files | Creation of GLF files]] instruction and outputs the result in the [[http://www.1000genomes.org/node/101 VCF]] format. The variant calling and ''de novo'' mutation detection are modelled jointly within families and can handle both nuclear and extended pedigrees without consanguinity loops. The input is a set of GLF files for each of family members and the relationships are specified through the .ped file.

* The evidence of variants and ''de novo'' mutations are assessed probabilistically. For a variant, the QUAL value is calculated as -10*log10(1-posterior(Variant | Data)) and for ''de novo'' mutation events a ''de novo'' quality (DQ) value is defined as log10(lk_denovo / lk_no_denovo) where lk_denovo and lk_no_denovo are the likelihoods of data allowing and disallowing ''de novo'' mutations respectively. Similarly, for each genotype, a genotype quality (GQ) value is defined as -10*log10(1-posterior(Genotype | Data)).

* Since unrelated individuals are kind of special case of families, unrelated individuals or a mixture of related and unrelated individuals can be handled.

* If some individuals in a family are not sequenced, this can be handled by setting the corresponding GLF file indices to zero for those family members who are not sequenced.

* NOTE: This version only works for autosomes. Variant calling for X, Y and MT is in the testing process and will be available in next version.

* See below for more details.

== Usage ==
A command without any input will display the basic usage

polymutt

The following parameters are in effect:
pedfile : (-pname)
datfile : (-dname)
glfIndexFile : (-gname)
posterior cutoff : 0.500 (-c99.999)

Additional Options
Map Quality Filter : --minMapQuality
Depth Filter : --minDepth, --maxDepth,
--minPercSampleWithData [0.00]
Scaled mutation rate : --theta [1.0e-03]
Prior of ts/tv ratio : --poly_tstv [2.00]
de novo mutation : --denovo, --rate_denovo [1.5e-08],
--tstv_denovo [2.00], --minLLR_denovo [1.00]
Optimization precision : --prec [1.0e-04]
Multiple threading : --nthreads [1]
Chromosomes to process : --chr2process []
Output : --vcf [variantCalls.vcf], --gl_off

An example command for variant calling looks like the following:
polymutt -p in.ped -d in.dat -g glfIndexFile --vcf out.vcf --nthreads 4

An example command for ''de novo'' mutation detection is as follows:
polymutt -p in.ped -d in.dat -g glfIndexFile --denovo --rate_denovo 1.5e-08 --min_denovo_LLR 1.0 --vcf out.denovo.vcf --nthreads 4

== Input files ==

Required input files are -p input.ped -d input.dat -g glfIndex files

* An example in.ped file looks like the following:
fam1 p1 0 0 1 1
fam1 p2 0 0 2 2
fam1 p3 p1 p2 1 3
fam2 p4 0 0 1 4
fam2 p5 0 0 2 5
fam2 p6 p4 p5 1 6
...

* An example in.dat file is like the following (for the 6th column above and in addition other traits/markers can be specified but will be ignored):
T GLF_Index

* An example glfIndex file is like the following and the numbers (except zeros) in the 6th column in the above in.ped file have to be present in the first column.
1 /home/me/sample1.glf
2 /home/me/sample2.glf
3 /home/me/sample3.glf
4 /home/me/sample4.glf
...

* If some of the members are not sequenced but are in the pedigree because of the relatedness with other members, the GLF_Index column (6th column) in the ped file should be set to zero
* For unrelated individuals, you can either (1) create a family for each unrelated individual as a founder or (2) put all unrelated individuals as founders in a single family.

== Other options ==

Some of command line options are explained below and others are self-explanatory.

-c : minimum cutoff of posterior probability to output a variant [''Default: 0.5'']
--theta : scaled mutation rate per site [''Default: 0.001'']
--tstv: prior of ts:tv ratio [''Default: 2.0'']
--nthreads : number of threads to run and it is recommended to use 4 threads for small number of input files [''Default: 1']

--denovo : a boolean flag to turn on ''de novo'' mutation detection. The following options take effect only when this flag is ON
--rate_denovo : mutation rate per haplotype per generation. [''Default: 1.5e-08'']
--tstv_denovo : the prior ts/tv ratio of ''de novo'' mutations. [''Default: 2.0'']
--minLLR_denovo : minimum value of log10 likelihood ratio of allowing vs. disallowing ''de novo'' mutations in the data to output [''Default: 1.0'']

--chr2process: the chromosome names to process. Default is empty and is to process all chromosomes in the input.
If multiple chromosomes are provided, they should be separated by comma, e.g. --chr2process 2,10 or --chr2process chr2,chr10

--gl_off: not to output genotype likelihood values for each individual. Default is to output 3 GLs for polymorphisms and 10 GLs for de novo mutations

== Output files ==
* The output file is a VCF file and the specification can be found [[http://www.1000genomes.org/node/101 here]]
* Since there is no standard to represent ''de novo'' mutations in the current VCF specification, actual genotypes (e.g. [ACGT]/[ACGT]) are output in the VCF file for ''de novo'' mutations.
* A summary about variant calling statistics is output to STDOUT and it may be redirected to a file for a record.

Summary of reference -- 9
Total Entry Count: 141213431
Total Base Cout: 120124735
Total '0' Base Count: 137
Non-Polymorphic Count: 655457
Transition Count: 6556
Transversion Count: 3127
Other Polymorphism Count: 0
Filter counts:
minMapQual 4550
minTotalDepth 1089
maxTotalDepth 736
Hard to call: 0
Skipped bases: 134

== Creation of GLF files ==

* The current version performs variant calling and ''de novo'' mutation detection from files in the genotype likelihood format (GLF). In future versions we plan to take [[http://samtools.sourceforge.net/ SAM/BAM]] files as input. See the following for instructions on how to create GLF files.
** Download a [[Samtools-hybrid|modified version of samtools]] ( [[https://github.com/statgen/samtools-0.1.7a-hybrid samtools-hybrid]] )
** Prepare the reference genome in fasta format and sequence alignments in [[http://samtools.sourceforge.net/ SAM/BAM]] format
** Generate BAQ adjusted GLF files using the following command
samtools-hybrid view -bh chr1.bam 1:0 | samtools-hybrid calmd -Abr - human.v37.fa 2> /dev/null | samtools-hybrid pileup - -g -f human.v37.fa > chr1.bam.glf
** If you want to clip overlapping reads, you can add clipOverlap command
samtools-hybrid view -bh chr1.bam 1:0 |samtools-hybrid calmd -Abr - human.v37.fa 2> /dev/null | bam clipOverlap --in -.bam --out -.ubam| samtools-hybrid pileup - -g -f human.v37.fa > chr1.bam.glf
* For other functionalities please refer to the [[http://samtools.sourceforge.net/ samtools]] website.

== Download ==
The latest version of source code v0.04 with test files can be [[Media:polymutt.0.04.tar.gz | downloaded]] here.

== Contact ==
For questions please contact the authors (Bingshan Li: [mailto:bingshan@umich.edu bingshan@umich.edu] or Goncalo Abecasis: [mailto:goncalo@umich.edu goncalo@umich.edu])

[[Category:Software]]

Polymutt

2012-10-05T16:20:15Z

Gjun: /* Creation of GLF files */

== Updates ==
The latest version of 0.04 is available for [[#Download | Download]].

== Introduction ==
* The program '''polymutt''' implemented a likelihood-based framework for calling '''single nucleotide variants''' and detecting '''''de novo''''' '''point mutation''' events in families for next-generation sequencing data. The program takes as input genotype likelihood format (GLF) files which can be generated following the [[#Creation of GLF files | Creation of GLF files]] instruction and outputs the result in the [[http://www.1000genomes.org/node/101 VCF]] format. The variant calling and ''de novo'' mutation detection are modelled jointly within families and can handle both nuclear and extended pedigrees without consanguinity loops. The input is a set of GLF files for each of family members and the relationships are specified through the .ped file.

* The evidence of variants and ''de novo'' mutations are assessed probabilistically. For a variant, the QUAL value is calculated as -10*log10(1-posterior(Variant | Data)) and for ''de novo'' mutation events a ''de novo'' quality (DQ) value is defined as log10(lk_denovo / lk_no_denovo) where lk_denovo and lk_no_denovo are the likelihoods of data allowing and disallowing ''de novo'' mutations respectively. Similarly, for each genotype, a genotype quality (GQ) value is defined as -10*log10(1-posterior(Genotype | Data)).

* Since unrelated individuals are kind of special case of families, unrelated individuals or a mixture of related and unrelated individuals can be handled.

* If some individuals in a family are not sequenced, this can be handled by setting the corresponding GLF file indices to zero for those family members who are not sequenced.

* NOTE: This version only works for autosomes. Variant calling for X, Y and MT is in the testing process and will be available in next version.

* See below for more details.

== Usage ==
A command without any input will display the basic usage

polymutt

The following parameters are in effect:
pedfile : (-pname)
datfile : (-dname)
glfIndexFile : (-gname)
posterior cutoff : 0.500 (-c99.999)

Additional Options
Map Quality Filter : --minMapQuality
Depth Filter : --minDepth, --maxDepth,
--minPercSampleWithData [0.00]
Scaled mutation rate : --theta [1.0e-03]
Prior of ts/tv ratio : --poly_tstv [2.00]
de novo mutation : --denovo, --rate_denovo [1.5e-08],
--tstv_denovo [2.00], --minLLR_denovo [1.00]
Optimization precision : --prec [1.0e-04]
Multiple threading : --nthreads [1]
Chromosomes to process : --chr2process []
Output : --vcf [variantCalls.vcf], --gl_off

An example command for variant calling looks like the following:
polymutt -p in.ped -d in.dat -g glfIndexFile --vcf out.vcf --nthreads 4

An example command for ''de novo'' mutation detection is as follows:
polymutt -p in.ped -d in.dat -g glfIndexFile --denovo --rate_denovo 1.5e-08 --min_denovo_LLR 1.0 --vcf out.denovo.vcf --nthreads 4

== Input files ==

Required input files are -p input.ped -d input.dat -g glfIndex files

* An example in.ped file looks like the following:
fam1 p1 0 0 1 1
fam1 p2 0 0 2 2
fam1 p3 p1 p2 1 3
fam2 p4 0 0 1 4
fam2 p5 0 0 2 5
fam2 p6 p4 p5 1 6
...

* An example in.dat file is like the following (for the 6th column above and in addition other traits/markers can be specified but will be ignored):
T GLF_Index

* An example glfIndex file is like the following and the numbers (except zeros) in the 6th column in the above in.ped file have to be present in the first column.
1 /home/me/sample1.glf
2 /home/me/sample2.glf
3 /home/me/sample3.glf
4 /home/me/sample4.glf
...

* If some of the members are not sequenced but are in the pedigree because of the relatedness with other members, the GLF_Index column (6th column) in the ped file should be set to zero
* For unrelated individuals, you can either (1) create a family for each unrelated individual as a founder or (2) put all unrelated individuals as founders in a single family.

== Other options ==

Some of command line options are explained below and others are self-explanatory.

-c : minimum cutoff of posterior probability to output a variant [''Default: 0.5'']
--theta : scaled mutation rate per site [''Default: 0.001'']
--tstv: prior of ts:tv ratio [''Default: 2.0'']
--nthreads : number of threads to run and it is recommended to use 4 threads for small number of input files [''Default: 1']

--denovo : a boolean flag to turn on ''de novo'' mutation detection. The following options take effect only when this flag is ON
--rate_denovo : mutation rate per haplotype per generation. [''Default: 1.5e-08'']
--tstv_denovo : the prior ts/tv ratio of ''de novo'' mutations. [''Default: 2.0'']
--minLLR_denovo : minimum value of log10 likelihood ratio of allowing vs. disallowing ''de novo'' mutations in the data to output [''Default: 1.0'']

--chr2process: the chromosome names to process. Default is empty and is to process all chromosomes in the input.
If multiple chromosomes are provided, they should be separated by comma, e.g. --chr2process 2,10 or --chr2process chr2,chr10

--gl_off: not to output genotype likelihood values for each individual. Default is to output 3 GLs for polymorphisms and 10 GLs for de novo mutations

== Output files ==
* The output file is a VCF file and the specification can be found [[http://www.1000genomes.org/node/101 here]]
* Since there is no standard to represent ''de novo'' mutations in the current VCF specification, actual genotypes (e.g. [ACGT]/[ACGT]) are output in the VCF file for ''de novo'' mutations.
* A summary about variant calling statistics is output to STDOUT and it may be redirected to a file for a record.

Summary of reference -- 9
Total Entry Count: 141213431
Total Base Cout: 120124735
Total '0' Base Count: 137
Non-Polymorphic Count: 655457
Transition Count: 6556
Transversion Count: 3127
Other Polymorphism Count: 0
Filter counts:
minMapQual 4550
minTotalDepth 1089
maxTotalDepth 736
Hard to call: 0
Skipped bases: 134

== Creation of GLF files ==

* The current version performs variant calling and ''de novo'' mutation detection from files in the genotype likelihood format (GLF). In future versions we plan to take [[http://samtools.sourceforge.net/ SAM/BAM]] files as input. See the following for instructions on how to create GLF files.
** Download a [[Samtools-hybrid|modified version of samtools]] ( [[https://github.com/statgen/samtools-0.1.7a-hybrid samtools-hybrid]] )
** Prepare the reference genome in fasta format and sequence alignments in [[http://samtools.sourceforge.net/ SAM/BAM]] format
** Generate BAQ adjusted GLF files using the following command
samtools-hybrid view -bh chr1.bam 1:0 | samtools-hybrid calmd -Abr - human.v37.fa 2> /dev/null | samtools-hybrid pileup - -g -f human.v37.fa > chr1.bam.glf
** If you want to clip overlapping reads, you can add clipOverlap command
samtools-hybrid view -bh chr1.bam 1:0 |bam clipOverlap --in -.bam --out -.ubam| samtools-hybrid calmd -Abr - human.v37.fa 2> /dev/null | samtools-hybrid pileup - -g -f human.v37.fa > chr1.bam.glf
* For other functionalities please refer to the [[http://samtools.sourceforge.net/ samtools]] website.

== Download ==
The latest version of source code v0.04 with test files can be [[Media:polymutt.0.04.tar.gz | downloaded]] here.

== Contact ==
For questions please contact the authors (Bingshan Li: [mailto:bingshan@umich.edu bingshan@umich.edu] or Goncalo Abecasis: [mailto:goncalo@umich.edu goncalo@umich.edu])

[[Category:Software]]

ContaminationDetection

2012-04-24T21:19:06Z

Gjun:

== Overview ==

DNA sample contamination (or sample swap) is quite common in sequencing studies. It is possible to detect DNA contamination from sequence reads of intensity data from SNP genotyping arrays by modeling the likelihood of sequence reads as a mixture of two samples and estimating the fraction of reads contributed by contaminating samples.

Depending on the types of available data, one can detect DNA sample contamination in the following ways.

# When a sequenced sample has also external array-based genotypes available
#* [[VerifyBamID]] software can estimate the levels of sample contamination and detect sample swap from alignment sequence reads and the external genotypes
#* The basic mathematical concept is described in [[Verifying_Sample_Identities_-_Implementation]]
# When a sequenced sample does not have external genotypes available
#* [[VerifyBamID]] software can still estimate sample contamination from aligned sequence reads and population minor allele frequency
#* The key idea of the method is to capture the excessive heterozygosity in the contaminated sample by modeling the sequence reads as mixture of independent samples.
# When a sample has array-based genotypes but not yet sequenced
#* [[VerifyIDintensity]] software can estimate the levels of DNA sample contamination from pre-computed intensity data using likelihood-based model.
#** The key idea is similar to that of [[VerifyBamID]] with sequence data alone
#** [[VerifyIDintensity]] models array intensity data instead of sequence reads
#** When a large number of samples are genotyped together, the 'multi-sample' option will estimate the intensity distribution for each marker across multiple samples. When only one a few samples are genotyped, 'per-sample' option will estimate the intensity distribution for each sample across all markers.
#* [[BAFRegress]] software can estimate levels of sample contamination and test the presence of contamination by regressing the allele frequency with respect to the B allele frequency (BAF).
#** This method is sensitive to estimate low levels of contamination.
#** The method provides a well calibrated Type I error to test the null hypothesis of no contamination from a DNA sample.

== History ==
* Initial description of method 1 can be found at [[Verifying_Sample_Identities_-_Implementation]] (last modified in April 29, 2010)

VerifyIDintensity

2012-04-24T21:16:20Z

Gjun: Created page with 'VerifyIDintensity '''verifyIDintensity''' is a software that detects and estimates sample contamination using intensity data from Illumina genotyping array…'

[[Category:Software|VerifyIDintensity]]

'''verifyIDintensity''' is a software that detects and estimates sample contamination using intensity data from Illumina genotyping arrays using a mixture model.

== Download verifyIDintensity ==

* [[File:VerifyIDintensity.tgz]]

== Build verifyID intensity ==

To build verifyIDintensity, run the following series of commands. You need [http://www.boost.org/ boost library] and [http://tclap.sourceforge.net/ tclap].
tar xzvf verifyIDintensity.tgz
make

== Basic Usage ==

verifyIDintensity [-t <float>] [-m <int>] -n <int> [-b <string>] [-s <string>] -i <string> [-v] [-p] [--] [--version] [-h]

== Options ==

-t <float>, --threshold <float>
Minimum allele frequency for likelihood estimation, default is 0.01

-m <int>, --marker <int>
(required) Number of markers

-n <int>, --number <int>
(required) Number of samples

-b <string>, --abf <string>
Allele frequency file (ABF), which is a plain text file with SNP_ID and Allele_B frequency.
SNP_IDs should be sorted in the same order as the intensity file

-s <string>, --stat <string>
Statistics file (created if not exist)

-i <string>, --in <string>
(required) Input pre-computed intensity (.adpc.bin) file

-v, --verbose
Turn on verbose mode

-p, --persample
Do per-sample analysis, default is per-marker analysis

--, --ignore_rest
Ignores the rest of the labeled arguments following this flag.

--version
Displays version information and exits.

-h, --help
Displays usage information and exits.

File:VerifyIDintensity.tgz

2012-04-24T21:12:43Z

Gjun:

ContaminationDetection

2012-04-24T20:48:21Z

Gjun: /* Overview */

== Overview ==

DNA sample contamination (or sample swap) is quite common in sequencing studies. It is possible to detect DNA contamination from sequence reads of intensity data from SNP genotyping arrays by modeling the likelihood of sequence reads as a mixture of two samples and estimating the fraction of reads contributed by contaminating samples.

Depending on the types of available data, one can detect DNA sample contamination in the following ways.

# When a sequenced sample has also external array-based genotypes available
#* [[VerifyBamID]] software can estimate the levels of sample contamination and detect sample swap from alignment sequence reads and the external genotypes
#* The basic mathematical concept is described in [[Verifying_Sample_Identities_-_Implementation]]
# When a sequenced sample does not have external genotypes available
#* [[VerifyBamID]] software can still estimate sample contamination from aligned sequence reads and population minor allele frequency
#* The key idea of the method is to capture the excessive heterozygosity in the contaminated sample by modeling the sequence reads as mixture of independent samples.
# When a sample has array-based genotypes but not yet sequenced
#* [[VerifyIDintensity]] software can estimate the levels of DNA sample contamination from aligned sequence reads and population allele frequency
#** The key idea is similar to that of [[VerifyBamID]] with sequence data alone
#** [[VerifyIDintensity]] models array intensity data instead of sequence reads
#** When a large number of samples are genotyped together, the 'multi-sample' option will estimate the intensity distribution for each marker across multiple samples. When only one a few samples are genotyped, 'per-sample' option will estimate the intensity distribution for each sample across all markers.
#* [[BAFRegress]] software can estimate levels of sample contamination and test the presence of contamination by regressing the allele frequency with respect to the B allele frequency (BAF).
#** This method is sensitive to estimate low levels of contamination.
#** The method provides a well calibrated Type I error to test the null hypothesis of no contamination from a DNA sample.

== History ==
* Initial description of method 1 can be found at [[Verifying_Sample_Identities_-_Implementation]] (last modified in April 29, 2010)