Variance adjusted UniFrac #20
Conversation
|
Once more, could you update this branch? |
|
@ElDeveloper all green |
| @@ -743,6 +743,48 @@ void test_generalized_unifrac() { | |||
| SUITE_END(); | |||
| } | |||
|
|
|||
| void test_vaw_unifrac_weighted_normalized() { | |||
There was a problem hiding this comment.
Minor stylistic note, I checked in vim and this test case has a mix of tabs and spaces. Up to you if you want to fix. 👍
sucpp/unifrac.cpp
Outdated
| @@ -322,6 +321,56 @@ void progressbar(float progress) { | |||
| std::cout.flush(); | |||
| } | |||
|
|
|||
| void initialize_embedded(double*& prop, const su::task_parameters* task_p) { | |||
| posix_memalign((void **)&prop, 32, sizeof(double) * task_p->n_samples * 2); | |||
There was a problem hiding this comment.
This line also has a tab instead of spaces ¯\_(ツ)_/¯
sucpp/unifrac.cpp
Outdated
| Method unifrac_method, | ||
| const su::task_parameters* task_p) { | ||
| for(unsigned int i = task_p->start; i < task_p->stop; i++){ | ||
| posix_memalign((void **)&dm_stripes[i], 32, sizeof(double) * task_p->n_samples); |
There was a problem hiding this comment.
Just in case, the interwebs seem to indicate that posix_memalign has a return value specific for out of memory, ENOMEM.
| break; | ||
| case generalized: | ||
| func = &su::_vaw_generalized_unifrac_task; | ||
| break; |
There was a problem hiding this comment.
Worth adding a default case where this returns an error and the program terminates?
There was a problem hiding this comment.
added. I'm not too concerned as the switch is on an enum so this can only arise if the switch receives an enum value that is not described from su::Method which I believe is a compile-time constraint. But it is defensive and could also arise from programmer error.
| * at the moment. basically, we can't assume the presence of avx2. | ||
| */ | ||
| for(unsigned int j = 0; j < task_p->n_samples / 4; j++) { | ||
| int k = j * 4; |
There was a problem hiding this comment.
k is used here and in a loop below, perhaps worth a different letter ;)
There was a problem hiding this comment.
Not sure I follow? k is instantiated here as it is the stride offset. Or is it because k is used in the in the cleanup loop as well? I also think it is appropriate to stick with k there as the variable should not leak scope, and the use of the variable are consistent.
There was a problem hiding this comment.
Fair enough, just thought it was confusing, I wasn't concerned with the leakage.
| dm_stripe_total = dm_stripes_total[stripe]; | ||
|
|
||
| for(unsigned int j = 0; j < task_p->n_samples / 4; j++) { | ||
| int k = j * 4; |
| double sum1 = (u1 + v1) / vaw; | ||
| double sub1 = fabs(u1 - v1) / vaw; | ||
| double sum_pow1 = pow(sum1, task_p->g_unifrac_alpha) * length; | ||
| dm_stripe[j] += sum_pow1 * (sub1 / sum1); |
There was a problem hiding this comment.
Could sum1 ever be zero i.e. u1 and v1 being zero?
There was a problem hiding this comment.
It's implicit and originally the zero check was done on sum1. The reason this should never happen is that, if sum1 is zero then vaw must be zero as well.
sucpp/unifrac_task.hpp
Outdated
| unsigned int tid; // thread ID | ||
| double g_unifrac_alpha; // generalized unifrac alpha | ||
| }; | ||
|
|
There was a problem hiding this comment.
This is going to need some documentation at some point ... 😟
| double sum1 = (u1 + v1) / vaw; | ||
| double sub1 = fabs(u1 - v1) / vaw; | ||
| double sum_pow1 = pow(sum1, task_p->g_unifrac_alpha) * length; | ||
| dm_stripe[j] += sum_pow1 * (sub1 / sum1); |
Dependent on #18 and #19.
Tests are not written for VA unweighted, weighted unnormalized and generalized unifrac. These need to be done by hand as there isn't a reference implementation.