13.6.13. Macro “modelerbuildSimpleGPEstimWithCov.C”

13.6.13.1. Objective

This macro is the one described in Prediction of a new data set, global approach, to create and use a simple gaussian process, whose training (utf_4D_train.dat) and testing (utf_4D_test.dat) database can both be found in the document folder of the Uranie installation (${URANIESYS}/share/uranie/docUMENTS). It uses the global approach, computing the input covariance matrix and translating it to the prediction covariance matrix.

13.6.13.2. Macro Uranie

    // Load observations
    TDataServer *tdsObs = new TDataServer("tdsObs", "observations");
    TDataServer *tdsEstim = new TDataServer("tdsEstim", "estimations");

    tdsObs->fileDataRead("utf_4D_train.dat");
        
    // Construct the GPBuilder
    TGPBuilder *gpb = new TGPBuilder(tdsObs,            // observations data 
                   "x1:x2:x3:x4",     // list of input variables
                   "y",               // output variable
                   "matern1/2");       // name of the correlation function
                      
    // Search for the optimal hyper-parameters
    gpb->findOptimalParameters("ML",            // optimisation criterion 
                 100,             // screening design size
                 "neldermead",    // optimisation algorithm
                 500);            // max. number of optimisation iterations
        
    // Construct the kriging model
    TKriging *krig = gpb->buildGP();
        
    // Display model information
    krig->printLog();
        
    // Load the data to estimate
    tdsEstim->fileDataRead("utf_4D_test.dat");

    // Reducing the database to 1000 first event (prediction cov matrix of a million value !)
    int nST=1000;
    tdsEstim->exportData("utf_4D_test_red.dat","",Form("tdsEstim__n__iter__<=%d",nST));
    tdsEstim->fileDataRead("utf_4D_test_red.dat",false,true); // Reload reduce sample    

    krig->estimateWithCov(tdsEstim, // data to estimate
                        "x1:x2:x3:x4",// list of the input variables
                        "yEstim:vEstim", // name given to the model's outputs
                        "y", //name of the true reference if validation database
                        ""); //options

    TCanvas *c2=NULL;
    // Residuals plots if true information provided
    if( tdsEstim->isAttribute("_Residuals_") ) 
    {
        c2 = new TCanvas("c2","c2",1200,800);
        c2->Divide(2,1);
        c2->cd(1);
        // Usual residual considering uncorrated input points
        tdsEstim->Draw("_Residuals_"); 
        c2->cd(2);
        // Corrected residuals, with prediction covariance matrix
        tdsEstim->Draw("_uncorrResiduals_"); 
    }

    // Retrieve all the prediction covariance coefficient
    tdsEstim->getTuple()->SetEstimate( nST * nST ); //allocate the correct size
    // Get a pointer to all values
    tdsEstim->getTuple()->Draw("_CovarianceMatrix_","","goff");
    double *cov=tdsEstim->getTuple()->GetV1();

    //Put these in a matrix nicely created
    TMatrixD Cov(nST,nST);
    Cov.Use(0,nST-1,0,nST-1,cov);

    //Print it if size is reasonnable
    if(nST<10)
        Cov.Print();

13.6.13.3. Graph

Figure 13.44 Graph of the macro “modelerbuildSimpleGPEstimWithCov.C”