src/jalview/schemes/HMMMatchScoreColourScheme.java

   1 package jalview.schemes;
   2
   3 import jalview.api.AlignViewportI;
   4 import jalview.datamodel.AnnotatedCollectionI;
   5 import jalview.datamodel.HiddenMarkovModel;
   6 import jalview.datamodel.SequenceI;
   7 import jalview.util.ColorUtils;
   8 import jalview.util.Comparison;
   9
  10 import java.awt.Color;
  11 import java.io.BufferedReader;
  12 import java.io.FileNotFoundException;
  13 import java.io.FileReader;
  14 import java.io.IOException;
  15 import java.util.ArrayList;
  16 import java.util.HashMap;
  17 import java.util.List;
  18 import java.util.Map;
  19
  20
  21
  22
  23 public class HMMMatchScoreColourScheme extends ResidueColourScheme
  24 {
  25
  26   private Map<Character, Map<Integer, Map<String, Double>>> probabilities;
  27
  28   private List<Integer> ranges;
  29
  30   private static double binSize;
  31
  32   public class MatchProbReader
  33   {
  34     private BufferedReader reader;
  35
  36     MatchProbReader() throws FileNotFoundException
  37     {
  38       reader = new BufferedReader(
  39               new FileReader("resources/ProbabilityOfMatch"));
  40     }
  41
  42     /*
  43     public Map<Character, Map<String, Double>>  getProbabilities() throws IOException
  44     {
  45       Map<Character, Map<String, Double>> probabilities = new HashMap<>();
  46
  47       String[] alphabet = reader.readLine().split("\\,");
  48       for (int i = 1; i < alphabet.length - 1; i++)
  49       {
  50         probabilities.put(alphabet[i].replaceAll("\\ ", "").charAt(0),
  51                 new HashMap<>());
  52       }
  53
  54       String line = reader.readLine();
  55       while(line != null)
  56       {
  57         String[] contents = line.split("\\,");
  58
  59         for(int i = 1; i < contents.length; i++)
  60         {
  61           probabilities.get(alphabet[i].replaceAll("\\ ", "").charAt(0))
  62                   .put(contents[0], Double
  63                           .valueOf(contents[i].replaceAll("\\ ", "")));
  64         }
  65         line = reader.readLine();
  66
  67       }
  68       reader.close();
  69       return probabilities;
  70     }
  71     */
  72
  73     public Map<Character, Map<Integer, Map<String, Double>>> getProbabilities()
  74             throws IOException
  75     {
  76
  77       Map<Character, Map<Integer, Map<String, Double>>> probabilities = new HashMap<>();
  78
  79       ranges = new ArrayList<>();
  80       ranges.add(0);
  81
  82       binSize = Double.valueOf((reader.readLine().replaceAll("\\ ", "")));
  83       String line = reader.readLine();
  84       char c = line.charAt(0);
  85
  86       while (line != null)
  87       {
  88         line = reader.readLine();
  89         while (line != null && line.split("\\,").length != 1)
  90         {
  91           String[] llrs = line.split("\\,");
  92           String[] counts = reader.readLine().split("\\,");
  93           int range = Integer.valueOf(llrs[0]);
  94
  95           if (!ranges.contains(range))
  96           {
  97             ranges.add(range);
  98           }
  99           if (!probabilities.containsKey(c))
 100           {
 101             probabilities.put(c, new HashMap<>());
 102           }
 103           probabilities.get(c).put(range, new HashMap<>());
 104
 105           for (int i = 1; i < llrs.length; i++)
 106           {
 107             probabilities.get(c).get(range).put(
 108                     llrs[i].replaceAll("\\ ", ""),
 109                     Double.valueOf(counts[i].replaceAll("\\ ", "")));
 110           }
 111
 112           line = reader.readLine();
 113         }
 114         if (line != null)
 115         {
 116           c = line.charAt(0);
 117         }
 118       }
 119
 120       return probabilities;
 121     }
 122
 123   }
 124
 125   private static final Color INSERTION_COLOUR = Color.white;
 126
 127   /*
 128    * the aligned HMM consensus sequence to use as reference for colouring
 129    */
 130   private SequenceI hmmSeq;
 131
 132   private HiddenMarkovModel hmm;
 133
 134   private Map<Character, Float> frequencies;
 135
 136   /**
 137    * Constructor given a list of Hidden Markov Model consensus sequences. The
 138    * first sequence provides the HMM profile from which we can read the emission
 139    * probabilities that determine the colour.
 140    *
 141    * @param hmmSeqs
 142    * @throws IOException
 143    */
 144   public HMMMatchScoreColourScheme(List<SequenceI> hmmSeqs)
 145   {
 146     hmmSeq = hmmSeqs.isEmpty() ? null : hmmSeqs.get(0);
 147     hmm = hmmSeq == null ? null : hmmSeq.getHMM();
 148
 149     try
 150     {
 151       MatchProbReader probabilityReader = new MatchProbReader();
 152       probabilities = probabilityReader.getProbabilities();
 153     } catch (IOException e)
 154     {
 155       System.out.println(e.getStackTrace());
 156     }
 157   }
 158
 159   /**
 160    * Default constructor (required by ColourSchemes.loadColourSchemes)
 161    */
 162   public HMMMatchScoreColourScheme()
 163   {
 164   }
 165
 166   @Override
 167   public Color findColour(char symbol, int column, SequenceI seq,
 168           String consensusResidue, float pid)
 169   {
 170     if (seq == null)
 171     {
 172       return null;
 173     }
 174     return findColour(symbol, column, seq.gapMap().length);
 175   }
 176
 177   // TODO change documentation
 178   /**
 179    * Returns the colour at a particular symbol at a column in the alignment:
 180    * <ul>
 181    * <li>white for a gap</li>
 182    * <li>red for an insertion</li>
 183    * <li>orange for negative information content</li>
 184    * <li>white to blue for increasing information content</li>
 185    * </ul>
 186    *
 187    * @param symbol
 188    * @param column
 189    * @return
 190    */
 191   private Color findColour(char symbol, int column, int length)
 192   {
 193     if (getHmm() == null || Comparison.isGap(symbol))
 194     {
 195       return Color.white;
 196     }
 197     if (Comparison.isGap(hmmSeq.getCharAt(column)))
 198     {
 199       return INSERTION_COLOUR;
 200     }
 201     if (Character.isLowerCase(symbol))
 202     {
 203       symbol = Character.toUpperCase(symbol);
 204     }
 205
 206     double prob = 0;
 207     if (hmm.getBackgroundFrequencies().containsKey(symbol))
 208     {
 209       int lengthBin = getLengthBin(length);
 210
 211       double llr = Math
 212               .log(getHmm().getMatchEmissionProbability(column, symbol)
 213                       / hmm.getBackgroundFrequencies().get(symbol));
 214
 215       if (!probabilities.get(symbol).get(lengthBin)
 216               .containsKey(format(llr)))
 217       {
 218         return Color.white;
 219       }
 220
 221
 222       prob = probabilities.get(symbol).get(lengthBin).get(format(llr));
 223     }
 224     else
 225     {
 226       return new Color(140, 140, 140);
 227     }
 228
 229     Color colour = Color.ORANGE;
 230     if (prob >= 0.5)
 231     {
 232
 233       colour = ColorUtils.getGraduatedColour((float) prob, 0.5f,
 234               Color.WHITE, 1f,
 235               Color.blue);
 236     }
 237     else
 238     {
 239       colour = ColorUtils.getGraduatedColour((float) prob, 0f, Color.red,
 240               0.5f, Color.WHITE);
 241     }
 242
 243     return colour;
 244   }
 245
 246   public static String format(Double d)
 247   {
 248     String formatArg = String.valueOf(binSize);
 249
 250     // if bin size, need format "%.n" where n is number of decimal places
 251     if (binSize < 1)
 252     {
 253       formatArg = "." + formatArg.split("\\.")[1].length();
 254     }
 255
 256     Double rounded = Math.round(d / binSize) * binSize;
 257     String formatted = String.format("%" + formatArg + "f", rounded);
 258
 259     // format sometimes returns a number rounded to 0 as -0
 260     // this ensures output will always be 0
 261     if (Double.valueOf(formatted) == 0)
 262     {
 263       formatted = "0";
 264     }
 265     return formatted;
 266
 267   }
 268
 269   /**
 270    * Answers the maximum possible value of information score (log ratio), for use
 271    * in scaling a graduated colour range
 272    *
 273    * @return
 274    */
 275   protected float getMaxInformationScore()
 276   {
 277     return 0f;
 278   }
 279
 280   /**
 281    * Answers a new colour scheme instance based on the HMM of the first sequence
 282    * in ac that has an HMM
 283    */
 284   @Override
 285   public ColourSchemeI getInstance(AlignViewportI viewport,
 286           AnnotatedCollectionI ac)
 287   {
 288     return newInstance(ac);
 289   }
 290
 291   /**
 292    * Constructor given a sequence collection
 293    *
 294    * @param ac
 295    */
 296   public HMMMatchScoreColourScheme(AnnotatedCollectionI ac)
 297   {
 298     this(ac.getHmmSequences());
 299   }
 300
 301   /**
 302    * Answers a new instance of the colour scheme for the given HMM
 303    *
 304    * @param ac
 305    * @return
 306    */
 307   protected HMMMatchScoreColourScheme newInstance(AnnotatedCollectionI ac)
 308   {
 309     return new HMMMatchScoreColourScheme(ac);
 310   }
 311
 312   @Override
 313   public boolean isSimple()
 314   {
 315     return false;
 316   }
 317
 318   /**
 319    * Answers true if the sequence collection has an HMM consensus sequence and
 320    * that the first HMM sequence contains background frequencies, else false
 321    */
 322   @Override
 323   public boolean isApplicableTo(AnnotatedCollectionI ac)
 324   {
 325     return !ac.getHmmSequences().isEmpty() && ac.getHmmSequences().get(0)
 326             .getHMM().getBackgroundFrequencies() != null;
 327   }
 328
 329   protected Map<Character, Float> getFrequencies()
 330   {
 331     return frequencies;
 332   }
 333
 334   protected void setFrequencies(Map<Character, Float> frequencies)
 335   {
 336     this.frequencies = frequencies;
 337   }
 338
 339   protected HiddenMarkovModel getHmm()
 340   {
 341     return hmm;
 342   }
 343
 344   protected SequenceI getHmmSequence()
 345   {
 346     return hmmSeq;
 347   }
 348
 349   @Override
 350   public String getSchemeName()
 351   {
 352     return JalviewColourScheme.HMMMatchScore.toString();
 353   }
 354
 355   private int getLengthBin(int l)
 356   {
 357     for (int i = 1; i < ranges.size(); i++)
 358     {
 359       if (l >= ranges.get(i - 1) && l < ranges.get(i))
 360       {
 361         return ranges.get(i);
 362       }
 363     }
 364     return -1;
 365   }
 366 }
 367
 368