Changeset 1331 for applications

Show
Ignore:
Timestamp:
04/14/11 19:24:43 (13 years ago)
Author:
sindj
Message:

Dodelan sampling sigma, zbyva doladit sampling alpha.JS

Location:
applications/robust
Files:
3 modified

Legend:

Unmodified
Added
Removed

  • applications/robust/main.cpp

    r1324 r1331  
    193193                                        //my_rarx->posterior->step_me(0); 
    194194 
    195                                         //my_rarx->posterior->sample_mat(1); 
     195                                        my_rarx->posterior->sample_mat(1); 
    196196 
    197197                                        cout << "MaxLik coords:" << my_rarx->posterior->minimal_vertex->get_coordinates() << endl; 

  • applications/robust/robustlib.cpp

    r1325 r1331  
    9797 
    9898                        emlig* current_emlig; 
    99                         simplex->gamma_parameters.clear(); 
    100                         simplex->gamma_sum = 0; 
     99                        simplex->clear_gammas(); 
    101100 
    102101                        if(this->my_emlig!=NULL) 
    … …  
    257256 
    258257                                double bracket = fact(as1_order-1)/pow(a_0-(*as_ref).first,sigma_order+1); 
    259                                  
    260                                 simplex->gamma_sum += bracket*multiplier; 
    261                                  
    262                                 multimap<double,double> map; 
    263                                 simplex->gamma_parameters.push_back(map); 
    264                                 simplex->gamma_parameters[0].insert(pair<double,double>(bracket*multiplier,a_0-(*as_ref).first)); 
     258                                                                 
     259                                simplex->insert_gamma(0,bracket*multiplier,a_0-(*as_ref).first);                                                                 
    265260 
    266261                                bracket *= fact(sigma_order); 
    … …  
    287282                                        } 
    288283 
    289                                         simplex->gamma_sum += local_bracket*multiplier; 
    290                                          
    291                                         if(simplex->gamma_parameters.size()<=k+1) 
    292                                         { 
    293                                                 multimap<double,double> loc_map; 
    294                                                 simplex->gamma_parameters.push_back(loc_map); 
    295                                         } 
    296  
    297                                         simplex->gamma_parameters[k+1].insert(pair<double,double>(local_bracket*multiplier,a_0-(*as_ref).first)); 
    298  
     284                                        simplex->insert_gamma(k+1,local_bracket*multiplier,a_0-(*as_ref).first); 
     285                                         
    299286                                        bracket += local_bracket*fact(sigma_order-k); 
    300287                                }                                
    … …  
    306293                        double correction_term_base = correction_term/pow(a_0,sigma_order+1); 
    307294 
    308                         simplex->gamma_sum += correction_term_base; 
    309                         simplex->gamma_parameters[0].insert(pair<double,double>(correction_term_base,a_0)); 
     295                        simplex->insert_gamma(0,correction_term_base,a_0);                       
    310296 
    311297                        correction_term = fact(sigma_order)*correction_term_base; 

  • applications/robust/robustlib.h

    r1325 r1331  
    6767class simplex 
    6868{ 
     69         
     70 
    6971public: 
    7072 
    … …  
    7375        double probability; 
    7476 
    75         vector<multimap<double,double>> gamma_parameters; 
    76  
    77         double gamma_sum; 
     77        vector<multimap<double,double>> positive_gamma_parameters; 
     78 
     79        vector<multimap<double,double>> negative_gamma_parameters; 
     80 
     81        double positive_gamma_sum; 
     82 
     83        double negative_gamma_sum; 
    7884         
    7985 
    … …  
    8894                this->vertices.insert(vertex); 
    8995                probability = 0; 
     96        } 
     97 
     98        void clear_gammas() 
     99        { 
     100                positive_gamma_parameters.clear(); 
     101                negative_gamma_parameters.clear(); 
     102                 
     103                 
     104                positive_gamma_sum = 0; 
     105                negative_gamma_sum = 0; 
     106        } 
     107 
     108        void insert_gamma(int order, double weight, double beta) 
     109        { 
     110                if(weight>=0) 
     111                { 
     112                        if(positive_gamma_parameters.size()<order+1) 
     113                        { 
     114                                multimap<double,double> map; 
     115                                positive_gamma_parameters.push_back(map); 
     116                        } 
     117 
     118                        positive_gamma_sum += weight; 
     119 
     120                        positive_gamma_parameters[order].insert(pair<double,double>(weight,beta));               
     121                } 
     122                else 
     123                { 
     124                        if(negative_gamma_parameters.size()<order+1) 
     125                        { 
     126                                multimap<double,double> map; 
     127                                negative_gamma_parameters.push_back(map); 
     128                        } 
     129 
     130                        negative_gamma_sum -= weight; 
     131 
     132                        negative_gamma_parameters[order].insert(pair<double,double>(-weight,beta)); 
     133                } 
    90134        } 
    91135}; 
    … …  
    17421786 
    17431787        mat sample_mat(int n) 
    1744         { 
    1745                  
    1746                  
     1788        {                
    17471789 
    17481790                /// \TODO tady je to spatne, tady nesmi byt conditions.size(), viz RARX.bayes() 
    … …  
    17961838                                // cout << &(*tri_ref) << endl; 
    17971839 
    1798                                 rnumber = randu(); 
    1799  
     1840                                bool have_sigma = false; 
    18001841                                double sigma = 0; 
    1801                                 double sum_g = 0; 
    1802                                 for(int i = 0;i<(*s_ref)->gamma_parameters.size();i++) 
    1803                                 { 
    1804                                         for(multimap<double,double>::iterator g_ref = (*s_ref)->gamma_parameters[i].begin();g_ref != (*s_ref)->gamma_parameters[i].end();g_ref++) 
    1805                                         { 
    1806                                                 sum_g += (*g_ref).first/(*s_ref)->gamma_sum; 
    1807  
    1808                                                 if(sum_g>rnumber) 
    1809                                                 { 
    1810                                                         itpp::Gamma_RNG* gamma = new itpp::Gamma_RNG(conditions.size()-number_of_parameters,(*g_ref).second); 
    1811                                                         sigma = (*gamma)(); 
     1842                                do 
     1843                                { 
     1844                                        rnumber = randu(); 
     1845                                         
     1846                                        double sum_g = 0; 
     1847                                        for(int i = 0;i<(*s_ref)->positive_gamma_parameters.size();i++) 
     1848                                        { 
     1849                                                for(multimap<double,double>::iterator g_ref = (*s_ref)->positive_gamma_parameters[i].begin();g_ref != (*s_ref)->positive_gamma_parameters[i].end();g_ref++) 
     1850                                                { 
     1851                                                        sum_g += (*g_ref).first/(*s_ref)->positive_gamma_sum; 
     1852 
     1853                                                        if(sum_g>rnumber) 
     1854                                                        { 
     1855                                                                itpp::Gamma_RNG* gamma = new itpp::Gamma_RNG(conditions.size()-number_of_parameters,1/(*g_ref).second); 
     1856                                                                sigma = 1/(*gamma)(); 
     1857                                                                break; 
     1858                                                        }                                                
     1859                                                } 
     1860 
     1861                                                if(sigma!=0) 
     1862                                                { 
    18121863                                                        break; 
    1813                                                 }                                                
    1814                                         } 
    1815  
    1816                                         if(sigma!=0) 
    1817                                         { 
    1818                                                 break; 
    1819                                         } 
    1820                                 } 
     1864                                                } 
     1865                                        } 
     1866 
     1867                                        rnumber = randu(); 
     1868 
     1869                                        double pg_sum = 0; 
     1870                                        for(vector<multimap<double,double>>::iterator v_ref = (*s_ref)->positive_gamma_parameters.begin();v_ref!=(*s_ref)->positive_gamma_parameters.end();v_ref++) 
     1871                                        { 
     1872                                                for(multimap<double,double>::iterator pg_ref = (*v_ref).begin();pg_ref!=(*v_ref).end();pg_ref++) 
     1873                                                { 
     1874                                                        pg_sum += exp(-(*pg_ref).second/sigma)*pow((*pg_ref).second/sigma,(int)conditions.size()-number_of_parameters-1)*(*pg_ref).second/fact(conditions.size()-number_of_parameters-1)*(*pg_ref).first; 
     1875                                                }                                        
     1876                                        } 
     1877 
     1878                                        double ng_sum = 0; 
     1879                                        for(vector<multimap<double,double>>::iterator v_ref = (*s_ref)->negative_gamma_parameters.begin();v_ref!=(*s_ref)->negative_gamma_parameters.end();v_ref++) 
     1880                                        { 
     1881                                                for(multimap<double,double>::iterator ng_ref = (*v_ref).begin();ng_ref!=(*v_ref).end();ng_ref++) 
     1882                                                { 
     1883                                                        ng_sum += exp(-(*ng_ref).second/sigma)*pow((*ng_ref).second/sigma,(int)conditions.size()-number_of_parameters-1)*(*ng_ref).second/fact(conditions.size()-number_of_parameters-1)*(*ng_ref).first; 
     1884                                                }                                        
     1885                                        } 
     1886                                         
     1887                                        if((pg_sum-ng_sum)/pg_sum>rnumber) 
     1888                                        { 
     1889                                                have_sigma = true; 
     1890                                        } 
     1891                                } 
     1892                                while(!have_sigma); 
    18211893 
    18221894                                int dimension = (*s_ref)->vertices.size()-1;