htswanalysis/src/SRLib/nuc.h

   1 #include <iostream>
   2
   3 #ifndef NUC_H
   4 #define NUC_H
   5
   6 class Nuc {
   7   protected:
   8     unsigned int n[4];
   9
  10     //background nucleotide probabilities
  11     //
  12     //can change this to a background model class later if needed
  13     static double qA;
  14     static double qC;
  15     static double qG;
  16     static double qT;
  17
  18     // pseudocount to avoid divide-by-zero errors
  19     static double pseudocount;
  20     static unsigned int strong_coverage_threshold;
  21
  22     // protected access directly to the data arrray
  23     unsigned int& operator[](size_t b);
  24
  25   public:
  26
  27     Nuc();
  28     Nuc(const Nuc& n);
  29     Nuc& operator=(const Nuc& n);
  30
  31     void add_nuc(char b);
  32
  33     //nth_nuc
  34     //=======
  35     //Useful to get the listing of all the nucleotides sequenced for this position.
  36     //Note that this function doesn't return the nucleotide in the order presented
  37     char nth_nuc(unsigned int i);
  38
  39     //N
  40     //====
  41     //returns the number of nucleotides at the position
  42     unsigned int N();
  43
  44     //at
  45     //==
  46     //allows read-only access to the dta array
  47     unsigned int at(size_t b) const;
  48     unsigned int A() const;
  49     unsigned int C() const;
  50     unsigned int G() const;
  51     unsigned int T() const;
  52
  53     //RE
  54     //==
  55     //returns the relative entropy of the nucleotide, taking into
  56     //account the background model and pseudocounts
  57     double RE();
  58
  59     //chastity
  60     //==
  61     //returns 2 x the ratio of the most prominent nucleotide to the first two nucleotides.
  62     //i.e. a homozygous position has chastity 1, and a heterozygous position has
  63     //     chastity ~0.5
  64     //note, the value is doubled to put it on the same scale and with roughly the same
  65     //interpretation as relative entropy (RE)
  66     double chastity();
  67
  68     //consensus:
  69     //=========
  70     //returns the consensus sequence, including IUPAC degenerate bases to represent heterozygosity
  71     char consensus();
  72 };
  73
  74 #endif
  75