src/jalview/schemes/HMMERAlignmentColourScheme.java

   1 package jalview.schemes;
   2
   3 import jalview.datamodel.AnnotatedCollectionI;
   4 import jalview.datamodel.HiddenMarkovModel;
   5 import jalview.datamodel.ResidueCount;
   6 import jalview.datamodel.SequenceCollectionI;
   7 import jalview.datamodel.SequenceI;
   8 import jalview.util.ColorUtils;
   9 import jalview.util.Comparison;
  10
  11 import java.awt.Color;
  12 import java.util.HashMap;
  13 import java.util.List;
  14 import java.util.Map;
  15
  16 public class HMMERAlignmentColourScheme extends ResidueColourScheme
  17 {
  18   Map<Character, Double> backgroundFrequencies = new HashMap<>();
  19
  20   Double maxLLR;
  21
  22   HiddenMarkovModel hmm;
  23
  24   boolean peptideSpecific;
  25
  26   boolean nucleotideSpecific;
  27
  28   public HMMERAlignmentColourScheme(HiddenMarkovModel markov)
  29   {
  30     hmm = markov;
  31   }
  32
  33   public HMMERAlignmentColourScheme()
  34   {
  35
  36   }
  37
  38   @Override
  39   public Color findColour(char symbol, int position, SequenceI seq,
  40           String consensusResidue, float pid)
  41   {
  42     if (hmm == null)
  43     {
  44       return Color.white;
  45     }
  46     return findColour(symbol, position);
  47   }
  48
  49   /**
  50    * Returns the colour at a particular symbol at a column in the alignment.
  51    *
  52    * @param symbol
  53    * @param position
  54    * @return Red for an insertion, white for a gap, orange for a negative
  55    *         information content, white to blue for increasing information
  56    *         content.
  57    */
  58   private Color findColour(char symbol, int position)
  59   {
  60
  61     if (Comparison.isGap(symbol))
  62     {
  63       return Color.white;
  64     }
  65     if (Character.isLowerCase(symbol))
  66     {
  67       symbol = Character.toUpperCase(symbol);
  68     }
  69     Double prob;
  70     prob = hmm.getMatchEmissionProbability(position, symbol);
  71     Double freq = backgroundFrequencies.get(symbol);
  72     if (freq == null)
  73     {
  74       return Color.white;
  75     }
  76     if (prob == 0)
  77     {
  78       return new Color(230, 0, 0);
  79     }
  80     Double value = Math.log(prob / freq);
  81     Color colour = null;
  82     if (value > 0)
  83     {
  84
  85       colour = ColorUtils.getGraduatedColour(value.floatValue(), 0,
  86               Color.WHITE, maxLLR.floatValue(), Color.blue);
  87     }
  88     else if (value < 0)
  89     {
  90       return Color.ORANGE;
  91
  92     }
  93     return colour;
  94   }
  95
  96   @Override
  97   public void alignmentChanged(AnnotatedCollectionI collection,
  98           Map<SequenceI, SequenceCollectionI> hiddenReps)
  99   {
 100     List<SequenceI> seqs = collection.getSequences();
 101     for (SequenceI seq : seqs)
 102     {
 103       if (seq.isHMMConsensusSequence())
 104       {
 105         hmm = seq.getHMM();
 106         break;
 107       }
 108     }
 109
 110     count(collection);
 111
 112   }
 113
 114   @Override
 115   public ColourSchemeI getInstance(AnnotatedCollectionI sg,
 116           Map<SequenceI, SequenceCollectionI> hiddenRepSequences)
 117   {
 118     HiddenMarkovModel markov = null;
 119     List<SequenceI> seqs = sg.getSequences();
 120     for (SequenceI seq : seqs)
 121     {
 122       if (seq.getHMM() != null)
 123       {
 124         markov = seq.getHMM();
 125         break;
 126       }
 127     }
 128
 129     count(sg);
 130
 131     HMMERAlignmentColourScheme colour = new HMMERAlignmentColourScheme(
 132             markov);
 133     return colour;
 134
 135   }
 136
 137   @Override
 138   public boolean isApplicableTo(AnnotatedCollectionI ac)
 139   {
 140     return true;
 141
 142   }
 143
 144   @Override
 145   public String getSchemeName()
 146   {
 147
 148     return JalviewColourScheme.HMMERA.name();
 149   }
 150
 151   @Override
 152   public boolean isSimple()
 153   {
 154     return false;
 155   }
 156
 157   public void count(AnnotatedCollectionI sg)
 158   {
 159     ResidueCount counts = new ResidueCount();
 160     for (SequenceI seq : sg.getSequences())
 161     {
 162       for (int i = 0; i < seq.getLength(); i++)
 163       {
 164         if (!Comparison.isGap(seq.getCharAt(i)))
 165         {
 166           counts.add(seq.getCharAt(i));
 167         }
 168       }
 169     }
 170
 171     int total = counts.getTotalCount();
 172
 173     for (char symbol : counts.getSymbolCounts().symbols)
 174     {
 175       double count = Double.valueOf(counts.getCount(symbol))
 176               / Double.valueOf(total);
 177       backgroundFrequencies.put(symbol, count);
 178     }
 179     maxLLR = Math.log(total);
 180   }
 181 }
 182