Skip to content

Commit

Permalink
windows!!! closes #2996
Browse files Browse the repository at this point in the history
  • Loading branch information
@petrelharp
petrelharp committed Sep 25, 2024
1 parent 9e34c22 commit 39174d2
Show file tree
Hide file tree
Showing 4 changed files with 136 additions and 85 deletions.
97 changes: 76 additions & 21 deletions c/tests/test_stats.c
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2001,10 +2001,12 @@ test_empty_genetic_relatedness_vector(void)
double *weights, *result;
tsk_size_t j;
tsk_size_t num_weights = 2;
double windows[] = { 0, 0 };

tsk_treeseq_from_text(
&ts, 1, single_tree_ex_nodes, "", NULL, NULL, NULL, NULL, NULL, 0);
num_samples = tsk_treeseq_get_num_samples(&ts);
windows[1] = tsk_treeseq_get_sequence_length(&ts);
weights = tsk_malloc(num_weights * num_samples * sizeof(double));
result = tsk_malloc(num_weights * num_samples * sizeof(double));
for (j = 0; j < num_samples; j++) {
Expand All @@ -2015,11 +2017,11 @@ test_empty_genetic_relatedness_vector(void)
}

ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, 0);
&ts, num_weights, weights, 1, windows, result, 0);
CU_ASSERT_EQUAL_FATAL(ret, 0);

ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, TSK_STAT_NONCENTRED);
&ts, num_weights, weights, 1, windows, result, TSK_STAT_NONCENTRED);
CU_ASSERT_EQUAL_FATAL(ret, 0);

tsk_treeseq_free(&ts);
Expand All @@ -2028,47 +2030,67 @@ test_empty_genetic_relatedness_vector(void)
}

static void
test_paper_ex_genetic_relatedness_vector(void)
verify_genetic_relatedness_vector(
tsk_treeseq_t *ts, tsk_size_t num_weights, tsk_size_t num_windows)
{
int ret;
tsk_treeseq_t ts;
tsk_size_t num_samples;
double *weights, *result;
tsk_size_t j;
tsk_size_t num_weights = 2;

tsk_treeseq_from_text(&ts, 10, paper_ex_nodes, paper_ex_edges, NULL, paper_ex_sites,
paper_ex_mutations, paper_ex_individuals, NULL, 0);
num_samples = tsk_treeseq_get_num_samples(&ts);
tsk_size_t j, k;
double *windows = tsk_malloc((num_windows + 1) * sizeof(*windows));
double L = tsk_treeseq_get_sequence_length(ts);

weights = tsk_malloc(num_weights * num_samples * sizeof(double));
result = tsk_malloc(num_weights * num_samples * sizeof(double));
for (j = 0; j < num_samples; j++) {
weights[j] = 1.0;
windows[0] = 0;
windows[num_windows] = L;
for (j = 1; j < num_windows; j++) {
windows[j] = ((double) j) * L / (double) num_windows;
}
num_samples = tsk_treeseq_get_num_samples(ts);

weights = tsk_malloc(num_weights * num_samples * sizeof(*weights));
result = tsk_malloc(num_windows * num_weights * num_samples * sizeof(*result));
for (j = 0; j < num_samples; j++) {
weights[j + num_samples] = (float) j;
for (k = 0; k < num_weights; k++) {
weights[j + k * num_samples] = 1.0 + (double) k;
}
}

ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, 0);
ts, num_weights, weights, num_windows, windows, result, 0);
CU_ASSERT_EQUAL_FATAL(ret, 0);

ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, TSK_STAT_NONCENTRED);
ts, num_weights, weights, num_windows, windows, result, TSK_STAT_NONCENTRED);
CU_ASSERT_EQUAL_FATAL(ret, 0);

tsk_set_debug_stream(_devnull);
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, TSK_DEBUG);
ts, num_weights, weights, num_windows, windows, result, TSK_DEBUG);
CU_ASSERT_EQUAL_FATAL(ret, 0);
tsk_set_debug_stream(stdout);

tsk_treeseq_free(&ts);
free(windows);
free(weights);
free(result);
}

static void
test_paper_ex_genetic_relatedness_vector(void)
{
tsk_treeseq_t ts;

tsk_treeseq_from_text(&ts, 10, paper_ex_nodes, paper_ex_edges, NULL, paper_ex_sites,
paper_ex_mutations, paper_ex_individuals, NULL, 0);

tsk_size_t j, k;
for (j = 1; j < 3; j++) {
for (k = 1; k < 3; k++) {
verify_genetic_relatedness_vector(&ts, j, k);
}
}
tsk_treeseq_free(&ts);
}

static void
test_paper_ex_genetic_relatedness_vector_errors(void)
{
Expand All @@ -2078,6 +2100,7 @@ test_paper_ex_genetic_relatedness_vector_errors(void)
double *weights, *result;
tsk_size_t j;
tsk_size_t num_weights = 2;
double windows[] = { 0, 0, 0 };

tsk_treeseq_from_text(&ts, 10, paper_ex_nodes, paper_ex_edges, NULL, paper_ex_sites,
paper_ex_mutations, paper_ex_individuals, NULL, 0);
Expand All @@ -2092,11 +2115,43 @@ test_paper_ex_genetic_relatedness_vector_errors(void)
weights[j + num_samples] = (float) j;
}

/* Window errors */
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 0, windows, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_NUM_WINDOWS);
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 0, NULL, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_NUM_WINDOWS);

ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 2, windows, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_WINDOWS);

windows[0] = -1;
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 2, windows, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_WINDOWS);

windows[0] = 10;
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 2, windows, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_WINDOWS);

windows[0] = 0;
windows[2] = 12;
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, 1, weights, 2, windows, result, TSK_STAT_BRANCH);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_BAD_WINDOWS);

/* unsupported mode errors */
windows[0] = 0.0;
windows[1] = 5.0;
windows[2] = 10.0;
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, TSK_STAT_SITE);
&ts, num_weights, weights, 2, windows, result, TSK_STAT_SITE);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_UNSUPPORTED_STAT_MODE);
ret = tsk_treeseq_genetic_relatedness_vector(
&ts, num_weights, weights, 0, NULL, result, TSK_STAT_NODE);
&ts, num_weights, weights, 2, windows, result, TSK_STAT_NODE);
CU_ASSERT_EQUAL_FATAL(ret, TSK_ERR_UNSUPPORTED_STAT_MODE);

tsk_treeseq_free(&ts);
Expand Down
Loading

0 comments on commit 39174d2

Please sign in to comment.