Provided by: qtltools_1.3.1+dfsg-4_amd64 bug

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)