Ordinary least-squares fit

1. (6 points) Ordinary least-squares fit by QR-decomposition

   1. Implement a function that makes a least-squares fit—using your QR-decomposition routines—of a given data-set, {x_i, y_i, δy_i}_i=1...n , with a linear combination F(x)≐∑ _k=1..m c_k f_k(x) of given functions f_k(x)|_k=1..m .

      The parameters to the function should be the data to fit, {x_i, y_i, δy_i}, and the set of functions, {f_k}, the linear combination of which should fit the data. The function must calculate the vector of the coefficients, {c_k}.

   2. Fit the following data,

      x  =     0.1    1.33    2.55    3.78       5    6.22    7.45    8.68     9.9
      y  =   -15.3    0.32    2.45    2.75    2.27    1.35   0.157   -1.23   -2.75
      dy =    1.04   0.594   0.983   0.998    1.11   0.398   0.535   0.968   0.478
      

      with the linear combination of functions {log(x), 1, x}.

   Hint: the set of fitting functions {f_k(x)} can be implemented, e.g., as

   • a function of two parameters (C,Fortran,C++,Python):

     double fitfunctions(int i, double x){return fi(x);}
     

     For example:

     #include<math.h> /* NAN is here */
     double funs(int i, double x){
        switch(i){
        case 0: return log(x); break;
        case 1: return 1.0;   break;
        case 2: return x;     break;
        default: {fprintf(stderr,"funs: wrong i:%d",i); return NAN;}
        }
     }
     

   • a vector of functions (C++11):

     std::vector<std::function<double(double)>> fitfunctions = {f1,...,fm};
     

   • a list of functions (Python):

     fitfunctions = [f1,...,fm]
     

2. (3 points) Uncertainties of the fitting coefficients

   • Modify you least-squares fitting function such that it also calculates the covariance matrix and the uncertainties of the fitting coefficients.

   • Make up some interesting fits and check that the reported uncertainties of the fitting coefficients are reasonable. For example,

     

3. (1 points) Ordinary least-squares solution by thin singular-value decomposition

   1. Implement a function that makes a singular value decomposition of a (tall) matrix A by diagonalising A^TA (with your implementation of the Jacobi eigenvalue algorithm), A^TA=VDV^T (where D is a diagonal matrix with eigenvalues of A^TA and V is the matrix of the corresponding eigenvectors) and then using A=UΣV^T where U=AVD^-1/2 and Σ=D^1/2.
   2. Implement a function that solves the linear least squares problem Ax=b using your implementation of singular-value decomposition.
   3. Implement a function that makes ordinary least-squares fit using your implementation of singular-value decomposition.