1LoadFunctionLibrary("GrabBag"); 2LoadFunctionLibrary("libv3/all-terms.bf"); 3 4/** @module parameters */ 5 6parameters.infinity = 1e10; 7 8/** 9 * Applies a namespace to parameter ids 10 * @name parameters.ApplyNameSpace 11 * @param {String} id 12 * @param {String} namespace 13 */ 14function parameters.ApplyNameSpace(id, namespace) { 15 if (Type(namespace) == "String") { 16 if (Abs(namespace) > 0) { 17 return namespace + "." + id; 18 } 19 } 20 return id; 21} 22 23/** 24 * @name parameters.UnconstrainParameterSet 25 * @param {LikelihoodFunction} lf - the likelihood function to operate on 26 * @param {Matrix} set - set of parameters to unconstrain 27 * @returns nothing 28 */ 29function parameters.UnconstrainParameterSet(lf, set) { 30 ExecuteCommands("GetString(parameters.UnconstrainParameterSet.info, `lf`, -1)"); 31 if (None == set) { 32 set = { 33 { 34 terms.parameters.global_constrained, terms.parameters.local_constrained 35 } 36 }; 37 } 38 for (parameters.UnconstrainParameterSet.s = 0; parameters.UnconstrainParameterSet.s < Columns(set); parameters.UnconstrainParameterSet.s += 1) { 39 parameters.UnconstrainParameterSet.m = parameters.UnconstrainParameterSet.info[set[parameters.UnconstrainParameterSet.s]]; 40 for (parameters.UnconstrainParameterSet.i = 0; parameters.UnconstrainParameterSet.i < Columns(parameters.UnconstrainParameterSet.m); parameters.UnconstrainParameterSet.i += 1) { 41 Eval(parameters.UnconstrainParameterSet.m[parameters.UnconstrainParameterSet.i] + "=" + Eval(parameters.UnconstrainParameterSet.m[parameters.UnconstrainParameterSet.i])); 42 } 43 } 44} 45 46/** 47 * @name parameters.DeclareGlobal 48 * @param {String} id 49 * @param {Matrix} cache 50 * @returns nothing 51 */ 52function parameters.DeclareGlobal(id, cache) { 53 if (Type(id) == "String") { 54 if (Abs(id)) { 55 if (Type(cache) == "AssociativeList") { 56 if (Abs(cache[id]) > 0) { 57 return; 58 } else { 59 cache[id] = 1; 60 } 61 } 62 ExecuteCommands("global `id` = 1;"); 63 } 64 } else { 65 if (Type(id) == "AssociativeList") { 66 parameters.DeclareGlobal.var_count = Abs(id); 67 parameters.DeclareGlobal.names = Columns(id); 68 for (parameters.DeclareGlobal.k = 0; parameters.DeclareGlobal.k < parameters.DeclareGlobal.var_count; parameters.DeclareGlobal.k += 1) { 69 parameters.DeclareGlobal(parameters.DeclareGlobal.names[parameters.DeclareGlobal.k], cache); 70 } 71 } else { 72 if (Type(id) == "Matrix") { 73 parameters.DeclareGlobal.var_count = Columns(id) * Rows(id); 74 for (parameters.DeclareGlobal.k = 0; parameters.DeclareGlobal.k < parameters.DeclareGlobal.var_count; parameters.DeclareGlobal.k += 1) { 75 parameters.DeclareGlobal(id[parameters.DeclareGlobal.k], cache); 76 } 77 } 78 } 79 } 80 81 return cache; 82} 83 84/** 85 * @name parameters.DeclareGlobalWithRanges 86 * @param {String} id variable id 87 * @param {Number} init initial value (could be None) 88 * @param {Number} lb lower bound (could be None) 89 * @param {Number} ub upper bound (could be None) 90 * @returns nothing 91 */ 92function parameters.DeclareGlobalWithRanges(id, init, lb, ub) { 93 if (Type (id) == "String") { 94 if (None != init) { 95 ExecuteCommands("global `id` = " + init); 96 } else { 97 ExecuteCommands("global `id`; "); 98 } 99 100 if (None != lb) { 101 ExecuteCommands ("`id` :> " + lb); 102 } 103 if (None != ub) { 104 ExecuteCommands ("`id` :< " + ub); 105 } 106 } else { 107 if (Type(id) == "AssociativeList") { 108 parameters.DeclareGlobalWithRanges.var_count = Abs(id); 109 parameters.DeclareGlobalWithRanges.names = Columns(id); 110 for (parameters.DeclareGlobalWithRanges.k = 0; parameters.DeclareGlobalWithRanges.k < parameters.DeclareGlobalWithRanges.var_count; parameters.DeclareGlobalWithRanges.k += 1) { 111 parameters.DeclareGlobalWithRanges(parameters.DeclareGlobalWithRanges.names[parameters.DeclareGlobalWithRanges.k], init, lb, ub); 112 } 113 } else { 114 if (Type(id) == "Matrix") { 115 parameters.DeclareGlobalWithRanges.var_count = Columns(id) * Rows(id); 116 for (parameters.DeclareGlobalWithRanges.k = 0; parameters.DeclareGlobalWithRanges.k < parameters.DeclareGlobalWithRanges.var_count; parameters.DeclareGlobalWithRanges.k += 1) { 117 parameters.DeclareGlobalWithRanges(id[parameters.DeclareGlobalWithRanges.k], init, lb, ub); 118 } 119 } 120 } 121 } 122} 123 124function parameters.DeclareCategory.helper (dict, key, default) { 125 if (dict / key) { 126 return dict[key]; 127 } 128 return default; 129} 130 131/** 132 * @name parameters.DeclareCategory 133 * @param {Dict} def category definition components 134 */ 135function parameters.DeclareCategory (def) { 136 ExecuteCommands ("category " + def[terms.id] + "= (" + 137 Join (",", 138 utility.Map ({"0": terms.category.bins, 139 "1": terms.category.weights, 140 "2": terms.category.represent, 141 "3": terms.category.PDF, 142 "4": terms.category.CDF, 143 "5": terms.lower_bound, 144 "6": terms.upper_bound, 145 "7": terms.category.dCDF, 146 "8": terms.category.HMM}, 147 "_value_", 148 'parameters.DeclareCategory.helper(def[terms.category.category_parameters], _value_, "")') 149 ) + ");"); 150 151} 152 153 154/** 155 * @name parameters.NormalizeRatio 156 * @param {Number} n 157 * @param {Number} d 158 * @returns n/d 159 */ 160function parameters.NormalizeRatio(n, d) { 161 if (d == 0) { 162 if (n == 0) { 163 return 1; 164 } else { 165 return parameters.infinity; 166 } 167 } 168 return n / d; 169} 170 171/** 172 * Sets value of passed parameter id 173 * @name parameters.SetValue 174 * @param {String} id - id of parameter to set value to 175 * @param {Number} value - value to set 176 * @returns nothing 177 */ 178function parameters.SetValue(id, value) { 179 Eval("`id` = " + value); 180} 181 182/** 183 * Sets value of passed parameter tree.branch.id 184 * @name parameters.SetLocalValue 185 * @param {String} tree - id of tree 186 * @param {String} branch - id of branch 187 * @param {String} id - id of parameter to set value to 188 * @param {Number} value - value to set 189 * @returns nothing 190 */ 191function parameters.SetLocalValue(tree, branch, id, value) { 192 Eval("`tree`.`branch`.`id` = " + value); 193} 194 195/** 196 * Sets value of passed parameter id 197 * @name parameters.SetValues 198 * @param {dict} desc -> {id : id, mle : value} 199 * @returns nothing 200 */ 201 202 203function parameters.SetValues(set) { 204 //console.log ("\n==== parameters.SetValues === \n"); 205 if (Type (set) == "AssociativeList") { 206 for (_key_, _value_; in; set) { 207 if (parameters.IsIndependent (_value_[terms.id])) { 208 //console.log (_value_[terms.id] + " (" + Eval (_value_[terms.id]) + ") =>" + _value_[terms.fit.MLE]); 209 parameters.SetValue (_value_[terms.id], _value_[terms.fit.MLE]); 210 } else { 211 //console.log (_value_[terms.id] + " is NOT independent"); 212 } 213 } 214 215 } 216} 217 218/** 219 * Ensures that the mean of parameters in a set is maintained 220 * @name parameters.ConstrainMeanOfSet 221 * @param {Dict/Matrix} set - list of variable ids 222 * @param {Dict/Matrix} weights - weights to apply 223 * @param {Number} mean - desired mean 224 * @param {String} namespace - desired mean 225 * @returns nothing 226 */ 227lfunction parameters.ConstrainMeanOfSet (set, weights, mean, namespace) { 228 if (Type (set) == "AssociativeList") { 229 unscaled = utility.Map (set, "_name_", "_name_ + '_scaler_variable'"); 230 constraint = utility.MapWithKey (unscaled, "_key_", "_name_", "_name_ + '*' + `&weights`[_key_]"); 231 } else { 232 if (Type (set) == "Matrix") { 233 unscaled = utility.Map (set, "_name_", "_name_ + '_scaler_variable'"); 234 constraint = utility.MapWithKey (unscaled, "_key_", "_name_", "_name_ + '*' + `&weights`[_key_[0]+_key_[1]]"); 235 } 236 else { 237 return; 238 } 239 } 240 241 scaler_variables = {}; 242 utility.ForEach (unscaled, "_name_", 'parameters.DeclareGlobal (_name_, null)'); 243 global_scaler = namespace + ".scaler_variable"; 244 parameters.SetConstraint (global_scaler, Join ("+", constraint), "global"); 245 utility.ForEach (set, "_name_", ' 246 `&scaler_variables`["Mean scaler variable for " + _name_] = _name_ + "_scaler_variable"; 247 parameters.SetValue (_name_ + "_scaler_variable", _name_); 248 parameters.SetConstraint (_name_, "(" + `&mean` + ")*" + _name_ + "_scaler_variable/`global_scaler`", ""); 249 '); 250 251 return {^'terms.global' : scaler_variables}; 252} 253 254/** 255 * Given a set of parameters [x1,x2,...] set x1 to a given value, and constrain x2 := x1, x3 := x1... 256 * @name parameters.ConstrainParameterSet 257 * @param {Dict} set - list of variable ids (as values) 258 * @param {Number/null} value - if is a Number, then set x1 to this value, 259 if null, set x1 to the mean of all values 260 * @returns {Dict} the list of constrained parameters 261 */ 262lfunction parameters.ConstrainParameterSet (set, value) { 263 if (Type (set) == "AssociativeList") { 264 if (utility.Array1D (set) > 1) { 265 _key_name = set["VALUEINDEXORDER"][0]; 266 if (None == value) { 267 value = + (utility.Map (set, "_id_", "Eval(_id_)")); 268 value = value / utility.Array1D (set); 269 } 270 parameters.SetValue (_key_name, value); 271 constraints = {}; 272 utility.ForEach (set, "_id_", ' 273 if (_id_ != `&_key_name`) { 274 parameters.SetConstraint (_id_, `&_key_name`, ""); 275 (`&constraints`) + _id_; 276 } 277 '); 278 279 return constraints; 280 } 281 } 282 283 return {}; 284 285} 286 287/** 288 * Given a set of parameters [x1,x2,...] set x_i := Eval (x_i) 289 * @name parameters.ConstrainParameterSet 290 * @param {Dict} set - list of variable ids (as values) 291 * @returns {Dict} the list of constrained parameters 292 */ 293lfunction parameters.FixParameterSet (set) { 294 if (Type (set) == "AssociativeList") { 295 if (utility.Array1D (set) > 1) { 296 constraints = {}; 297 utility.ForEach (set, "_id_", ' 298 parameters.SetConstraint (_id_, "" + Eval (_id_), ""); 299 (`&constraints`) + _id_; 300 '); 301 302 return constraints; 303 } 304 } 305 306 return {}; 307} 308 309 310 311/** 312 * Returns mean of values 313 * @name parameters.Mean 314 * @param {Matrix} values - values to return mean of 315 * @param {Matrix} weights - weights to multiply values by 316 * @param {Number} d - does nothing 317 * @returns {Number} mean 318 */ 319lfunction parameters.Mean(values, weights, d) { 320 m = 0; 321 d = Rows(values) * Columns(values); 322 for (i = 0; i < d; i += 1) { 323 m += Eval(values[i]) * Eval(weights[i]); 324 } 325 return m; 326} 327 328/** 329 * Quotes the argument 330 * @name parameters.Quote 331 * @param {String} arg - string to be quoted 332 * @returns {String} string in quotes 333 */ 334function parameters.Quote(arg) { 335 if (Type(arg) == "String") { 336 return "\"" + (arg && 2) + "\""; 337 } 338 return arg; 339} 340 341/** 342 * @name parameters.AppendMultiplicativeTerm 343 * @param {String} expression - the matrix to modify 344 * @param {String} term - the multiplier to append 345 * @returns {String} (expression) * (term) 346 */ 347lfunction parameters.AppendMultiplicativeTerm (expression, term) { 348 if (Type (expression) == "String") { 349 if (Abs (expression)) { 350 return "(" + expression + ")*(" + term + ")"; 351 } 352 return term; 353 } 354 return expression; 355} 356 357/** 358 * @name parameters.AddMultiplicativeTerm 359 * @param {Matrix} matrix - matrix to scale 360 * @param {Number} term - scalar to multiply matrix by 361 * @param {Number} do_empties - if element matrix is empty, fill with term 362 * @returns {Matrix} New matrix 363 */ 364lfunction parameters.AddMultiplicativeTerm(matrix, term, do_empties) { 365 366 if (Abs(term) > 0) { 367 __N = Rows(matrix); 368 369 for (__r = 0; __r < __N; __r += 1) { 370 for (__c = 0; __c < __N; __c += 1) { 371 if (__r != __c) { 372 if (Abs(matrix[__r][__c])) { 373 matrix[__r][__c] = "(" + matrix[__r][__c] + ")*(" + term + ")"; 374 } else { 375 if (do_empties) { 376 matrix[__r][__c] = term; 377 } 378 } 379 } 380 } 381 } 382 } 383 384 return matrix; 385} 386 387/** 388 * @name parameters.StringMatrixToFormulas 389 * @param {String} id - matrix to scale 390 * @param {Matrix} matrix - if element matrix is empty, fill with term 391 * @returns nothing 392 */ 393function parameters.StringMatrixToFormulas(id, matrix) { 394 __N = Rows(matrix); 395 __M = Columns(matrix); 396 397 ExecuteCommands("`id` = {__N,__M}"); 398 399 __expr := "`id`[__r][__c] := " + matrix[__r][__c]; 400 401 for (__r = 0; __r < __N; __r += 1) { 402 for (__c = 0; __c < __M; __c += 1) { 403 404 if (Abs(matrix[__r][__c])) { 405 ExecuteCommands(__expr); 406 } 407 } 408 } 409} 410 411/** 412 * @name parameters.GenerateAttributedNames 413 * @param {String} prefix 414 * @param {Dictionary} attributes 415 * @param {String} delimiter 416 */ 417function parameters.GenerateAttributedNames(prefix, attributes, delimiter) { 418 if (delimiter == None) { 419 delimiter = "_"; 420 } 421 parameters.generate_names.holder = {}; 422 for (parameters.generate_names.k = 0; parameters.generate_names.k < Columns(attributes); parameters.generate_names.k += 1) { 423 parameters.generate_names.holder + (prefix + delimiter + attributes[parameters.generate_names.k]); 424 } 425 return parameters.generate_names.holder; 426} 427 428/** 429 * @name parameters.GenerateSequentialNames 430 * @param {String} prefix 431 * @param {Number} count 432 * @param {String} delimiter 433 * @returns {Matrix} 1 x <count> row vector of generated names 434 */ 435 436lfunction parameters.GenerateSequentialNames(prefix, count, delimiter) { 437 if (delimiter == None) { 438 delimiter = "_"; 439 } 440 holder = {}; 441 for (k = 0; k < count; k += 1) { 442 holder + (prefix + delimiter + k); 443 } 444 return holder; 445} 446 447/** 448 * @name parameters.GetRange 449 * @param id 450 * @returns variable range 451 */ 452lfunction parameters.GetRange(id) { 453 454 if (Type(id) == "String") { 455 GetInformation (range, ^id); 456 return { 457 ^"terms.lower_bound" : range[1], 458 ^"terms.upper_bound" : range[2] 459 }; 460 } 461 io.ReportAnExecutionError ("An invalid combination of parameters was passed to parameters.GetRange. ID = " + id); 462 return None; 463} 464 465 466/** 467 * @name parameters.SetRange 468 * @param id 469 * @param ranges 470 * @returns nothing 471 */ 472function parameters.SetRange(id, ranges) { 473 if (Type(id) == "String") { 474 if (Abs(id)) { 475 if (Type(ranges) == "AssociativeList") { 476 if (Abs(ranges[terms.lower_bound])) { 477 ExecuteCommands("`id` :> " + ranges[terms.lower_bound]); 478 } 479 if (Abs(ranges[terms.upper_bound])) { 480 ExecuteCommands("`id` :< " + ranges[terms.upper_bound]); 481 } 482 return 0; 483 } 484 } 485 } else { 486 if (Type(id) == "AssociativeList") { 487 parameters.SetRange.var_count = Abs(id); 488 for (parameters.SetRange.k = 0; parameters.SetRange.k < parameters.SetRange.var_count; parameters.SetRange.k += 1) { 489 parameters.SetRange(id[parameters.SetRange.k], ranges); 490 } 491 return 0; 492 } 493 } 494 io.ReportAnExecutionError ("An invalid combination of parameters was passed to parameters.SetRange. ID = " + id + ", range = " + ranges); 495 496} 497 498 499/** 500 * Check if parameter is independent 501 * @name parameters.IsIndependent 502 * @param parameter - id of parameter to check 503 * @returns {Bool} TRUE if independent, FALSE otherwise 504 */ 505function parameters.IsIndependent(parameter) { 506 507 SetParameter(HBL_EXECUTION_ERROR_HANDLING,1,0); 508 GetString(info, ^parameter, -1); 509 SetParameter(HBL_EXECUTION_ERROR_HANDLING,0,0); 510 if (Abs (LAST_HBL_EXECUTION_ERROR)) { // variable does not exist 511 return TRUE; 512 } 513 if (Type(info) == "AssociativeList") { 514 return (utility.CheckKey(info, "Local", "Matrix") && utility.CheckKey(info, "Global", "Matrix")) == FALSE; 515 } 516 return TRUE; 517} 518 519lfunction parameters.GetConstraint(parameter) { 520 GetString(info, ^parameter, -2); 521 return info; 522} 523 524/** 525 * sets constraint on parameter 526 * @name parameters.SetConstraint 527 * @param {String} or {AssociativeList} id - id(s) of parameter(s) to set constraint on 528 * @param {Number} value - the constraint to set on the parameter 529 * @param {String} global_tag - the global namespace of the parameter 530 * @returns nothing 531 */ 532function parameters.SetConstraint(id, value, global_tag) { 533 if (Type(id) == "String") { 534 if (Abs(id)) { 535 ExecuteCommands("`global_tag` `id` := " + value); 536 } 537 } else { 538 if (Type(id) == "AssociativeList" && Type(value) == "AssociativeList") { 539 parameters.SetConstraint.var_count = Abs(id); 540 for (parameters.SetConstraint.k = 0; parameters.SetConstraint.k < parameters.SetConstraint.var_count; parameters.SetConstraint.k += 1) { 541 parameters.SetConstraint(id[parameters.SetConstraint.k], 542 value[parameters.SetConstraint.k], 543 global_tag); 544 } 545 } 546 } 547} 548 549/** 550 * constrain x to be x := C * (current value of x) 551 * @name parameters.SetProprtionalConstraint 552 * @param {String} id of parameter(s) to set constraint on 553 * @param {String} scaler variable - the ID of the 'C' scaler variable above; could also be an expression 554 * @returns nothing 555 */ 556function parameters.SetProprtionalConstraint(id, scaler) { 557 if (Type(id) == "String") { 558 if (Abs(id)) { 559 //console.log ("`id` => " + Eval (id)); 560 ExecuteCommands("`id` := (`scaler`)*" + Eval (id)); 561 } 562 } 563} 564 565/** 566 * constraint set of parameters 567 * @name parameters.ConstrainSets 568 * @param {AssociativeList} set1 - 569 * @param {AssociativeList} set2 - 570 * @returns nothing 571 */ 572function parameters.ConstrainSets(set1, set2) { 573 parameters.ConstrainSets.tags = Rows(set1); 574 for (parameters.ConstrainSets.k = 0; parameters.ConstrainSets.k < Abs(set1); parameters.ConstrainSets.k += 1) { 575 576 if (Type(set2[parameters.ConstrainSets.tags[parameters.ConstrainSets.k]]) == "String") { 577 ExecuteCommands(set2[parameters.ConstrainSets.tags[parameters.ConstrainSets.k]] + ":=" + 578 set1[parameters.ConstrainSets.tags[parameters.ConstrainSets.k]]); 579 } 580 } 581} 582 583/** 584 * Removes a constraint from a parameter 585 * @name parameters.RemoveConstraint 586 * @param {String} id - id of parameter to remove constraint from 587 * @returns nothing 588 */ 589function parameters.RemoveConstraint(id) { 590 if (Type(id) == "String") { 591 if (Abs(id)) { 592 Eval("`id` = " + Eval(id)); 593 } 594 } else { 595 if (Type(id) == "AssociativeList") { 596 return parameters.RemoveConstraint(Columns(id)); 597 } 598 if (Type(id) == "Matrix") { 599 parameters.RemoveConstraint.var_count = Columns(id) * Rows(id); 600 for (parameters.RemoveConstraint.k = 0; parameters.RemoveConstraint.k < parameters.RemoveConstraint.var_count; parameters.RemoveConstraint.k += 1) { 601 parameters.RemoveConstraint(id[parameters.RemoveConstraint.k]); 602 } 603 } 604 } 605} 606 607 608/** 609 * Copies definitions from target to source 610 * @name parameters.helper.copy_definitions 611 * @param {Dictionary} target - the target dictionary 612 * @param {Dictionary} source - the source element to copy to target 613 * @returns nothing 614 */ 615function parameters.helper.copy_definitions(target, source) { 616 parameters.helper.copy_definitions.key_iterator = {1,2}; 617 parameters.helper.copy_definitions.key_iterator [0] = terms.local; 618 parameters.helper.copy_definitions.key_iterator [1] = terms.global; 619 620 for (parameters.helper.copy_definitions.i = 0; parameters.helper.copy_definitions.i < Columns(parameters.helper.copy_definitions.key_iterator); parameters.helper.copy_definitions.i += 1) { 621 parameters.helper.copy_definitions.key = parameters.helper.copy_definitions.key_iterator[parameters.helper.copy_definitions.i]; 622 if (Type(source[parameters.helper.copy_definitions.key]) == "AssociativeList") { 623 utility.EnsureKey (target, parameters.helper.copy_definitions.key); 624 target[parameters.helper.copy_definitions.key] * source[parameters.helper.copy_definitions.key]; 625 } 626 } 627 628 if (utility.Has (source, terms.category, "AssociativeList")) { 629 utility.EnsureKey (target, terms.category); 630 (target[terms.category])[(source[terms.category])[terms.id]] = (source[terms.category])[terms.description]; 631 } 632 633 return target; 634} 635 636/** 637 * @name pparameters.SetStickBreakingDistribution 638 * @param {AssociativeList} parameters 639 * @param {Matrix} values 640 * @returns nothing 641 */ 642 643lfunction parameters.SetStickBreakingDistribution (parameters, values) { 644 rate_count = Rows (values); 645 left_over = 1; 646 647 for (i = 0; i < rate_count; i += 1) { 648 parameters.SetValue ((parameters["rates"])[i], values[i][0]); 649 if (i < rate_count - 1) { 650 break_here = values[i][1] / left_over; 651 parameters.SetValue ((parameters["weights"])[i], break_here); 652 left_over = left_over * (1-break_here); 653 } 654 } 655} 656 657/** 658 * @name pparameters.GetStickBreakingDistribution 659 * @param {AssociativeList} parameters 660 * @returns {Matrix} computed distribution (Nx2) 661 */ 662 663lfunction parameters.GetStickBreakingDistribution (parameters) { 664 rate_count = utility.Array1D (parameters["rates"]); 665 distribution = {rate_count, 2}; 666 667 current_weight = 1; 668 669 for (i = 0; i < rate_count; i += 1) { 670 distribution [i][0] = Eval ((parameters["rates"])[i]); 671 if (i < rate_count - 1) { 672 distribution [i][1] = current_weight * Eval ((parameters["weights"])[i]); 673 current_weight = current_weight * (1-Eval ((parameters["weights"])[i])); 674 } else { 675 distribution [i][1] = current_weight; 676 } 677 } 678 return distribution; 679} 680 681 682/** 683 * @name parameters.helper.stick_breaking 684 * @param {AssociativeList} parameters 685 * @param {Matrix} initial_values 686 * @returns weights 687 */ 688lfunction parameters.helper.stick_breaking(parameters, initial_values) { 689 left_over = "1*"; // handle the case of ONE component 690 weights = {}; 691 accumulator = 1; 692 693 694 for (k = 0; k < Abs(parameters); k += 1) { 695 if (None != initial_values) { 696 vid = parameters[k]; ^ (vid) = initial_values[k] / accumulator; 697 accumulator = accumulator * (1 - ^ (vid)); 698 } 699 weights[k] = left_over + parameters[k]; 700 left_over += "(1-" + parameters[k] + ")*"; 701 } 702 703 weights[k] = left_over[0][Abs(left_over) - 2]; 704 return weights; 705} 706 707/** 708 * Prints matrix to screen 709 * @name parameters.helper.dump_matrix 710 * @param {Matrix} matrix 711 * @returns nothing 712 */ 713lfunction parameters.helper.dump_matrix(matrix) { 714 for (i = 0; i < Rows( ^ matrix); i += 1) { 715 for (j = 0; j < Columns( ^ matrix); j += 1) { 716 ExecuteCommands("GetString (cell, `matrix`, i, j)"); 717 fprintf(stdout, "`matrix`[", i, "][", j, "] := ", cell, "\n"); 718 } 719 } 720 return None; 721} 722 723/** 724 * Sets tree lengths to initial values 725 * @name parameters.helper.tree_lengths_to_initial_values 726 * @param dict - a [0 to N-1] dictrionary of tree objects 727 * @param type - codon or nucleotide 728 * @returns {Dictionary} dictionary of initial branch lengths 729 */ 730lfunction parameters.helper.tree_lengths_to_initial_values(dict, type) { 731 732 components = utility.Array1D(dict); 733 734 if (type == "codon") { // TODO 735 factor = 1; 736 } else { 737 factor = 1; 738 } 739 740 result = {}; 741 742 for (i = 0; i < components; i += 1) { 743 this_component = {}; 744 745 746 utility.ForEachPair((dict[i])[ utility.getGlobalValue("terms.branch_length")], "_branch_name_", "_branch_length_", 747 " 748 `&this_component`[_branch_name_] = {utility.getGlobalValue('terms.fit.MLE') : `&factor`*_branch_length_} 749 "); 750 result[i] = this_component; 751 752 } 753 754 return { utility.getGlobalValue("terms.branch_length"): result 755 }; 756} 757 758/** 759 * Profiles likelihood function based on covariance precision level 760 * @name parameters.GetProfileCI 761 * @param {String} id - covariance parameter id 762 * @param {LikelihoodFunction} lf - likelihood function to profile 763 * @param {Number} - Covariance precision level 764 * @returns {Dictionary} a dictionary containing profiling information 765 */ 766function parameters.GetProfileCI(id, lf, level) { 767 768 utility.ToggleEnvVariable("COVARIANCE_PRECISION", level); 769 utility.ToggleEnvVariable("COVARIANCE_PARAMETER", id); 770 CovarianceMatrix(parameters.GetProfileCI.mx, * lf); 771 utility.ToggleEnvVariable("COVARIANCE_PRECISION", None); 772 utility.ToggleEnvVariable("COVARIANCE_PARAMETER", None); 773 774 //TODO: used to be "`terms.lower_bound`". 775 return { 776 terms.lower_bound: parameters.GetProfileCI.mx[0], 777 terms.fit.MLE: parameters.GetProfileCI.mx[1], 778 terms.upper_bound: parameters.GetProfileCI.mx[2] 779 }; 780} 781 782/** 783 * Geneate an HBL string needed to define a parameter 784 * @name parameters.ExportParameterDefinition 785 * @param {String} id - the name of the parameter to export; 786 * @returns {String} the string for an HBL definition of the parameter 787 * e.g. "global x := 2.3; x :> 1; x :< 3"; note that the definition will 788 * be NOT recursive, so if x depends on y and z, then y and z need to be 789 * exported separately 790 */ 791 792lfunction parameters.ExportParameterDefinition (id) { 793 GetString (parameter_definition, ^id, -3); 794 return parameter_definition; 795} 796 797/** 798 * Copy local parameters from a node to 'template' variables 799 * @name parameters.SetLocalModelParameters 800 * @param {Dict} model - model description 801 * @param {String} tree id 802 * @param {String} node id 803 804 * e.g. 805 * set alpha = Tree.Node.alpha 806 * set beta = Tree.Node.beta 807 */ 808 809lfunction parameters.SetLocalModelParameters (model, tree, node) { 810 node_name = tree + "." + node + "."; 811 utility.ForEach ((model[utility.getGlobalValue("terms.parameters")])[utility.getGlobalValue("terms.local")], "_parameter_", ' 812 ^_parameter_ = ^(`&node_name` + _parameter_); 813 '); 814 815} 816 817/** 818 * Set category variables to their mean values for branch length calculations 819 * @name parameters.SetLocalModelParameters 820 * @param {Dict} model - model description 821*/ 822 823lfunction parameters.SetCategoryVariables (model) { 824 825 cat_vars = (model[utility.getGlobalValue("terms.parameters")])[utility.getGlobalValue("terms.category")]; 826 cat_var_count = utility.Array1D (cat_vars); 827 cat_var_substitutions = {}; 828 829 if (cat_var_count) { 830 cat_vars = Rows (cat_vars); 831 for (i = 0; i < cat_var_count; i+=1) { 832 cv = rate_variation.compute_mean (cat_vars[i]); 833 cat_var_substitutions [cat_vars[i]] = cv; 834 parameters.SetValue (cat_vars[i], cv); 835 } 836 } 837 838 return cat_var_substitutions; 839} 840 841/** 842 * Given a set of strings, convert them to valid IDs that don't clash following renames 843 * @name parameters.ValidateID 844 * @param {Array/Dict} ids - the ids to validate (should be unique), if dict, should be i -> ID 845 * @returns {Dictionary} a dictionary containing old -> new id map 846*/ 847 848lfunction parameters.ValidateIDs (ids) { 849 result = {}; 850 N = utility.Array1D (ids); 851 for (i = 0; i < N; i+=1) { 852 try_name = ids [i] && 7; 853 if (result / try_name ) { 854 for (k = 2; ; k+=1) { 855 tagged = try_name + "_" + k; 856 if (!(result/tagged)) { 857 try_name = tagged; 858 break; 859 } 860 } 861 } 862 result [ids[i]] = try_name; 863 } 864 return result; 865} 866 867