r1013: changed to a more accurate similarity est

lh3 · lh3 · commit 507d39af15c3 · 2021-04-08T13:57:53.000-04:00
Based on DOI:10.1101/2021.01.15.426881. One minimap2 reviewer suggested the
right formula to me but I thought the difference would be insignificant. I was
wrong.
diff --git a/esterr.c b/esterr.c
@@ -59,6 +59,6 @@ void mm_est_err(const mm_idx_t *mi, int qlen, int n_regs, mm_reg1_t *regs, const
 		n_tot = en - st + 1;
 		if (r->qs > avg_k && r->rs > avg_k) ++n_tot;
 		if (qlen - r->qs > avg_k && l_ref - r->re > avg_k) ++n_tot;
-		r->div = logf((float)n_tot / n_match) / avg_k;
+		r->div = n_match >= n_tot? 0.0f : (float)(1.0 - pow((double)n_match / n_tot, 1.0 / avg_k));
 	}
 }
diff --git a/main.c b/main.c
@@ -7,7 +7,7 @@
 #include "mmpriv.h"
 #include "ketopt.h"
 
-#define MM_VERSION "2.17-r1012-dirty"
+#define MM_VERSION "2.17-r1013-dirty"
 
 #ifdef __linux__
 #include <sys/resource.h>

Original file line number	Diff line number	Diff line change
`@@ -59,6 +59,6 @@ void mm_est_err(const mm_idx_t mi, int qlen, int n_regs, mm_reg1_t regs, const`
`59`	`59`	`n_tot = en - st + 1;`
`60`	`60`	`if (r->qs > avg_k && r->rs > avg_k) ++n_tot;`
`61`	`61`	`if (qlen - r->qs > avg_k && l_ref - r->re > avg_k) ++n_tot;`
`62`		`- r->div = logf((float)n_tot / n_match) / avg_k;`
	`62`	`+ r->div = n_match >= n_tot? 0.0f : (float)(1.0 - pow((double)n_match / n_tot, 1.0 / avg_k));`
`63`	`63`	`}`
`64`	`64`	`}`