binaries/src/mafft/core/blosum.c

   1 #define DEFAULTGOP_B -1530
   2 #define DEFAULTGEP_B   -00
   3 #define DEFAULTOFS_B  -123   /* +10 -- -50  teido ka ? */
   4
   5 void BLOSUMmtx( int n, double **matrix, double *freq, char *amino, char *amino_grp )
   6 {
   7         /*
   8         char locaminod[26] = "GASTPLIMVDNEQFYWKRHCXXX.-U";
   9         */
  10 //      char locaminod[] = "ARNDCQEGHILKMFPSTWYVBZX.-U";
  11         char locaminod[] = "ARNDCQEGHILKMFPSTWYVBZX.-J";
  12         char locgrpd[] =
  13         {
  14                 0, 3, 2, 2, 5, 2, 2, 0, 3, 1, 1, 3, 1, 4, 0, 0, 0, 4, 4, 1, 2, 2,
  15                 6, 6, 6, 1,
  16         };
  17         double freqd[20] =
  18         {
  19             0.077,
  20             0.051,
  21             0.043,
  22             0.052,
  23             0.020,
  24             0.041,
  25             0.062,
  26             0.074,
  27             0.023,
  28             0.052,
  29             0.091,
  30             0.059,
  31             0.024,
  32             0.040,
  33             0.051,
  34             0.069,
  35             0.059,
  36             0.014,
  37             0.032,
  38             0.066,
  39         };
  40
  41         double tmpmtx30[] =
  42         {
  43     4,
  44    -1,     8,
  45     0,    -2,     8,
  46     0,    -1,     1,     9,
  47    -3,    -2,    -1,    -3,    17,
  48     1,     3,    -1,    -1,    -2,     8,
  49     0,    -1,    -1,     1,     1,     2,     6,
  50     0,    -2,     0,    -1,    -4,    -2,    -2,     8,
  51    -2,    -1,    -1,    -2,    -5,     0,     0,    -3,    14,
  52     0,    -3,     0,    -4,    -2,    -2,    -3,    -1,    -2,     6,
  53    -1,    -2,    -2,    -1,     0,    -2,    -1,    -2,    -1,     2,     4,
  54     0,     1,     0,     0,    -3,     0,     2,    -1,    -2,    -2,    -2,     4,
  55     1,     0,     0,    -3,    -2,    -1,    -1,    -2,     2,     1,     2,     2,     6,
  56    -2,    -1,    -1,    -5,    -3,    -3,    -4,    -3,    -3,     0,     2,    -1,    -2,    10,
  57    -1,    -1,    -3,    -1,    -3,     0,     1,    -1,     1,    -3,    -3,     1,    -4,    -4,    11,
  58     1,    -1,     0,     0,    -2,    -1,     0,     0,    -1,    -1,    -2,     0,    -2,    -1,    -1,     4,
  59     1,    -3,     1,    -1,    -2,     0,    -2,    -2,    -2,     0,     0,    -1,     0,    -2,     0,     2,     5,
  60    -5,     0,    -7,    -4,    -2,    -1,    -1,     1,    -5,    -3,    -2,    -2,    -3,     1,    -3,    -3,    -5,    20,
  61    -4,     0,    -4,    -1,    -6,    -1,    -2,    -3,     0,    -1,     3,    -1,    -1,     3,    -2,    -2,    -1,     5,     9,
  62     1,    -1,    -2,    -2,    -2,    -3,    -3,    -3,    -3,     4,     1,    -2,     0,     1,    -4,    -1,     1,    -3,     1,     5,
  63     0,    -2,     4,     5,    -2,    -1,     0,     0,    -2,    -2,    -1,     0,    -2,    -3,    -2,     0,     0,    -5,    -3,    -2,     5,
  64     0,     0,    -1,     0,     0,     4,     5,    -2,     0,    -3,    -1,     1,    -1,    -4,     0,    -1,    -1,    -1,    -2,    -3,     0,     4,
  65     0,    -1,     0,    -1,    -2,     0,    -1,    -1,    -1,     0,     0,     0,     0,    -1,    -1,     0,     0,    -2,    -1,     0,    -1,     0,    -1,
  66         };
  67
  68         double tmpmtx45[] =
  69         {
  70       5,
  71      -2,      7,
  72      -1,      0,      6,
  73      -2,     -1,      2,      7,
  74      -1,     -3,     -2,     -3,     12,
  75      -1,      1,      0,      0,     -3,      6,
  76      -1,      0,      0,      2,     -3,      2,      6,
  77       0,     -2,      0,     -1,     -3,     -2,     -2,      7,
  78      -2,      0,      1,      0,     -3,      1,      0,     -2,     10,
  79      -1,     -3,     -2,     -4,     -3,     -2,     -3,     -4,     -3,      5,
  80      -1,     -2,     -3,     -3,     -2,     -2,     -2,     -3,     -2,      2,      5,
  81      -1,      3,      0,      0,     -3,      1,      1,     -2,     -1,     -3,     -3,      5,
  82      -1,     -1,     -2,     -3,     -2,      0,     -2,     -2,      0,      2,      2,     -1,      6,
  83      -2,     -2,     -2,     -4,     -2,     -4,     -3,     -3,     -2,      0,      1,     -3,      0,      8,
  84      -1,     -2,     -2,     -1,     -4,     -1,      0,     -2,     -2,     -2,     -3,     -1,     -2,     -3,      9,
  85       1,     -1,      1,      0,     -1,      0,      0,      0,     -1,     -2,     -3,     -1,     -2,     -2,     -1,      4,
  86       0,     -1,      0,     -1,     -1,     -1,     -1,     -2,     -2,     -1,     -1,     -1,     -1,     -1,     -1,      2,      5,
  87      -2,     -2,     -4,     -4,     -5,     -2,     -3,     -2,     -3,     -2,     -2,     -2,     -2,      1,     -3,     -4,     -3,     15,
  88      -2,     -1,     -2,     -2,     -3,     -1,     -2,     -3,      2,      0,      0,     -1,      0,      3,     -3,     -2,     -1,      3,      8,
  89       0,     -2,     -3,     -3,     -1,     -3,     -3,     -3,     -3,      3,      1,     -2,      1,      0,     -3,     -1,      0,     -3,     -1,      5,
  90         };
  91     double tmpmtx50[] =
  92     {
  93        5,
  94       -2,      7,
  95       -1,     -1,      7,
  96       -2,     -2,      2,      8,
  97       -1,     -4,     -2,     -4,     13,
  98       -1,      1,      0,      0,     -3,      7,
  99       -1,      0,      0,      2,     -3,      2,      6,
 100        0,     -3,      0,     -1,     -3,     -2,     -3,      8,
 101       -2,      0,      1,     -1,     -3,      1,      0,     -2,     10,
 102       -1,     -4,     -3,     -4,     -2,     -3,     -4,     -4,     -4,      5,
 103       -2,     -3,     -4,     -4,     -2,     -2,     -3,     -4,     -3,      2,      5,
 104       -1,      3,      0,     -1,     -3,      2,      1,     -2,      0,     -3,     -3,      6,
 105       -1,     -2,     -2,     -4,     -2,      0,     -2,     -3,     -1,      2,      3,     -2,      7,
 106       -3,     -3,     -4,     -5,     -2,     -4,     -3,     -4,     -1,      0,      1,     -4,      0,      8,
 107       -1,     -3,     -2,     -1,     -4,     -1,     -1,     -2,     -2,     -3,     -4,     -1,     -3,     -4,     10,
 108        1,     -1,      1,      0,     -1,      0,     -1,      0,     -1,     -3,     -3,      0,     -2,     -3,     -1,      5,
 109        0,     -1,      0,     -1,     -1,     -1,     -1,     -2,     -2,     -1,     -1,     -1,     -1,     -2,     -1,      2,      5,
 110       -3,     -3,     -4,     -5,     -5,     -1,     -3,     -3,     -3,     -3,     -2,     -3,     -1,      1,     -4,     -4,     -3,     15,
 111       -2,     -1,     -2,     -3,     -3,     -1,     -2,     -3,      2,     -1,     -1,     -2,      0,      4,     -3,     -2,     -2,      2,      8,
 112        0,     -3,     -3,     -4,     -1,     -3,     -3,     -4,     -4,      4,      1,     -3,      1,     -1,     -3,     -2,      0,     -3,     -1,      5,
 113     };
 114         double tmpmtx62[] =
 115         {
 116       6,
 117      -2,      8,
 118      -2,     -1,      8,
 119      -3,     -2,      2,      9,
 120      -1,     -5,     -4,     -5,     13,
 121      -1,      1,      0,      0,     -4,      8,
 122      -1,      0,      0,      2,     -5,      3,      7,
 123       0,     -3,     -1,     -2,     -4,     -3,     -3,      8,
 124      -2,      0,      1,     -2,     -4,      1,      0,     -3,     11,
 125      -2,     -4,     -5,     -5,     -2,     -4,     -5,     -6,     -5,      6,
 126      -2,     -3,     -5,     -5,     -2,     -3,     -4,     -5,     -4,      2,      6,
 127      -1,      3,      0,     -1,     -5,      2,      1,     -2,     -1,     -4,     -4,      7,
 128      -1,     -2,     -3,     -5,     -2,     -1,     -3,     -4,     -2,      2,      3,     -2,      8,
 129      -3,     -4,     -4,     -5,     -4,     -5,     -5,     -5,     -2,      0,      1,     -5,      0,      9,
 130      -1,     -3,     -3,     -2,     -4,     -2,     -2,     -3,     -3,     -4,     -4,     -2,     -4,     -5,     11,
 131       2,     -1,      1,      0,     -1,      0,      0,      0,     -1,     -4,     -4,      0,     -2,     -4,     -1,      6,
 132       0,     -2,      0,     -2,     -1,     -1,     -1,     -2,     -3,     -1,     -2,     -1,     -1,     -3,     -2,      2,      7,
 133      -4,     -4,     -6,     -6,     -3,     -3,     -4,     -4,     -4,     -4,     -2,     -4,     -2,      1,     -5,     -4,     -4,     16,
 134      -3,     -3,     -3,     -5,     -4,     -2,     -3,     -5,      3,     -2,     -2,     -3,     -1,      4,     -4,     -3,     -2,      3,     10,
 135       0,     -4,     -4,     -5,     -1,     -3,     -4,     -5,     -5,      4,      1,     -3,      1,     -1,     -4,     -2,      0,     -4,     -2,      6,
 136         };
 137         double tmpmtx80[] =
 138         {
 139       7,
 140      -3,      9,
 141      -3,     -1,      9,
 142      -3,     -3,      2,     10,
 143      -1,     -6,     -5,     -7,     13,
 144      -2,      1,      0,     -1,     -5,      9,
 145      -2,     -1,     -1,      2,     -7,      3,      8,
 146       0,     -4,     -1,     -3,     -6,     -4,     -4,      9,
 147      -3,      0,      1,     -2,     -7,      1,      0,     -4,     12,
 148      -3,     -5,     -6,     -7,     -2,     -5,     -6,     -7,     -6,      7,
 149      -3,     -4,     -6,     -7,     -3,     -4,     -6,     -7,     -5,      2,      6,
 150      -1,      3,      0,     -2,     -6,      2,      1,     -3,     -1,     -5,     -4,      8,
 151      -2,     -3,     -4,     -6,     -3,     -1,     -4,     -5,     -4,      2,      3,     -3,      9,
 152      -4,     -5,     -6,     -6,     -4,     -5,     -6,     -6,     -2,     -1,      0,     -5,      0,     10,
 153      -1,     -3,     -4,     -3,     -6,     -3,     -2,     -5,     -4,     -5,     -5,     -2,     -4,     -6,     12,
 154       2,     -2,      1,     -1,     -2,     -1,     -1,     -1,     -2,     -4,     -4,     -1,     -3,     -4,     -2,      7,
 155       0,     -2,      0,     -2,     -2,     -1,     -2,     -3,     -3,     -2,     -3,     -1,     -1,     -4,     -3,      2,      8,
 156      -5,     -5,     -7,     -8,     -5,     -4,     -6,     -6,     -4,     -5,     -4,     -6,     -3,      0,     -7,     -6,     -5,     16,
 157      -4,     -4,     -4,     -6,     -5,     -3,     -5,     -6,      3,     -3,     -2,     -4,     -3,      4,     -6,     -3,     -3,      3,     11,
 158      -1,     -4,     -5,     -6,     -2,     -4,     -4,     -6,     -5,      4,      1,     -4,      1,     -2,     -4,     -3,      0,     -5,     -3,      7,
 159         };
 160         double tmpmtx0[] =
 161         {
 162      2.4,
 163     -0.6,    4.7,
 164     -0.3,    0.3,    3.8,
 165     -0.3,   -0.3,    2.2,    4.7,
 166      0.5,   -2.2,   -1.8,   -3.2,   11.5,
 167     -0.2,    1.5,    0.7,    0.9,   -2.4,    2.7,
 168      0.0,    0.4,    0.9,    2.7,   -3.0,    1.7,    3.6,
 169      0.5,   -1.0,    0.4,    0.1,   -2.0,   -1.0,   -0.8,    6.6,
 170     -0.8,    0.6,    1.2,    0.4,   -1.3,    1.2,    0.4,   -1.4,    6.0,
 171     -0.8,   -2.4,   -2.8,   -3.8,   -1.1,   -1.9,   -2.7,   -4.5,   -2.2,    4.0,
 172     -1.2,   -2.2,   -3.0,   -4.0,   -1.5,   -1.6,   -2.8,   -4.4,   -1.9,    2.8,    4.0,
 173     -0.4,    2.7,    0.8,    0.5,   -2.8,    1.5,    1.2,   -1.1,    0.6,   -2.1,   -2.1,    3.2,
 174     -0.7,   -1.7,   -2.2,   -3.0,   -0.9,   -1.0,   -2.0,   -3.5,   -1.3,    2.5,    2.8,   -1.4,    4.3,
 175     -2.3,   -3.2,   -3.1,   -4.5,   -0.8,   -2.6,   -3.9,   -5.2,   -0.1,    1.0,    2.0,   -3.3,    1.6,    7.0,
 176      0.3,   -0.9,   -0.9,   -0.7,   -3.1,   -0.2,   -0.5,   -1.6,   -1.1,   -2.6,   -2.3,   -0.6,   -2.4,   -3.8,    7.6,
 177      1.1,   -0.2,    0.9,    0.5,    0.1,    0.2,    0.2,    0.4,   -0.2,   -1.8,   -2.1,    0.1,   -1.4,   -2.8,    0.4,    2.2,
 178      0.6,   -0.2,    0.5,    0.0,   -0.5,    0.0,   -0.1,   -1.1,   -0.3,   -0.6,   -1.3,    0.1,   -0.6,   -2.2,    0.1,    1.5,    2.5,
 179     -3.6,   -1.6,   -3.6,   -5.2,   -1.0,   -2.7,   -4.3,   -4.0,   -0.8,   -1.8,   -0.7,   -3.5,   -1.0,    3.6,   -5.0,   -3.3,   -3.5,   14.2,
 180     -2.2,   -1.8,   -1.4,   -2.8,   -0.5,   -1.7,   -2.7,   -4.0,    2.2,   -0.7,    0.0,   -2.1,   -0.2,    5.1,   -3.1,   -1.9,   -1.9,    4.1,    7.8,
 181      0.1,   -2.0,   -2.2,   -2.9,    0.0,   -1.5,   -1.9,   -3.3,   -2.0,    3.1,    1.8,   -1.7,    1.6,    0.1,   -1.8,   -1.0,    0.0,   -2.6,   -1.1,    3.4,
 182         };
 183
 184         int i, j, count;
 185         double av;
 186         double *tmpmtx;
 187
 188         if( n == 30 ) tmpmtx = tmpmtx30;
 189         else if( n == 45 ) tmpmtx = tmpmtx45;
 190         else if( n == 50 ) tmpmtx = tmpmtx50;
 191         else if( n == 62 ) tmpmtx = tmpmtx62;
 192         else if( n == 80 ) tmpmtx = tmpmtx80;
 193         else if( n == 0 ) tmpmtx = tmpmtx0;
 194         else if( n == -1 ) tmpmtx = loadaamtx();
 195         else
 196         {
 197                 fprintf( stderr, "blosum %d ?\n", n );
 198                 exit( 1 );
 199         }
 200
 201         count = 0;
 202         for( i=0; i<20; i++ )
 203         {
 204                 for( j=0; j<=i; j++ )
 205                 {
 206                         matrix[i][j] = matrix[j][i] = (double)tmpmtx[count++];
 207                 }
 208         }
 209         if( n == -1 && tmpmtx[400] != -1.0 )
 210         {
 211                 for( i=0; i<20; i++ ) freq[i] = tmpmtx[400+i];
 212                 av = 0.0;
 213                 for( i=0; i<20; i++ ) av += freq[i];
 214                 for( i=0; i<20; i++ ) freq[i] /= av;
 215         }
 216         else
 217                 for( i=0; i<20; i++ ) freq[i] = freqd[i];
 218
 219 #if 0
 220         av = 0.0;
 221         for( i=0; i<20; i++ )
 222                 av += matrix[i][i];
 223         av /= 20;
 224         fprintf( stdout, "av = %f\n", av );
 225
 226         for( i=0; i<20; i++ ) for( j=0; j<20; j++ )
 227                 matrix[i][j] /= av;
 228
 229         av = wav = 0;
 230         count = 0;
 231         wcount = 0.0;
 232         tmptmp = 0.0;
 233         for( i=0; i<20; i++ )
 234         {
 235                 fprintf( stdout, "freq[%d] = %f\n", i, freq[i] );
 236                 tmptmp += freq[i];
 237                 for( j=0; j<20; j++ )
 238                 {
 239                         av += matrix[i][j];
 240                         wav += freq[i] * freq[j] * matrix[i][j];
 241                         count++;
 242                         wcount += freq[i] * freq[j];
 243                 }
 244         }
 245
 246         av /= count;
 247         wav /= wcount;
 248         fprintf( stdout, "av = %f\n", av );
 249         fprintf( stdout, "wav = %f\n", wav );
 250         fprintf( stdout, "wcount = %f\n", wcount );
 251         fprintf( stdout, "tmptmp = %f\n", tmptmp );
 252
 253         for( i=0; i<20; i++ )
 254         {
 255                 for( j=0; j<=i; j++ )
 256                 {
 257                         fprintf( stderr, "## %d-%d, %f\n", i, j, matrix[i][j] );
 258                 }
 259         }
 260
 261         exit( 1 );
 262 #endif
 263
 264     for( i=0; i<26; i++ ) amino[i] = locaminod[i];
 265     for( i=0; i<26; i++ ) amino_grp[(int)amino[i]] = locgrpd[i];
 266 }