2323 */
2424
2525#include "testlib.h"
26+ #include <math.h>
2627#include <tskit/stats.h>
2728
2829#include <unistd.h>
@@ -602,7 +603,7 @@ verify_pair_coalescence_rates(tsk_treeseq_t *ts)
602603 epochs [B ] = DBL_MAX ;
603604 ret = tsk_treeseq_pair_coalescence_rates (ts , P , sample_set_sizes , sample_sets , I ,
604605 index_tuples , T , breakpoints , B , node_bin_map , epochs , 0 , C );
605- CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS );
606+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS_END );
606607 epochs [B ] = INFINITY ;
607608
608609 node_bin_map [0 ] = (tsk_id_t ) B ;
@@ -868,6 +869,84 @@ verify_one_way_stat_func_errors(tsk_treeseq_t *ts, one_way_sample_stat_method *m
868869 CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_WINDOWS );
869870}
870871
872+ // Temporary definition for time_windows in tsk_treeseq_allele_frequency_spectrum
873+ typedef int one_way_sample_stat_method_tw (const tsk_treeseq_t * self ,
874+ tsk_size_t num_sample_sets , const tsk_size_t * sample_set_sizes ,
875+ const tsk_id_t * sample_sets , tsk_size_t num_windows , const double * windows ,
876+ tsk_size_t num_time_windows , const double * time_windows , tsk_flags_t options ,
877+ double * result );
878+
879+ // Temporary duplicate for time-windows-having methods
880+ static void
881+ verify_one_way_stat_func_errors_tw (
882+ tsk_treeseq_t * ts , one_way_sample_stat_method_tw * method )
883+ {
884+ int ret ;
885+ tsk_id_t num_nodes = (tsk_id_t ) tsk_treeseq_get_num_nodes (ts );
886+ tsk_id_t samples [] = { 0 , 1 , 2 , 3 };
887+ tsk_size_t sample_set_sizes = 4 ;
888+ double windows [] = { 0 , 0 , 0 };
889+ double time_windows [] = { -1 , 0.5 , INFINITY };
890+ double result ;
891+
892+ ret = method (ts , 0 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
893+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_INSUFFICIENT_SAMPLE_SETS );
894+
895+ samples [0 ] = TSK_NULL ;
896+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
897+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_NODE_OUT_OF_BOUNDS );
898+ samples [0 ] = -10 ;
899+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
900+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_NODE_OUT_OF_BOUNDS );
901+ samples [0 ] = num_nodes ;
902+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
903+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_NODE_OUT_OF_BOUNDS );
904+ samples [0 ] = num_nodes + 1 ;
905+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
906+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_NODE_OUT_OF_BOUNDS );
907+
908+ samples [0 ] = num_nodes - 1 ;
909+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
910+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_SAMPLES );
911+
912+ samples [0 ] = 1 ;
913+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
914+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_DUPLICATE_SAMPLE );
915+
916+ samples [0 ] = 0 ;
917+ sample_set_sizes = 0 ;
918+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , & result );
919+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_EMPTY_SAMPLE_SET );
920+
921+ sample_set_sizes = 4 ;
922+ /* Window errors */
923+ ret = method (ts , 1 , & sample_set_sizes , samples , 0 , windows , 0 , NULL , 0 , & result );
924+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_NUM_WINDOWS );
925+
926+ ret = method (ts , 1 , & sample_set_sizes , samples , 2 , windows , 0 , NULL , 0 , & result );
927+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_WINDOWS );
928+
929+ /* Time window errors */
930+ ret = method (
931+ ts , 1 , & sample_set_sizes , samples , 0 , NULL , 0 , time_windows , 0 , & result );
932+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS_DIM );
933+
934+ ret = method (
935+ ts , 1 , & sample_set_sizes , samples , 0 , NULL , 2 , time_windows , 0 , & result );
936+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS );
937+
938+ time_windows [0 ] = 0.1 ;
939+ ret = method (
940+ ts , 1 , & sample_set_sizes , samples , 0 , NULL , 2 , time_windows , 0 , & result );
941+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS );
942+
943+ time_windows [0 ] = 0 ;
944+ time_windows [1 ] = 0 ;
945+ ret = method (
946+ ts , 1 , & sample_set_sizes , samples , 0 , NULL , 2 , time_windows , 0 , & result );
947+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_BAD_TIME_WINDOWS );
948+ }
949+
871950static void
872951verify_two_way_stat_func_errors (
873952 tsk_treeseq_t * ts , general_sample_stat_method * method , tsk_flags_t options )
@@ -1180,6 +1259,7 @@ verify_afs(tsk_treeseq_t *ts)
11801259 int ret ;
11811260 tsk_size_t n = tsk_treeseq_get_num_samples (ts );
11821261 tsk_size_t sample_set_sizes [2 ];
1262+ double time_windows [] = { 0 , 1 };
11831263 const tsk_id_t * samples = tsk_treeseq_get_samples (ts );
11841264 double * result = tsk_malloc (n * n * sizeof (* result ));
11851265
@@ -1188,23 +1268,28 @@ verify_afs(tsk_treeseq_t *ts)
11881268 sample_set_sizes [0 ] = n - 2 ;
11891269 sample_set_sizes [1 ] = 2 ;
11901270 ret = tsk_treeseq_allele_frequency_spectrum (
1191- ts , 2 , sample_set_sizes , samples , 0 , NULL , 0 , result );
1271+ ts , 2 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , result );
11921272 CU_ASSERT_EQUAL_FATAL (ret , 0 );
11931273
11941274 ret = tsk_treeseq_allele_frequency_spectrum (
1195- ts , 2 , sample_set_sizes , samples , 0 , NULL , TSK_STAT_POLARISED , result );
1275+ ts , 2 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , TSK_STAT_POLARISED , result );
11961276 CU_ASSERT_EQUAL_FATAL (ret , 0 );
11971277
11981278 ret = tsk_treeseq_allele_frequency_spectrum (ts , 2 , sample_set_sizes , samples , 0 ,
1199- NULL , TSK_STAT_POLARISED | TSK_STAT_SPAN_NORMALISE , result );
1279+ NULL , 0 , NULL , TSK_STAT_POLARISED | TSK_STAT_SPAN_NORMALISE , result );
12001280 CU_ASSERT_EQUAL_FATAL (ret , 0 );
12011281
12021282 ret = tsk_treeseq_allele_frequency_spectrum (ts , 2 , sample_set_sizes , samples , 0 ,
1203- NULL , TSK_STAT_BRANCH | TSK_STAT_POLARISED | TSK_STAT_SPAN_NORMALISE , result );
1283+ NULL , 0 , NULL , TSK_STAT_BRANCH | TSK_STAT_POLARISED | TSK_STAT_SPAN_NORMALISE ,
1284+ result );
12041285 CU_ASSERT_EQUAL_FATAL (ret , 0 );
12051286
12061287 ret = tsk_treeseq_allele_frequency_spectrum (ts , 2 , sample_set_sizes , samples , 0 ,
1207- NULL , TSK_STAT_BRANCH | TSK_STAT_SPAN_NORMALISE , result );
1288+ NULL , 0 , NULL , TSK_STAT_BRANCH | TSK_STAT_SPAN_NORMALISE , result );
1289+ CU_ASSERT_EQUAL_FATAL (ret , 0 );
1290+
1291+ ret = tsk_treeseq_allele_frequency_spectrum (ts , 2 , sample_set_sizes , samples , 0 ,
1292+ NULL , 1 , time_windows , TSK_STAT_BRANCH | TSK_STAT_SPAN_NORMALISE , result );
12081293 CU_ASSERT_EQUAL_FATAL (ret , 0 );
12091294
12101295 free (result );
@@ -2398,21 +2483,26 @@ test_paper_ex_afs_errors(void)
23982483 tsk_size_t sample_set_sizes [] = { 2 , 2 };
23992484 tsk_id_t samples [] = { 0 , 1 , 2 , 3 };
24002485 double result [10 ]; /* not thinking too hard about the actual value needed */
2486+ double time_windows [] = { 0 , 1 };
24012487 int ret ;
24022488
24032489 tsk_treeseq_from_text (& ts , 10 , paper_ex_nodes , paper_ex_edges , NULL , paper_ex_sites ,
24042490 paper_ex_mutations , paper_ex_individuals , NULL , 0 );
24052491
2406- verify_one_way_stat_func_errors (& ts , tsk_treeseq_allele_frequency_spectrum );
2492+ verify_one_way_stat_func_errors_tw (& ts , tsk_treeseq_allele_frequency_spectrum );
24072493
24082494 ret = tsk_treeseq_allele_frequency_spectrum (
2409- & ts , 2 , sample_set_sizes , samples , 0 , NULL , TSK_STAT_NODE , result );
2495+ & ts , 2 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , TSK_STAT_NODE , result );
24102496 CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_UNSUPPORTED_STAT_MODE );
24112497
24122498 ret = tsk_treeseq_allele_frequency_spectrum (& ts , 2 , sample_set_sizes , samples , 0 ,
2413- NULL , TSK_STAT_BRANCH | TSK_STAT_SITE , result );
2499+ NULL , 0 , NULL , TSK_STAT_BRANCH | TSK_STAT_SITE , result );
24142500 CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_MULTIPLE_STAT_MODES );
24152501
2502+ ret = tsk_treeseq_allele_frequency_spectrum (& ts , 2 , sample_set_sizes , samples , 0 ,
2503+ NULL , 1 , time_windows , TSK_STAT_SITE , result );
2504+ CU_ASSERT_EQUAL_FATAL (ret , TSK_ERR_UNSUPPORTED_STAT_MODE );
2505+
24162506 tsk_treeseq_free (& ts );
24172507}
24182508
@@ -2430,14 +2520,14 @@ test_paper_ex_afs(void)
24302520 /* we have two singletons and one tripleton */
24312521
24322522 ret = tsk_treeseq_allele_frequency_spectrum (
2433- & ts , 1 , sample_set_sizes , samples , 0 , NULL , 0 , result );
2523+ & ts , 1 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , 0 , result );
24342524 CU_ASSERT_EQUAL_FATAL (ret , 0 );
24352525 CU_ASSERT_EQUAL_FATAL (result [0 ], 0 );
24362526 CU_ASSERT_EQUAL_FATAL (result [1 ], 3.0 );
24372527 CU_ASSERT_EQUAL_FATAL (result [2 ], 0 );
24382528
24392529 ret = tsk_treeseq_allele_frequency_spectrum (
2440- & ts , 1 , sample_set_sizes , samples , 0 , NULL , TSK_STAT_POLARISED , result );
2530+ & ts , 1 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , TSK_STAT_POLARISED , result );
24412531 CU_ASSERT_EQUAL_FATAL (ret , 0 );
24422532 CU_ASSERT_EQUAL_FATAL (result [0 ], 0 );
24432533 CU_ASSERT_EQUAL_FATAL (result [1 ], 2.0 );
@@ -2465,6 +2555,39 @@ test_paper_ex_divergence_matrix(void)
24652555 tsk_treeseq_free (& ts );
24662556}
24672557
2558+ static void
2559+ test_unary_ex_afs (void )
2560+ {
2561+ tsk_treeseq_t ts ;
2562+ tsk_id_t samples [] = { 0 , 2 , 3 };
2563+ tsk_size_t sample_set_sizes [] = { 3 , 0 };
2564+ double result [25 ];
2565+ int ret ;
2566+
2567+ tsk_treeseq_from_text (& ts , 100 , unary_ex_nodes , unary_ex_edges , NULL , unary_ex_sites ,
2568+ unary_ex_mutations , NULL , NULL , 0 );
2569+ /* we have a singleton and a doubleton */
2570+
2571+ ret = tsk_treeseq_allele_frequency_spectrum (
2572+ & ts , 1 , sample_set_sizes , samples , 0 , NULL , 0 , NULL , TSK_STAT_POLARISED , result );
2573+ CU_ASSERT_EQUAL_FATAL (ret , 0 );
2574+ CU_ASSERT_EQUAL_FATAL (result [0 ], 0 );
2575+ CU_ASSERT_EQUAL_FATAL (result [1 ], 1.0 );
2576+ CU_ASSERT_EQUAL_FATAL (result [2 ], 1.0 );
2577+ CU_ASSERT_EQUAL_FATAL (result [3 ], 0.0 );
2578+
2579+ ret = tsk_treeseq_allele_frequency_spectrum (& ts , 1 , sample_set_sizes , samples , 0 ,
2580+ NULL , 0 , NULL , TSK_STAT_BRANCH | TSK_STAT_POLARISED , result );
2581+ CU_ASSERT_EQUAL_FATAL (ret , 0 );
2582+ CU_ASSERT_TRUE_FATAL (result [0 ] > 0 );
2583+ CU_ASSERT_TRUE_FATAL (result [1 ] > 0 );
2584+ CU_ASSERT_TRUE_FATAL (result [2 ] > 0 );
2585+ CU_ASSERT_EQUAL_FATAL (result [3 ], 0.0 );
2586+
2587+ verify_afs (& ts );
2588+ tsk_treeseq_free (& ts );
2589+ }
2590+
24682591static void
24692592test_nonbinary_ex_ld (void )
24702593{
@@ -4000,6 +4123,7 @@ main(int argc, char **argv)
40004123 { "test_paper_ex_afs" , test_paper_ex_afs },
40014124 { "test_paper_ex_divergence_matrix" , test_paper_ex_divergence_matrix },
40024125
4126+ { "test_unary_ex_afs" , test_unary_ex_afs },
40034127 { "test_nonbinary_ex_ld" , test_nonbinary_ex_ld },
40044128 { "test_nonbinary_ex_mean_descendants" , test_nonbinary_ex_mean_descendants },
40054129 { "test_nonbinary_ex_genealogical_nearest_neighbours" ,
0 commit comments