Provided by: qtltools_1.3.1+dfsg-2build2_amd64 
NAME
QTLtools ase - Measure ASE from RNA-seq
SYNOPSIS
QTLtools ase --bam [sample.bam|sample.sam|sample.cram] --vcf [in.vcf|in.bcf|in.vcf.gz]
--sample sample_name_in_vcf --mapq integer --out output_file_prefix [OPTIONS]
DESCRIPTION
This mode measures allele specific expression (ASE) from RNAseq for transcribed
heterozygous SNPs as detailed in <https://www.ncbi.nlm.nih.gov/pmc/articles/PMC3918453/>.
In brief, the reference allele mapping bias is calculated for each of the 12 REF/ALT pairs
separately provided that there are at least --sites number of REF/ALT sites that have a
minimum --cov-bias number of reads overlapping. For REF/ALT pairs that fail these
criteria, reference allele mapping bias is calculated from all sites. Reference allele
mapping bias is the number of all reads across many sites that contained the reference
allele for a given REF/ALT pair over the total number of reads overlapping. --subsample
controls which percentile of the highest covered sites are subsampled, which is necessary
so that the reference allele mapping bias is not estimated mostly from very high coverage
sites, but from sites covering all the genome. The ASE p-value for each site is then
calculated as a two-tailed binomial test checking if the observed number of reference
allele reads is significantly different from random, given the total number of reads, and
the probability of observing a reference allele, which is the reference allele mapping
bias for a certain REF/ALT pair.
The defaults for options should work well for most RNAseq experiments. It is NOT
advisable to decrease the --baseq below 10, --cov and --cov-bias below 8, and setting
--subsample to 1. It is NOT recommended to use the following options regarding filtering
--keep-bans-for-bias, --keep-discordant-for-bias, --filter-duplicates, --ignore-
orientation, and --legacy-options. Also refrain from using --auto-flip, but rather create
correct VCF/BCF files.
We highly recommended using a BCF file rather than a VCF due to performance benefits,
--fasta to provide the genome sequence used, --blacklist to filter low mappability
regions, and --gtf to annotate the SNPs. If you are using unfiltered imputed genotypes
then consider using --imp-qual and --geno-prob. For trouble shooting purposes you can use
--filtered, which will list why certain variants are filtered from the analysis. If you
are I/O bound you may try using --group-by for better performance. If your blacklist file
contains many overlapping or contiguous regions you can decrease the memory usage with
--merge-on-the-fly.
OPTIONS
--vcf, -v [in.vcf|in.bcf|in.vcf.gz]
Genotypes in VCF/BCF format sorted by chromosome and then position. Can contain
multiple samples. REQUIRED and RECOMMENDED to use a BCF file for performance.
--bam, -b [in.bam|in.sam|in.cram]
Sequence data in BAM/SAM/CRAM format sorted by chromosome and then position. One
sample per BAM file. REQUIRED.
--ind, -i sample_name
Sample name in the VCF corresponding to the BAM file. REQUIRED.
--out, -o output
Output prefix. This will generate output.ase and output.ref_bias files. Or if you
give output.gz or output.bz2 then output.ase.gz etc. REQUIRED.
--mapq, -q integer
Minimum mapping quality for a read or read pair to be considered. Set this to only
include uniquely mapped reads. REQUIRED.
--fasta, -f genome_sequence.fa
Genome sequence in FASTA format. Make sure this the sequence for the correct
genome build. RECOMMENDED.
--blacklist, -B poor_mappability_regions.bed
Poor mappability regions in BED format. ASE estimates will be unreliable in
regions with low mappability. RECOMMENDED.
--merge-on-the-fly, -t
Merges the blacklisted regions on the fly. This can reduce memory usage if there
are many overlapping or contiguous regions in the blacklist.
--gtf, -g gene_annotation.gtf
Gene annotations in GTF format. These can be obtained from
<https://www.gencodegenes.org/>. RECOMMENDED.
--filtered, -l filename
File to output filtered variants. RECOMMENDED for troubleshooting especially if
--suppress-warnings.
--reg, -r chr:start-end
Genomic region to be processed. E.g. chr4:12334456-16334456, or chr5
--fix-chr, -F
Attempt to match chromosome names to the BAM file by adding or removing chr to
chromosome names. Does not apply to --{include,exclude}-positions options. These
should be in the VCF chromosome names.
--fix-id, -R
Convert missing VCF variant IDs to chr_pos_refalt.
--auto-flip, -x
Attempt to fix reference allele mismatches. Requires a fasta file for the
reference sequence by --fasta. NOT RECOMMENDED.
--print-stats, -P
Print out stats for the filtered reads for ASE sites. This will be slower and if
there are sites with more reads than --max-depth, then will potentially give
different results.
--suppress-warnings, -k
Suppress the warnings about individual variants.
--illumina13, -j
Base quality is in the Illumina-1.3+ encoding.
--group-by, -G integer
Group variants separated by this much into batches. This allows you not to stream
the whole BAM file and may improve running time.
--max-depth, -d integer
Pileup max-depth DEFAULT=1000000. Set to more reasonable value if you experience
memory or performance issues. This is set to a high value by default since with
RNA-seq you can have very high coverage sites. Set to 0 for max.
--baseq, -Q integer
Minimum phred quality for a base to be considered DEFAULT=10.
--pvalue, -p float
Binomial p-value threshold for ASE output DEFAULT=1.0.
--cov, -c integer
Minimum coverage for a genotype to be considered in ASE analysis DEFAULT=16.
--cov-bias, -C integer
Minimum coverage for a genotype to be considered in reference allele mapping bias
analysis DEFAULT=10.
--sites, -s integer
Minimum number of sites to calculate a reference allele mapping bias from for a
specific REF/ALT pair DEFAULT=200. The reference allele mapping bias for pairs
with less than this many sites will be calculated from all sites.
--imp-qual-id, -I string
The VCF INFO field ID of the imputation score in the VCF DEFAULT=INFO.
--geno-prob-id, -L string
The VCF FORMAT field ID of the genotype posterior probabilities for RR/RA/AA in the
VCF DEFAULT=GP.
--imp-qual, -W float
Minimum imputation score for a variant to be considered DEFAULT=0.0.
--geno-prob, -V float
Minimum posterior probability for a genotype to be considered DEFAULT=0.0.
--subsample, -S float
Randomly subsample sites that have greater coverage than this percentile of all the
sites in reference allele mapping bias calculations DEFAULT=0.75. Set to 1 to turn
off which is NOT RECOMMENDED. Set to 0 to subsample all sites to --cov-bias.
--both-alleles-seen, -a
Require both alleles to be observed in RNA-seq reads for a site for ASE
calculations.
--keep-bans-for-bias, -A
DON'T require both alleles to be observed in RNAseq reads for a site for reference
allele mapping bias calculations. NOT RECOMMENDED.
--keep-discordant-for-bias, -E
If given, sites with more discordant alleles than REF or ALT alleles will be
included in the reference allele mapping bias bias calculations. NOT RECOMMENDED.
--filter-indel-reads, -D
Remove reads that contain indels.
--keep-failed-qc, -e
Keep fastq reads that fail sequencing QC (as indicated by the sequencer).
--keep-orphan-reads, -O
Keep paired end reads where one of mates is unmapped.
--check-proper-pairing, -y
If provided only properly paired reads according to the aligner will be considered.
--ignore-orientation, -X
If NOT provided only mate pairs where both mates are on the same chromosome and
where the first mate is on the +ve strand and the second is on the -ve strand will
be considered. NOT RECOMMENDED.
--filter-duplicates, -u
Remove reads designated as duplicate by the aligner. NOT RECOMMENDED.
--filter-supp, -m
Remove supplementary (non-linear) alignments.
--legacy-options, -J
Replicate legacy options used. NOT RECOMMENDED.
OUTPUT FILES
.ase
This file is the main ASE results file with the following columns. Columns after the
23rd column are only printed if --print-stats is provided.
1 INDIVIDUAL The sample id
2 RSID The SNP ID from the VCF file
3 CHR Chromosome of the SNP
4 POS Position of the SNP
5 ALLELES The SNP's alleles
6 BOTH_ALLELES_SEEN Whether or not both of the SNP's alleles were seen in the
RNAseq reads
7 MIN_ALLELE_RATIO The minor allele ration among RNAseq reads
8 REF_COUNT Number of reference alleles in reads
9 NONREF_COUNT Number of alternative alleles in reads
10 TOTAL_COUNT Number of read overlapping the SNP
11 WEIGHTED_REF_COUNT REF_COUNT adjusted for ref mapping bias
12 WEIGHTED_NONREF_COUNT NONREF_COUNT adjusted for ref mapping bias
13 WRC_MINUS_WNC WEIGHTED_REF_COUNT - WEIGHTED_NONREF_COUNT
14 ALLELES_SEEN Alleles observed in RNAseq reads
15 REF_ALLELE Reference allele
16 ALT_ALLELE Alternative allele
17 OTHER_COUNT Number of discordant reads (not REF or ALT)
18 EXPECTED_DISCORDANT Expected number of discordant alleles. Calculated by adding
up all the base error probabilities of RNAseq read positions
overlapping the SNP
19 DISCORDANT_PVAL P-value for observed number of discordant reads being more
than expected. One can Bonferroni correct these p-values,
and exclude the significant ones from downstream analyses.
20 REF_RATIO Reference allele mapping bias
21 PVALUE ASE p-value
22 CONCERN If there were any concerns potentially rendering this SNP
unusable, they will be coded here. See below OUTPUT FILE
CODES for how to decode
23 EXON_INFO Exons that overlap this SNP will be listed here if a GTF
file is provided. Exon names, which are formed by
concatenating gene id, transcript id, exon position, and
gene name delimited with a colon, are separated with
semicolons
24 SECONDARY The number of secondary alignments that were filtered out
25 SUPPLEMENTARY If --filter-supp is provided, then the number of
supplementary alignments that were filtered out
26 FAIL_MAPQ The number of alignments that were filtered out due to low
mapping quality
27 FAILQC If --keep-failed-qc is not provided, then the number of
filtered reads that failed qc according to the sequencer
28 DUPLICATE If --filter-duplicates is provided, then the number of
filtered reads that were labeled as duplicate by the aligner
29 MATE_UNMAPPED If --keep-orphan-reads is not provided, then the number of
reads that were filtered since one mate was unmapped
30 WRONG_ORIENTATION If --ignore-orientation is not provided, then the number of
reads that were filtered since they were in the wrong
orientation
31 NOT_PROPER_PAIR If --check-proper-pairing is provided, then the number of
filtered reads that were not properly paired according to
the aligner
32 SKIPPED The number of alignments that were filtered out since they
did not have actual bases overlapping the SNP
33 FAIL_BASEQ The number of filtered reads where the base overlapping the
SNP had a base quality less than --baseq
34 INDEL If --filter-indel-reads is provided, then the number of
filtered reads that contained indels
35 DEPTH The number alignments that overlap with the SNP position
.ref_bias
Details the reference allele mapping bias results, and has the following columns:
1 INDIVIDUAL The sample id
2 ALLELES The reference alternative allele pair
3 REF_ALL Total number of reference alleles observed across all sites
4 NONREF_ALL Total number of alternative alleles observed across all sites
5 SITES Number of sites for this REF/ALT pair that pass the thresholds
6 SUBSAMPLED_SITES Number of sites that were subsampled to SUBSAMPLED_TO since they
had too high a coverage
7 SUBSAMPLED_TO The coverage SUBSAMPLED_SITES were subsampled to
8 PERC Reference allele mapping bias
9 SOURCE Whether the reference allele mapping bias was calculated from
allele specific sites of all sites
--filtered filename
This file lists the variants that were omitted from the analysis and why. The first
column is the coded omission reason, followed by the variant ID. The codes and their
meaning are listed in the following section.
OUTPUT FILE CODES
.ase file codes for the CONCERN column
You probably want to exclude SNPs with RM, NRA, MDTA, or MDTR concerns from your analyses
as these likely have wrong genotypes. SNPs with other concern may be OK
RM Reference allele in the VCF file mismatches the reference sequence. Requires
--fasta
DP Multiple variants observed at the same position in the VCF file
NRA No reference or alternative allele observed in the RNAseq reads
MDTA More discordant (not reference or alternative) alleles than alternative alleles in
the RNAseq reads
MDTR More discordant (not reference or alternative) alleles than reference alleles in
the RNAseq reads
BANS Both alleles of the SNP were not observed in RNAseq reads
LMAR Among the RNAseq reads the minor allele ratio was less than 2%
PD Pileup depth (--max-depth) was potentially exceeded around this SNP. This may
prevent some overlapping reads from being counted. We recommend rerunning with a
higher --max-depth.
--filtered file codes
VMA Multi-allelic variant
VU Variant position or ID is excluded by the user
VCNIB Variant chromosome is not in the BAM file
VB Variant is in a blacklisted region
VI Variant is an indel
VMRA Variant is missing either the reference or the alternative allele
VNIF Variant is not in the reference genome. Only if --fasta was provided
VS Reference and alternative alleles were swapped (not excluded). Only if
--auto-flip and --fasta were provided
VF Reference and alternative alleles' strand was flipped (not excluded). Only if
--auto-flip and --fasta were provided
VWR Reference allele does not match the reference sequence. Only if --auto-flip and
--fasta were provided
VBI Variant failed the imputation quality filter
VMGT Variant with missing GT field or variant is not diploid
VMG Missing genotype
VH Homozygous variant
VBG Variant failed the genotype probability filter
VD Duplicate variant position
VMI Variant with a missing ID (not excluded). Will be renamed if --fix-id is
provided
VDK Heterozygous variant has another variant with the same position (not excluded)
BC Variant did not have enough coverage for reference allele mapping bias
calculations
BBANS Both alleles were not observed for reference allele mapping bias calculations
BMDTRA Variant had more discordant alleles than reference or alternative alleles thus
was excluded from reference allele mapping bias calculations
AC Variant did not have enough coverage ASE calculations
ABANS Both alleles were not observed for ASE calculations
EXAMPLES
o ASE analysis of a sample mapped with STAR:
QTLtools ase --vcf multi_sample.bcf --bam sample1.bam --ind sample1 --mapq 255 --out
sample1 --filtered sample1.filtered.gz --gtf gencode.v19.annotation.gtf.gz --blacklist
poor_mappability_regions.bed --fasta hg19.fa --fix-chr --fix-id
SEE ALSO
QTLtools(1)
QTLtools website: <https://qtltools.github.io/qtltools>
BUGS
Please submit bugs to <https://github.com/qtltools/qtltools>
CITATION
Delaneau, O., Ongen, H., Brown, A. et al. A complete tool set for molecular QTL discovery
and analysis. Nat Commun 8, 15452 (2017). <https://doi.org/10.1038/ncomms15452>
AUTHORS
Halit Ongen (halitongen@gmail.com), Olivier Delaneau (olivier.delaneau@gmail.com)