Changeset 1412 for applications/robust

Timestamp:

12/03/11 20:12:14 (14 years ago)

Author:

sindj

Message:

Main upravena pro experimenty na burzovnich datech. JS

Files:

• applications/robust/main.cpp (modified) (11 diffs)

applications/robust/main.cpp

@@ -16 +16 @@
 #include "ddeml.h"
 #include "stdio.h"
-#include <itpp/itoptim.h>
+#include <dirent.h>
+//#include <itpp/itoptim.h>
 
 //#include "DDEClient.h"
@@ -28 +29 @@
 //const int utility_constant = 5;
 
-const int max_model_order = 2;
-const double apriorno     = 0.01;
-const int max_window_size = 40;
-const int utility_order   = 19;
-const int prediction_time = 30;
+const int max_model_order         = 2;
+const double apriorno             = 0.01;
+const int max_window_size         = 40;
+const int utility_order           = 25;
+const int prediction_time         = 30;
+const double min_utility_argument = 0.001;
+const double max_investment       = 10.0;
+const char* commodity                 = "CL\\";
 
 /*
@@ -69 +73 @@
 }
 */
+vector<char*> listFiles(char* dir){
+        vector<char*> files;
+                DIR *pDIR;
+        struct dirent *entry;
+        if( pDIR=opendir(dir) ){
+                while(entry = readdir(pDIR)){
+                        if( strcmp(entry->d_name, ".") != 0 && strcmp(entry->d_name, "..") != 0 )
+                                                {                                                       
+                                                        char *file_string = new char[strlen(entry->d_name) + 1];
+                                                        strcpy(file_string, entry->d_name);                                                                                             
+                                                        files.push_back(file_string);
+                                                }
+                }
+                closedir(pDIR);
+        }
+                return files;
+}
 
 double valueCRRAUtility(const double &position, const vec &samples, const int order)
@@ -77 +98 @@
         {
                 double sample = samples.get(i);
-                value += sample/pow(position*sample+1,order+1);
+                if((position*sample+1)>min_utility_argument)
+                {
+                        value += sample/pow(position*sample+1,order+1);
+                }
+                else
+                {
+                        value += (min_utility_argument-1)/position/pow(min_utility_argument,order+1);
+                }
         }       
         
@@ -90 +118 @@
         {
                 double sample = samples.get(i);
-                value += (-(order+1)*pow(sample,2))/pow(position*sample+1,order+2);
+                if((position*sample+1)>min_utility_argument)
+                {
+                        value += (-(order+1)*pow(sample,2))/pow(position*sample+1,order+2);
+                }
         }
 
@@ -129 +160 @@
         }
 }
-        
-        
-
-        
-
-
 
 class model
@@ -242 +267 @@
         }
 
+        ~model()
+        {
+                if(my_rarx!=NULL)
+                {
+                        delete my_rarx;
+                }
+                else
+                {
+                        delete my_arx;
+                }
+        }
+
         void data_update(int time) 
         {
@@ -309 +346 @@
                                 gau_sample.ins(gau_sample.size(),randn());
 
-                                sample_prediction.ins(0,gau_sample*samples.get_col(t));                         
+                                vec param_sample = samples.get_col(t);
+                                param_sample.set(param_sample.size()-1,sqrt(param_sample[param_sample.size()-1]));
+                                sample_prediction.ins(0,gau_sample*param_sample);                               
                         }
 
@@ -365 +404 @@
 };
 
