X-Git-Url: https://pintos-os.org/cgi-bin/gitweb.cgi?a=blobdiff_plain;f=src%2Fmath%2Fcovariance-matrix.c;h=f929a370cf5470456ba0c33018a7f292d2e4e3a4;hb=b734a970ee8bfda57e72398d4f424fb16f2ea80a;hp=85998dfe49d80e1e20081f2bb503ce15bfc230e9;hpb=685407ad62911e5edb1ec093a01ec9e46563af44;p=pspp diff --git a/src/math/covariance-matrix.c b/src/math/covariance-matrix.c index 85998dfe49..f929a370cf 100644 --- a/src/math/covariance-matrix.c +++ b/src/math/covariance-matrix.c @@ -39,80 +39,96 @@ void covariance_matrix_destroy (struct design_matrix *x) } /* - Update the covariance matrix with the new entries, assuming that V1 - is categorical and V2 is numeric. + Update the covariance matrix with the new entries, assuming that ROW + corresponds to a categorical variable and V2 is numeric. */ static void covariance_update_categorical_numeric (struct design_matrix *cov, double mean, - double weight, double ssize, const struct variable *v1, - const struct variable *v2, const union value *val1, const union value *val2) + size_t row, + const struct variable *v2, double x, const union value *val2) { - double x; - size_t i; size_t col; - size_t row; + double tmp; - assert (var_is_alpha (v1)); assert (var_is_numeric (v2)); - row = design_matrix_var_to_column (cov, v1); col = design_matrix_var_to_column (cov, v2); - for (i = 0; i < cat_get_n_categories (v1); i++) + assert (val2 != NULL); + tmp = gsl_matrix_get (cov->m, row, col); + gsl_matrix_set (cov->m, row, col, (val2->f - mean) * x + tmp); + gsl_matrix_set (cov->m, col, row, (val2->f - mean) * x + tmp); +} +static void +column_iterate (struct design_matrix *cov, const struct variable *v, + double ssize, double x, const union value *val1, size_t row) +{ + size_t col; + size_t i; + double y; + double tmp; + const union value *tmp_val; + + col = design_matrix_var_to_column (cov, v); + for (i = 0; i < cat_get_n_categories (v) - 1; i++) { - row += i; - x = -1.0 * cat_get_n_categories (v1) / ssize; - if (i == cat_value_find (v1, val1)) + col += i; + y = -1.0 * cat_get_category_count (i, v) / ssize; + tmp_val = cat_subscript_to_value (i, v); + if (compare_values (tmp_val, val1, v)) { - x += 1.0; + y += -1.0; } - assert (val2 != NULL); - gsl_matrix_set (cov->m, row, col, (val2->f - mean) * x * weight); + tmp = gsl_matrix_get (cov->m, row, col); + gsl_matrix_set (cov->m, row, col, x * y + tmp); + gsl_matrix_set (cov->m, col, row, x * y + tmp); } } - /* - Call this function in the first data pass. The central moments are + Call this function in the second data pass. The central moments are MEAN1 and MEAN2. Any categorical variables should already have their values summarized in in its OBS_VALS element. */ -void covariance_pass_one (struct design_matrix *cov, double mean1, double mean2, - double weight, double ssize, const struct variable *v1, +void covariance_pass_two (struct design_matrix *cov, double mean1, double mean2, + double ssize, const struct variable *v1, const struct variable *v2, const union value *val1, const union value *val2) { size_t row; size_t col; size_t i; double x; - double y; + const union value *tmp_val; if (var_is_alpha (v1)) { - if (var_is_numeric (v2)) - { - covariance_update_categorical_numeric (cov, mean2, weight, ssize, v1, - v2, val1, val2); - } - else + row = design_matrix_var_to_column (cov, v1); + for (i = 0; i < cat_get_n_categories (v1) - 1; i++) { - row = design_matrix_var_to_column (cov, v1); - col = design_matrix_var_to_column (cov, v2); - for (i = 0; i < cat_get_n_categories (v2); i++) + row += i; + x = -1.0 * cat_get_category_count (i, v1) / ssize; + tmp_val = cat_subscript_to_value (i, v1); + if (compare_values (tmp_val, val1, v1)) { - col += i; - y = -1.0 * cat_get_n_categories (v2) / ssize; - if (i == cat_value_find (v2, val2)) - { - y += 1.0; - } - gsl_matrix_set (cov->m, row, col, x * y * weight); - gsl_matrix_set (cov->m, col, row, x * y * weight); + x += 1.0; + } + if (var_is_numeric (v2)) + { + covariance_update_categorical_numeric (cov, mean2, row, + v2, x, val2); + } + else + { + column_iterate (cov, v1, ssize, x, val1, row); + column_iterate (cov, v2, ssize, x, val2, row); } } } else if (var_is_alpha (v2)) { - covariance_update_categorical_numeric (cov, mean1, weight, ssize, v2, - v1, val2, val1); + /* + Reverse the orders of V1, V2, etc. and put ourselves back + in the previous IF scope. + */ + covariance_pass_two (cov, mean2, mean1, ssize, v2, v1, val2, val1); } else { @@ -121,7 +137,8 @@ void covariance_pass_one (struct design_matrix *cov, double mean1, double mean2, */ row = design_matrix_var_to_column (cov, v1); col = design_matrix_var_to_column (cov, v2); - x = (val1->f - mean1) * (val2->f - mean2) * weight; + x = (val1->f - mean1) * (val2->f - mean2); + x += gsl_matrix_get (cov->m, col, row); gsl_matrix_set (cov->m, row, col, x); gsl_matrix_set (cov->m, col, row, x); }