--- /dev/null
+#include "muscle.h"\r
+#include "msa.h"\r
+\r
+static int Blosum62[23][23] =\r
+ {\r
+// A B C D E F G H I K L M N P Q R S T V W X Y Z \r
+ +4, -2, +0, -2, -1, -2, +0, -2, -1, -1, -1, -1, -2, -1, -1, -1, +1, +0, +0, -3, -1, -2, -1, // A\r
+ -2, +6, -3, +6, +2, -3, -1, -1, -3, -1, -4, -3, +1, -1, +0, -2, +0, -1, -3, -4, -1, -3, +2, // B\r
+ +0, -3, +9, -3, -4, -2, -3, -3, -1, -3, -1, -1, -3, -3, -3, -3, -1, -1, -1, -2, -1, -2, -4, // C\r
+ -2, +6, -3, +6, +2, -3, -1, -1, -3, -1, -4, -3, +1, -1, +0, -2, +0, -1, -3, -4, -1, -3, +2, // D\r
+ -1, +2, -4, +2, +5, -3, -2, +0, -3, +1, -3, -2, +0, -1, +2, +0, +0, -1, -2, -3, -1, -2, +5, // E\r
+ \r
+ -2, -3, -2, -3, -3, +6, -3, -1, +0, -3, +0, +0, -3, -4, -3, -3, -2, -2, -1, +1, -1, +3, -3, // F\r
+ +0, -1, -3, -1, -2, -3, +6, -2, -4, -2, -4, -3, +0, -2, -2, -2, +0, -2, -3, -2, -1, -3, -2, // G\r
+ -2, -1, -3, -1, +0, -1, -2, +8, -3, -1, -3, -2, +1, -2, +0, +0, -1, -2, -3, -2, -1, +2, +0, // H\r
+ -1, -3, -1, -3, -3, +0, -4, -3, +4, -3, +2, +1, -3, -3, -3, -3, -2, -1, +3, -3, -1, -1, -3, // I\r
+ -1, -1, -3, -1, +1, -3, -2, -1, -3, +5, -2, -1, +0, -1, +1, +2, +0, -1, -2, -3, -1, -2, +1, // K\r
+ \r
+ -1, -4, -1, -4, -3, +0, -4, -3, +2, -2, +4, +2, -3, -3, -2, -2, -2, -1, +1, -2, -1, -1, -3, // L\r
+ -1, -3, -1, -3, -2, +0, -3, -2, +1, -1, +2, +5, -2, -2, +0, -1, -1, -1, +1, -1, -1, -1, -2, // M\r
+ -2, +1, -3, +1, +0, -3, +0, +1, -3, +0, -3, -2, +6, -2, +0, +0, +1, +0, -3, -4, -1, -2, +0, // N\r
+ -1, -1, -3, -1, -1, -4, -2, -2, -3, -1, -3, -2, -2, +7, -1, -2, -1, -1, -2, -4, -1, -3, -1, // P\r
+ -1, +0, -3, +0, +2, -3, -2, +0, -3, +1, -2, +0, +0, -1, +5, +1, +0, -1, -2, -2, -1, -1, +2, // Q\r
+ \r
+ -1, -2, -3, -2, +0, -3, -2, +0, -3, +2, -2, -1, +0, -2, +1, +5, -1, -1, -3, -3, -1, -2, +0, // R\r
+ +1, +0, -1, +0, +0, -2, +0, -1, -2, +0, -2, -1, +1, -1, +0, -1, +4, +1, -2, -3, -1, -2, +0, // S\r
+ +0, -1, -1, -1, -1, -2, -2, -2, -1, -1, -1, -1, +0, -1, -1, -1, +1, +5, +0, -2, -1, -2, -1, // T\r
+ +0, -3, -1, -3, -2, -1, -3, -3, +3, -2, +1, +1, -3, -2, -2, -3, -2, +0, +4, -3, -1, -1, -2, // V\r
+ -3, -4, -2, -4, -3, +1, -2, -2, -3, -3, -2, -1, -4, -4, -2, -3, -3, -2, -3,+11, -1, +2, -3, // W\r
+ \r
+ -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, -1, // X\r
+ -2, -3, -2, -3, -2, +3, -3, +2, -1, -2, -1, -1, -2, -3, -1, -2, -2, -2, -1, +2, -1, +7, -2, // Y\r
+ -1, +2, -4, +2, +5, -3, -2, +0, -3, +1, -3, -2, +0, -1, +2, +0, +0, -1, -2, -3, -1, -2, +5, // Z\r
+ };\r
+\r
+static int toi_tab[26] =\r
+ {\r
+ 0, // A\r
+ 1, // B\r
+ 2, // C\r
+ 3, // D\r
+ 4, // E\r
+ 5, // F\r
+ 6, // G\r
+ 7, // H\r
+ 8, // I\r
+ -1, // J\r
+ 9, // K\r
+ 10, // L\r
+ 11, // M\r
+ 12, // N\r
+ -1, // O\r
+ 13, // P\r
+ 14, // Q\r
+ 15, // R\r
+ 16, // S\r
+ 17, // T\r
+ 17, // U\r
+ 18, // V\r
+ 19, // W\r
+ 20, // X\r
+ 21, // Y\r
+ 22, // Z\r
+ };\r
+\r
+static int toi(char c)\r
+ {\r
+ c = toupper(c);\r
+ return toi_tab[c - 'A'];\r
+ }\r
+\r
+static int BlosumScore(char c1, char c2)\r
+ {\r
+ int i1 = toi(c1);\r
+ int i2 = toi(c2);\r
+ return Blosum62[i1][i2];\r
+ }\r
+\r
+/***\r
+Consider a column with 5 As and 3 Bs.\r
+There are:\r
+ 5x4 pairs of As.\r
+ 3x2 pairs of Bs.\r
+ 5x3x2 AB pairs\r
+ 8x7 = 5x4 + 3x2 + 5x3x2 pairs of letters\r
+***/\r
+static double BlosumScoreCol(const MSA &a, unsigned uColIndex)\r
+ {\r
+ int iCounts[23];\r
+ memset(iCounts, 0, sizeof(iCounts));\r
+ const unsigned uSeqCount = a.GetSeqCount();\r
+ unsigned uCharCount = 0;\r
+ for (unsigned uSeqIndex = 0; uSeqIndex < uSeqCount; ++uSeqIndex)\r
+ {\r
+ char c = a.GetChar(uSeqIndex, uColIndex);\r
+ if (IsGapChar(c))\r
+ continue;\r
+ int iChar = toi(c);\r
+ ++iCounts[iChar];\r
+ ++uCharCount;\r
+ }\r
+ if (uCharCount < 2)\r
+ return -9;\r
+ int iTotalScore = 0;\r
+ for (int i1 = 0; i1 < 23; ++i1)\r
+ {\r
+ int iCounts1 = iCounts[i1];\r
+ iTotalScore += iCounts1*(iCounts1 - 1)*Blosum62[i1][i1];\r
+ for (int i2 = i1 + 1; i2 < 23; ++i2)\r
+ iTotalScore += iCounts[i2]*iCounts1*2*Blosum62[i1][i2];\r
+ }\r
+ int iPairCount = uCharCount*(uCharCount - 1);\r
+ return (double) iTotalScore / (double) iPairCount;\r
+ }\r
+\r
+/***\r
+Consider a column with 5 As and 3 Bs.\r
+A residue of type Q scores:\r
+ 5xAQ + 3xBQ\r
+***/\r
+static void AssignColorsCol(const MSA &a, unsigned uColIndex, int **Colors)\r
+ {\r
+ int iCounts[23];\r
+ memset(iCounts, 0, sizeof(iCounts));\r
+ const unsigned uSeqCount = a.GetSeqCount();\r
+ unsigned uCharCount = 0;\r
+ for (unsigned uSeqIndex = 0; uSeqIndex < uSeqCount; ++uSeqIndex)\r
+ {\r
+ char c = a.GetChar(uSeqIndex, uColIndex);\r
+ if (IsGapChar(c))\r
+ continue;\r
+ int iChar = toi(c);\r
+ ++iCounts[iChar];\r
+ ++uCharCount;\r
+ }\r
+ int iMostConservedType = -1;\r
+ int iMostConservedCount = -1;\r
+ for (unsigned i = 0; i < 23; ++i)\r
+ {\r
+ if (iCounts[i] > iMostConservedCount)\r
+ {\r
+ iMostConservedType = i;\r
+ iMostConservedCount = iCounts[i];\r
+ }\r
+ }\r
+\r
+ double dColScore = BlosumScoreCol(a, uColIndex);\r
+ int c;\r
+ if (dColScore >= 3.0)\r
+ c = 3;\r
+ //else if (dColScore >= 1.0)\r
+ // c = 2;\r
+ else if (dColScore >= 0.2)\r
+ c = 1;\r
+ else\r
+ c = 0;\r
+\r
+ int Color[23];\r
+ for (unsigned uLetter = 0; uLetter < 23; ++uLetter)\r
+ {\r
+ double dScore = Blosum62[uLetter][iMostConservedType];\r
+ if (dScore >= dColScore)\r
+ Color[uLetter] = c;\r
+ else\r
+ Color[uLetter] = 0;\r
+ }\r
+\r
+ for (unsigned uSeqIndex = 0; uSeqIndex < uSeqCount; ++uSeqIndex)\r
+ {\r
+ char c = a.GetChar(uSeqIndex, uColIndex);\r
+ if (IsGapChar(c))\r
+ {\r
+ Colors[uSeqIndex][uColIndex] = 0;\r
+ continue;\r
+ }\r
+ int iLetter = toi(c);\r
+ if (iLetter >= 0 && iLetter < 23)\r
+ Colors[uSeqIndex][uColIndex] = Color[iLetter];\r
+ else\r
+ Colors[uSeqIndex][uColIndex] = 0;\r
+ }\r
+ }\r
+\r
+void AssignColors(const MSA &a, int **Colors)\r
+ {\r
+ const unsigned uColCount = a.GetColCount();\r
+ for (unsigned uColIndex = 0; uColIndex < uColCount; ++uColIndex)\r
+ AssignColorsCol(a, uColIndex, Colors);\r
+ }\r
git://source.jalview.org / jabaws.git / blobdiff