-
-
-
+// ****************************************************
+//                           MAIN MAIN MAIN MAIN MAIN
+// ****************************************************
 int main ( int argc, char* argv[] ) 
 {
@@ -375 +414 @@
     // can be compared to the classical setup.
         
-        itpp::Laplace_RNG LapRNG = Laplace_RNG();       
-        
-        vector<vector<string>> strings;
-
-        char* file_string =  "C:\\CD2"; // "C:\\dataADClosePercDiff"; //  
-
-        char dfstring[80];
-        strcpy(dfstring,file_string);
-        strcat(dfstring,".txt");
-        
-        
-        mat data_matrix;
-        ifstream myfile(dfstring);
-        if (myfile.is_open())
-        {               
-                string line;
-                while(getline(myfile,line))
-                {                       
-                        vec data_vector;
-                        while(line.find(',') != string::npos) //zmenil som ciarku za medzeru
-                        {
-                                //line.erase(0,1); // toto som sem pridal
-                                int loc2 = line.find('\n');
-                                int loc  = line.find(',');
-                                data_vector.ins(data_vector.size(),atof(line.substr(0,loc).c_str()));                           
-                                line.erase(0,loc+1);                                    
-                        }
-
-                        data_matrix.ins_row(data_matrix.rows(),data_vector);
-                }               
-
-                myfile.close(); 
-        }
-        else
-        {
-                cout << "Can't open data file!" << endl;
-        }
-        
-        //konec nacitavania dat
-        set<set<pair<int,int>>> model_types = model::possible_models_recurse(max_model_order,data_matrix.rows()); //volanie funkce kde robi kombinace modelov
-                                                                                                //to priradime do model_types, data_matrix.row urcuje pocet kanalov dat
-        vector<model*> models;
-        for(set<set<pair<int,int>>>::iterator model_type = model_types.begin();model_type!=model_types.end();model_type++)
-        {// prechadza rozne typy kanalov, a poctu regresorov
-                for(int window_size = max_window_size-1;window_size < max_window_size;window_size++)
-                {
-                        models.push_back(new model((*model_type),true,true,window_size,0,&data_matrix));   // to su len konstruktory, len inicializujeme rozne typy
-                        models.push_back(new model((*model_type),false,true,window_size,0,&data_matrix));
-                        models.push_back(new model((*model_type),true,false,window_size,0,&data_matrix));
-                        models.push_back(new model((*model_type),false,false,window_size,0,&data_matrix));              
-                }
-
-                //set<pair<int,int>> empty_list;
-                //models.push_back(new model(empty_list,false,true,100,0,&data_matrix));
-        }
-        
-        mat result_lognc;
-        // mat result_preds;
-
-        ofstream myfilew;
-        char fstring[80];                                       
-        strcpy(fstring,file_string);
-        //strcat(fstring,"lognc.txt");
-        strcat(fstring,"preds.txt");    
-
-        for(int time = max_model_order;time<data_matrix.cols();time++) //time<data_matrix.cols() 
-        {       
-                cout << "Steps: " << time << endl;
+        itpp::Laplace_RNG LapRNG = Laplace_RNG();
+
+        char* folder_string  =  "C:\\RobustExperiments\\"; // "C:\\dataADClosePercDiff"; //
+        char* data_folder    = "data\\";
+        char* results_folder = "results\\";     
+
+        char dfstring[150];
+        strcpy(dfstring,folder_string);
+        strcat(dfstring,data_folder);
+        strcat(dfstring,commodity);
+        vector<char*> files = listFiles(dfstring);
+        
+        for(int contract=0;contract<files.size();contract++)
+        {
+                char *cdf_str = new char[strlen(dfstring) + 1];
+                strcpy(cdf_str, dfstring);              
+                strcat(cdf_str,files[contract]);        
+        
+                mat data_matrix;
+                ifstream myfile(cdf_str);
+                if (myfile.is_open())
+                {               
+                        string line;
+                        while(getline(myfile,line))
+                        {                       
+                                vec data_vector;
+                                while(line.find(',') != string::npos) //zmenil som ciarku za medzeru
+                                {
+                                        //line.erase(0,1); // toto som sem pridal
+                                        int loc2 = line.find('\n');
+                                        int loc  = line.find(',');
+                                        data_vector.ins(data_vector.size(),atof(line.substr(0,loc).c_str()));                           
+                                        line.erase(0,loc+1);                                    
+                                }
+
+                                data_matrix.ins_row(data_matrix.rows(),data_vector);
+                        }               
+
+                        myfile.close(); 
+                }
+                else
+                {
+                        cout << "Can't open data file!" << endl;
+                }
+        
+                //konec nacitavania dat
+                set<set<pair<int,int>>> model_types = model::possible_models_recurse(max_model_order,data_matrix.rows()); //volanie funkce kde robi kombinace modelov
+                                                                                                        //to priradime do model_types, data_matrix.row urcuje pocet kanalov dat
+                vector<model*> models;
+                for(set<set<pair<int,int>>>::iterator model_type = model_types.begin();model_type!=model_types.end();model_type++)
+                {// prechadza rozne typy kanalov, a poctu regresorov
+                        for(int window_size = max_window_size-1;window_size < max_window_size;window_size++)
+                        {
+                                models.push_back(new model((*model_type),true,true,window_size,0,&data_matrix));   // to su len konstruktory, len inicializujeme rozne typy
+                                models.push_back(new model((*model_type),false,true,window_size,0,&data_matrix));
+                                models.push_back(new model((*model_type),true,false,window_size,0,&data_matrix));
+                                models.push_back(new model((*model_type),false,false,window_size,0,&data_matrix));              
+                        }
+
+                        //set<pair<int,int>> empty_list;
+                        //models.push_back(new model(empty_list,false,true,100,0,&data_matrix));
+                }
+        
+                mat result_lognc;
+                // mat result_preds;
+
+                ofstream myfilew;
+                char rfstring[150];                                     
+                strcpy(rfstring,folder_string);
+                strcat(rfstring,results_folder);
+                strcat(rfstring,commodity);
+                strcat(rfstring,files[contract]);
+
+                /*
+                char predstring[150];                                   
+                strcpy(predstring,folder_string);
+                strcat(predstring,results_folder);
+                strcat(predstring,commodity);
+                strcat(predstring,"PRED");
+                strcat(predstring,files[contract]);
+                */
+
+                for(int time = max_model_order;time<data_matrix.cols();time++) //time<data_matrix.cols() 
+                {       
+                        cout << "Steps: " << time << endl;
                 
-                /*
-                if(time==max_window_size)
-                {
-                        exit(1);
-                }*/
+                        if(time == max_model_order)
+                        { 
+                                
+                                myfilew.open(rfstring,ios::app);
+                                for(int i = 0;i<models.size();i++)
+                                {
+                                        for(set<pair<int,int>>::iterator ar_ref = models[i]->ar_components.begin();ar_ref != models[i]->ar_components.end();ar_ref++)
+                                        {
+                                                myfilew << (*ar_ref).second << (*ar_ref).first;                                                 
+                                        }
+
+                                        myfilew << ".";
+
+                                        if(models[i]->my_arx == NULL)
+                                        {
+                                                myfilew << "31";
+                                        }
+                                        else
+                                        {
+                                                myfilew << "30";
+                                        }
+                                        
+                                        if(models[i]->has_constant)
+                                        {
+                                                myfilew << "61";
+                                        }
+                                        else
+                                        {
+                                                myfilew << "60";
+                                        }
+
+                                        myfilew << ",999,999,999,";
+                                }
+
+                                myfilew << "888" << endl;
+                                myfilew.close();
+                        }
                 
-                //pocet stlpcov data_matrix je pocet casovych krokov
-                vec cur_res_lognc;
-                // vec preds;
-                vector<string> nazvy;
-                for(vector<model*>::iterator model_ref = models.begin();model_ref!=models.end();model_ref++)
-                {//posuvam s apo models, co je pole modelov urobene o cyklus vyssie. Teda som v case time a robim to tam pre vsetky typy modelov, kombinace regresorov
-                        (*model_ref)->data_update(time); //pozret sa preco je toto tu nutne
-
-                        //cout << "Updated." << endl;
-                        //if (time = max_model_order) nazvy.push_back(models.model_ref]);// ako by som mohol dostat nazov modelu?
-
-                        if((*model_ref)->my_rarx!=NULL) //vklada normalizacni faktor do cur_res_lognc
-                        {
-                                //cout << "Maxlik vertex:" << (*model_ref)->my_rarx->posterior->minimal_vertex->get_coordinates() << endl;
-                                cur_res_lognc.ins(cur_res_lognc.size(),(*model_ref)->my_rarx->posterior->_ll());                                
-                        }
-                        else
-                        {
-                                double cur_lognc = (*model_ref)->my_arx->posterior().lognc();
-                                double cur_ll = cur_lognc-(*model_ref)->previous_lognc;
+                        //pocet stlpcov data_matrix je pocet casovych krokov
+                        double cur_loglikelihood;
+                        // vec preds;
+                        int  prev_samples_nr;
+                        bool previous_switch = true;
+                        for(vector<model*>::iterator model_ref = models.begin();model_ref!=models.end();model_ref++)
+                        {//posuvam s apo models, co je pole modelov urobene o cyklus vyssie. Teda som v case time a robim to tam pre vsetky typy modelov, kombinace regresorov
+                                (*model_ref)->data_update(time); //pozret sa preco je toto tu nutne
+
+                                //cout << "Updated." << endl;                           
+
+                                if((*model_ref)->my_rarx!=NULL) //vklada normalizacni faktor do cur_res_lognc
+                                {
+                                        //cout << "Maxlik vertex(robust):" << (*model_ref)->my_rarx->posterior->minimal_vertex->get_coordinates() << endl;
+                                        cur_loglikelihood = (*model_ref)->my_rarx->posterior->_ll();                            
+                                }
+                                else
+                                {
+                                        //cout << "Mean(classical):" << (*model_ref)->my_arx->posterior().mean() << endl;
+                                        double cur_lognc = (*model_ref)->my_arx->posterior().lognc();
+                                        cur_loglikelihood = cur_lognc-(*model_ref)->previous_lognc;
+                                                                
+                                        (*model_ref)->previous_lognc = cur_lognc;                               
+                                }
                                 
                                 /*
-                                if(time<(*model_ref)->window_size)
-                                {
-                                        cur_ll-=1.83787706640935;
-                                }*/
-                                
-                                (*model_ref)->my_arx->_ll();
-
-                                cur_res_lognc.ins(cur_res_lognc.size(),cur_ll);
-                                
-                                (*model_ref)->previous_lognc = cur_lognc;                               
-                        }
-                        
-                        /*
                         if(time == max_window_size-1)
                         {
@@ -541 +620 @@
                                 myfilew.close();                                
                         }
-                        */              
-
-                        if(time>prediction_time)
-                        {                       
-                                // PREDICTIONS
-                                pair<vec,vec> predictions = (*model_ref)->predict(5000,time,&LapRNG);
-
-                                /*
-                                cout << predictions.first << endl << endl << predictions.second << endl << "*************************************" ;
-                                pause(5);
-                                */      
-
-                                double optimalInvestment = newtonRaphson(0,0.00001,predictions.second,utility_order);
-
-                                /*
-                                vec utilityValues;
-                                for(int j=0;j<1000;j++)
-                                {
-                                        utilityValues.ins(utilityValues.length(),valueCRRAUtility(-0.5+0.001*j, predictions.second, utility_order));
-                                }*/
-
-                                double avg_prediction = (predictions.first*predictions.second)/(predictions.first*ones(predictions.first.size()));
-
-                                (*model_ref)->predictions.ins((*model_ref)->predictions.size(),avg_prediction);
+                        */                              
+
+                                if(time>prediction_time)
+                                {                       
+                                        int samples_nr;                                 
+                                        if(previous_switch)
+                                        {
+                                                samples_nr = 10000;
+                                        }
+                                        else
+                                        {
+                                                samples_nr = prev_samples_nr;
+                                        }
+                                        
+                                        
+                                        // PREDICTIONS
+                                        pair<vec,vec> predictions = (*model_ref)->predict(samples_nr,time,&LapRNG);
+
+                                        /*
+                                        myfilew.open(predstring,ios::app);
+                                        for(int i=0;i<10000;i++)
+                                        {
+                                                if(i<predictions.second.size())
+                                                {
+                                                        myfilew << predictions.second.get(i) << ",";
+                                                }
+                                                else
+                                                {
+                                                        myfilew << ",";
+                                                }
+                                        }
+                                        myfilew << endl;
+                                        myfilew.close();
+                                        */
+
+                                        if(previous_switch)
+                                        {
+                                                prev_samples_nr = predictions.second.size();
+                                                samples_nr = prev_samples_nr;
+                                        }
+
+                                        previous_switch = !previous_switch;
+
+                                        double optimalInvestment = newtonRaphson(0,0.00001,predictions.second,utility_order);
+
+                                        if(abs(optimalInvestment)>max_investment)
+                                        {
+                                                optimalInvestment = max_investment*sign(optimalInvestment);
+                                        }
+
+
+                                        /*
+                                        vec utilityValues;
+                                        for(int j=0;j<1000;j++)
+                                        {
+                                                utilityValues.ins(utilityValues.length(),valueCRRAUtility(-0.5+0.001*j, predictions.second, utility_order));
+                                        }*/
+
+                                        double avg_prediction = (predictions.first*predictions.second)/(predictions.first*ones(predictions.first.size()));
+
+                                        (*model_ref)->predictions.ins((*model_ref)->predictions.size(),avg_prediction);
                                                 
-                                myfilew.open(fstring,ios::app);
-                                
-                                /*
-                                for(int j=0;j<utilityValues.length();j++)
-                                {
-                                        myfilew << utilityValues.get(j) << ",";
-                                }
-                                myfilew << endl;
-                                */
-
-                                myfilew << avg_prediction << "," << optimalInvestment << ","; 
+                                        myfilew.open(rfstring,ios::app);
+                                
+                                        /*
+                                        for(int j=0;j<utilityValues.length();j++)
+                                        {
+                                                myfilew << utilityValues.get(j) << ",";
+                                        }
+                                        myfilew << endl;
+                                        */
+
+                                        myfilew << avg_prediction << "," << optimalInvestment << "," << samples_nr << "," << cur_loglikelihood << ","; 
+                                        myfilew.close();
+                                }                       
+                        }
+
+                        if(time>prediction_time&&(time+1)<data_matrix.cols())
+                        {
+                                // REAL PRICE           
+                                myfilew.open(rfstring,ios::app);
+                                myfilew << data_matrix.get(0,time+1) << endl;
                                 myfilew.close();
-                        }
-                        
-
-                        //preds.ins(0,data_matrix.get(0,time+1));
-                }
-
-                if(time>prediction_time)
-                {
-                        // REAL PRICE           
-                        myfilew.open(fstring,ios::app);
-                        myfilew << data_matrix.get(0,time+1) << endl;
-                        myfilew.close();
-                }
-                
-                
-
-                result_lognc.ins_col(result_lognc.cols(),cur_res_lognc);
-                //result_preds.ins_col(result_preds.cols(),preds);
-
-                // cout << "Updated." << endl;  
-                                                        
-                /*
-                myfilew.open(fstring,ios::app);
-                
-                // myfile << my_rarx->posterior->minimal_vertex->get_coordinates()[0];
-                
-                if(time == max_model_order)
-                { 
-                        for(int i = 0;i<cur_res_lognc.size();i++)
-                        {
-                                for(set<pair<int,int>>::iterator ar_ref = models[i]->ar_components.begin();ar_ref != models[i]->ar_components.end();ar_ref++)
-                                {
-                                        myfilew << (*ar_ref).second << (*ar_ref).first;                                                 
-                                }
-
-                                myfilew << ".";
-
-                                if(models[i]->my_arx == NULL)
-                                {
-                                        myfilew << "1";
-                                }
-                                else
-                                {
-                                        myfilew << "0";
-                                }
-                                        
-                                if(models[i]->has_constant)
-                                {
-                                        myfilew << "1";
-                                }
-                                else
-                                {
-                                        myfilew << "0";
-                                }
-
-                                myfilew << ",";
-                        }
-
-                        myfilew << endl;
-                }
-                
-                
-                for(int i = 0;i<cur_res_lognc.size();i++)
-                {
-                   myfilew << cur_res_lognc[i] << ' ';//zmenil som ciarku ze medzeru
-                }
-
-                myfilew << endl;                
-                myfilew.close();        
-                */
-}
-
-        // EXPERIMENT: One step ahead price prediction. Comparison of classical and robust model using optimal trading
-    //             with maximization of logarithm of one-step ahead wealth.
-
-        
-                
-                /*
-                cout << "One experiment finished." << endl;
-
-                ofstream myfile;
-                myfile.open("c:\\robust_ar1.txt",ios::app);
-                myfile << endl;
-                myfile.close();
-
-                myfile.open("c:\\robust_ar2.txt",ios::app);
-                myfile << endl;
-                myfile.close();*/
-        
-
-        //emlig* emlig1 = new emlig(emlig_size);
-
-        //emlig1->step_me(0);
-        //emlig* emlig2 = new emlig(emlig_size);
-        
-        /*
-        emlig1->set_correction_factors(4);
-
-        for(int j = 0;j<emlig1->correction_factors.size();j++)
-        {
-                for(set<my_ivec>::iterator vec_ref = emlig1->correction_factors[j].begin();vec_ref!=emlig1->correction_factors[j].end();vec_ref++)
-                {
-                        cout << j << "    ";
-                        
-                        for(int i=0;i<(*vec_ref).size();i++)
-                        {
-                                cout << (*vec_ref)[i];
-                        }
-
-                        cout << endl;
-                }
-        }*/
-        
-        /*
-    vec condition5 = "1.0 1.0 1.01";//"-0.3 1.7 1.5";
-
-        emlig1->add_condition(condition5);
-        //emlig1->step_me(0);
-
-
-        vec condition1a = "-1.0 1.02 0.5";
-        //vec condition1b = "1.0 1.0 1.01";
-        emlig1->add_condition(condition1a);
-        //emlig2->add_condition(condition1b);
-
-        vec condition2a = "-0.3 1.7 1.5";
-        //vec condition2b = "-1.0 1.0 1.0";
-        emlig1->add_condition(condition2a);
-        //emlig2->add_condition(condition2b);
-
-        vec condition3a = "0.5 -1.01 1.0";
-        //vec condition3b = "0.5 -1.01 1.0";
-
-        emlig1->add_condition(condition3a);
-        //emlig2->add_condition(condition3b);   
-
-        vec condition4a = "-0.5 -1.0 1.0";
-        //vec condition4b = "-0.5 -1.0 1.0";    
-
-        emlig1->add_condition(condition4a);
-        //cout << "************************************************" << endl;
-        //emlig2->add_condition(condition4b);
-        //cout << "************************************************" << endl;
-        
-        //cout << emlig1->minimal_vertex->get_coordinates();
-        
-        //emlig1->remove_condition(condition3a);
-        //emlig1->step_me(0);
-        //emlig1->remove_condition(condition2a);
-        //emlig1->remove_condition(condition1a);
-        //emlig1->remove_condition(condition5);
-        
-
-        //emlig1->step_me(0);
-        //emlig2->step_me(0);
-        
-
-        // DA SE POUZIT PRO VYPIS DO SOUBORU
-        // emlig1->step_me(0);
-
-        //emlig1->remove_condition(condition1);
-        
-        
-
-        
-        
-        /*
-        for(int i = 0;i<100;i++)
-        {
-                cout << endl << "Step:" << i << endl;           
-
-                double condition[emlig_size+1];         
-
-                for(int k = 0;k<=emlig_size;k++)
-                {
-                        condition[k] = (rand()-RAND_MAX/2)/1000.0;              
-                }
-                        
-
-                vec* condition_vec = new vec(condition,emlig_size+1);
-                emlig1->add_condition(*condition_vec);
-
-                /*
-                for(polyhedron* toprow_ref = emlig1->statistic.rows[emlig_size]; toprow_ref != emlig1->statistic.end_poly; toprow_ref = toprow_ref->next_poly)
-                {
-                        cout << ((toprow*)toprow_ref)->probability << endl;
-                }
-                */
-                /*
-                cout << emlig1->statistic_rowsize(emlig_size) << endl << endl;
-        
-                /*
-                if(i-emlig1->number_of_parameters >= 0)
-                {
-                        pause(30);
-                }
-                */
-
-                // emlig1->step_me(i);
-                
-                /*
-                vector<int> sizevector;
-                for(int s = 0;s<=emlig1->number_of_parameters;s++)
-                {
-                        sizevector.push_back(emlig1->statistic_rowsize(s));
-                }
-                */
-        //}
-    
-
-
-        
-        /*
-        emlig1->step_me(1);
-
-        vec condition = "2.0 0.0 1.0";  
-
-        emlig1->add_condition(condition);
-
-        vector<int> sizevector;
-        for(int s = 0;s<=emlig1->number_of_parameters;s++)
-        {
-                sizevector.push_back(emlig1->statistic_rowsize(s));
-        }
-
-        emlig1->step_me(2);
-
-        condition = "2.0 1.0 0.0";
-
-        emlig1->add_condition(condition);
-
-        sizevector.clear();
-        for(int s = 0;s<=emlig1->number_of_parameters;s++)
-        {
-                sizevector.push_back(emlig1->statistic_rowsize(s));
-        }
-        */
+                        }       
+                }
+
+                for(vector<model*>::reverse_iterator model_ref = models.rbegin();model_ref!=models.rend();model_ref++)
+                {
+                        delete *model_ref;
+                }
+        }
 
         return 0;