Как использовать реализацию GSL коэффициента корреляции Пирсона?

Question

У меня есть два вектора поплавков x и y, и я хочу вычислить коэффициенты корреляции Pearson. Поскольку я должен делать это на большом количестве данных (например, 10 миллионов разных векторов x и 20 тысяч разных векторов y), я использую C++, а точнее функцию gsl_stats_correlation от GSL.

Вот мой C++ код:

#include <iostream> 
#include <vector> 
using namespace std; 

#include <gsl/gsl_vector.h> 
#include <gsl/gsl_statistics.h> 

int main (int argc, char ** argv) 
{ 
    vector<double> x, y; 
    size_t n = 5; 
    x.push_back(1.0); y.push_back(1.0); 
    x.push_back(3.1); y.push_back(3.2); 
    x.push_back(2.0); y.push_back(1.9); 
    x.push_back(5.0); y.push_back(4.9); 
    x.push_back(2.0); y.push_back(2.1); 
    for(size_t i=0; i<n; ++i) 
     printf ("x[%ld]=%.1f y[%ld]=%.1f\n", i, x[i], i, y[i]); 

    gsl_vector_const_view gsl_x = gsl_vector_const_view_array(&x[0], x.size()); 
    gsl_vector_const_view gsl_y = gsl_vector_const_view_array(&y[0], y.size()); 

    double pearson = gsl_stats_correlation((double*) gsl_x.vector.data, sizeof(double), 
              (double*) gsl_y.vector.data, sizeof(double), 
              n); 
    printf ("Pearson correlation = %f\n", pearson); 

    return 0; 
} 

Он успешно компилируется (НКУ -Wall -g pearson.cpp -lstdC++ -lgsl -lgslcblas -o пирсон), но когда я запускаю его здесь есть выход:

x[0]=1.0 y[0]=1.0 
x[1]=3.1 y[1]=3.2 
x[2]=2.0 y[2]=1.9 
x[3]=5.0 y[3]=4.9 
x[4]=2.0 y[4]=2.1 
Pearson correlation = 1.000000 

Хотя очевидно, что результаты не должны быть точно 1, как показано с помощью следующего кода R:

x <- c(1.0,3.1,2.0,5.0,2.0); y <-c(1.0,3.2,1.9,4.9,2.1) 
cor(x, y, method="pearson") # 0.99798 

Что мне не хватает?

Answer 1

Изменить строку:

double pearson = gsl_stats_correlation((double*) gsl_x.vector.data, sizeof(double), 
              (double*) gsl_y.vector.data, sizeof(double), 
              n); 

к:

double pearson = gsl_stats_correlation((double*) gsl_x.vector.data, 1, 
              (double*) gsl_y.vector.data, 1, 
              n); 

или, если вы хотите, чтобы избежать повторения "магическое число" 1:

const size_t stride = 1; 
    double pearson = gsl_stats_correlation((double*) gsl_x.vector.data, stride, 
              (double*) gsl_y.vector.data, stride, 
              n); 

gsl_stats_correlation предполагает double, второй и четвертый аргументы - количество удвоений для «шагания», поэтому, давая ему sizeof(double) было прыгнуто на sizeof(double)*sizeof(double) байт.