1/********* HIVm MODELS (within/between) OF PROTEIN EVOLUTION ************/ 2LoadFunctionLibrary("../empirical.bf"); 3LoadFunctionLibrary("../../protein.bf"); 4LoadFunctionLibrary("../../parameters.bf"); 5LoadFunctionLibrary("../../frequencies.bf"); 6LoadFunctionLibrary("../../../UtilityFunctions.bf"); 7LoadFunctionLibrary("../../../all-terms.bf"); 8 9 10 11/**************************************** HIVBm functions *************************************/ 12 13 14/** 15 * @name models.protein.HIVBm.ModelDescription 16 * @description Create the baseline schema (dictionary) for the HIVBm model of protein evolution 17 * @returns {Dictionary} model description 18 * @param {String} type 19 */ 20function models.protein.HIVBm.ModelDescription(type) { 21 models.protein.HIVBm.ModelDescription.model_definition = models.protein.empirical.ModelDescription(type); 22 models.protein.HIVBm.ModelDescription.model_definition [terms.model.empirical_rates] = models.protein.HIVBm.Rij; 23 models.protein.HIVBm.ModelDescription.model_definition [terms.model.frequency_estimator] = "models.protein.HIVBm.frequencies"; 24 return models.protein.HIVBm.ModelDescription.model_definition; 25} 26 27/** 28 * @name models.protein.HIVBmF.ModelDescription 29 * @description Create the baseline schema (dictionary) for the HIVBm+F model of protein evolution 30 * @returns {Dictionary} model description 31 * @param {String} type 32 */ 33function models.protein.HIVBmF.ModelDescription(type) { 34 models.protein.HIVBmF.ModelDescription.model_definition = models.protein.HIVBm.ModelDescription(type); 35 models.protein.HIVBmF.ModelDescription.model_definition [terms.model.frequency_estimator] = "frequencies.empirical.protein"; 36 models.protein.HIVBmF.ModelDescription.model_definition [terms.model.efv_estimate_name] = utility.getGlobalValue("terms.frequencies._20x1"); 37 return models.protein.HIVBmF.ModelDescription.model_definition; 38} 39 40/** 41 * @name models.protein.HIVBmML.ModelDescription 42 * @description Create the baseline schema (dictionary) for the HIVBm+ML model of protein evolution 43 * @returns {Dictionary} model description 44 * @param {String} type 45 */ 46function models.protein.HIVBmML.ModelDescription(type) { 47 models.protein.HIVBmML.ModelDescription.model_definition = models.protein.HIVBm.ModelDescription(type); 48 models.protein.HIVBmML.ModelDescription.model_definition [terms.model.frequency_estimator] = "frequencies.ML.protein"; 49 models.protein.HIVBmML.ModelDescription.model_definition [terms.model.efv_estimate_name] = utility.getGlobalValue("terms.frequencies.MLE"); 50 return models.protein.HIVBmML.ModelDescription.model_definition; 51} 52 53 54 55 56models.protein.HIVBm.Rij = { 57 "A":{ 58 "C":0.123758, 59 "D":1.45504, 60 "E":1.48135, 61 "F":0.0141269, 62 "G":2.13536, 63 "H":0.0847613, 64 "I":0.005, 65 "K":0.005, 66 "L":0.215256, 67 "M":0.0186643, 68 "N":0.005, 69 "P":2.12217, 70 "Q":0.0551128, 71 "R":0.307507, 72 "S":2.46633, 73 "T":15.9183, 74 "V":7.61428, 75 "W":0.005, 76 "Y":0.005 77 }, 78 "C":{ 79 "D":0.005, 80 "E":0.005, 81 "F":9.29815, 82 "G":0.897871, 83 "H":0.240073, 84 "I":0.005, 85 "K":0.005, 86 "L":0.129777, 87 "M":0.005, 88 "N":0.08606419999999999, 89 "P":0.005, 90 "Q":0.005, 91 "R":0.351721, 92 "S":4.69314, 93 "T":0.739969, 94 "V":0.420027, 95 "W":2.63277, 96 "Y":7.57932 97 }, 98 "D":{ 99 "E":10.5872, 100 "F":0.005, 101 "G":2.83806, 102 "H":1.9169, 103 "I":0.0176792, 104 "K":0.005, 105 "L":0.008760479999999999, 106 "M":0.005, 107 "N":17.6612, 108 "P":0.0342658, 109 "Q":0.005, 110 "R":0.005, 111 "S":0.52823, 112 "T":0.274724, 113 "V":1.04793, 114 "W":0.005, 115 "Y":0.6746529999999999 116 }, 117 "E":{ 118 "F":0.005, 119 "G":3.92775, 120 "H":0.11974, 121 "I":0.00609079, 122 "K":4.61482, 123 "L":0.005, 124 "M":0.175789, 125 "N":0.07926329999999999, 126 "P":0.0120226, 127 "Q":2.5602, 128 "R":0.0749218, 129 "S":0.005, 130 "T":0.289774, 131 "V":1.02847, 132 "W":0.005, 133 "Y":0.07926329999999999 134 }, 135 "F":{ 136 "G":0.291561, 137 "H":0.145558, 138 "I":3.39836, 139 "K":0.0342658, 140 "L":8.524839999999999, 141 "M":0.188025, 142 "N":0.005, 143 "P":0.005, 144 "Q":0.005, 145 "R":0.005, 146 "S":0.956472, 147 "T":0.0141269, 148 "V":0.723274, 149 "W":0.8293430000000001, 150 "Y":15.34 151 }, 152 "G":{ 153 "H":0.005, 154 "I":0.005, 155 "K":0.521705, 156 "L":0.005, 157 "M":0.005, 158 "N":0.323401, 159 "P":0.005, 160 "Q":0.0619137, 161 "R":3.65345, 162 "S":4.38041, 163 "T":0.369615, 164 "V":0.953155, 165 "W":1.21674, 166 "Y":0.005 167 }, 168 "H":{ 169 "I":0.103111, 170 "K":0.005, 171 "L":1.74171, 172 "M":0.005, 173 "N":7.64585, 174 "P":2.45318, 175 "Q":7.05545, 176 "R":9.04044, 177 "S":0.382747, 178 "T":0.7115939999999999, 179 "V":0.005, 180 "W":0.06951789999999999, 181 "Y":18.6943 182 }, 183 "I":{ 184 "K":0.322319, 185 "L":5.95879, 186 "M":11.2065, 187 "N":0.680565, 188 "P":0.0410593, 189 "Q":0.005, 190 "R":0.677289, 191 "S":1.21803, 192 "T":8.612170000000001, 193 "V":17.7389, 194 "W":0.005, 195 "Y":0.148168 196 }, 197 "K":{ 198 "L":0.0814995, 199 "M":1.28246, 200 "N":7.90443, 201 "P":0.0313862, 202 "Q":6.54737, 203 "R":20.45, 204 "S":0.504111, 205 "T":4.67142, 206 "V":0.265829, 207 "W":0.005, 208 "Y":0.005 209 }, 210 "L":{ 211 "M":5.31961, 212 "N":0.005, 213 "P":2.07757, 214 "Q":1.49456, 215 "R":0.701427, 216 "S":0.927656, 217 "T":0.0437673, 218 "V":1.41036, 219 "W":0.748843, 220 "Y":0.111986 221 }, 222 "M":{ 223 "N":0.005, 224 "P":0.005, 225 "Q":0.303676, 226 "R":2.51394, 227 "S":0.005, 228 "T":4.94026, 229 "V":6.8532, 230 "W":0.089078, 231 "Y":0.005 232 }, 233 "N":{ 234 "P":0.00739578, 235 "Q":0.672052, 236 "R":0.295543, 237 "S":13.1447, 238 "T":6.88667, 239 "V":0.026656, 240 "W":0.005, 241 "Y":1.76417 242 }, 243 "P":{ 244 "Q":4.47211, 245 "R":1.28355, 246 "S":5.37762, 247 "T":2.01417, 248 "V":0.005, 249 "W":0.0444506, 250 "Y":0.0304381 251 }, 252 "Q":{ 253 "R":3.4215, 254 "S":0.116311, 255 "T":0.243589, 256 "V":0.0209153, 257 "W":0.026656, 258 "Y":0.113033 259 }, 260 "R":{ 261 "S":3.4791, 262 "T":2.86868, 263 "V":0.0812454, 264 "W":0.9913380000000001, 265 "Y":0.00991826 266 }, 267 "S":{ 268 "T":8.93107, 269 "V":0.0749218, 270 "W":0.0248728, 271 "Y":0.648024 272 }, 273 "T":{ 274 "V":0.709226, 275 "W":0.005, 276 "Y":0.105652 277 }, 278 "V":{ 279 "W":0.005, 280 "Y":0.0410593 281 }, 282 "W":{ 283 "Y":1.28022 284 } 285}; 286 287 288lfunction models.protein.HIVBm.frequencies (model, namespace, datafilter) { 289 model[utility.getGlobalValue("terms.efv_estimate")] = 290 {{ 0.060490222} 291 { 0.020075899} 292 { 0.042109048} 293 { 0.071567447} 294 { 0.028809447} 295 { 0.072308239} 296 { 0.022293943} 297 { 0.069730629} 298 { 0.056968211} 299 { 0.098851122} 300 { 0.019768318} 301 { 0.044127815} 302 { 0.046025282} 303 { 0.053606488} 304 { 0.066039665} 305 { 0.05060433} 306 { 0.053636813} 307 { 0.061625237} 308 { 0.033011601} 309 { 0.028350243}}; 310 model[utility.getGlobalValue("terms.model.efv_estimate_name")] = utility.getGlobalValue("terms.frequencies.predefined"); 311 (model[utility.getGlobalValue("terms.parameters")])[utility.getGlobalValue("terms.model.empirical")] = 0; 312 return model; 313} 314 315 316 317 318 319 320 321 322 323 324 325 326/**************************************** HIVWm functions *************************************/ 327 328 329/** 330 * @name models.protein.HIVWm.ModelDescription 331 * @description Create the baseline schema (dictionary) for the HIVWm model of protein evolution 332 * @returns {Dictionary} model description 333 * @param {String} type 334 */ 335function models.protein.HIVWm.ModelDescription(type) { 336 models.protein.HIVWm.ModelDescription.model_definition = models.protein.empirical.ModelDescription(type); 337 models.protein.HIVWm.ModelDescription.model_definition [terms.model.empirical_rates] = models.protein.HIVWm.Rij; 338 models.protein.HIVWm.ModelDescription.model_definition [terms.model.frequency_estimator] = "models.protein.HIVWm.frequencies"; 339 return models.protein.HIVWm.ModelDescription.model_definition; 340} 341 342/** 343 * @name models.protein.HIVWmF.ModelDescription 344 * @description Create the baseline schema (dictionary) for the HIVWm+F model of protein evolution 345 * @returns {Dictionary} model description 346 * @param {String} type 347 */ 348function models.protein.HIVWmF.ModelDescription(type) { 349 models.protein.HIVWmF.ModelDescription.model_definition = models.protein.HIVWm.ModelDescription(type); 350 models.protein.HIVWmF.ModelDescription.model_definition [terms.model.frequency_estimator] = "frequencies.empirical.protein"; 351 models.protein.HIVWmF.ModelDescription.model_definition [terms.model.efv_estimate_name] = utility.getGlobalValue("terms.frequencies._20x1"); 352 return models.protein.HIVWmF.ModelDescription.model_definition; 353} 354 355/** 356 * @name models.protein.HIVWmML.ModelDescription 357 * @description Create the baseline schema (dictionary) for the HIVWm+ML model of protein evolution 358 * @returns {Dictionary} model description 359 * @param {String} type 360 */ 361function models.protein.HIVWmML.ModelDescription(type) { 362 models.protein.HIVWmML.ModelDescription.model_definition = models.protein.HIVWm.ModelDescription(type); 363 models.protein.HIVWmML.ModelDescription.model_definition [terms.model.frequency_estimator] = "frequencies.ML.protein"; 364 models.protein.HIVWmML.ModelDescription.model_definition [terms.model.efv_estimate_name] = utility.getGlobalValue("terms.frequencies.MLE"); 365 return models.protein.HIVWmML.ModelDescription.model_definition; 366} 367 368 369 370 371 372 373models.protein.HIVWm.Rij = { 374 "A":{ 375 "C":0.167653, 376 "D":4.43521, 377 "E":5.56325, 378 "F":0.597923, 379 "G":1.8685, 380 "H":0.005, 381 "I":0.005, 382 "K":0.592784, 383 "L":0.16024, 384 "M":0.005, 385 "N":0.617509, 386 "P":1.00981, 387 "Q":0.005, 388 "R":0.0744808, 389 "S":8.594200000000001, 390 "T":24.1422, 391 "V":24.8094, 392 "W":0.005, 393 "Y":0.005 394 }, 395 "C":{ 396 "D":0.005, 397 "E":0.005, 398 "F":0.362959, 399 "G":0.0489798, 400 "H":0.005, 401 "I":0.005, 402 "K":0.005, 403 "L":0.005, 404 "M":0.005, 405 "N":0.0604932, 406 "P":0.005, 407 "Q":0.005, 408 "R":2.86364, 409 "S":1.12195, 410 "T":0.005, 411 "V":0.005, 412 "W":5.49894, 413 "Y":8.34835 414 }, 415 "D":{ 416 "E":12.1233, 417 "F":0.005, 418 "G":10.3969, 419 "H":2.31779, 420 "I":0.145124, 421 "K":0.894313, 422 "L":0.005, 423 "M":0.005, 424 "N":29.4087, 425 "P":0.005, 426 "Q":0.005, 427 "R":0.0674539, 428 "S":0.427881, 429 "T":0.630395, 430 "V":2.91786, 431 "W":0.005, 432 "Y":2.28154 433 }, 434 "E":{ 435 "F":0.005, 436 "G":14.7801, 437 "H":0.005, 438 "I":0.0390512, 439 "K":23.9626, 440 "L":0.129839, 441 "M":0.005, 442 "N":0.201526, 443 "P":0.005, 444 "Q":3.20656, 445 "R":0.0251632, 446 "S":0.005, 447 "T":0.458743, 448 "V":2.19952, 449 "W":0.005, 450 "Y":0.005 451 }, 452 "F":{ 453 "G":0.005, 454 "H":0.005, 455 "I":1.48288, 456 "K":0.005, 457 "L":7.48781, 458 "M":0.005, 459 "N":0.005, 460 "P":0.0342252, 461 "Q":0.005, 462 "R":0.005, 463 "S":4.27939, 464 "T":0.114512, 465 "V":2.28, 466 "W":0.005, 467 "Y":4.12728 468 }, 469 "G":{ 470 "H":0.005, 471 "I":0.005, 472 "K":0.279425, 473 "L":0.0489798, 474 "M":0.0489798, 475 "N":0.0604932, 476 "P":0.005, 477 "Q":0.0604932, 478 "R":13.4379, 479 "S":6.27966, 480 "T":0.0489798, 481 "V":2.79622, 482 "W":2.8258, 483 "Y":0.005 484 }, 485 "H":{ 486 "I":0.005, 487 "K":0.22406, 488 "L":1.76382, 489 "M":0.005, 490 "N":8.59876, 491 "P":13.9444, 492 "Q":18.5465, 493 "R":6.84405, 494 "S":0.7251570000000001, 495 "T":0.95956, 496 "V":0.827479, 497 "W":0.005, 498 "Y":47.4889 499 }, 500 "I":{ 501 "K":0.817481, 502 "L":9.102460000000001, 503 "M":17.3064, 504 "N":0.987028, 505 "P":0.005, 506 "Q":0.0342252, 507 "R":1.34069, 508 "S":0.7400910000000001, 509 "T":9.36345, 510 "V":24.8231, 511 "W":0.005, 512 "Y":0.114512 513 }, 514 "K":{ 515 "L":0.005, 516 "M":4.09564, 517 "N":10.6655, 518 "P":0.111928, 519 "Q":13.0705, 520 "R":39.8897, 521 "S":0.005, 522 "T":4.04802, 523 "V":0.128065, 524 "W":0.005, 525 "Y":0.005 526 }, 527 "L":{ 528 "M":11.3839, 529 "N":0.005, 530 "P":9.83095, 531 "Q":2.89048, 532 "R":0.586757, 533 "S":6.14396, 534 "T":0.005, 535 "V":2.95344, 536 "W":1.37031, 537 "Y":0.005 538 }, 539 "M":{ 540 "N":0.201526, 541 "P":0.005, 542 "Q":0.005, 543 "R":3.28652, 544 "S":0.392575, 545 "T":7.41313, 546 "V":14.7683, 547 "W":0.005, 548 "Y":0.579198 549 }, 550 "N":{ 551 "P":0.344848, 552 "Q":0.342068, 553 "R":0.16024, 554 "S":14.5699, 555 "T":4.54206, 556 "V":0.0744808, 557 "W":0.005, 558 "Y":5.06475 559 }, 560 "P":{ 561 "Q":3.04502, 562 "R":0.404723, 563 "S":14.249, 564 "T":4.33701, 565 "V":0.005, 566 "W":0.005, 567 "Y":0.005 568 }, 569 "Q":{ 570 "R":10.6746, 571 "S":0.16024, 572 "T":0.203091, 573 "V":0.005, 574 "W":0.0443298, 575 "Y":0.005 576 }, 577 "R":{ 578 "S":8.350239999999999, 579 "T":0.928203, 580 "V":0.279425, 581 "W":5.96564, 582 "Y":0.005 583 }, 584 "S":{ 585 "T":6.34079, 586 "V":0.862637, 587 "W":1.10156, 588 "Y":0.933142 589 }, 590 "T":{ 591 "V":0.005, 592 "W":0.005, 593 "Y":0.490608 594 }, 595 "V":{ 596 "W":0.005, 597 "Y":1.35482 598 }, 599 "W":{ 600 "Y":0.005 601 } 602}; 603 604 605 606lfunction models.protein.HIVWm.frequencies (model, namespace, datafilter) { 607 model[utility.getGlobalValue("terms.efv_estimate")] = 608 {{0.0377494} 609 {0.0240105} 610 {0.0342034} 611 {0.0618606} 612 {0.0422741} 613 {0.0838496} 614 {0.0156076} 615 {0.0983641} 616 {0.0641682} 617 {0.0577867} 618 {0.0158419} 619 {0.0891129} 620 {0.0458601} 621 {0.0437824} 622 {0.057321} 623 {0.0550846} 624 {0.0813774} 625 {0.0515639} 626 {0.019597} 627 {0.0205847}}; 628 model[utility.getGlobalValue("terms.model.efv_estimate_name")] = utility.getGlobalValue("terms.frequencies.predefined"); 629 (model[utility.getGlobalValue("terms.parameters")])[utility.getGlobalValue("terms.model.empirical")] = 0; 630 return model; 631}