Changeset 1383
- Timestamp:
- 08/18/11 16:02:53 (13 years ago)
- Files:
-
- 7 modified
-
applications/robust/main.cpp (modified) (15 diffs)
-
applications/robust/robustlib.cpp (modified) (6 diffs)
-
applications/robust/robustlib.h (modified) (5 diffs)
-
library/bdm/base/itpp/CMakeLists.txt (modified) (1 diff)
-
library/bdm/estim/arx.h (modified) (1 diff)
-
library/bdm/stat/exp_family.h (modified) (1 diff)
-
library/tests/tutorial/arx_simple.cpp (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
applications/robust/main.cpp
r1379 r1383 27 27 //const int utility_constant = 5; 28 28 29 const int max_model_order = 2;30 const double apriorno =0.005;29 const int max_model_order = 1; 30 const double apriorno = 0.01; 31 31 32 32 /* … … 80 80 int predicted_channel; 81 81 mat* data_matrix; 82 vec predictions; 82 83 83 84 model(set<pair<int,int>> ar_components, //funkcie treidz model-konstruktor … … 116 117 { 117 118 V0 = apriorno * eye(ar_components.size()+2); //aj tu konst 118 V0(0,0) = 1;119 V0(0,0) = 0; 119 120 my_arx->set_constant(true); 120 121 } … … 123 124 124 125 V0 = apriorno * eye(ar_components.size()+1);//menit konstantu 125 V0(0,0) = 1;126 V0(0,0) = 0; 126 127 my_arx->set_constant(false); 127 128 … … 143 144 } 144 145 if(my_rarx!=NULL) 145 { //pre robus r priradi az tu do data_vector aj rpedikciu146 { //pre robust priradi az tu do data_vector aj predikciu 146 147 data_vector.ins(0,(*data_matrix).get(predicted_channel,time)); 147 148 my_rarx->bayes(data_vector); … … 188 189 else 189 190 { 190 mat samples = my_arx->posterior().sample_mat(sample_size); 191 mat samples = my_arx->posterior().sample_mat(sample_size); 191 192 192 193 vec sample_prediction; … … 200 201 } 201 202 202 gau_sample.ins( 0,randn());203 gau_sample.ins(gau_sample.size(),randn()); 203 204 204 205 sample_prediction.ins(0,gau_sample*samples.get_col(t)); … … 262 263 263 264 264 int main ( int argc, char* argv[] ) { 265 266 /* 267 DWORD Id; 268 HANDLE hThrd = CreateThread( NULL, 0, (LPTHREAD_START_ROUTINE)ThrdFunc, (LPVOID)1, 0, &Id); 269 270 if ( !hThrd ) 271 { 272 cout<<"Error Creating Threads,,,,.exiting"<<endl; 273 return -1; 274 } 275 Sleep ( 100 ); 276 277 278 char szApp[] = "MT4"; 279 char szTopic[] = "QUOTE"; 280 char szItem1[] = "EURUSD"; 281 282 //DDE Initialization 283 DWORD idInst=0; 284 UINT iReturn; 285 iReturn = DdeInitialize(&idInst, (PFNCALLBACK)DdeCallback, 286 APPCLASS_STANDARD | APPCMD_CLIENTONLY, 0 ); 287 if (iReturn!=DMLERR_NO_ERROR) 288 { 289 printf("DDE Initialization Failed: 0x%04x\n", iReturn); 290 Sleep(1500); 291 return 0; 292 } 293 294 //DDE Connect to Server using given AppName and topic. 295 HSZ hszApp, hszTopic; 296 HCONV hConv; 297 hszApp = DdeCreateStringHandle(idInst, szApp, 0); 298 hszTopic = DdeCreateStringHandle(idInst, szTopic, 0); 299 hConv = DdeConnect(idInst, hszApp, hszTopic, NULL); 300 //DdeFreeStringHandle(idInst, hszApp); 301 //DdeFreeStringHandle(idInst, hszTopic); 302 if (hConv == NULL) 303 { 304 printf("DDE Connection Failed.\n"); 305 Sleep(1500); DdeUninitialize(idInst); 306 return 0; 307 } 308 309 //Execute commands/requests specific to the DDE Server. 310 311 DDERequest(idInst, hConv, szItem1); 312 313 while(1) 314 { 315 MSG msg; 316 BOOL MsgReturn = GetMessage ( &msg , NULL , 0 , 0 ); 317 318 if(MsgReturn) 319 { 320 TranslateMessage(&msg); 321 DispatchMessage(&msg); 322 } 323 } 324 325 //DDE Disconnect and Uninitialize. 326 DdeDisconnect(hConv); 327 DdeUninitialize(idInst); 328 */ 329 330 331 332 /* 333 // EXPERIMENT: 100 AR model generated time series of length of 30 from y_t=0.95*y_(t-1)+0.05*y_(t-2)+0.2*e_t, 334 // where e_t is normally, student(4) and cauchy distributed are tested using robust AR model, to obtain the 335 // variance of location parameter estimators and compare it to the classical setup. 336 vector<vector<vector<string>>> string_lists; 337 string_lists.push_back(vector<vector<string>>()); 338 string_lists.push_back(vector<vector<string>>()); 339 string_lists.push_back(vector<vector<string>>()); 340 341 char* file_strings[3] = {"c:\\ar_normal.txt", "c:\\ar_student.txt", "c:\\ar_cauchy.txt"}; 342 343 344 for(int i = 0;i<3;i++) 345 { 346 ifstream myfile(file_strings[i]); 347 if (myfile.is_open()) 348 { 349 while ( myfile.good() ) 350 { 351 string line; 352 getline(myfile,line); 353 354 vector<string> parsed_line; 355 while(line.find(',') != string::npos) 356 { 357 int loc = line.find(','); 358 parsed_line.push_back(line.substr(0,loc)); 359 line.erase(0,loc+1); 360 } 361 362 string_lists[i].push_back(parsed_line); 363 } 364 myfile.close(); 365 } 366 } 367 368 for(int j = 0;j<string_lists.size();j++) 369 { 370 371 for(int i = 0;i<string_lists[j].size()-1;i++) 372 { 373 vector<vec> conditions; 374 //emlig* emliga = new emlig(2); 375 RARX* my_rarx = new RARX(2,30); 376 377 for(int k = 1;k<string_lists[j][i].size();k++) 378 { 379 vec condition; 380 //condition.ins(0,1); 381 condition.ins(0,string_lists[j][i][k]); 382 conditions.push_back(condition); 383 384 //cout << "orig:" << condition << endl; 385 386 if(conditions.size()>1) 387 { 388 conditions[k-2].ins(0,string_lists[j][i][k]); 389 390 } 391 392 if(conditions.size()>2) 393 { 394 conditions[k-3].ins(0,string_lists[j][i][k]); 395 396 //cout << "modi:" << conditions[k-3] << endl; 397 398 my_rarx->bayes(conditions[k-3]); 399 400 401 //if(k>5) 402 //{ 403 // cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl; 404 //} 405 406 } 407 408 } 409 410 //emliga->step_me(0); 411 /* 412 ofstream myfile; 413 myfile.open("c:\\robust_ar1.txt",ios::app); 414 myfile << my_rarx->minimal_vertex->get_coordinates()[0] << ";"; 415 myfile.close(); 416 417 myfile.open("c:\\robust_ar2.txt",ios::app); 418 myfile << emliga->minimal_vertex->get_coordinates()[1] << ";"; 419 myfile.close(); 420 421 422 cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl; 423 cout << "Step: " << i << endl; 424 } 425 426 cout << "One experiment finished." << endl; 427 428 ofstream myfile; 429 myfile.open("c:\\robust_ar1.txt",ios::app); 430 myfile << endl; 431 myfile.close(); 432 433 myfile.open("c:\\robust_ar2.txt",ios::app); 434 myfile << endl; 435 myfile.close(); 436 }*/ 265 int main ( int argc, char* argv[] ) 266 { 437 267 438 268 // EXPERIMENT: A moving window estimation and prediction of RARX is tested on data generated from … … 444 274 vector<vector<string>> strings; 445 275 446 char* file_string = "c:\\ar_normal_single"; // "c:\\dataTYClosePercDiff"; // 276 char* file_string = "c:\\ar_normal_single"; // "c:\\dataTYClosePercDiff"; // 447 277 448 278 char dfstring[80]; … … 461 291 while(line.find(',') != string::npos) //zmenil som ciarku za medzeru 462 292 { 463 line.erase(0,1); // toto som sem pridal293 //line.erase(0,1); // toto som sem pridal 464 294 int loc2 = line.find('\n'); 465 295 int loc = line.find(','); … … 484 314 for(set<set<pair<int,int>>>::iterator model_type = model_types.begin();model_type!=model_types.end();model_type++) 485 315 {// prechadza rozne typy kanalov, a poctu regresorov 486 for(int window_size = 15;window_size < 16;window_size++)487 { 488 models.push_back(new model((*model_type),true,true,window_size,0,&data_matrix)); // to su len konstruktory, len inicializujeme rozne typy489 models.push_back(new model((*model_type),false,true,window_size,0,&data_matrix));316 for(int window_size = 50;window_size < 51;window_size++) 317 { 318 //models.push_back(new model((*model_type),true,true,window_size,0,&data_matrix)); // to su len konstruktory, len inicializujeme rozne typy 319 //models.push_back(new model((*model_type),false,true,window_size,0,&data_matrix)); 490 320 models.push_back(new model((*model_type),true,false,window_size,0,&data_matrix)); 491 models.push_back(new model((*model_type),false,false,window_size,0,&data_matrix));321 //models.push_back(new model((*model_type),false,false,window_size,0,&data_matrix)); 492 322 } 493 323 … … 498 328 mat result_lognc; 499 329 // mat result_preds; 330 331 ofstream myfilew; 332 char fstring[80]; 333 strcpy(fstring,file_string); 334 strcat(fstring,"lognc.txt"); 335 //strcat(fstring,"preds.txt"); 500 336 501 337 for(int time = max_model_order;time<data_matrix.cols();time++) //time<data_matrix.cols() … … 513 349 if((*model_ref)->my_rarx!=NULL) //vklada normalizacnz faktor do cur_res_lognc 514 350 { 515 cur_res_lognc.ins(cur_res_lognc.size(),(*model_ref)->my_rarx->posterior->log_nc); 351 cur_res_lognc.ins(cur_res_lognc.size(),(*model_ref)->my_rarx->posterior->log_nc); 516 352 } 517 353 else … … 520 356 } 521 357 522 // pair<vec,vec> predictions = (*model_ref)->predict(200,time,&LapRNG); 523 524 // preds.ins(preds.size(),(predictions.first*predictions.second)/(predictions.first*ones(predictions.first.size()))); 525 // preds.ins(0,data_matrix.get(0,time+1)); 526 527 } 358 pair<vec,vec> predictions = (*model_ref)->predict(20,time,&LapRNG); 359 360 cout << predictions.first << endl << predictions.second << endl; 361 362 double avg_prediction = (predictions.first*predictions.second)/(predictions.first*ones(predictions.first.size())); 363 364 (*model_ref)->predictions.ins((*model_ref)->predictions.size(),avg_prediction); 365 366 /* 367 myfilew.open(fstring,ios::app); 368 myfilew << avg_prediction << ","; 369 myfilew.close(); 370 */ 371 372 //preds.ins(0,data_matrix.get(0,time+1)); 373 } 374 375 /* 376 myfilew.open(fstring,ios::app); 377 myfilew << data_matrix.get(0,time+1) << endl; 378 myfilew.close(); 379 */ 528 380 529 381 result_lognc.ins_col(result_lognc.cols(),cur_res_lognc); … … 531 383 532 384 // cout << "Updated." << endl; 533 534 /* 535 vector<vec> conditions; 536 //emlig* emliga = new emlig(2); 537 RARX* my_rarx = new RARX(2,10,false); 538 539 540 mat V0 = 0.0001 * eye ( 3 ); 541 ARX* my_arx = new ARX(0.85); 542 my_arx->set_statistics ( 1, V0 ); //nu is default (set to have finite moments) 543 my_arx->set_constant ( false ); 544 my_arx->validate(); 545 546 547 for(int k = 1;k<strings[j].size();k++) 548 { 549 vec condition; 550 //condition.ins(0,1); 551 condition.ins(0,strings[j][k]); 552 conditions.push_back(condition); 553 554 //cout << "orig:" << condition << endl; 555 556 if(conditions.size()>1) 557 { 558 conditions[k-2].ins(0,strings[j][k]); 385 386 387 388 389 390 myfilew.open(fstring,ios::app); 391 392 // myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0]; 393 394 if(time == max_model_order) 395 { 396 for(int i = 0;i<cur_res_lognc.size();i++) 397 { 398 for(set<pair<int,int>>::iterator ar_ref = models[i]->ar_components.begin();ar_ref != models[i]->ar_components.end();ar_ref++) 399 { 400 myfilew << (*ar_ref).second << (*ar_ref).first; 401 } 402 403 myfilew << "."; 404 405 if(models[i]->my_arx == NULL) 406 { 407 myfilew << "1"; 408 } 409 else 410 { 411 myfilew << "0"; 412 } 559 413 560 } 561 562 if(conditions.size()>2) 563 { 564 conditions[k-3].ins(0,strings[j][k]); 565 566 // cout << "Condition:" << conditions[k-3] << endl; 567 568 my_rarx->bayes(conditions[k-3]); 569 //my_rarx->posterior->step_me(1); 570 571 vec cond_vec; 572 cond_vec.ins(0,conditions[k-3][0]); 573 574 my_arx->bayes(cond_vec,conditions[k-3].right(2)); 575 576 /* 577 if(k>8) 578 { 579 //my_rarx->posterior->step_me(0); 580 581 //mat samples = my_rarx->posterior->sample_mat(10); 582 583 pair<vec,mat> imp_samples = my_rarx->posterior->importance_sample(1000); 584 585 //cout << imp_samples.first << endl; 586 587 vec sample_prediction; 588 vec averaged_params = zeros(imp_samples.second.rows()); 589 for(int t = 0;t<imp_samples.first.size();t++) 590 { 591 vec lap_sample = conditions[k-3].left(2); 592 //lap_sample.ins(lap_sample.size(),1.0); 593 594 lap_sample.ins(0,LapRNG()); 595 596 sample_prediction.ins(0,lap_sample*imp_samples.second.get_col(t)); 597 598 averaged_params += imp_samples.first[t]*imp_samples.second.get_col(t); 599 } 600 601 averaged_params = averaged_params*(1/(imp_samples.first*ones(imp_samples.first.size()))); 602 603 // cout << "Averaged estimated parameters: " << averaged_params << endl; 604 605 vec sample_pow = sample_prediction; 606 607 // cout << sample_prediction << endl; 608 vec poly_coefs; 609 double prediction; 610 bool stop_iteration = false; 611 int en = 1; 612 do 613 { 614 double poly_coef = imp_samples.first*sample_pow/(imp_samples.first*ones(imp_samples.first.size())); 615 616 if(en==1) 617 { 618 prediction = poly_coef; 619 } 620 621 poly_coef = poly_coef*en*fact(utility_constant-2+en)/fact(utility_constant-2); 622 623 if(abs(poly_coef)>numeric_limits<double>::epsilon()) 624 { 625 sample_pow = elem_mult(sample_pow,sample_prediction); 626 poly_coefs.ins(0,pow(-1.0,en+1)*poly_coef); 627 } 628 else 629 { 630 stop_iteration = true; 631 } 632 633 en++; 634 635 if(en>20) 636 { 637 stop_iteration = true; 638 } 639 } 640 while(!stop_iteration); 641 642 /* 643 ofstream myfile_coef; 644 645 myfile_coef.open("c:\\coefs.txt",ios::app); 646 647 for(int t = 0;t<poly_coefs.size();t++) 648 { 649 myfile_coef << poly_coefs[t] << ","; 650 } 651 652 myfile_coef << endl; 653 myfile_coef.close(); 654 */ 655 656 //cout << "Coefficients: " << poly_coefs << endl; 657 658 /* 659 vec bas_coef = vec("1.0 2.0 -8.0"); 660 cout << "Coefs: " << bas_coef << endl; 661 cvec actions2 = roots(bas_coef); 662 cout << "Roots: " << actions2 << endl; 663 */ 664 665 /* 666 667 cvec actions = roots(poly_coefs); 668 669 670 bool is_max = false; 671 for(int t = 0;t<actions.size();t++) 672 { 673 if(actions[t].imag() == 0) 674 { 675 double second_derivative = 0; 676 for(int q = 1;q<poly_coefs.size();q++) 677 { 678 second_derivative+=poly_coefs[q]*pow(actions[t].real(),q-1)*q; 679 } 680 681 if(second_derivative<0) 682 { 683 cout << "Action:" << actions[t].real() << endl; 684 685 is_max = true; 686 } 687 } 688 } 689 690 if(!is_max) 691 { 692 cout << "No maximum." << endl; 693 } 694 695 // cout << "MaxLik coords:" << my_rarx->posterior->minimal_vertex->get_coordinates() << endl; 696 697 /* 698 double prediction = 0; 699 for(int s = 1;s<samples.rows();s++) 700 { 701 702 double avg_parameter = imp_samples.get_row(s)*ones(samples.cols())/samples.cols(); 703 704 prediction += avg_parameter*conditions[k-3][s-1]; 705 706 707 708 /* 709 ofstream myfile; 710 char fstring[80]; 711 strcpy(fstring,file_strings[j]); 712 713 char es[5]; 714 strcat(fstring,itoa(s,es,10)); 715 716 strcat(fstring,"_res.txt"); 717 718 719 myfile.open(fstring,ios::app); 720 721 //myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0]; 722 myfile << avg_parameter; 723 724 if(k!=strings[j].size()-1) 725 { 726 myfile << ","; 727 } 728 else 729 { 730 myfile << endl; 731 } 732 myfile.close(); 733 */ 734 735 736 //} 737 738 // cout << "Prediction: "<< prediction << endl; 739 /* 740 enorm<ldmat>* pred_mat = my_arx->epredictor(conditions[k-3].left(2)); 741 double prediction2 = pred_mat->mean()[0]; 742 */ 743 744 745 ofstream myfile; 746 char fstring[80]; 747 strcpy(fstring,file_string); 748 749 strcat(fstring,"lognc.txt"); 750 751 myfile.open(fstring,ios::app); 752 753 // myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0]; 754 755 if(time == max_model_order) 756 { 757 for(int i = 0;i<cur_res_lognc.size();i++) 758 { 759 for(set<pair<int,int>>::iterator ar_ref = models[i]->ar_components.begin();ar_ref != models[i]->ar_components.end();ar_ref++) 760 { 761 myfile << (*ar_ref).second << (*ar_ref).first; 762 } 763 764 myfile << "."; 765 766 if(models[i]->my_arx == NULL) 767 { 768 myfile << "1"; 769 } 770 else 771 { 772 myfile << "0"; 773 } 774 775 if(models[i]->has_constant) 776 { 777 myfile << "1"; 778 } 779 else 780 { 781 myfile << "0"; 782 } 783 784 myfile << ","; 785 } 786 787 myfile << endl; 788 } 789 790 for(int i = 0;i<cur_res_lognc.size();i++) 791 { 792 myfile << cur_res_lognc[i] << ' ';//zmenil som ciarku ze medzeru 793 } 794 795 myfile << endl; 796 797 myfile.close(); 798 } 799 /* 800 myfile.open(f2string,ios::app); 801 myfile << prediction2; 802 803 if(k!=strings[j].size()-1) 804 { 805 myfile << ","; 806 } 807 else 808 { 809 myfile << endl; 810 } 811 myfile.close(); 812 //*//* 813 814 } 815 } */ 816 817 //emliga->step_me(0); 818 /* 819 ofstream myfile; 820 myfile.open("c:\\robust_ar1.txt",ios::app); 821 myfile << my_rarx->minimal_vertex->get_coordinates()[0] << ";"; 822 myfile.close(); 823 824 myfile.open("c:\\robust_ar2.txt",ios::app); 825 myfile << emliga->minimal_vertex->get_coordinates()[1] << ";"; 826 myfile.close(); 827 828 829 cout << "MaxLik coords:" << emliga->minimal_vertex->get_coordinates() << endl; 830 cout << "Step: " << i << endl;*/ 831 //} 832 833 834 //} 414 if(models[i]->has_constant) 415 { 416 myfilew << "1"; 417 } 418 else 419 { 420 myfilew << "0"; 421 } 422 423 myfilew << ","; 424 } 425 426 myfilew << endl; 427 } 428 429 for(int i = 0;i<cur_res_lognc.size();i++) 430 { 431 myfilew << cur_res_lognc[i] << ' ';//zmenil som ciarku ze medzeru 432 } 433 434 myfilew << endl; 435 436 myfilew.close(); 437 438 } 439 835 440 836 441 -
applications/robust/robustlib.cpp
r1379 r1383 143 143 extended_coords.ins(0,-1); 144 144 145 cout << "Ext. coords:" << extended_coords << "Condition sum:" << as_toprow->condition_sum << endl; 146 145 147 double a = extended_coords*as_toprow->condition_sum; 146 148 if(a<current_emlig->min_ll) … … 152 154 //cout << "a0:" << a_0 << " a0 coords:" << base_vertex->get_coordinates() << " am:" << a_m << " am coords:" << (*vert_ref)->get_coordinates() << endl; 153 155 154 //cout << "Absolute coords:(V" << row_count << ")" << (*vert_ref)->get_coordinates() << endl;156 cout << "Absolute coords:(V" << row_count << ")" << (*vert_ref)->get_coordinates() << endl; 155 157 //cout << "Relative coords:(V" << row_count << ")" << relative_coords << endl; 156 158 … … 209 211 210 212 if(k!=0) 211 { 213 { 214 double value = 0; 212 215 ivec control_vec = ivec(); 213 216 control_vec.ins(0,my_emlig->number_of_parameters-a_order+1); … … 223 226 } 224 227 225 double value = bracket_factor*factor_multiplier*k_multiplier; 226 227 simplex->insert_gamma(k,value,(*a_ref).first); 228 gamma_facs[k] += value; 228 value += bracket_factor*factor_multiplier*k_multiplier; 229 229 } 230 231 simplex->insert_gamma(k,value,(*a_ref).first); 232 gamma_facs[k] += value; 230 233 } 231 234 else … … 234 237 235 238 simplex->insert_gamma(0,value,(*a_ref).first); 236 gamma_facs[ k] += value;239 gamma_facs[0] += value; 237 240 } 238 241 } … … 254 257 255 258 simplex->probability = int_value; 256 //cout << "Probability:" << int_value << endl;259 cout << "Probability:" << int_value << endl; 257 260 return int_value; 258 261 } -
applications/robust/robustlib.h
r1379 r1383 1077 1077 this->number_of_parameters = number_of_parameters; 1078 1078 1079 condition_order = number_of_parameters+ 3;1079 condition_order = number_of_parameters+2; 1080 1080 1081 1081 create_statistic(number_of_parameters, soft_prior_parameter); … … 2056 2056 //sigma = 1/(*gamma)(); 2057 2057 2058 GamRNG.setup(conditions.size()-number_of_parameters ,(*g_ref).second);2058 GamRNG.setup(conditions.size()-number_of_parameters+3,(*g_ref).second); 2059 2059 2060 2060 sigma = 1/GamRNG(); … … 2463 2463 2464 2464 double exponent = extended_coords*condition_and_simplex.first; 2465 double sample_prob = 1/condition_and_simplex.second->probability/pow(probability_and_sigma.second,(int)conditions.size()-number_of_parameters )*exp((-1)/probability_and_sigma.second*exponent);2465 double sample_prob = 1/condition_and_simplex.second->probability/pow(probability_and_sigma.second,(int)conditions.size()-number_of_parameters+3)*exp((-1)/probability_and_sigma.second*exponent); 2466 2466 sample_prob *= probability_and_sigma.first; 2467 2467 … … 2820 2820 this->has_constant = has_constant; 2821 2821 2822 posterior = new emlig(number_of_parameters,0.1 );2822 posterior = new emlig(number_of_parameters,0.141); 2823 2823 2824 2824 this->window_size = window_size; … … 2845 2845 //posterior->step_me(0); 2846 2846 2847 cout << "Current condition:" << yt << endl; 2848 2847 2849 /// \TODO tohle je spatne, tady musi byt jiny vypocet poctu podminek, kdyby nejaka byla multiplicitni, tak tohle bude spatne 2848 2850 if(conditions.size()>window_size && window_size!=0) 2849 { 2851 { 2850 2852 posterior->add_and_remove_condition(yt,conditions.front()); 2851 2853 conditions.pop_front(); -
library/bdm/base/itpp/CMakeLists.txt
r812 r1383 13 13 SET(itpp_math base/math/elem_math.cpp base/math/error.cpp base/math/integration.cpp base/math/misc.cpp) 14 14 15 SET(itpp_signal signal/fastica.cpp signal/filter.cpp signal/filter_design.cpp signal/freq_filt.cpp signal/poly.cpp signal/resampling.cpp signal/source.cpp signal/transforms.cpp signal/window.cpp) 16 17 SET(itpp_stat stat/misc_stat.cpp stat/mog_diag.cpp stat/mog_diag_em.cpp stat/mog_diag_kmeans.cpp stat/mog_generic.cpp) 18 15 19 # add BDMLIB compile flag 16 20 #ADD_DEFINITIONS(-DNDEBUG) 17 21 18 22 # Normal BDM library 19 add_library (itpp STATIC ${itpp_base} ${itpp_algebra} ${itpp_bessel} ${itpp_math} )23 add_library (itpp STATIC ${itpp_base} ${itpp_algebra} ${itpp_bessel} ${itpp_math} ${itpp_signal} ${itpp_stat}) -
library/bdm/estim/arx.h
r1359 r1383 311 311 312 312 void bayes ( const vec &val, const vec &cond ) { 313 // fill window 314 Y.append_col(val); 315 Cond.append_col(cond); 316 if (Y.cols()>win_length){ 317 // shift the buffer 318 Y=Y.get_cols(1,Y.cols()-1); 319 Cond=Cond.get_cols(1,Cond.cols()-1); 313 314 if(cond.size()>0) 315 { 316 // fill window 317 Y.append_col(val); 318 Cond.append_col(cond); 319 if (Y.cols()>win_length){ 320 // shift the buffer 321 Y=Y.get_cols(1,Y.cols()-1); 322 Cond=Cond.get_cols(1,Cond.cols()-1); 323 } 324 325 est._V()=V0; 326 est._nu()=nu0; 327 for ( int t = 0; t < Y.cols(); t++ ) { 328 ARX::bayes ( Y.get_col ( t ), Cond.get_col ( t ) ); 329 } 330 } 331 else 332 { 333 Y.append_col(val); 334 335 if (Y.cols()>win_length){ 336 // shift the buffer 337 Y=Y.get_cols(1,Y.cols()-1); 338 } 339 340 est._V()=V0; 341 est._nu()=nu0; 342 343 for ( int t = 0; t < Y.cols(); t++ ) { 344 ARX::bayes (Y.get_col ( t )); 345 } 320 346 } 321 347 322 est._V()=V0; 323 est._nu()=nu0; 324 for ( int t = 0; t < Y.cols(); t++ ) { 325 ARX::bayes ( Y.get_col ( t ), Cond.get_col ( t ) ); 326 } 327 328 } 348 } 349 329 350 void from_setting(const Setting &set){ 330 351 ARX::from_setting(set); -
library/bdm/stat/exp_family.h
r1196 r1383 28 28 //! Global Gamma_RNG 29 29 extern Gamma_RNG GamRNG; 30 30 31 31 32 /*! -
library/tests/tutorial/arx_simple.cpp
r1323 r1383 1 1 #include "estim/arx.h" 2 #include <vector> 3 #include <iostream> 4 #include <fstream> 2 5 using namespace bdm; 3 6 4 7 // estimation of AR(0) model 5 int main() { 8 int main(){ 9 //data 10 vector<vector<vector<string>>> string_lists; 11 string_lists.push_back(vector<vector<string>>()); 12 string_lists.push_back(vector<vector<string>>()); 13 string_lists.push_back(vector<vector<string>>()); 14 15 char* file_strings[3] = {"c:\\ar_normal.txt", "c:\\ar_student.txt", "c:\\ar_cauchy.txt"}; 16 17 18 for(int i = 0;i<3;i++) 19 { 20 ifstream myfile(file_strings[i]); 21 if (myfile.is_open()) 22 { 23 while ( myfile.good() ) 24 { 25 string line; 26 getline(myfile,line); 27 28 vector<string> parsed_line; 29 while(line.find(',') != string::npos) 30 { 31 int loc = line.find(','); 32 parsed_line.push_back(line.substr(0,loc)); 33 line.erase(0,loc+1); 34 } 35 36 string_lists[i].push_back(parsed_line); 37 } 38 myfile.close(); 39 } 40 } 41 6 42 //prior 7 mat V0 = 0.00001 * eye ( 3 ); 8 V0 ( 0, 0 ) = 0.1; // 9 ARX Ar; 10 Ar.set_statistics ( 1, V0 ); //nu is default (set to have finite moments) 11 Ar.set_constant ( false ); 12 Ar.validate(); 13 // forgetting is default: 1.0 14 Ar.bayes_batch ( randn ( 1, 3400 ), randn(2,3400) ); 43 mat V0 = 0.0001 * eye ( 3 ); 44 //V0 ( 0, 0 ) = 0.1; // 15 45 16 cout << "Expected value of Theta is: " << Ar.posterior().mean() << endl; 17 cout << "NC of posterior: " << Ar.posterior().lognc() << endl; 46 for(int j = 0;j<string_lists.size();j++) 47 { 48 49 for(int i = 0;i<string_lists[j].size();i++) 50 { 51 ARX Ar; 52 Ar.set_statistics ( 1, V0 ); //nu is default (set to have finite moments) 53 Ar.set_constant ( false ); 54 Ar.validate(); 55 // forgetting is default: 1.0 56 57 vector<vec> conditions; 58 for(int k = 1;k<string_lists[j][i].size();k++) 59 { 60 vec condition; 61 condition.ins(0,string_lists[j][i][k]); 62 conditions.push_back(condition); 63 64 if(conditions.size()>1) 65 { 66 conditions[k-2].ins(0,string_lists[j][i][k]); 67 68 } 69 70 if(conditions.size()>2) 71 { 72 conditions[k-3].ins(0,string_lists[j][i][k]); 73 74 //cout << conditions[k-3] << endl;// << conditions[k-3].left(1) << conditions[k-3].right(2); 75 76 Ar.bayes(conditions[k-3].left(1),conditions[k-3].right(2)); 77 } 78 } 79 80 ofstream myfile; 81 myfile.open("c:\\classic_ar1.txt",ios::app); 82 myfile << Ar.posterior().mean()[0] << ";"; 83 myfile.close(); 84 85 myfile.open("c:\\classic_ar2.txt",ios::app); 86 myfile << Ar.posterior().mean()[1] << ";"; 87 myfile.close(); 88 89 90 91 } 92 93 ofstream myfile; 94 myfile.open("c:\\classic_ar1.txt",ios::app); 95 myfile << endl; 96 myfile.close(); 97 98 myfile.open("c:\\classic_ar2.txt",ios::app); 99 myfile << endl; 100 myfile.close(); 101 } 18 102 } 103
