SciRuby/rb-gsl

QR_lssolve() seems to insist on only solving square problems, unlike gsl_linalg_QR_lssolve()

crazzell opened this issue · 2 comments

I am trying to use gsl (1.16.0.6) for solving over-determined least-squares problems, and a simple test case is not working for me.

require("gsl")
include GSL

m = GSL::Matrix::alloc([0.18, 0.60, 0.57], [0.24, 0.99, 0.58],
            [0.14, 0.30, 0.97], [0.51,  0.19, 0.85], [0.34, 0.91, 0.18])

B = GSL::Vector[1, 2, 3, 4, 5]
qr, tau = m.QR_decomp
x, r = qr.QR_lssolve(tau,B)

Whereas, if I write the same problem in plain C, it works find:

#include <stdio.h>
#include <gsl/gsl_linalg.h>

 int
 main (void)
 {
   double a_data[] = {0.18, 0.60, 0.57, 0.24, 0.99, 0.58,
                0.14, 0.30, 0.97, 0.51,  0.19, 0.85, 0.34, 0.91, 0.18};

   double b_data[] = {1.0,2.0,3.0,4.0,5.0};

   gsl_matrix_view m = gsl_matrix_view_array (a_data, 5, 3);

   gsl_vector_view b = gsl_vector_view_array (b_data, 5);

   gsl_vector *x = gsl_vector_alloc (3); // size equal to n
   gsl_vector *residual = gsl_vector_alloc (5); // size equal to m
   gsl_vector *tau = gsl_vector_alloc (3); //size equal to min(m,n)
   gsl_linalg_QR_decomp (&m.matrix, tau); // 
   gsl_linalg_QR_lssolve(&m.matrix, tau, &b.vector, x, residual);

   printf ("x = \n");
   gsl_vector_fprintf (stdout, x, "%g");
   gsl_vector_free (x);
   gsl_vector_free (tau);
   gsl_vector_free (residual);
   return 0;
}

x =
8.06835
0.884443
0.231883

The only way I can get the ruby code to work is by making the problem square, but that defeats the overall purpose of a least-squares fit in my use-case. Of course, it is just possible that I am messing up the row and column dimensions in the Ruby code in a way that I can't see for looking, but I simply can't find another permutation that should work. I am guessing the ruby-wrapper is mistakenly assuming that over-determined problems are always user-error, whereas this is precisely where a LS fit is most useful.

Thanks! Charles.

Forgot to add the error message from Ruby with the above script:

../testls.rb:9:in 'QR_lssolve': Ruby/GSL error code 19, matrix size must match solution size (file qr.c, line >193), matrix/vector sizes are not conformant (GSL::ERROR::EBADLEN)
from ../testls.rb:9:in '

'

Apologies for the false alarm, I found that the following method works:

require("gsl")
include GSL

m = GSL::Matrix::alloc([0.18, 0.60, 0.57], [0.24, 0.99, 0.58],
                [0.14, 0.30, 0.97], [0.51,  0.19, 0.85], [0.34, 0.91, 0.18])

B = GSL::Vector[1, 2, 3, 4, 5]
qr, tau = m.QR_decomp
x = GSL::Vector.alloc(3)
r = GSL::Vector.alloc(5)
qr.QR_lssolve(tau,B,x,r)

p x

I was not expecting to have to put the outputs {x,r} as arguments to the function, so this took me a while to find.