gsl_matrix_set (out, i, j, x);
}
}
-
+
gsl_matrix_free (in);
return out;
struct covariance
{
+ /* True if the covariances are centerered. (ie Real covariances) */
+ bool centered;
+
/* The variables for which the covariance matrix is to be calculated. */
size_t n_vars;
const struct variable *const *vars;
double *cm;
int n_cm;
- /* 1 for single pass algorithm;
+ /* 1 for single pass algorithm;
2 for double pass algorithm
*/
short passes;
/*
0 : No pass has been made
1 : First pass has been started
- 2 : Second pass has been
-
+ 2 : Second pass has been
+
IE: How many passes have been (partially) made. */
short state;
*/
struct covariance *
covariance_1pass_create (size_t n_vars, const struct variable *const *vars,
- const struct variable *weight, enum mv_class exclude)
+ const struct variable *weight, enum mv_class exclude,
+ bool centered)
{
size_t i;
struct covariance *cov = xzalloc (sizeof *cov);
+ cov->centered = centered;
cov->passes = 1;
cov->state = 0;
cov->pass_one_first_case_seen = cov->pass_two_first_case_seen = false;
-
+
cov->vars = vars;
cov->wv = weight;
cov->dim = n_vars;
cov->moments = xmalloc (sizeof *cov->moments * n_MOMENTS);
-
+
for (i = 0; i < n_MOMENTS; ++i)
cov->moments[i] = gsl_matrix_calloc (n_vars, n_vars);
struct covariance *
covariance_2pass_create (size_t n_vars, const struct variable *const *vars,
struct categoricals *cats,
- const struct variable *wv, enum mv_class exclude)
+ const struct variable *wv, enum mv_class exclude,
+ bool centered)
{
size_t i;
struct covariance *cov = xmalloc (sizeof *cov);
+ cov->centered = centered;
cov->passes = 2;
cov->state = 0;
cov->pass_one_first_case_seen = cov->pass_two_first_case_seen = false;
-
+
cov->vars = vars;
cov->wv = wv;
cov->dim = n_vars;
cov->moments = xmalloc (sizeof *cov->moments * n_MOMENTS);
-
+
for (i = 0; i < n_MOMENTS; ++i)
cov->moments[i] = gsl_matrix_calloc (n_vars, n_vars);
return cov;
}
-/* Return an integer, which can be used to index
+/* Return an integer, which can be used to index
into COV->cm, to obtain the I, J th element
of the covariance matrix. If COV->cm does not
contain that element, then a negative value
int as;
const int n2j = cov->dim - 2 - j;
const int nj = cov->dim - 2 ;
-
+
assert (i >= 0);
assert (j < cov->dim);
if (j >= cov->dim - 1)
return -1;
- if ( i <= j)
+ if ( i <= j)
return -1 ;
as = nj * (nj + 1) ;
- as -= n2j * (n2j + 1) ;
+ as -= n2j * (n2j + 1) ;
as /= 2;
return i - 1 + as;
/*
- Returns true iff the variable corresponding to the Ith element of the covariance matrix
+ Returns true iff the variable corresponding to the Ith element of the covariance matrix
has a missing value for case C
*/
static bool
is_missing (const struct covariance *cov, int i, const struct ccase *c)
{
const struct variable *var = i < cov->n_vars ?
- cov->vars[i] :
+ cov->vars[i] :
categoricals_get_interaction_by_subscript (cov->categoricals, i - cov->n_vars)->vars[0];
const union value *val = case_data (c, var);
categoricals_done (cov->categoricals);
cov->dim = cov->n_vars;
-
+
if (cov->categoricals)
cov->dim += categoricals_df_total (cov->categoricals);
*x += s;
}
- ss =
+ ss =
(v1 - gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
- *
+ *
(v2 - gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
* weight
;
}
-/*
+/*
Allocate and return a gsl_matrix containing the covariances of the
data.
*/
}
}
- /* Centre the moments */
- for ( j = 0 ; j < cov->dim - 1; ++j)
+ if (cov->centered)
{
- for (i = j + 1 ; i < cov->dim; ++i)
+ /* Centre the moments */
+ for ( j = 0 ; j < cov->dim - 1; ++j)
{
- double *x = &cov->cm [cm_idx (cov, i, j)];
-
- *x /= gsl_matrix_get (cov->moments[0], i, j);
+ for (i = j + 1 ; i < cov->dim; ++i)
+ {
+ double *x = &cov->cm [cm_idx (cov, i, j)];
- *x -=
- gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
- *
- gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i);
+ *x /= gsl_matrix_get (cov->moments[0], i, j);
+
+ *x -=
+ gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
+ *
+ gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i);
+ }
}
}
switch (cov->passes)
{
case 1:
- return covariance_calculate_single_pass (cov);
+ return covariance_calculate_single_pass (cov);
break;
case 2:
- return covariance_calculate_double_pass (cov);
+ return covariance_calculate_double_pass (cov);
break;
default:
NOT_REACHED ();
{
size_t i, j;
- for (i = 0 ; i < cov->dim; ++i)
+ if (cov->centered)
{
- for (j = 0 ; j < cov->dim; ++j)
+ for (i = 0 ; i < cov->dim; ++i)
{
- double *x = gsl_matrix_ptr (cov->moments[MOMENT_VARIANCE], i, j);
- *x -= pow2 (gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
- / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+ for (j = 0 ; j < cov->dim; ++j)
+ {
+ double *x = gsl_matrix_ptr (cov->moments[MOMENT_VARIANCE], i, j);
+ *x -= pow2 (gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j))
+ / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+ }
}
- }
- for ( j = 0 ; j < cov->dim - 1; ++j)
- {
- for (i = j + 1 ; i < cov->dim; ++i)
+
+ for ( j = 0 ; j < cov->dim - 1; ++j)
{
- double *x = &cov->cm [cm_idx (cov, i, j)];
-
- *x -=
- gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
- *
- gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i)
- / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+ for (i = j + 1 ; i < cov->dim; ++i)
+ {
+ double *x = &cov->cm [cm_idx (cov, i, j)];
+
+ *x -=
+ gsl_matrix_get (cov->moments[MOMENT_MEAN], i, j)
+ *
+ gsl_matrix_get (cov->moments[MOMENT_MEAN], j, i)
+ / gsl_matrix_get (cov->moments[MOMENT_NONE], i, j);
+ }
}
}
-
+
return cm_to_gsl (cov);
}
switch (cov->passes)
{
case 1:
- cov->unnormalised = covariance_calculate_single_pass_unnormalized (cov);
+ cov->unnormalised = covariance_calculate_single_pass_unnormalized (cov);
break;
case 2:
- cov->unnormalised = covariance_calculate_double_pass_unnormalized (cov);
+ cov->unnormalised = covariance_calculate_double_pass_unnormalized (cov);
break;
default:
NOT_REACHED ();
tab_title (t, "Covariance Encoding");
- tab_box (t,
+ tab_box (t,
TAL_2, TAL_2, 0, 0,
0, 0, tab_nc (t) - 1, tab_nr (t) - 1);