// Example of using:
// nagcpp::correg::lars (g02ma)
#include <iomanip>
#include <iostream>
#include <math.h>
#include <vector>
#include "g02/nagcpp_g02ma.hpp"
// a simple sample matrix class
#include "include/nag_my_matrix.hpp"
int main(void) {
std::cout << "nagcpp::correg::lars Example" << std::endl;
nagcpp::types::f77_integer mtype = 1;
// clang-format off
std::vector<double> rd = {
10.28, 1.77, 9.69, 15.58, 8.23, 10.44,
9.08, 8.99, 11.53, 6.57, 15.89, 12.58,
17.98, 13.10, 1.04, 10.45, 10.12, 16.68,
14.82, 13.79, 12.23, 7.00, 8.14, 7.79,
17.53, 9.41, 6.24, 3.75, 13.12, 17.08,
7.78, 10.38, 9.83, 2.58, 10.13, 4.25,
11.95, 21.71, 8.83, 11.00, 12.59, 10.52,
14.60, 10.09, -2.70, 9.89, 14.67, 6.49,
3.63, 9.07, 12.59, 14.09, 9.06, 8.19,
6.35, 9.79, 9.40, 12.79, 8.38, 16.79,
4.66, 3.55, 16.82, 13.83, 21.39, 13.88,
8.32, 14.04, 17.17, 7.93, 7.39, -1.09,
10.86, 13.68, 5.75, 10.44, 10.36, 10.06,
4.76, 4.92, 17.83, 2.90, 7.58, 11.97,
5.05, 10.41, 9.89, 9.04, 7.90, 13.12,
5.41, 9.32, 5.27, 15.53, 5.06, 19.84,
9.77, 2.37, 9.54, 20.23, 9.33, 8.82,
14.28, 4.34, 14.23, 14.95, 18.16, 11.03,
10.17, 6.80, 3.17, 8.57, 16.07, 15.93,
5.39, 2.67, 6.37, 13.56, 10.68, 7.35
};
// clang-format on
MyMatrix<double> d(20, 6, rd);
std::vector<double> y = {-46.47, -35.80, -129.22, -42.44, -73.51,
-26.61, -63.90, -76.73, -32.64, -83.29,
-16.31, -5.82, -47.75, 18.38, -54.71,
-55.62, -45.28, -22.76, -104.32, -55.94};
nagcpp::types::f77_integer ip, nstep;
MyMatrix<double> b, fitsum;
// we are going to be using the default values for the optional
// arguments, however we want to print out any warning messages
nagcpp::correg::OptionalG02MA opt;
try {
nagcpp::correg::lars(mtype, d, nullptr, y, ip, nstep, b, fitsum, opt);
} catch (nagcpp::error_handler::Exception &e) {
std::cout << e.msg << std::endl;
return 1;
}
if (opt.fail.warning_thrown) {
std::cout << opt.fail.msg << std::endl;
}
std::cout.precision(3);
std::cout.setf(std::ios::fixed);
// display the parameter estimates
std::cout << " Step " << std::string(((ip > 2) ? (ip - 2) : 0) * 5, ' ')
<< "Parameter Estimate" << std::endl;
std::cout << std::string(5 + ip * 10, '-') << std::endl;
for (int k = 0; k < nstep; k++) {
std::cout << std::setw(4) << k + 1;
for (int i = 0; i < ip; i++) {
std::cout << std::setw(10) << b(i, k);
}
std::cout << std::endl;
}
std::cout << std::endl;
std::cout << " alpha: " << std::setw(10) << fitsum(0, nstep) << std::endl;
std::cout << std::endl;
std::cout << " Step Sum RSS df Cp Ck "
<< "Step Size" << std::endl;
std::cout << std::string(64, '-') << std::endl;
for (int k = 0; k < nstep; k++) {
std::cout << std::setw(4) << k + 1 << std::setw(10) << fitsum(0, k)
<< std::setw(10) << fitsum(1, k) << std::setw(7)
<< static_cast<int>(floor(fitsum(2, k) + 0.5)) << std::setw(10)
<< fitsum(3, k) << std::setw(10) << fitsum(4, k) << std::setw(10)
<< fitsum(5, k) << std::endl;
}
std::cout << std::endl;
std::cout << "sigma^2: " << fitsum(4, nstep) << std::endl;
return 0;
}